Dear Michael,
I am using DESeq2 for a factorial design. I have Celltype and condition as colData. I have used the formula ~ Celltype + Cellypte:condition as suggested by Steve Lianoglou elsewhere and then I look at the single comparisons in only one Celltype

Everything works fine and I have good adj p values, genes make sene...

The problem arises when I want to make an MA plot of the 2 conditions for one cell type. I call the function as DESeq2:: since I known there are some. Here is my code:

> ddsFull

class: DESeqDataSet
dim: 54668 18
exptData(0):
assays(1): counts
rownames(54668): 5S_rRNA 61E3.4 ... snosnR66 yR211F11.2
rowData metadata column names(0):
colnames(18): r3_1394_A r3_2435_A ... r3_1395_A r3_1395_B
colData names(2): condition Celltype

> design(ddsFull)<-formula(~ Celltype + Celltype:condition)

> ddsFull.fact<-DESeq(ddsFull)

estimating size factors
estimating dispersions
gene-wise dispersion estimates
mean-dispersion relationship
final dispersion estimates
fitting model and testing
18 rows did not converge in beta, labelled in mcols(object)$betaConv. Use larger maxit argument with nbinomWaldTest

> resultsNames(ddsFull.fact)
[1] "Intercept"              "Celltype_R6_vs_R3"      "CelltypeR3.conditionKD" "CelltypeR6.conditionKD" "CelltypeR3.conditionOE" "CelltypeR6.conditionOE"

Now, the problem:

If I use plotMA as for the vignette:

> DESeq2::plotMA(ddsFull.fact)

I obtain a right graph but for the last value in the design formula:


But if I try to plot for my contrasts of interest:

> DESeq2::plotMA(ddsFull.fact, lfcColname="CelltypeR3.conditionKD")


It looks as if the red points are spread over with no correlation to fold change, even with red points over the 0 fold change.

I did not find the pvalueColname variable anymore on plotMA, I guess is the same as log2FC.

If I create the results objects by:

> resFull.KD<-results(ddsFull.fact,name="CelltypeR3.conditionKD")
> resFull.OE<-results(ddsFull.fact,name="CelltypeR3.conditionOE")

and try to make my own MA plot by:

> col=ifelse(resFull.KD$padj>=0.1,"black","red3");plot(resFull.KD$baseMean,resFull.KD$log2FoldChange,log="x",col=col)


It looks pretty convincing with the red points being those external to the cloud of maximal fold changes an high baseMean as it should

> sessionInfo()
R version 3.0.2 (2013-09-25)
Platform: x86_64-apple-darwin10.8.0 (64-bit)

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] grid      parallel  stats     graphics  grDevices utils     datasets  methods   base

other attached packages:
 [1] parathyroidSE_1.0.2     pasilla_0.2.16          DESeq_1.13.3            lattice_0.20-23         locfit_1.5-9.1          DEXSeq_1.7.22           BiocInstaller_1.12.0
 [8] gplots_2.11.3           MASS_7.3-29             KernSmooth_2.23-10      caTools_1.14            gdata_2.13.2            gtools_3.1.0            AnnotationDbi_1.23.27
[15] Biobase_2.21.7          RSQLite_0.11.4          DBI_0.2-7               edgeR_3.3.8             limma_3.17.26           DESeq2_1.2.0            RcppArmadillo_0.3.910.0
[22] Rcpp_0.10.4             GenomicRanges_1.13.51   XVector_0.1.4           IRanges_1.19.38         BiocGenerics_0.7.5

loaded via a namespace (and not attached):
 [1] annotate_1.39.0          biomaRt_2.17.3           Biostrings_2.29.19       biovizBase_1.9.4         bitops_1.0-6             BSgenome_1.29.1          cluster_1.14.4
 [8] colorspace_1.2-3         dichromat_2.0-0          digest_0.6.3             evaluate_0.5             formatR_0.9              genefilter_1.43.0        geneplotter_1.39.05
[15] GenomicFeatures_1.13.47  ggplot2_0.9.3.1          GO.db_2.10.1             gridSVG_1.3-1            gtable_0.1.2             Gviz_1.5.15              Hmisc_3.12-2
[22] httpuv_1.1.0             hwriter_1.3              interactiveDisplay_1.0.0 knitr_1.5                labeling_0.2             latticeExtra_0.6-26      munsell_0.4.2
[29] plyr_1.8                 proto_0.3-10             RColorBrewer_1.0-5       RCurl_1.95-4.1           reshape2_1.2.2           RJSONIO_1.0-3            rpart_4.1-3
[36] Rsamtools_1.13.48        rtracklayer_1.21.12      scales_0.2.3             shiny_0.7.0              splines_3.0.2            statmod_1.4.18           stats4_3.0.2
[43] stringr_0.6.2            survival_2.37-4          tools_3.0.2              XML_3.95-0.2             xtable_1.7-1             zlibbioc_1.7.0


The same code had been used with older versions of DESeq2 (<1.13) with good results, the ones that did not contain the 'contrast' parameter in results(). Now with the newest I have this issue.
What can be the trouble?
Thanks



Jose M. Garcia Manteiga PhD

Data analyst in Functional Genomics
Center for Translational Genomics and Bioinformatics
DIBIT2-A3 Room 21
San Raffaele Scientific Institute
Via Olgettina 58
20132 Milano
Italy

Office: +39 02 26439114







	[[alternative HTML version deleted]]

