| subsetByChr<-function(methylRawList.obj,my.chr="chr21"){ | |
| sub=lapply(methylRawList.obj, function(x) x[x$chr==my.chr] ) | |
| new("methylRawList",sub,treatment=methylRawList.obj@treatment) | |
| } | |
| my.chrs= unique(methylRawList.obj[[1]]$chr) |
| .readTableFast<-function(filename,header=T,skip=0,sep="") | |
| { | |
| tab5rows <- read.table(filename, header = header,skip=skip,sep=sep, | |
| nrows = 100,stringsAsFactors=FALSE) | |
| classes <- sapply(tab5rows, class) | |
| classes[classes=="logical"]="character" | |
| return( read.table(filename, header = header,skip=skip,sep=sep, | |
| colClasses = classes) ) | |
| } |
| # returns a matrix of mean methylation values for groups defined | |
| # by treatment vector in the methylBase obj | |
| meanMethGroup<-function(methylBase.obj,weighted=TRUE ){ | |
| data=getData(methylBase.obj) # get data frame part of the object | |
| treatment=methylBase.obj@treatment # get the treatment vector from the object | |
| # create the the empty resulting matrix | |
| result=matrix(0,ncol=length(unique(methylBase.obj@treatment)) ,nrow=nrow(methylBase.obj) ) | |
| colnames(result)=paste('group',unique(treatment),sep='') # column names are from treatmet vector |
| require(data.table) | |
| .structureAMPoutput<-function(data) | |
| { | |
| strand=rep("+",nrow(data)) | |
| strand[data[,4]=="R"]="-" | |
| numCs=round(data[,5]*data[,6]/100) | |
| numTs=round(data[,5]*data[,7]/100) | |
| data.frame(chr=data[,2],start=data[,3],end=data[,3] | |
| ,strand=strand,coverage=data[,5],numCs=numCs,numTs=numTs) |
| #!/bin/bash | |
| #$ -N run_bowtie | |
| #$ -cwd | |
| #$ -pe smp 2 | |
| #$ -l h_vmem=6G | |
| infile=/data/bioinfo/genome_data/sample.ce10.fastq | |
| outfile=~/my.alignment.sam | |
| btindex=/data/bioinfo/genome_data/Caenorhabditis_elegans/UCSC/ce10/Sequence/BowtieIndex/genome |
The git command-line utility has plenty of inconsistencies http://steveko.wordpress.com/2012/02/24/10-things-i-hate-about-git/
A GUI like http://sourcetreeapp.com is often helpful, but staying on the command line usually quicker. This is a list of the commands I use most frequently, listed by funcional category:
git status list which (unstaged) files have changed
This document shows how to use packrat to create reproducible R projects. packrat allows reproducuiblty of R code in an OS independent manner. First we need to install the packrat package.
install.packages("packrat")
setwd("~/mdc_desktop_projects/packrat_trial3")
now we can initialize the packrat project. This will create local directories for packages and .Rprofile file. Whenever we run R in this directory, we will be using locally installed R packages and even if we copy this directory to some other computer (even with different OS) we will still have the same setup.
| y=readRDS("data/othersmRNA.pca.table.rds") | |
| p1=nPlot(PC2 ~ PC1, data = y,group="cellDisease", type = "scatterChart", | |
| filter="disease") | |
| p1$xAxis(axisLabel = 'PC1') | |
| p1$yAxis(axisLabel = 'PC2') | |
| p1$chart(color = y$color) # color for each dot | |
| p1$chart(sizeRange = c(100,100)) # for dot sizes |
| download.file("https://dl.dropboxusercontent.com/u/1373164/H1.chr21.chr22.rds",destfile="H1.chr21.chr22.rds",method="curl") | |
| mbw=readRDS("H1.chr21.chr22.rds") | |
| # it finds the optimal number of componets as 6 | |
| res=methSeg(mbw,diagnostic.plot=TRUE,maxInt=100,minSeg=10) | |
| # however the BIC stabilizes after 4, we can also try 4 componets | |
| res=methSeg(mbw,diagnostic.plot=TRUE,maxInt=100,minSeg=10,G=1:4) | |
| # get segments to BED file |
| .Rproj.user | |
| .Rhistory | |
| .RData | |
| *.Rproj | |
| *.html |