Question

Contrasts & coefficients in limma

0

Entering edit mode

prathap1708 • 0

@prathap1708-14714

Last seen 13 months ago

Duke NUS medical school, Singapore

Hi, I have a doubt regarding contrasts. In my case, i understood,the 1st coefficient is comparison of MHO with MHNW, 2nd coefficient is comparison of MAO with MHNW and 3rd coefficient is comparison of MAO with MHW . Please let me know, if I understood wrong. Please find my code below.

design.adipose = model.matrix(~0+AdiposeCPMCountGT1_min20percentsamples$samples$group.1+AdiposeCPMCountGT1_min20percentsamples$samples$Gender)

# To rename column names
colnames(design.adipose)[1] = "MAO"
colnames(design.adipose)[2] = "MHNW"
colnames(design.adipose)[3] = "MHO"
colnames(design.adipose)[4] = "Gender"
# To do comparisons between groups
contr.matrix.Adipose = makeContrasts(MHOvsMHNW = MHO-MHNW,MAOvsMHNW = MAO-MHNW,MAOvsMHO = MAO-MHO,levels = colnames(design.adipose)[1:4])
# Voom function v = voom(AdiposeDGE,design.adipose,lib.size = colSums(AdiposeDGE$counts)*AdiposeDGE$samples$norm.factors,normalize.method="quantile",plot = TRUE)

design and contrast matrix • 2.8k views

ADD COMMENT • link updated 7.0 years ago by Gordon Smyth 52k • written 7.0 years ago by prathap1708 • 0

score 0 · Answer 1 · 2018-01-09

0

Entering edit mode

Gordon Smyth 52k

@gordon-smyth

Last seen 4 hours ago

WEHI, Melbourne, Australia

Have you looked at the limma User's Guide, for example Sections 9.1-9.3? You have to fit a linear model and then use contrasts.fit before your contrast matrix will have any effect.

At the moment, you don't yet have any coefficients and you haven't yet made any comparisons.

To make the comparisons you want, code might go like this:

design.adipose <- model.matrix(~0+group.1+Gender,
                  data=AdiposeDGE$samples)
colnames(design.adipose) <- c("MAO","MHNW","MHO","Gender")
contr.matrix.Adipose = makeContrasts(
       MHOvsMHNW = MHO-MHNW,
       MAOvsMHNW = MAO-MHNW,
       MAOvsMHO = MAO-MHO,
       levels = design.adipose)
v <- voom(AdiposeDGE,design.adipose, normalize.method="quantile", plot = TRUE)
vfit <- lmFit(v, design.adipose)
vfitAdipose <- contrasts.fit(vfit, contr.matrix.Adipose)
efit <- eBayes(vfitAdipose)
voom.Results_MHO.vs.MHNW <- topTable(efit,coef="MHOvsMHNW", n=Inf)

and so on.

ADD COMMENT • link 7.0 years ago Gordon Smyth 52k

0

Entering edit mode

I did linear modeling fitting and then the contrast was done.I didn't post the last few lines of code in the previous question. Sorry for that. Please see the remaining code here. Thanks for your spontaneous reply.

# Usual limma pipeline
vfit <- lmFit(v, design.adipose)
vfitAdipose <- contrasts.fit(vfit, contrasts=contr.matrix.Adipose)
efit = eBayes(vfitAdipose)
plotSA(efit)
summary(decideTests(efit))

## Access Results
# Examining individual genes from top to bottom
voom.Results_MHO.vs.MHNW = toptable(efit,coef = 1,adjust="BH",number = nrow(Adipose.Dataset))
voom.Results_MAO.vs.MHNW = toptable(efit,coef = 2,adjust="BH",number = nrow(Adipose.Dataset))
voom.Results_MAO.vs.MHO = toptable(efit,coef = 3,adjust="BH",number = nrow(Adipose.Dataset))

ADD REPLY • link 7.0 years ago prathap1708 • 0

0

Entering edit mode

You should use topTable instead of toptable.

Also, you seem to be using two different DGEList data objects, one called AdiposeDGE and the other called AdiposeCPMCountGT1_min20percentsamples. I can't check whether these two datasets match up in terms of sample information. You obviously need to be consistent.

Otherwise your code looks correct.

Your code is longer than it needs to be in some respects, but it should work ok.

ADD REPLY • link 7.0 years ago Gordon Smyth 52k

0

Entering edit mode

Thanks prof. The two DGE list objects are same but we removed genes which were expressed very low in the 2nd DGE object.

ADD REPLY • link 7.0 years ago prathap1708 • 0