How to apply two tailed paired t-test in limma for DEGs selection?

Question

[prithvi.mastermind]

I want to apply paired t-test in limma for DEGs selection. I have 46 pairs of samples (tumor and normal each) making a total of 92 samples.

The processed matrix file can be found at https://od.lk/s/ODdfMjgyMTA2MDNf/Data%20for%20Limma.txt

The phenotype file can be found at https://od.lk/s/ODdfMjgyMTA2MDFf/LUSC_Phenotype.txt

###### Load necessary packages######
library(NOISeq)
library(edgeR)
library(DESeq2)
library(ggplot2)
library(gplots)
library(RColorBrewer)
library(limma)
library(sva)
library(biomaRt)

 ###### Load processed matrix file and phenotyp files######

LIMS <- read.delim(file.choose(), row.names=1)
PHENO1 <- read.delim(file = file.choose(), as.is = T, row.names = 1)

###### Create design matrix######
    trt <- factor(rep(1:2, 46))
    pairs <- factor(rep(1:46, each = 2))
    design.gse74706 = model.matrix(~trt+pairs)

###### LIMMA for DEGs screening######
data.fit.gse74706 = lmFit(LIMS, design.gse74706)
fit2 <- eBayes(data.fit.gse74706)
top.table<- topTable(fit2, n = Inf)
length(which(top.table$P.Value < 0.05))

But at the end it states that 47 DEGs at p-value < 0.05. Is the code right?
However, when I apply unpaired t-test I get 2809 DEGs at p-value < 0.05. This code I used for unpaired t-test.

    ###### Create design matrix######
    groups.gse74706 = PHENO1$Type
    f.gse74706 = factor(groups.gse74706, levels=c("Tumor", "Normal"))
    design.gse74706 = model.matrix(~0+f.gse74706)
    colnames(design.gse74706) = c("Tumor", "Normal")
    statusCol <- as.numeric(factor(PHENO1$Type)) + 1

<h6>LIMMA for DEGs screening</h6>

data.fit.gse74706 = lmFit(LIMS, design.gse74706)
contrast.matrix.gse74706 = makeContrasts(Tumor-Normal, levels=colnames(coef(data.fit.gse74706)))
data.fit.con.gse74706 = contrasts.fit(data.fit.gse74706, contrast.matrix.gse74706)
data.fit.con.gse74706 = eBayes(data.fit.con.gse74706)
top.table <- topTable(data.fit.con.gse74706, sort.by = "P", n = Inf)
length(which(top.table$P.Value < 0.05))

[ATpoint]

What is the question? I am not sure your paired t-test is correct. Be sure to set a contrast there as well to make sure you are really extracting the treatment comparison. I am not sure what the default coef=NULL in topTable extracts without a contrast or contrast.fit.

Answer 1

[Gordon Smyth]

In your paired t-test code, you need to add the argument coef=2 to the topTable call:

top.table<- topTable(fit2, n = Inf, coef = 2)

You should also be using the adj.P.Val column rather than the unadjusted p-values.

See the limma User's Guide for more detail.

The links you give to the data matrix and simple information do not work. The files appear to have been deleted. I wonder what sort of processed data file you have. Your code seems to omit normalization and filtering steps that would be part of most analysis pipelines.

[ATpoint]

Gordon Smyth Out of curiosity, if coef=NULL as in OPs code, it internally does coef <- 1:ncol(fit), is that an ANOVA-like test for any differences based on any covariate or how would one interpret this?

[Gordon Smyth]

If coef=NULL, then topTable will return ANOVA-like F-tests for any differences, which in this case means between pairs or between treatments. That is equivalent to setting coef=2:ncol(fit), i.e., all coefficients except for the intercept.

[prithvi.mastermind]

Thank you sir for your assistance and much needed valuable suggestions.