Tutorials › High-Dimensional Genomic Data Analysis › Design Matrices For Differential Expression

Count-Based Differential Expression for RNA-seq · Tutorial 427 of 795

Design Matrices For Differential Expression

Build two-group and multi-group RNA-seq design matrices, interpret their coefficients alongside sample-specific offsets, and check that the coding makes the intended group effects identifiable.

Intermediate 11 min read

What You'll Learn

  • Arrange design-matrix rows in the same sample order as the count matrix columns.
  • Construct two-group and multi-group matrices using reference and group-means coding.
  • Distinguish coefficient contributions from the full linear predictor, which also includes the offset.
  • Check whether a proposed design is full rank and its coefficients are identifiable.
  • Diagnose common coding, sample-alignment, and missing-group problems.

From Sample Groups to a Design Matrix

A differential-expression model is fitted separately for each gene, but the sample-level experiment determines the predictors used in every one of those gene-specific models. The design matrix is the bookkeeping device that makes this connection explicit: each row describes one sample's predictor values, and each column represents a model term. A design matrix is not a table of counts, and it does not contain the response values themselves.

For \(m\) genes and \(n\) samples, let \(Y_{ij}\) be the count for gene \(i\) in sample \(j\). Write \(X\) for the \(n\times p\) design matrix, with sample \(j\)'s row denoted by \(\mathbf{x}_j^\mathsf{T}\). In the log-linked count model established in the earlier material, the expected count \(\mu_{ij}\) is represented as

$$ \log(\mu_{ij})=\log(s_j)+\mathbf{x}_j^\mathsf{T}\boldsymbol{\beta}_i. $$

Here \(s_j\) is the sample-specific size factor, and \(\boldsymbol{\beta}_i\) is the \(p\)-component coefficient vector for gene \(i\). The same \(X\) is used for every gene, while the coefficients vary from gene to gene. The term \(\mathbf{x}_j^\mathsf{T}\boldsymbol{\beta}_i\) is the sample's design contribution to the linear predictor. The full linear predictor also includes the offset \(\log(s_j)\). Keeping those parts distinct is important when interpreting any coding scheme.

Definition: A design matrix \(X\) has one row per sample and one column per model term. Its entry \(x_{jk}\) records the value of predictor \(k\) for sample \(j\). The row's product with the gene-specific coefficient vector gives that sample's design contribution, \(\mathbf{x}_j^\mathsf{T}\boldsymbol{\beta}_i\); adding the offset gives the full linear predictor.

The model structure described here assumes that each sample has been assigned the intended predictor values, that the row order of \(X\) matches the sample order of the count matrix, and that the specified log-mean relationship is appropriate for the comparison. For a simple group-only design, group categories must be mutually exclusive and correctly recorded. The design must also contain enough independent columns to identify its coefficients. Constructing \(X\) does not itself verify the count distribution, estimate coefficients, or establish that samples are independent; those are separate modeling questions.

Two Groups: Reference Coding

Suppose samples belong to a control group or a treatment group. A common coding uses an intercept and a treatment indicator, with control as the reference group. The indicator equals zero for control and one for treatment.

$$ X= \begin{pmatrix} 1&0\\ 1&0\\ 1&1\\ 1&1 \end{pmatrix}, \qquad \boldsymbol{\beta}_i= \begin{pmatrix} \beta_{i0}\\ \beta_{i1} \end{pmatrix}. $$

The first two rows represent control samples, and the last two represent treatment samples. Multiplication by the coefficient vector gives a design contribution of \(\beta_{i0}\) for a control sample and \(\beta_{i0}+\beta_{i1}\) for a treatment sample. Thus \(\beta_{i0}\) is the reference group's design contribution, and \(\beta_{i1}\) is the treatment indicator's contribution relative to that reference. For sample \(j\), the full linear predictor is \(\log(s_j)+\beta_{i0}\) in control and \(\log(s_j)+\beta_{i0}+\beta_{i1}\) in treatment.

This distinction also clarifies how the offset works: samples in the same group share a design contribution under this model, but their offsets can differ. The offset is added for each sample; it is not another group coefficient.

Worked Example: A Two-Group Matrix and Its Predictions

Consider four hypothetical samples, listed in the same order as the columns of the gene-by-sample count matrix: C1, C2, T1, T2. Their groups are control, control, treatment, treatment. The size factors are \(1.00, 1.20, 0.80,\) and \(1.10\), respectively. For illustration, suppose a particular gene has coefficients \(\beta_{i0}=\log(20)\) and \(\beta_{i1}=\log(1.5)\). These are hypothetical values, not estimates from a real study.

The intercept-and-indicator design is

