rm(list=ls())
source("robust-fdr.R")
setwd("~... ... CNV/pvalues,qvalues/RTFET")
infiles = list.files(, pattern="txt")
for( i in 1:length(infiles)){
tb = read.table(infiles[i], sep="\t", header=F)
tmp1 = tb[tb$V4<1.0, ];
if( length(tmp1[,1]) >2 ) {
hist(tmp1$V2)
ret = robust.fdr(tmp1$V4, sides=1, discrete=T )
tmp1$q1 = ret$q
tb$q1 = tmp1$q1[match(tb$V1, tmp1$V1)]
} else {
tb$q1 = NA;
}
summary(tb)
tb[tb$q1<0.15 & (!is.na(tb$q1)), ]
outfile = paste("output_robustFDR/", infiles[i], ".csv", sep='')
write.csv(tb, outfile,, row.names=F, quote=F)
}
No comments:
Post a Comment