[BioC] limma design matrix and contrast specification for time course PCR data with two conditions and variable replicate number

[Rao,Xiayu]

Hi, Linus

There is a section called Time course Experiments (9.6) in the Limma user's guide. You can find the design similar to yours there as well as other designs. 

Thanks,
Xiayu

-----Original Message-----
From: bioconductor-bounces at r-project.org [mailto:bioconductor-bounces at r-project.org] On Behalf Of Linus Schumacher
Sent: Monday, August 11, 2014 7:38 PM
To: bioconductor at r-project.org
Subject: [BioC] limma design matrix and contrast specification for time course PCR data with two conditions and variable replicate number

Dear list,

I'm trying to find differential expression between two conditions, as well as within a condition over time. The data are RTqPCR, log10 relative expression of 96 genes (by 55), under 2 conditions, at 8 unevenly spaced timepoints with some missing data (NA) and variable number of replicates (2-3):

nReplicates <- c(3,3,3,3,3,3,2,3, 2,3,2,3,2,2,3,3) replicates <- rep(c(1:16), nReplicates) timePoints <- rep(rep(c(2,4,8,16,30,45,60,90),2), nReplicates) cond <- factor(rep(rep(0:1, each=8), nReplicates))

I've managed to run the analysis, but I'm unfamiliar with the syntax for specifying the design matrix and would greatly appreciate if someone could confirm the design matrix is set up correctly. I couldn't find an intuitive explanation, so I have been working of examples from related cases.

design <- model.matrix(~cond+timePoints:cond)
   (Intercept) cond1 cond0:timePoints cond1:timePoints
1            1     0                2                0
2            1     0                2                0
...
23           1     0               90                0
24           1     1                0                2
...
43           1     1                0               90
attr(,"assign")
[1] 0 1 2 2
attr(,"contrasts")
attr(,"contrasts")$cond
[1] "contr.treatment"

Now to determine which gene are differentially expressed between the two conditions,

dupfit <- duplicateCorrelation(data,design,block=replicates)
fit <- lmFit(data,design,block=replicates,correlation=dupfit$consensus)
fit <- eBayes(fit)
topTable(fit,coef="cond1")

(I have checked dupfit$consensus and is positive, if small (0.11))

And to determine difference in slope of genes under two conditions, I tried

cont.VEGF <- makeContrasts(contrasts="cond0:timePoints-cond1:timePoints",levels=design)
fit2 <- contrasts.fit(fit, cont.VEGF)
fit2 <- ebayes(fit2)
topTable(fit2)

but I seem to have the wrong contrast specification, as I got the following error:

Error in makeContrasts(contrasts = "cond0:timePoints-cond1:timePoints",  :
  The levels must by syntactically valid names in R, see help(make.names).  Non-valid names: cond0:timePoints,cond1:timePoints

Is it the semi-colon? If so, how do I change this in the design matrix specification (where the semicolon comes from), while still getting the same model fitting?

Apologies if this question is already answered elsewhere, I could only find parts of it in the list archives and the documentation.

With kind regards,

Linus Schumacher


	[[alternative HTML version deleted]]