Limma Paired Analysis - Force intercept to 0
0
Entering edit mode
@moltenlight-24367
Last seen 3 days ago

I have a question regarding paired analysis in the Limma package. I am analyzing Illumina methylation arrays of paired tumor samples from two different sites. I came up with the following script to analyze the data, which I put together from old scripts I found in our lab. But there is something that bothers me, I do not understand why the design matrix has an intercept that is forced to 0. From a lecture on linear models I remembered that forcing the intercept to 0 is bad practice in many cases, but I'm not sure if I'm missing the point here? The results are substantially different between r~0+Location+Pairs and r~Location+Pairs.

# this is the factor of interest
Location <- factor(pd$Sample_Group, levels = c("M","P")) # the individual effect, each pair is a subject Pairs <- factor(pd$Pair)

# design matrix
design <- model.matrix(~0+Location+Pairs, data=pd)
colnames(design) <- c(levels(Tissuetype),levels(Pairs)[-1])

# fit the linear model
fit <- lmFit(Beta_To_M(myCombatLimma), design)

# create a contrast matrix for specific comparisons
contMatrix <- makeContrasts(P-M, levels=design)

# fit the contrasts, estimate sample-variance, compute moderated t-statistics
fit2 <- contrasts.fit(fit, contMatrix)
fit2 <- eBayes(fit2)

# pairwise analysis
summary(decideTests(fit2,
lfc = 0, # cutoff based on beta values
p.value = 0.05,
method = "separate"))

limma MethylationArray • 53 views
2
Entering edit mode
@gordon-smyth
Last seen 34 minutes ago
WEHI, Melbourne, Australia

The intercept isn't being forced to zero, it is just a reparametrization. You will get exactly the same results regardless of whether you specify 0+ or not if you make the corresponding change to the contrast. See Section 9.2 of the limma User's Guide.