Answers.Body: <p>It is actually simple. The first two conditions are nested with the Pt's starting with m and the other two conditions are nested with those PT's starting with c. That means you cannot run this design. The only solution is to subset the data into two chunks and analyze separately.</p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2022-07-26T12:47:42.553Z
Answers.LastActivityDate: 2022-07-26T12:47:42.553Z
Answers.UserDisplayName: charcool88
Answers.Body: <p>The <a href="https://www.bioconductor.org/packages/devel/bioc/vignettes/DESeq2/inst/doc/DESeq2.html" rel="nofollow noreferrer">DESeq2 documentation</a> explains this error as follows:</p>
  <blockquote>
  <p>While most experimental designs run easily using design formula, some design formulas can cause problems and result in the DESeq function returning an error with the text: “the model matrix is not full rank, so the model cannot be fit as specified.” There are two main reasons for this problem: either one or more columns in the model matrix are linear combinations of other columns, or there are levels of factors or combinations of levels of multiple factors which are missing samples.</p>
  </blockquote>
  <p>In your case, I think the error is popping up because there are too many defined values in the &quot;Pts&quot; variable. DESeq2 doesn't require you to separately define replicates as a variable; it assumes that different data columns are taken from different biological replicates.</p>
  <p>Getting the models right is tricky, and needs a deeper discussion about what the goal of analysis should be. My guess is that the following should work, but only if you're not interested in comparing expression results from different patients:</p>
  <pre><code>coldata$Pts &lt;- c(
    &quot;m&quot;, &quot;m&quot;, &quot;m&quot;,
    &quot;m&quot;, &quot;m&quot;, &quot;m&quot;,
    &quot;c&quot;, &quot;c&quot;, &quot;c&quot;,
    &quot;c&quot;, &quot;c&quot;, &quot;c&quot;
  )
  coldata$class &lt;- c(
    &quot;Treated&quot;, &quot;Treated&quot;, &quot;Treated&quot;,
    &quot;Distal&quot;,  &quot;Distal&quot;,  &quot;Distal&quot;,
    &quot;CLEFT&quot;,   &quot;CLEFT&quot;,   &quot;CLEFT&quot;,
    &quot;CRIGHT&quot;,  &quot;CRIGHT&quot;,  &quot;CRIGHT&quot;
  )
  </code></pre>
  <p>This may not be what you want to do; other ways to deal with this are explained in the DESeq2 documentation on the Bioconductor website:</p>
  <p><a href="https://www.bioconductor.org/packages/devel/bioc/vignettes/DESeq2/inst/doc/DESeq2.html#model-matrix-not-full-rank" rel="nofollow noreferrer">https://www.bioconductor.org/packages/devel/bioc/vignettes/DESeq2/inst/doc/DESeq2.html#model-matrix-not-full-rank</a></p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2022-07-26T22:11:24.543Z
Answers.LastActivityDate: 2022-07-26T22:11:24.543Z
Answers.UserDisplayName: gringer
Answers.Body: <p>check if <code>Pts</code> and <code>class</code> are of <code>factor</code> class or not by using
  <code>is.factor(coldata$Pts)</code> and <code>is.factor(coldata$class)</code></p>
  <p>if it is not then convert it into <code>factor</code> by</p>
  <p><code>coldata$Pts = factor(coldata$Pts)</code></p>
  <p><code>coldata$class = factor(coldata$class)</code></p>
  <p>Then try to design your matrix (A design matrix is a matrix of values of the grouping variable).
  I think this would work as it works with</p>
  <p><code>design &lt;- model.matrix(~0+coldata$Pts+coldata$class)</code> in <code>limma</code> package.</p>
  <p>You can also check in <a href="https://master.bioconductor.org/packages/release/workflows/vignettes/rnaseqGene/inst/doc/rnaseqGene.html#the-deseqdataset-object-sample-information-and-the-design-formula" rel="nofollow noreferrer">this RNA-seq analysis workflow</a> of <code>DESeq2</code> where they also used nested samples for their design matrix.</p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2022-08-14T20:32:22.837Z
