<p>Your table is not very clear. If I understood it correctly, there are three groups of t-cells, SPtetramer+, Ly49- and Ly49+ and within it there is a MOGSP vs MOG comparison. I would set colData like this:</p>
<pre><code>                   name             Ly treatment
1  SPtetramer+CD8+TCell SPtetramerplus     MOGSP
2  SPtetramer+CD8+TCell SPtetramerplus     MOGSP
3  SPtetramer+CD8+TCell SPtetramerplus       MOG
4  SPtetramer+CD8+TCell SPtetramerplus       MOG
5  SPtetramer+CD8+TCell SPtetramerplus       MOG
6        Ly49-CD8+TCell      Ly49minus     MOGSP
7        Ly49+CD8+TCell       Ly49plus     MOGSP
8        Ly49-CD8+TCell      Ly49minus     MOGSP
9        Ly49+CD8+TCell       Ly49plus     MOGSP
10       Ly49-CD8+TCell      Ly49minus     MOGSP
11       Ly49+CD8+TCell       Ly49plus     MOGSP
</code></pre>
<p>I generate an example dataset:</p>
<pre><code>counts = matrix(rnbinom(11000,mu=100,size=1),ncol=11)
</code></pre>
<p>And in the design you specify these two factors, meaning gene expression can be affected by these two conditions:</p>
<pre><code>dds = DESeqDataSetFromMatrix(counts,colData=colData,~Ly +treatment)
dds = DESeq(dds)
</code></pre>
<p>You can compare them like this:</p>
<pre><code>head(results(dds,c(&quot;Ly&quot;,&quot;Ly49plus&quot;,&quot;Ly49minus&quot;)))
log2 fold change (MLE): Ly Ly49plus vs Ly49minus 
Wald test p-value: Ly Ly49plus vs Ly49minus 
DataFrame with 6 rows and 6 columns
          baseMean     log2FoldChange            lfcSE               stat
         &lt;numeric&gt;          &lt;numeric&gt;        &lt;numeric&gt;          &lt;numeric&gt;
1 91.4458114779915  0.712987165797688 1.28291188260395  0.555756927241586
2 81.7945234843954 -0.418266530066077 1.14378290088875 -0.365686993345568
3 94.4351681904941 -0.775592496878288 1.32965355129901 -0.583304196886903
4 107.484335203298  -1.36694885973712 1.38595911743023 -0.986283680770933
5 61.2205696982841 -0.760219864564472 1.16779888688769 -0.650985262188887
6 95.5908168437137  0.497427543096335 1.12768919469662   0.44110340458672
             pvalue              padj
          &lt;numeric&gt;         &lt;numeric&gt;
1 0.578377034355333 0.977537706056193
2    0.714598652706 0.992321035598301
3 0.559688537757548 0.977537706056193
4 0.323993925580938 0.954772990635094
5 0.515056000552832 0.977537706056193
6 0.659138138867294 0.992321035598301
</code></pre>
<p>To compare groups within another factor:</p>
<pre><code>results(dds,c(&quot;treatment&quot;,&quot;MOGSP&quot;,&quot;MOG&quot;))
log2 fold change (MLE): treatment MOGSP vs MOG 
Wald test p-value: treatment MOGSP vs MOG 
DataFrame with 6 rows and 6 columns
          baseMean     log2FoldChange            lfcSE               stat
         &lt;numeric&gt;          &lt;numeric&gt;        &lt;numeric&gt;          &lt;numeric&gt;
1 91.4458114779915 -0.875905107403943 1.43007567177694 -0.612488642867122
2 81.7945234843954  0.595588650280833 1.28109080878515  0.464907441530728
3 94.4351681904941 0.0255578400583986 1.48592401497315 0.0171999643325371
4 107.484335203298  -1.46652057431159 1.55285829334496 -0.944400774105796
5 61.2205696982841  -1.27442596494316 1.30489822606005  -0.97664778715433
6 95.5908168437137  -2.35402796630456 1.26608624796739  -1.85929510733078
              pvalue              padj
           &lt;numeric&gt;         &lt;numeric&gt;
1  0.540214509839105 0.944881849859735
2  0.641997741235223 0.973651054534338
3  0.986277090644156 0.994944323499595
4  0.344964886056295  0.87421063696937
5  0.328743552427222                NA
6 0.0629853201285643 0.525889480109579
</code></pre>
<p>The colData i used:</p>
<pre><code>structure(list(name = structure(c(3L, 3L, 3L, 3L, 3L, 1L, 2L, 
1L, 2L, 1L, 2L), .Label = c(&quot;Ly49-CD8+TCell&quot;, &quot;Ly49+CD8+TCell&quot;, 
&quot;SPtetramer+CD8+TCell&quot;), class = &quot;factor&quot;), Ly = c(&quot;SPtetramerplus&quot;, 
&quot;SPtetramerplus&quot;, &quot;SPtetramerplus&quot;, &quot;SPtetramerplus&quot;, &quot;SPtetramerplus&quot;, 
&quot;Ly49minus&quot;, &quot;Ly49plus&quot;, &quot;Ly49minus&quot;, &quot;Ly49plus&quot;, &quot;Ly49minus&quot;, 
&quot;Ly49plus&quot;), treatment = structure(c(2L, 2L, 1L, 1L, 1L, 2L, 
2L, 2L, 2L, 2L, 2L), .Label = c(&quot;MOG&quot;, &quot;MOGSP&quot;), class = &quot;factor&quot;)), class = &quot;data.frame&quot;, row.names = c(NA, 
-11L))
</code></pre>

