library(edgeR) files <- dir(pattern="*\\.tab$") RG <- readDGE(files, skip=3) table<-RG$counts; head(table); nn<-colnames(table) nn<-sub(".ReadsPerGene.out","",nn) colnames(table)<-nn colnames(RG$counts)<-nn group <- factor(rep(c("T0","T1","T6"),2)) counts<-RG$counts > dim(counts)# 48526 6 nonzero<-counts[rowSums(counts)>0,]; dim(nonzero)# 30772 6 y <- DGEList(counts=nonzero, group=group) keep <- rowSums(cpm(y)>1) >= 4; #por lo menos el 66% de las muestras tengan 1 cpm dim (y[keep,])# 13467 6 y <- y[keep, , keep.lib.sizes=FALSE] y <- calcNormFactors(y) y <- estimateDisp(y) y <- estimateTrendedDisp(y) y <- estimateTagwiseDisp(y) png("bcv.png") plotBCV(y) dev.off() png("smear.png") plotSmear(y) dev.off() geneNames<-read.table("Hsa.ID.names.txt")