[Bioc-devel] renameModelMatrixColumns mishap + patch for DESeq2

Steve Lianoglou lianoglou.steve at gene.com
Mon Jul 15 17:45:47 CEST 2013


Thanks, Mike!

On Mon, Jul 15, 2013 at 8:04 AM, Michael Love
<michaelisaiahlove at gmail.com> wrote:
> hi Steve,
>
> I have fixed this bug in version 1.0.18 and 1.1.22, and I added a unit
> test for such a design formula in the devel version.
>
> best,
>
> Mike
>
> On Wed, Jul 10, 2013 at 8:52 AM, Michael Love
> <michaelisaiahlove at gmail.com> wrote:
>> Hi Steve,
>>
>> Thanks for pointing out this bug in the renaming code. I will fix this as
>> soon as possible.
>>
>> Best,
>>
>> Mike
>>
>> On Jul 10, 2013 1:00 AM, "Steve Lianoglou" <lianoglou.steve at gene.com> wrote:
>>>
>>> Greetings DESeq2ers,
>>>
>>> I've been playing around with different ways to encode a 2x2 factorial
>>> design I'm working with: 4 cell types, 4 treatments and exploring
>>> what's cooking in there.
>>>
>>> I thought that the nested interaction formula used in the limma user's
>>> guide (section 8.5.3) would be an easy way to explore some obvious
>>> questions since differential expression from different treatments
>>> within each cell type (among other things) shake out quite cleanly
>>> from the results of running the wald test since these end up as
>>> columns in the design matrix, eg:
>>>
>>> ~ cell + cell:treatment
>>>
>>> as long as the control treatment is the first level of the treatment
>>> factor, we got columns for all cell:treatment effects, eg.
>>> cellA:treatment1, cellA:treatment2, ..., cellD:treatment4)
>>>
>>> Yay!
>>>
>>> ... almost ...
>>>
>>> The problem is that the "renaming column" stuff in the fitNbinomGLMs
>>> function seems to assume that we'll have columns for all main effects
>>> in the design matrix, so the following line triggers an error:
>>>
>>> modelMatrixNames[match(convertNames$from, modelMatrixNames)] <-
>>> convertNames$to
>>>
>>> since you will find elements like "treatment2, ..., treatment4" among
>>> the elements in `covertNames$from`. The call to match then returns NA
>>> for these, and *boom*.
>>>
>>> An easy fix would be to add the following line after the call to
>>> `renameModelMatrixColumns`:
>>>
>>> convertNames <- subset(convertNames, from %in% modelMatrixNames)
>>>
>>> So the "if (renameCols)" block (line 1066 in core.R) now becomes:
>>>
>>>   if (renameCols) {
>>>     names(modelMatrixNames) <- modelMatrixNames
>>>     convertNames <- renameModelMatrixColumns(modelMatrixNames,
>>>
>>> as.data.frame(colData(object)),
>>>                                              modelFormula)
>>>     convertNames <- subset(convertNames, from %in% modelMatrixNames)
>>>     modelMatrixNames[match(convertNames$from, modelMatrixNames)] <-
>>> convertNames$to
>>>   }
>>>
>>>
>>>
>>> --
>>> Steve Lianoglou
>>> Computational Biologist
>>> Bioinformatics and Computational Biology
>>> Genentech
>
> _______________________________________________
> Bioc-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/bioc-devel



--
Steve Lianoglou
Computational Biologist
Bioinformatics and Computational Biology
Genentech



More information about the Bioc-devel mailing list