I choose microarray data GSE75693 of 30 patients with stable kidney
transplantation and 15 with BKVN to identify differentially expressed genes
(DEGs). I performed this in GEO2R and find R script there and Runs R script
Successfully on R studio as well. The R script is :
# Differential expression analysis with limma
library(Biobase)
library(GEOquery)
library(limma)
# load series and platform data from GEO
gset <- getGEO("GSE75693", GSEMatrix =TRUE, AnnotGPL=TRUE)if
(length(gset) > 1) idx <- grep("GPL570", attr(gset,
"names")) else idx
<- 1
gset <- gset[[idx]]
# make proper column names to match toptable
fvarLabels(gset) <- make.names(fvarLabels(gset))
# group names for all samples
gsms <-
paste0("000000000000000000000000000000XXXXXXXXXXXXXXX11111",
"1111111111XXXXXXXXXXXXXXXXXXX")
sml <- c()for (i in 1:nchar(gsms)) { sml[i] <- substr(gsms,i,i) }
# eliminate samples marked as "X"
sel <- which(sml != "X")
sml <- sml[sel]
gset <- gset[ ,sel]
# log2 transform
exprs(gset) <- log2(exprs(gset))
# set up the data and proceed with analysis
sml <- paste("G", sml, sep="") # set group names
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=1250)
tT <- subset(tT,
select=c("ID","adj.P.Val","P.Value","t","B","logFC","Gene.symbol","Gene.title"))
DEGs = subset(tT, P.Value < 0.01 & abs(logFC) > 2)
After running this no genes are found plz help me
[[alternative HTML version deleted]]
https://www.bioconductor.org/help/ On August 1, 2020 4:01:08 AM PDT, Anas Jamshed <anasjamshed1994 at gmail.com> wrote:>I choose microarray data GSE75693 of 30 patients with stable kidney >transplantation and 15 with BKVN to identify differentially expressed >genes >(DEGs). I performed this in GEO2R and find R script there and Runs R >script >Successfully on R studio as well. The R script is : > > # Differential expression analysis with limma > >library(Biobase) >library(GEOquery) >library(limma) ># load series and platform data from GEO > >gset <- getGEO("GSE75693", GSEMatrix =TRUE, AnnotGPL=TRUE)if >(length(gset) > 1) idx <- grep("GPL570", attr(gset, "names")) else idx ><- 1 >gset <- gset[[idx]] ># make proper column names to match toptable >fvarLabels(gset) <- make.names(fvarLabels(gset)) ># group names for all samples >gsms <- paste0("000000000000000000000000000000XXXXXXXXXXXXXXX11111", > "1111111111XXXXXXXXXXXXXXXXXXX") >sml <- c()for (i in 1:nchar(gsms)) { sml[i] <- substr(gsms,i,i) } ># eliminate samples marked as "X" >sel <- which(sml != "X") >sml <- sml[sel] >gset <- gset[ ,sel] ># log2 transform >exprs(gset) <- log2(exprs(gset)) ># set up the data and proceed with analysis >sml <- paste("G", sml, sep="") # set group names >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=1250) > >tT <- subset(tT, >select=c("ID","adj.P.Val","P.Value","t","B","logFC","Gene.symbol","Gene.title")) >DEGs = subset(tT, P.Value < 0.01 & abs(logFC) > 2) > >After running this no genes are found plz help me > > [[alternative HTML version deleted]] > >______________________________________________ >R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see >https://stat.ethz.ch/mailman/listinfo/r-help >PLEASE do read the posting guide >http://www.R-project.org/posting-guide.html >and provide commented, minimal, self-contained, reproducible code.-- Sent from my phone. Please excuse my brevity.
As with the previous post, I agree that Bioconductor will be a better place to ask this question. As a quick thought you also might try to adjust the p-value in the last line: DEGs = subset(tT, P.Value < 0.01 & abs(logFC) > 2). You could play around with the P.Value, 0.01 is pretty low, you could try 0.05 and maybe abs(logFC) > 1. But, first you should try to print out tT with something like write.table(tT, file = TopTable.txt, sep = "\t"). This will write out tT to a tab-delimited text file (in the directory that you are working in) that you can import into Excel and then inspect the logFC and p-values for the top 1250 genes. Matthew On 8/1/20 1:13 PM, Jeff Newmiller wrote:> External Email - Use Caution > > https://www.bioconductor.org/help/ > > On August 1, 2020 4:01:08 AM PDT, Anas Jamshed <anasjamshed1994 at gmail.com> wrote: >> I choose microarray data GSE75693 of 30 patients with stable kidney >> transplantation and 15 with BKVN to identify differentially expressed >> genes >> (DEGs). I performed this in GEO2R and find R script there and Runs R >> script >> Successfully on R studio as well. The R script is : >> >> # Differential expression analysis with limma >> >> library(Biobase) >> library(GEOquery) >> library(limma) >> # load series and platform data from GEO >> >> gset <- getGEO("GSE75693", GSEMatrix =TRUE, AnnotGPL=TRUE)if >> (length(gset) > 1) idx <- grep("GPL570", attr(gset, "names")) else idx >> <- 1 >> gset <- gset[[idx]] >> # make proper column names to match toptable >> fvarLabels(gset) <- make.names(fvarLabels(gset)) >> # group names for all samples >> gsms <- paste0("000000000000000000000000000000XXXXXXXXXXXXXXX11111", >> "1111111111XXXXXXXXXXXXXXXXXXX") >> sml <- c()for (i in 1:nchar(gsms)) { sml[i] <- substr(gsms,i,i) } >> # eliminate samples marked as "X" >> sel <- which(sml != "X") >> sml <- sml[sel] >> gset <- gset[ ,sel] >> # log2 transform >> exprs(gset) <- log2(exprs(gset)) >> # set up the data and proceed with analysis >> sml <- paste("G", sml, sep="") # set group names >> 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=1250) >> >> tT <- subset(tT, >> select=c("ID","adj.P.Val","P.Value","t","B","logFC","Gene.symbol","Gene.title")) >> DEGs = subset(tT, P.Value < 0.01 & abs(logFC) > 2) >> >> After running this no genes are found plz help me >> >> [[alternative HTML version deleted]] >> >> ______________________________________________ >> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see >> https://stat.ethz.ch/mailman/listinfo/r-help >> PLEASE do read the posting guide >> http://www.R-project.org/posting-guide.html >> and provide commented, minimal, self-contained, reproducible code.[[alternative HTML version deleted]]