Question

Weird MA-plot (DESeq2)

0

Entering edit mode

jgarces • 0

@jgarces-20286

Last seen 2.7 years ago

Spain

Hi,

I'm running a RNAseq experiment comparing two different conditions on the same subject... my count matrix is some like that (but much bigger):

         341_MO_r150 341_SP_r150 345_MO_r129 345_SP_r129 350_MO_r178 350_SP_r178
A1BG               0           0           2           0          50           0
A1BG-AS1           0           0           0          11           3           0
A1CF               0           0           0           0          37           0
A2M                1           0           0           0           0           0
A2M-AS1            1           0           0           0           9           0
A2ML1              0           0           1           0           6           0

After creating a DESeq2 object in that way...

> dds <- DESeqDataSetFromMatrix(countData = as.matrix(cnt), colData = mde, design = ~ condition+ subject)
> dds <- dds[rowSums(counts(dds)) >= 10,]
> dds <- DESeq(dds, parallel = T)
> res <- results(dds)
> plotMA(res, ylim = c(-2,2), main = "res")

... I'm getting this (very) weird MA-plot: enter image description here

Any idea what I'm missing, please? Thanks a lot!

deseq2 RNAseq • 1.7k views

ADD COMMENT • link 4.6 years ago jgarces • 0

score 0 · Answer 1 · 2020-04-17

0

Entering edit mode

Michael Love 43k

@mikelove

Last seen 3 days ago

United States

Note that this is not the condition effect, results(), without being told what contrast to perform by the user, will use the last coefficient in the design (see vignette).

ADD COMMENT • link 4.6 years ago Michael Love 43k

0

Entering edit mode

Thanks for so quick response, Michael. Although I didn't indicate a specific contrast, the DE is correctly comparing my two conditions:

log2 fold change (MLE): sample SP vs MO 
Wald test p-value: sample SP vs MO 
DataFrame with 36432 rows and 6 columns
                                baseMean      log2FoldChange              lfcSE               stat             pvalue               padj
                               <numeric>           <numeric>          <numeric>          <numeric>          <numeric>          <numeric>
A1BG                    1.74363594715959  -0.579334448248012  0.828054617068813 -0.699633135672581  0.484156443520277  0.657718235467392
A1BG-AS1                 3.7262462242163  -0.704298317437928  0.581062162809431  -1.21208772919003  0.225478800381989  0.351558642610188
A1CF                   0.436922029696208   0.817919276432277   2.94859379266782  0.277392999492223  0.781478356401314                 NA
A2M                    0.681894561684853     0.3756186780853   2.94778300705695  0.127424127619324  0.898604720812493                 NA
A2M-AS1                0.862822450286881   0.459158270253347   2.94682082534243  0.155814790741474  0.876179011699376                 NA

ADD REPLY • link 4.6 years ago jgarces • 0

0

Entering edit mode

Ok, I just was using the code you pasted which is different than the code you used (~condition+subject above).

Otherwise, the plot above may reflect some confounding. What are the library sizes across the two conditions?

ADD REPLY • link 4.6 years ago Michael Love 43k

0

Entering edit mode

In theory, library sizes are equal since both conditions were simultaneously processed.

ADD REPLY • link 4.6 years ago jgarces • 0

0

Entering edit mode

I'd look into possible confounding of counts with your condition of interest, unless you have reason to believe it is not artifactual that low count genes should always favor SP. Take a look and perhaps consult a bioinformatician who could help diagnose technical problems.

ADD REPLY • link 4.6 years ago Michael Love 43k

0

Entering edit mode

Thanks, Michael. I've plotted the number of reads and the CT values before library preparation and there's an effect on the CT values... do you think that could explain this?

enter image description here

ADD REPLY • link 4.6 years ago jgarces • 0

0

Entering edit mode

I'm not going to be able to help diagnose the problem here - I'm just too busy and can only provide software support right now.