[R] RNA Seq Analysis in R

Matthew McCormack mccorm@ck @end|ng |rom mo|b|o@mgh@h@rv@rd@edu
Sat Aug 1 21:52:19 CEST 2020


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 using 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 using 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]]



More information about the R-help mailing list