topTable in limma

Question

[AHW]

I am using limma to perform differential expression of illumina microarray gene expression data. I get the topTable by using the following line of code

top <- topTable(fit2, coef = i, n = Inf)

                logFC   AveExpr             t     P.Value adj.P.Val         B
ILMN_1748476 -3.302386e-01  9.429974 -3.088904e+00 0.005235213 0.9999924 -3.701710
ILMN_1813544 -2.985103e-01  7.751847 -2.834129e+00 0.009476781 0.9999924 -3.831167
ILMN_1679438 -2.431550e-01  7.305140 -2.737038e+00 0.011831301 0.9999924 -3.880349
ILMN_1690546 -3.134771e-01  9.486155 -2.728182e+00 0.012071672 0.9999924 -3.884826
ILMN_1676423 -3.119174e-01  8.284954 -2.658629e+00 0.014126178 0.9999924 -3.919914
ILMN_1703555 -2.914465e-01  7.673945 -2.590389e+00 0.016458251 0.9999924 -3.954196
ILMN_1691418 -2.362332e-01  7.486428 -2.576331e+00 0.016981505 0.9999924 -3.961237
ILMN_1832155 -2.510850e-01  7.035359 -2.575989e+00 0.016994431 0.9999924 -3.961408
ILMN_1692223  9.603060e-01  8.908288  2.536541e+00 0.018548104 0.9999924 -3.981123
ILMN_1806056  9.548615e-01  9.079619  2.476310e+00 0.021177783 0.9999924 -4.011091
ILMN_1689294 -2.109505e-01  6.794665 -2.469909e+00 0.021476760 0.9999924 -4.014266
ILMN_1679045 -4.135768e-01  9.117017 -2.463320e+00 0.021788585 0.9999924 -4.017531
ILMN_1687857  4.310552e-01  8.288689  2.457058e+00 0.022088852 0.9999924 -4.020632
ILMN_1658149 -2.335111e-01  9.369083 -2.440079e+00 0.022922402 0.9999924 -4.029031
ILMN_1796038 -2.484446e-01 11.072939 -2.422479e+00 0.023817087 0.9999924 -4.037720
ILMN_1702858 -3.202939e-01  7.977537 -2.421800e+00 0.023852240 0.9999924 -4.038055
ILMN_1777849 -3.336440e-01 12.242453 -2.410675e+00 0.024435156 0.9999924 -4.043538
ILMN_1723035  7.375359e-01  7.704460  2.407119e+00 0.024624210 0.9999924 -4.045289
ILMN_1746666 -2.647887e-01  7.165868 -2.390870e+00 0.025505513 0.9999924 -4.053281
ILMN_1748915  3.498114e-01 10.570429  2.381106e+00 0.026048981 0.9999924 -4.058075

But when I try to filter the topTable based on p-value and fold-change I get no genes. The following code I used to filter the topTable. The parameter of p-value and fold-change are used to get all the genes whose p-value and fold-change is less than or equal to the supplied parameters.

topTable(fit2, coef = i, n = Inf, p.value = 1, lfc = 10, adjust.method = FDRMethod)

With this I get no genes. What I am doing wrong.

[cpad0112]

lfc=10 means FC=2^10=1024. You are looking for genes that changed by 1024 times (between case vs control) and p-value filter 1 is not advisable. A reasonable one would be 0.05. However, you already have adjusted pvalue. Print the ranges of your FC and your p-value and adjusted p value. You can filter within those ranges. It is not advisable to filter on FC and p-value at the same time as per limma toptable documentation.

Although topTable enables users to set p-value and lfc cutoffs simultaneously, this is not generally recommended.  from topTable function in limma.

[AHW]

I don't understand how lfc = 10 is equal to 1024 (Is it given in the documentation). If what you say is right also then it should select some genes as you can see there are many genes with lfc < 1024 . I also understand that p-value > 0.054 is not recommended but I just tried by setting different parameters of fold-change and p-value and all i got is no genes.

[Devon Ryan]

lfc means log2(fold-change), so unless you have genes going from no expression to moderate expression you won't have any of those. You're setting a minimum amount with this, the lfc must be at least 10 for genes to be included.

[AHW]

Thank you. I think lfc = 10 is the cutoff value, so all genes having less or equal should be reported. I checked with different lfc values like 0.1 but no genes reported. On the basis of above pasted data what should be my lfc value to filter?.

[Devon Ryan]

BTW, lfc=10 will return genes with (absolute) lfc of at least 10. This is mentioned in the documentation and fits how people want to filter things.

[AHW]

Thanks Devon Ryan, you are right.

[Devon Ryan]

If you filter by (adjusted) p-values then don't filter by fold-change.

[cpad0112]

lfc is a minimum cutoff, not a maximum, for filtering DEGs i.e all DEGs should satisfy minimum of lfc=10 (between groups, in OP). Refer to page 51 of https://www.bioconductor.org/packages/devel/bioc/manuals/limma/man/limma.pdf.

lfc=10 means log2FC=10. FC would be 2^10 (1024). If lfc=1, then FC =2.

To understand the lfc conversion, (from logFC values to FC value): R: How to convert log2FC (Fold Change) obtained by limma's topTable() function to FC.

Copy/pasted from Limma manual: For example, one might choose lfc=log2(1.5) to restrict to 50% changes or lfc=1 for 2-fold changes.

In OP, code is (copy/pasted):

topTable(fit2, coef = i, n = Inf, p.value = 1, lfc = 10, adjust.method = FDRMethod)

• p.value=1. This is not advisable in general practice. At least I haven't come across such a cut off in DEGs
• lfc = 10 - too high cutoff
• adjust.method = FDRMethod - toptable has adjusted p-values (you can see it in OP). It is preferable to use adjusted p-values than p-values them selves.

However, Limma manual advises not to filter DEGs by both FC and p-value at the same time.