$$ X= \begin{pmatrix} 1&0\\ 1&0\\ 1&1\\ 1&1 \end{pmatrix}. $$

For C2, the design contribution is \(1\cdot\log(20)+0\cdot\log(1.5)=\log(20)\). The full linear predictor includes its offset: \(\log(1.20)+\log(20)=\log(24)\). Therefore its expected count is \(\mu_{i,\mathrm{C2}}=\exp(\log(24))=24\) counts.

For T2, the design contribution is \(\log(20)+\log(1.5)=\log(30)\). Including its size-factor offset gives \(\log(1.10)+\log(30)=\log(33)\), so its expected count is \(33\) counts. The two samples have different expected counts partly because their offsets differ, even though both have the same treatment design contribution.

The calculation follows directly from exponentiating the full linear predictor: \(\mu_{ij}=s_j\exp(\mathbf{x}_j^\mathsf{T}\boldsymbol{\beta}_i)\). In this example, the control design contribution corresponds to \(20\) after exponentiation, and the treatment design contribution corresponds to \(30\); the sample-specific factor \(s_j\) is then applied.

More Than Two Groups

For three or more groups, the same principles apply. With groups A, B, and C, reference coding can use A as the reference and include indicators for B and C. A sample from A has row \((1,0,0)\); a sample from B has row \((1,1,0)\); a sample from C has row \((1,0,1)\). The coefficients are gene-specific: \(\beta_{i0}\) is the A design contribution, while the B and C indicators add their respective contributions relative to A. Each sample's offset still has to be added to obtain its full linear predictor.

An alternative is group-means coding, also called cell-means coding: use one indicator column for each group and omit the intercept. A sample in group A has row \((1,0,0)\), one in B has row \((0,1,0)\), and one in C has row \((0,0,1)\). The coefficients then give the separate group-specific design contributions. They are not, by themselves, full sample-specific linear predictors, because each sample still has its own offset.

Worked Example: Three Groups With Reference Coding

Suppose five hypothetical samples, in count-matrix column order, belong to A, B, A, C, and B. Their size factors are \(1.00, 0.90, 1.10, 1.30,\) and \(1.00\). Using A as the reference gives the following rows:

$$ X= \begin{pmatrix} 1&0&0\\ 1&1&0\\ 1&0&0\\ 1&0&1\\ 1&1&0 \end{pmatrix}. $$

For a hypothetical gene, let \(\beta_{i0}=\log(10)\), \(\beta_{iB}=\log(1.4)\), and \(\beta_{iC}=\log(0.7)\). In group B, the design contribution is \(\log(10)+\log(1.4)=\log(14)\). For the first B sample, whose size factor is \(0.90\), the full linear predictor is \(\log(0.90)+\log(14)=\log(12.6)\), so its expected count is \(12.6\).

For the C sample, the design contribution is \(\log(10)+\log(0.7)=\log(7)\). Its size factor is \(1.30\), so its full linear predictor is \(\log(1.30)+\log(7)=\log(9.1)\), and its expected count is \(9.1\). The design contributions identify the group terms; the different size factors scale the expected counts for their particular samples.

With group-means coding, the same sample order would instead produce rows \((1,0,0)\), \((0,1,0)\), \((1,0,0)\), \((0,0,1)\), and \((0,1,0)\). The coefficient vector would contain a separate design contribution for each group. The two codings describe the same group structure when their columns are properly specified; they label the coefficients differently.

Why Identifiability and Sample Order Matter

A design matrix must have full column rank for all its coefficients to be separately identifiable. In practical terms, no column can be written as an exact linear combination of the other columns. If columns are redundant, different coefficient vectors can produce the same fitted design contributions, so the data cannot distinguish those coefficients.

For example, with two groups, including an intercept and indicators for both control and treatment creates three columns when only two independent group patterns exist. The intercept column equals the sum of the two group-indicator columns:

$$ \begin{pmatrix}1\\1\\1\\1\end{pmatrix} = \begin{pmatrix}1\\1\\0\\0\end{pmatrix} + \begin{pmatrix}0\\0\\1\\1\end{pmatrix}. $$

This is an exact dependence, so the three coefficients cannot all be uniquely determined. Use an intercept and omit one group's indicator, or omit the intercept and include indicators for both groups. These are two ways to avoid that particular redundancy.

Correct sample order is just as essential. A valid matrix with rows in the wrong order assigns predictors to the wrong count columns. The resulting computation may be numerically possible, yet it models a different and incorrectly labelled experiment. Always verify sample identifiers and row order rather than relying on visual similarity.

Worked Example: Detecting Redundant Columns