CC BY-SA 4.0
2020-09-16T05:27:46.427Z
so it can used for multiple variable or cell type comparisons instead of many pairwise comparison which i do normally
kcm
CC BY-SA 4.0
2020-09-16T16:15:23.957Z
it depends on the design but in general yes. In this example, you can see it's necessary because the groups have some overlap
StupidWolf
CC BY-SA 4.0
2021-01-21T06:36:38.923Z
conceptual doubt so if i have to perform multiple pairwise comparison do i need to set reference level for every comparison or i can simply use "results(dds,c("treatment","MOGSP","MOG"))" or the contrast function?
kcm
CC BY-SA 4.0
2021-01-22T01:09:10.543Z
using ```results(dds,c("treatment","MOGSP","MOG"))``` takes care of everything. The reference level of the factor does not matter here
StupidWolf
CC BY-SA 4.0
2021-01-22T11:12:18.613Z
will try this..
kcm
CC BY-SA 4.0
2020-09-13T18:34:53.213Z
2020-09-13T18:34:53.213Z
StupidWolf
<p><em>This question was also asked on <a href="https://www.biostars.org/p/461091/" rel="nofollow noreferrer">BioStars</a></em></p>
<p>I am trying to recreate <a href="https://www.nature.com/articles/s41586-019-1467-x/figures/6" rel="nofollow noreferrer">this heatmap</a>. How would I compare the four variables: Ly49+/- and MOG/MOGSP to get a single heatmap? Thank you for helping me through this.</p>
<p>I need to compare both Ly49+ versus Ly49− and MOG versus MOG plus SP. Column data in at the bottom. The function I used for DESeq2 is,</p>
<pre class="lang-r prettyprint-override"><code>dds = DESeqDataSetFromMatrix(countData=countData,
                                colData=colData,
                                design= ~ cell_treatment)
</code></pre>
<p>I am not sure if the design is right. For the results, how should I write the contrasts? I am very confused and would really appreciate any help!</p>
<pre class="lang-r prettyprint-override"><code>DataFrame with 11 rows and 2 columns
cell treatment
&lt;character&gt; &lt;character&gt;
MOGSP_2 SPtetramer+CD8+TCell MOGSP
MOGSP_3 SPtetramer+CD8+TCell MOGSP
MOG_3 SPtetramer+CD8+TCell MOG
MOG_4 SPtetramer+CD8+TCell MOG
MOG_5 SPtetramer+CD8+TCell MOG
Ly49N_1 Ly49-CD8+TCell MOGSP
Ly49P_1 Ly49+CD8+TCell MOGSP
Ly49N_2 Ly49-CD8+TCell MOGSP
Ly49P_2 Ly49+CD8+TCell MOGSP
Ly49N_3 Ly49-CD8+TCell MOGSP
Ly49P_3 Ly49+CD8+TCell MOGSP
</code></pre>
<p>Thank you!</p>

<rna-seq><deseq2><design>
Setting Contrast for DESeq2 results
# Setting Contrast for DESeq2 results

(Asked by: username_0 on 2020-09-13)

This question was also asked on BioStars
I am trying to recreate this heatmap. How would I compare the four variables: Ly49+/- and MOG/MOGSP to get a single heatmap? Thank you for helping me through this.
I need to compare both Ly49+ versus Ly49− and MOG versus MOG plus SP. Column data in at the bottom. The function I used for DESeq2 is,
dds = DESeqDataSetFromMatrix(countData=countData,
                                colData=colData,
                                design= ~ cell_treatment)

I am not sure if the design is right. For the results, how should I write the contrasts? I am very confused and would really appreciate any help!
DataFrame with 11 rows and 2 columns
cell treatment
 
MOGSP_2 SPtetramer+CD8+TCell MOGSP
MOGSP_3 SPtetramer+CD8+TCell MOGSP
MOG_3 SPtetramer+CD8+TCell MOG
MOG_4 SPtetramer+CD8+TCell MOG
MOG_5 SPtetramer+CD8+TCell MOG
Ly49N_1 Ly49-CD8+TCell MOGSP
Ly49P_1 Ly49+CD8+TCell MOGSP
Ly49N_2 Ly49-CD8+TCell MOGSP
Ly49P_2 Ly49+CD8+TCell MOGSP
Ly49N_3 Ly49-CD8+TCell MOGSP
Ly49P_3 Ly49+CD8+TCell MOGSP

Thank you!

**Tags:** rna-seq, deseq2, design

---

## 1 answers

### Answer by username_1 on 2020-09-13. Score: 3 | Accepted answer

