Skip to content

Instantly share code, notes, and snippets.

View al2na's full-sized avatar
💭
fishing 🎣

Altuna Akalin al2na

💭
fishing 🎣
View GitHub Profile
@al2na
al2na / mread.methylKit.R
Last active August 29, 2015 14:08
faster read for methylKit using fread from data.table package
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)
@al2na
al2na / meanMethGroup.R
Created September 27, 2014 17:34
get mean methylation per group from a methylBase object
# 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
@al2na
al2na / ideoDMC.R
Last active June 7, 2023 16:02
ideogram for differential methylation values
#' function for making ideogram for differential methylation values
#' requires methylKit, ggbio and GenomicRanges
#'
#' @example
#' library(BSgenome)
#' library("BSgenome.Hsapiens.UCSC.hg18")
#' chr.len = seqlengths(Hsapiens) # get chromosome lengths
#' # remove X,Y,M and random chromosomes
#' chr.len = chr.len[grep("_|M|X|Y", names(chr.len), invert = T)]
#'
@al2na
al2na / readTableFast.R
Last active August 29, 2015 13:56
read.table faster with known column classes and data.table::fread. I wonder which one is
.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) )
}
@al2na
al2na / uniteByChr.R
Created February 6, 2014 22:19
methylKit unite chromosome by chromosome
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)
@al2na
al2na / getOE.R
Last active January 4, 2016 02:59
get O/E ration and CpG ratio and GC content for given GRanges object and BSgenome. Returns a matrix of O/E ratio, CpGcontent and GC content
require(GenomicRanges)
require(Biostrings)
#' get OE ratio and GC content for a given set of DNAstrings
getOE.strset<-function(str.set)
{
di.mat=dinucleotideFrequency( str.set )
a.mat =alphabetFrequency( str.set ,baseOnly=TRUE )
exp=(a.mat[,2]*a.mat[,3])/width(str.set )
@al2na
al2na / callPrimer3.R
Last active December 11, 2021 19:48
calling primer3 from R
#' call primer3 for a given set of DNAstringSet object
#'
#' @param seq DNAstring object, one DNA string for the given amplicon
#' @param size_range default: '151-500'
#' @param Tm melting temprature parameters default:c(55,57,58)
#' @param name name of the amplicon in chr_start_end format
#' @param primer3 primer3 location
#' @param therme.param thermodynamic parameters folder
#' @param settings text file for parameters
#' @author Altuna Akalin modified Arnaud Krebs' original function
@al2na
al2na / tileMethylCounts2.R
Last active December 30, 2015 21:49
a version of tileMethylCounts where the chromosome names obtained from GRanges out of unique seqnames, so there will be no chr names that do not have coverage. see the discussion here: https://groups.google.com/forum/#!topic/methylkit_discussion/75y6dbrbCWI
setGeneric("tileMethylCounts2",
function(object,win.size=1000,step.size=1000,cov.bases=0)
standardGeneric("tileMethylCounts2") )
setMethod("tileMethylCounts2", signature(object="methylRaw"),
function(object,win.size,step.size,cov.bases){
g.meth =as(object,"GRanges")
chrs =as.character(unique(seqnames(g.meth)))
widths =seqlengths(g.meth)
@al2na
al2na / getCMethMatrix.R
Created November 11, 2013 15:13
getting a matrix methylation status matrix via QuasR
require(data.table)
require(QuasR)
# arguments:
# proj: qProject object
# range: GRanges object with ONE!!!! range
# samp: sample.name
#
getCMethMatrix<-function(proj,range,samp){
Cs=qMeth(proj, query=range,mode="allC",reportLevel="alignment")
@al2na
al2na / summarizeScoresOnRegions.R
Last active March 23, 2020 22:20
getting average score of a GRanges over another GRanges object. summarizeScores(loc.gr,score.gr,score.col)
.libPaths("/work2/gschub/altuna/RlibsDev")
require(data.table)
summarizeScores<-function(loc.gr,score.gr,score.col){
ov <- as.matrix(findOverlaps(loc.gr,score.gr))
dt=data.table(id=ov[,1],score=values(score.gr)[ov[,2],which(names(values(score.gr))==score.col)])
dt=dt[,list(av=mean(score)),by=id]