Suppose a proposed two-group design has an intercept, a control indicator, and a treatment indicator, with samples ordered C1, C2, T1, T2. Its rows are \((1,1,0)\), \((1,1,0)\), \((1,0,1)\), and \((1,0,1)\). Let its three coefficients for a gene be \(a,b,c\). A control row gives design contribution \(a+b\), while a treatment row gives \(a+c\).

For any number \(t\), changing the coefficients to \((a+t,b-t,c-t)\) leaves both contributions unchanged: \(a+t+b-t=a+b\), and \(a+t+c-t=a+c\). Thus infinitely many coefficient vectors give the same design contributions. The issue is not a small sample size or an unusual count; it is exact column redundancy.

Removing the intercept and retaining the two group indicators gives rows \((1,0)\) for control and \((0,1)\) for treatment. The resulting coefficients \(b\) and \(c\) directly label the two group design contributions. Alternatively, retaining the intercept and only the treatment indicator makes control the reference. Either choice avoids the redundant three-column design.

Choosing and Checking a Coding

Reference coding is convenient when one group is a natural comparison baseline, because coefficients are expressed relative to that group. Group-means coding can be easier to read when the main goal is to represent each group's design contribution directly. Neither coding changes the sample-level group assignments, counts, or offsets; it changes the coordinate system used for the coefficients. A coefficient's label and interpretation must therefore be read from the actual columns of \(X\), not inferred from software defaults.

FeatureReference codingGroup-means coding
Columns for \(K\) groupsIntercept plus \(K-1\) group indicators\(K\) group indicators and no intercept
Coefficient labelsReference-group contribution and additions for other groupsOne design contribution per group
Full linear predictorOffset plus the row's design contributionOffset plus the row's design contribution
Main cautionState the reference group explicitlyDo not also include an intercept with all group indicators

A practical check is to write down the group represented by each row pattern and multiply that row by a symbolic coefficient vector. This simple calculation reveals what each coefficient contributes and whether the intercept is present. Also check that every intended group occurs in the data: a group with no samples cannot have a separately supported group contribution, and its indicator column contains only zeros.

Several failure signs are worth investigating before fitting a model:

  • Unexpected row patterns. A sample may have been assigned the wrong group, or the sample order of the design may not match the count matrix.
  • Redundant columns. An intercept plus indicators for every group is a standard source of exact dependence.
  • An all-zero group column. This can indicate a missing group or an unused category. Its coefficient is not estimable from these samples.
  • Coefficient descriptions that omit offsets. A group coefficient or group design contribution is not the full sample linear predictor; add \(\log(s_j)\) for that sample.
  • Overinterpreting coding. Re-labelling the columns does not change the observed data or make a group comparison causal.

A design matrix is a compact representation of the planned comparison, not a guarantee that the comparison is scientifically appropriate. A group-only design does not encode other sources of sample variation. If the experiment's structure requires additional terms, the matrix must represent them; that extension is a separate modeling decision. For the two-group and multi-group designs considered here, the essential checks are clear group labels, aligned sample order, a stated coding scheme, and identifiable columns.

Key takeaway: Build one design row per sample, interpret every coefficient through its column, and distinguish the design contribution \(\mathbf{x}_j^\mathsf{T}\boldsymbol{\beta}_i\) from the full linear predictor \(\log(s_j)+\mathbf{x}_j^\mathsf{T}\boldsymbol{\beta}_i\). Check rank and sample alignment before interpreting model output.

Check Your Understanding

Answer each question, then compare your reasoning with the solutions.

  1. In a two-group reference-coded design with control as the reference, what are the rows for a control sample and a treatment sample?
  2. For a treatment sample with row \((1,1)\), coefficients \((\log(8),\log(2))\), and size factor \(1.25\), calculate its design contribution, full linear predictor, and expected count.
  3. Why is an intercept together with indicators for all three groups redundant?
  4. In group-means coding for groups A, B, and C, what row represents a sample in group B, and what must be added to its design contribution to obtain the full linear predictor?
  5. What is the problem if the rows of a correctly constructed design matrix are in a different order from the count-matrix sample columns?

Solutions. (1) The control row is \((1,0)\), and the treatment row is \((1,1)\). (2) The design contribution is \(\log(8)+\log(2)=\log(16)\). The full linear predictor is \(\log(1.25)+\log(16)=\log(20)\), so the expected count is \(\exp(\log(20))=20\) counts. (3) The intercept column is the sum of the three group-indicator columns, so the columns are linearly dependent and the coefficients are not separately identifiable. (4) The B row is \((0,1,0)\); add that sample's offset \(\log(s_j)\) to the design contribution. (5) The predictor values would be assigned to the wrong samples, so counts would be modeled using incorrect group labels.

Clinical Biostats