Why my limma DEG result in top Table all without significant meaning,showing logFC <1.5 or logFC>-1.5 ?

Question

[BioLite]

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

[Kevin Blighe]

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

[BioLite]

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!

[Kevin Blighe]

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)

[BioLite]

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!

[BioLite]

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)

##loading data
targets<-read.targets("targets.txt")
raw_data<-read.AgilentFE(targets = targets,makePLOT = FALSE)

##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)

[Devon Ryan]

Does Result indicate you have any DE genes?

[BioLite]

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.

[Devon Ryan]

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.

[Kevin Blighe]

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.

[Fabio Marroni]

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?

[BioLite]

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!

[WouterDeCoster]

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:

[Kevin Blighe]

In addition, cross-posted: https://support.bioconductor.org/p/112875/

[BioLite]

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

Answer 1

[Gordon Smyth]

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/

[BioLite]

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!