Answers.LastActivityDate: 2022-08-15T05:34:57.340Z
Answers.LastEditDate: 2022-08-15T05:34:57.340Z
Answers.UserDisplayName: Nitesh Shriwash
Answers.Body: <p>Your original was,</p>
  <pre><code>ddsMat &lt;- DESeqDataSetFromMatrix(countData = as.matrix(matrix), # colData = coldata, design = ~ class+Pts)
  </code></pre>
  <p>What happens when you use</p>
  <pre><code>ddsMat &lt;- DESeqDataSetFromMatrix(countData = as.matrix(matrix), colData = coldata, design = ~ class+Pts)
  </code></pre>
  <p>Thus <code>#</code> is removed. The suspicion is <code>#</code> is causing the code to fall over. <code>#</code>  is used for comments in <code>R</code> so the presence of this within a line of code will cause the following code <code>colData = coldata, design = ~ class+Pts</code> to not be interpreted because it can't be seen.</p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2022-07-25T20:31:44.410Z
Answers.Text: THe #was only added by mistake when I wrote the question. In the original code the colData is not commented
Answers.UserDisplayName: GuidoL
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2022-07-25T18:40:44.237Z
Answers.LastActivityDate: 2022-07-25T18:47:31.870Z
Answers.LastEditDate: 2022-07-25T18:47:31.870Z
Answers.UserDisplayName: M__
Body: <p>I'm running a DEG analysis with Deseq2 by specifing the following design that includes 12 samples from 6 pts in 4 different conditions:</p>
  <pre><code>coldata &lt;- data.frame(&quot;sample&quot;=c(
    &quot;T58&quot;, &quot;T86&quot;, &quot;T87&quot;,
    &quot;D58&quot;, &quot;D86&quot;, &quot;D87&quot;,
    &quot;C62&quot;, &quot;C80&quot;, &quot;C84&quot;,
    &quot;C62&quot;, &quot;C80&quot;, &quot;C84&quot;
  ))
  coldata$Pts &lt;- c(
    &quot;m58&quot;, &quot;m86&quot;, &quot;m87&quot;,
    &quot;m58&quot;, &quot;m86&quot;, &quot;m87&quot;,
    &quot;c62&quot;, &quot;c80&quot;, &quot;c84&quot;,
    &quot;c62&quot;, &quot;c80&quot;, &quot;c84&quot;
  )
  coldata$class &lt;- c(
    &quot;Treated&quot;, &quot;Treated&quot;, &quot;Treated&quot;,
    &quot;Distal&quot;, &quot;Distal&quot;, &quot;Distal&quot;,
    &quot;CLEFT&quot;, &quot;CLEFT&quot;, &quot;CLEFT&quot;,
    &quot;CRIGHT&quot;, &quot;CRIGHT&quot;, &quot;CRIGHT&quot;
  )
  </code></pre>
  <p>I would model the differences among the conditions by considering that the same individuals are nested in different conditions</p>
  <p>I tried the following code</p>
  <pre><code>ddsMat &lt;- DESeqDataSetFromMatrix(
    countData = as.matrix(matrix),
    colData = coldata,
    design = ~ class+Pts
  )
  </code></pre>
  <p>but this produces the famous error</p>
  <blockquote>
  <p>&quot;the model matrix is not full rank...&quot;</p>
  </blockquote>
  <p>If I reduce the previous dataset to</p>
  <pre><code>coldata2 &lt;- as.data.frame(c(
    &quot;T58&quot;, &quot;T86&quot;, &quot;T87&quot;,
    &quot;D58&quot;, &quot;D86&quot;, &quot;D87&quot;
  ))
  coldata2$Pts &lt;- c(
    &quot;m58&quot;, &quot;m86&quot;, &quot;m87&quot;,
    &quot;m58&quot;, &quot;m86&quot;, &quot;m87&quot;
  )
  coldata2$class &lt;- c(
    &quot;Treated&quot;, &quot;Treated&quot;, &quot;Treated&quot;,
    &quot;Distal&quot;, &quot;Distal&quot;, &quot;Distal&quot;
  )
  </code></pre>
  <p>then everything goes well and I'm able to create a <code>DESeq</code> object with the following design:</p>
  <pre><code>ddsMat &lt;- DESeqDataSetFromMatrix(
    countData = as.matrix(matrix[,1:6]),
    colData = coldata2,
    design = ~ class+Pts
  )
  </code></pre>
  <p>Please can you provide to me some suggestion on how to solve the problem?</p>
  <p>Thank you</p>
  
Tags: <rna-seq><deseq2>
Title: nested samples in design DEseq2
ThreadText: # nested samples in design DEseq2
  
  (Asked by: username_0 on 2022-07-25)
  
  I'm running a DEG analysis with Deseq2 by specifing the following design that includes 12 samples from 6 pts in 4 different conditions:
  
  ```
  coldata <- data.frame("sample"=c(
    "T58", "T86", "T87",
    "D58", "D86", "D87",
    "C62", "C80", "C84",
    "C62", "C80", "C84"
  ))
  coldata$Pts <- c(
    "m58", "m86", "m87",
    "m58", "m86", "m87",
    "c62", "c80", "c84",
    "c62", "c80", "c84"
  )
  coldata$class <- c(
    "Treated", "Treated", "Treated",
    "Distal", "Distal", "Distal",
    "CLEFT", "CLEFT", "CLEFT",
    "CRIGHT", "CRIGHT", "CRIGHT"
  )
  
  ```
  
  I would model the differences among the conditions by considering that the same individuals are nested in different conditions
  
  I tried the following code
  
  ```
  ddsMat <- DESeqDataSetFromMatrix(
    countData = as.matrix(matrix),
    colData = coldata,
    design = ~ class+Pts
  )
  
  ```
  
  but this produces the famous error
  
  > "the model matrix is not full rank..."
  
  If I reduce the previous dataset to
  
  ```
  coldata2 <- as.data.frame(c(
    "T58", "T86", "T87",
    "D58", "D86", "D87"
  ))
  coldata2$Pts <- c(
    "m58", "m86", "m87",
    "m58", "m86", "m87"
  )
  coldata2$class <- c(
    "Treated", "Treated", "Treated",
    "Distal", "Distal", "Distal"
  )
  
  ```
  
  then everything goes well and I'm able to create a `DESeq` object with the following design:
  
  ```
  ddsMat <- DESeqDataSetFromMatrix(
    countData = as.matrix(matrix[,1:6]),
    colData = coldata2,
    design = ~ class+Pts
  )
  
  ```
  
  Please can you provide to me some suggestion on how to solve the problem?
  
  Thank you
  
  **Tags:** rna-seq, deseq2
  
  
  ---
  
  ## 4 answers
  
  ### Answer by username_1 on 2022-07-26. Score: 2
  
  It is actually simple. The first two conditions are nested with the Pt's starting with m and the other two conditions are nested with those PT's starting with c. That means you cannot run this design. The only solution is to subset the data into two chunks and analyze separately.
  
  
  ---
  
  ### Answer by username_2 on 2022-07-26. Score: 1
  
  The [DESeq2 documentation](https://www.bioconductor.org/packages/devel/bioc/vignettes/DESeq2/inst/doc/DESeq2.html) explains this error as follows:
  
  > While most experimental designs run easily using design formula, some design formulas can cause problems and result in the DESeq function returning an error with the text: “the model matrix is not full rank, so the model cannot be fit as specified.” There are two main reasons for this problem: either one or more columns in the model matrix are linear combinations of other columns, or there are levels of factors or combinations of levels of multiple factors which are missing samples.
  
  In your case, I think the error is popping up because there are too many defined values in the "Pts" variable. DESeq2 doesn't require you to separately define replicates as a variable; it assumes that different data columns are taken from different biological replicates.
  
  Getting the models right is tricky, and needs a deeper discussion about what the goal of analysis should be. My guess is that the following should work, but only if you're not interested in comparing expression results from different patients:
  
  ```
  coldata$Pts <- c(
    "m", "m", "m",
    "m", "m", "m",
    "c", "c", "c",
    "c", "c", "c"
  )
  coldata$class <- c(
    "Treated", "Treated", "Treated",
    "Distal",  "Distal",  "Distal",
    "CLEFT",   "CLEFT",   "CLEFT",
    "CRIGHT",  "CRIGHT",  "CRIGHT"
  )
  
  ```
  
  This may not be what you want to do; other ways to deal with this are explained in the DESeq2 documentation on the Bioconductor website:
  
  <https://www.bioconductor.org/packages/devel/bioc/vignettes/DESeq2/inst/doc/DESeq2.html#model-matrix-not-full-rank>
  
  
  ---
  
  ### Answer by username_3 on 2022-08-14. Score: 1
  
  check if `Pts` and `class` are of `factor` class or not by using
  `is.factor(coldata$Pts)` and `is.factor(coldata$class)`
  
  if it is not then convert it into `factor` by
  
  `coldata$Pts = factor(coldata$Pts)`
  
  `coldata$class = factor(coldata$class)`
  
  Then try to design your matrix (A design matrix is a matrix of values of the grouping variable).
  I think this would work as it works with
  
  `design <- model.matrix(~0+coldata$Pts+coldata$class)` in `limma` package.
  
  You can also check in [this RNA-seq analysis workflow](https://master.bioconductor.org/packages/release/workflows/vignettes/rnaseqGene/inst/doc/rnaseqGene.html#the-deseqdataset-object-sample-information-and-the-design-formula) of `DESeq2` where they also used nested samples for their design matrix.
  
  
  ---
  
  ### Answer by username_4 on 2022-07-25. Score: 0
  
  Your original was,
  
  ```
  ddsMat <- DESeqDataSetFromMatrix(countData = as.matrix(matrix), # colData = coldata, design = ~ class+Pts)
  
  ```
  
  What happens when you use
  
  ```
  ddsMat <- DESeqDataSetFromMatrix(countData = as.matrix(matrix), colData = coldata, design = ~ class+Pts)
  
  ```
  
  Thus `#` is removed. The suspicion is `#` is causing the code to fall over. `#` is used for comments in `R` so the presence of this within a line of code will cause the following code `colData = coldata, design = ~ class+Pts` to not be interpreted because it can't be seen.
  
  > **Comments:**
  > **username_0**: THe #was only added by mistake when I wrote the question. In the original code the colData is not commented  
  
  ---
  
  