Why my limma DEG result in top Table all without significant meaning,showing logFC <1.5 or logFC>-1.5 ?
1
0
Entering edit mode
4.2 years ago
BioLite ▴ 20

Hi, all I'm a newbie, and I'm dealing with Agilent single channel microarray when I filtering DEGs using limma, I find all my results showing logFC>-1.5 or logFC<1.5, without significant expression. I have been searching reason for a long time, but I still don't know what my faults are, I really need some help, thanks! My code as following:

treatment<-targets$Treatment Level<-c("neurologically healthy control","Parkinson disease") Group<-factor(treatment,levels = Level) Design<-model.matrix(~0+Group) colnames(Design) <- c("Control","PD") fita<-lmFit(exprs(eSet),Design) cont.matrix <- makeContrasts(Diff=PD-Control, levels=Design) fitb<- contrasts.fit(fita, cont.matrix) fitc<- eBayes(fitb) ##result Result<-topTable(fitc,adjust="BH",number = Inf,p.value = 0.05) dt_result<-decideTests(fitc,lfc =1.5,p.value=0.05,method="separate") summary(dt_result)range (Result$logFC)
[1] -1.471461  1.192352
limma deg toptable ;logfc • 3.1k views
1
Entering edit mode

One issue is here:

dt_result<-decideTests(fitc,lfc = log2(1.5),p.value=0.05,method="separate")

The parameter, lfc, expects a value that is already on the log (base 2) scale. So, the value that you are passing as the cut-off is actually:

log2(1.5) = 0.5849625

You need to use:

dt_result <- decideTests(fitc, lfc = 1.5, p.value=0.05, method="separate")

Kevin

0
Entering edit mode

Yes Kevin, thank for your answer. But my key question still alive. Once I set parameter lfc=1.5, summary(dt_result) showing Diff Down 0 NotSig 38822 Up 0 That's because of my topTable results with no changes, with nothing significant expression. I really don't know how to get rid of this situation. Could you help me, Kevin? Very thanks!

0
Entering edit mode

You are in good hands, now that Devon has added a comment. Please follow up with both him and Fabio (Fabio is another great guy)

0
Entering edit mode

Kevin, thank you very much! Devon and Fabio both experienced than me, I will learn modestly from them. Thanks again, wish you have a good day!

0
Entering edit mode

Here is my complete code, really hope to get some help.

rm(list = ls())
setwd("E:/Practice/Aglient/PD")
library(Agi4x44PreProcess)
library(hgug4112a.db)
library(limma)

##Background Correction and Normalization between arrays
data_NORM<-BGandNorm(raw_data,
BGmethod="half",
NORMmethod="quantile",
foreground="MeanSignal",
background="BGMedianSignal",
offset=50,makePLOTpre=FALSE,
makePLOTpost=FALSE)

##filtering probe
ddFILT<-filter.probes(data_NORM,
control=TRUE,
wellaboveBG=TRUE,
isfound=TRUE,
wellaboveNEG=TRUE,
sat=TRUE,
PopnOL=TRUE,
NonUnifOL=TRUE,
nas=TRUE,
limWellAbove=75,
limISF=75,
limNEG=75,
limSAT=75,
limPopnOL=75,
limNonUnifOL=75,
limNAS=100,
makePLOT=TRUE,
annotation.package="hgug4112a.db",
flag.counts=TRUE,targets)

##summarizing
data_filter<-summarize.probe(ddFILT, makePLOT=FALSE, targets)
dim(data_filter)

##create an expressionSet
eSet<-build.eset(data_filter, targets, makePLOT=FALSE,annotation.package="hgug4112a.db")

##DEG filter
treatment<-targets\$Treatment
Level<-c("neurologically healthy control","Parkinson disease")
Group<-factor(treatment,levels = Level)
Design<-model.matrix(~0+Group)
colnames(Design) <- c("Control","PD")
fita<-lmFit(exprs(eSet),Design)
cont.matrix <- makeContrasts(Diff=PD-Control, levels=Design)
fitb<- contrasts.fit(fita, cont.matrix)
fitc<- eBayes(fitb)

##DEG result
Result<-topTable(fitc,adjust="BH",number = Inf,p.value = 0.05,sort.by = "p")
dt_result<-decideTests(fitc,lfc = 1.5,p.value=0.05,method="separate")
summary(dt_result)
0
Entering edit mode

Does Result indicate you have any DE genes?

0
Entering edit mode

Nothing. I just checked my code again, finding some wrongs. So sad, I always make mistakes. Let me correct it and then discuss with you.

Devon, thanks for your caring! Have a good day!

BTW, my data from this article,https://journals.plos.org/plosgenetics/article?id=10.1371/journal.pgen.1002794. I am trying to replicate the Agilent part results (the first part) of this article. As a newbie, It's apparently I failed completely.

Best wish for you again.

0
Entering edit mode

If there really are DE genes in that dataset (from the paper I guess one would expect so) then it's likely that there's either an outlier sample that needs excluding or some of the sample labels are swapped. Make a PCA plot and do some clustering and see if that indicates either of those.

0
Entering edit mode

So sad, I always make mistakes.

Do not worry - we all make mistakes. At least you are honest about it, whereas others would prefer to lie about their mistakes. You are a great person.

0
Entering edit mode

How many replicates do you have per condition?

Do you see very high variation between replicates?

Did you try doing the analysis with a toy data set (probably provided by the package developer) to see if using that you got significant results?

Are there additional software packages you may try to see if the result is due to some problem in the analysis with this specific package?

0
Entering edit mode

Yes, Fabio. Q1, 27 PD, and 26 control samples;

Q2, my analysis data from this published article,https://journals.plos.org/plosgenetics/article?id=10.1371/journal.pgen.1002794. The bad samples have been removed by the authors, I get the cleaned data from ArrayExpress,https://www.ebi.ac.uk/arrayexpress/experiments/E-MTAB-812/

Q3, not yet.

Q4, I made this analysis using limma at first and got the same results.

I reviewed my practice article again, thinking there may be some faults in my limma lmFit object and model design. I will correct them firstly.

Thanks for your precious time. Fabio, wish you have a great day!

0
Entering edit mode

I added code markup to your post for increased readability. You can do this by selecting the text and clicking the 101010 button. When you compose or edit a post that button is in your toolbar, see image below:

0
Entering edit mode

0
Entering edit mode

Yes, thanks for your sweet tip, I will notice my code layout next time.

1
Entering edit mode
4.2 years ago
Gordon Smyth ★ 5.1k

The Agi4x44PreProcess package was removed from Bioconductor many years ago. You must be using very old versions of R and Bioconductor. I would certainly not advise you to preprocess the Agilent microarray data in the complicated way that you have. Why not use more up-to-date software and a pipeline as recommended in the limma User's Guide?

I analysed this data myself and had no difficulty finding DE genes. See my answer on the Bioconductor Support site https://support.bioconductor.org/p/112875/

0
Entering edit mode

Dear Gordon,

Actually, I worked with limma at first. However, I couldn't understand the probe filter part in the user guide.

Following your tips, I will analysis with limma again.

Thanks!