demo_geo2r.r

demo_geo2r.r

gset <- getGEO("GSE7442", GSEMatrix =TRUE)
if (length(gset) > 1) idx <- grep("GPL5058", attr(gset, "names")) else idx <- 1
gset <- gset[[idx]]

# group names for all samples
sml <- c("X","X","X","X","X","G0","X","X","X","X","G0","X","G0","G1","G1","G0","X","X","X","X","X","X","X","X","X","X","X","X","X","X","X","X","X","G1","X","G0","X","X","X","G1","X","X","G1");

# eliminate samples marked as "X"
sel <- which(sml != "X")
sml <- sml[sel]
gset <- gset[ ,sel]

# set up the data and proceed with analysis
fl <- as.factor(sml)
gset$description <- fl
design <- model.matrix(~ description + 0, gset)
colnames(design) <- levels(fl)
fit <- lmFit(gset, design)
cont.matrix <- makeContrasts(G1-G0, levels=design)
fit2 <- contrasts.fit(fit, cont.matrix)
fit2 <- eBayes(fit2, 0.01)
tT <- topTable(fit2, adjust="fdr", sort.by="B", number=250)

tT <- subset(tT, select=c("ID","adj.P.Val","P.Value","t","B","logFC","SPOT_ID","GB_LIST"))
write.table(tT, file=stdout(), row.names=F, sep="\t")