Your table is not very clear. If I understood it correctly, there are three groups of t-cells, SPtetramer+, Ly49- and Ly49+ and within it there is a MOGSP vs MOG comparison. I would set colData like this:
                   name             Ly treatment
1  SPtetramer+CD8+TCell SPtetramerplus     MOGSP
2  SPtetramer+CD8+TCell SPtetramerplus     MOGSP
3  SPtetramer+CD8+TCell SPtetramerplus       MOG
4  SPtetramer+CD8+TCell SPtetramerplus       MOG
5  SPtetramer+CD8+TCell SPtetramerplus       MOG
6        Ly49-CD8+TCell      Ly49minus     MOGSP
7        Ly49+CD8+TCell       Ly49plus     MOGSP
8        Ly49-CD8+TCell      Ly49minus     MOGSP
9        Ly49+CD8+TCell       Ly49plus     MOGSP
10       Ly49-CD8+TCell      Ly49minus     MOGSP
11       Ly49+CD8+TCell       Ly49plus     MOGSP

I generate an example dataset:
counts = matrix(rnbinom(11000,mu=100,size=1),ncol=11)

And in the design you specify these two factors, meaning gene expression can be affected by these two conditions:
dds = DESeqDataSetFromMatrix(counts,colData=colData,~Ly +treatment)
dds = DESeq(dds)

You can compare them like this:
head(results(dds,c("Ly","Ly49plus","Ly49minus")))
log2 fold change (MLE): Ly Ly49plus vs Ly49minus 
Wald test p-value: Ly Ly49plus vs Ly49minus 
DataFrame with 6 rows and 6 columns
          baseMean     log2FoldChange            lfcSE               stat
                                     
1 91.4458114779915  0.712987165797688 1.28291188260395  0.555756927241586
2 81.7945234843954 -0.418266530066077 1.14378290088875 -0.365686993345568
3 94.4351681904941 -0.775592496878288 1.32965355129901 -0.583304196886903
4 107.484335203298  -1.36694885973712 1.38595911743023 -0.986283680770933
5 61.2205696982841 -0.760219864564472 1.16779888688769 -0.650985262188887
6 95.5908168437137  0.497427543096335 1.12768919469662   0.44110340458672
             pvalue              padj
                   
1 0.578377034355333 0.977537706056193
2    0.714598652706 0.992321035598301
3 0.559688537757548 0.977537706056193
4 0.323993925580938 0.954772990635094
5 0.515056000552832 0.977537706056193
6 0.659138138867294 0.992321035598301

To compare groups within another factor:
results(dds,c("treatment","MOGSP","MOG"))
log2 fold change (MLE): treatment MOGSP vs MOG 
Wald test p-value: treatment MOGSP vs MOG 
DataFrame with 6 rows and 6 columns
          baseMean     log2FoldChange            lfcSE               stat
                                     
1 91.4458114779915 -0.875905107403943 1.43007567177694 -0.612488642867122
2 81.7945234843954  0.595588650280833 1.28109080878515  0.464907441530728
3 94.4351681904941 0.0255578400583986 1.48592401497315 0.0171999643325371
4 107.484335203298  -1.46652057431159 1.55285829334496 -0.944400774105796
5 61.2205696982841  -1.27442596494316 1.30489822606005  -0.97664778715433
6 95.5908168437137  -2.35402796630456 1.26608624796739  -1.85929510733078
              pvalue              padj
                    
1  0.540214509839105 0.944881849859735
2  0.641997741235223 0.973651054534338
3  0.986277090644156 0.994944323499595
4  0.344964886056295  0.87421063696937
5  0.328743552427222                NA
6 0.0629853201285643 0.525889480109579

The colData i used:
structure(list(name = structure(c(3L, 3L, 3L, 3L, 3L, 1L, 2L, 
1L, 2L, 1L, 2L), .Label = c("Ly49-CD8+TCell", "Ly49+CD8+TCell", 
"SPtetramer+CD8+TCell"), class = "factor"), Ly = c("SPtetramerplus", 
"SPtetramerplus", "SPtetramerplus", "SPtetramerplus", "SPtetramerplus", 
"Ly49minus", "Ly49plus", "Ly49minus", "Ly49plus", "Ly49minus", 
"Ly49plus"), treatment = structure(c(2L, 2L, 1L, 1L, 1L, 2L, 
2L, 2L, 2L, 2L, 2L), .Label = c("MOG", "MOGSP"), class = "factor")), class = "data.frame", row.names = c(NA, 
-11L))

> **Comments:**
> **username_2**: so it can used for multiple variable or cell type comparisons instead of many pairwise comparison which i do normally  
> **username_1**: it depends on the design but in general yes. In this example, you can see it's necessary because the groups have some overlap  
> **username_2**: conceptual doubt so if i have to perform multiple pairwise comparison do i need to set reference level for every comparison or i can simply use "results(dds,c("treatment","MOGSP","MOG"))" or the contrast function?  
> **username_1**: using ```results(dds,c("treatment","MOGSP","MOG"))``` takes care of everything. The reference level of the factor does not matter here  
> **username_2**: will try this..  

---

