Correct use of arrayWeights alongside duplicateCorrelation for microarrays in Limma?
1
1
Entering edit mode
willj ▴ 30
@willj-8763
Last seen 6.6 years ago
France

I'm trying to use arrayWeights for the experimental design described below which features technical replicates and the duplicateCorrelation function. I would greatly appreciate it if someone could confirm that what i'm doing is sensible.

I want to determine differential expression for mutant v control cell lines. The study features 5 cell lines (2 control and 3 mutant) with 4 'technical' (cell culture) replicates per line, as shown:

sample    reps    groups
Line1Rep1    1    Control
Line1Rep2    1    Control
Line1Rep3    1    Control
Line1Rep4    1    Control
Line2Rep1    2    Control
Line2Rep2    2    Control
Line2Rep3    2    Control
Line2Rep4    2    Control
Line3Rep1    3    Mutant
Line3Rep2    3    Mutant
Line3Rep3    3    Mutant
Line3Rep4    3    Mutant
Line4Rep1    4    Mutant
Line4Rep2    4    Mutant
Line4Rep3    4    Mutant
Line4Rep4    4    Mutant
Line5Rep1    5    Mutant
Line5Rep2    5    Mutant
Line5Rep3    5    Mutant
Line5Rep4    5    Mutant

 

I created a matrix as follows:

groups <- factor(rep(c("Control","Mutant"), times = c(8,12)))
reps <- factor(rep(c(1,2,3,4,5), each=4))
design <- model.matrix(~0+groups)

 

My normalized data are called res.filt2

colnames(design) <- levels(groups)
aw <- arrayWeights(res.filt2,design)
corfit <- duplicateCorrelation(res.filt2, design, block = reps)
fit <- lmFit(res.filt2, design, block = reps, cor = corfit$consensus, weights=aw)
contrasts <- makeContrasts(Mutant-Control, levels=design)
contr.fit <- eBayes(contrasts.fit(fit,contrasts))

If I understand well, the code above calculates arrayWeights, taking into account both the reps and the groups. What I'm curious about is how does lmFit then use these arrayWeights in combination with the blocks of technical replicates? Is it correct to combine array weights along with the block parameter in lmFit in the way that I've done here?

Also, I notice that the duplicateCorrelation function has a "weights" parameter but I have not used it because (if I've understood correctly) it is for individual spot weights and not for array weights.

Many thanks for any help.

limma duplicatecorrelation arrayweights • 1.4k views
ADD COMMENT
1
Entering edit mode
@ryan-c-thompson-5618
Last seen 8 months ago
Scripps Research, La Jolla, CA

I think that you actually should pass the weights to duplicateCorrelation, because the documentation says "If smaller than ‘object’ then it will be filled out the same size." Presumably every spot on an array gets the weight of that array if you don't have additional information about the relative weights of each spot/probe.

The main question I would have is whether the array weight calculation would theoretically have to take account of the duplicate correlation and vice versa. If both are true, then you might have to iterate them back and forth until convergence or something like that. It would a similar problem to the one solved by voomWithQualityWeights (which is for RNA-seq data), where the authors determined that alternating voom and arrayWeights just twice produced sufficiently stable results (i.e. close enough to iterating until convergence). In this case, I don't see a way to pass the duplicate correlation to arrayWeights, but that might just mean it hasn't been implemented.

ADD COMMENT
0
Entering edit mode

Dear Ryan, thanks for your advice.

I've repeated the analysis, this time using:

corfit <- duplicateCorrelation(res.filt2, design, block = reps, weights=aw)

The total number of differentially expressed genes is now more than doubled, but includes all of the DE genes identified in the previous analysis. Passing the weights to duplicateCorrelation seemed to preserve approximately the same ranking (when the probes are ranked by p-value) but makes the p-values a bit stronger.

 

ADD REPLY

Login before adding your answer.

Traffic: 870 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6