Tutorials › Statistical Methods for Single-Cell RNA Sequencing › Condition-by-Time Interactions

Differential Expression With Pseudobulk and Mixed Models · Tutorial 493 of 700

Condition-by-Time Interactions

Learn to specify, test, and interpret condition-by-time interactions in repeated-measures pseudobulk analyses, including how to distinguish slope differences from time-specific contrasts.

Advanced 11 min read

What You'll Learn

  • Write a negative-binomial pseudobulk model with condition, time, and their interaction.
  • Distinguish a difference in linear time slopes from a joint test of categorical-time interactions.
  • Interpret interaction coefficients as changes in condition contrasts over time.
  • Calculate a Wald test and a joint interaction test from fitted estimates.
  • Recognize design and model problems that can make interaction results unreliable.

Why Test a Condition-by-Time Interaction?

In a repeated-measures experiment, a condition difference at one time point does not by itself tell us whether expression changes differently over time between conditions. The relevant question is whether the condition contrast changes as time advances. A condition-by-time interaction represents that change: it tests whether the expression trajectories differ between conditions, rather than merely whether their average levels differ.

For pseudobulk data, the observations are gene counts aggregated within biological samples, such as donor-by-condition-by-time combinations. The model must account for repeated measurements from the same donor. The earlier material on modeling time trends provides the basic time terms; here the focus is on adding and testing the interaction terms that compare trajectories.

Definition: A condition-by-time interaction is a model term allowing the association between condition and expected expression to vary with time. With two conditions and numeric time, it is commonly represented by the product of a condition indicator and time. With categorical time, it is represented by product terms between condition and each non-reference time indicator.

Model Specification and Interpretation

Let \(Y_{igct}\) be the pseudobulk count for gene \(g\), donor \(i\), condition \(c\), and time \(t\). For a negative-binomial model, the count has conditional mean \(\mu_{igct}\) and variance \(\mu_{igct}+\phi_g\mu_{igct}^2\), where \(\phi_g\) is the gene-specific dispersion. A library-size offset accounts for differences in sequencing depth. A repeated-measures model can include a donor random intercept \(b_i\), with \(b_i\sim N(0,\tau_g^2)\), or another suitable covariance structure.

For two conditions and numeric time, let \(C=0\) denote the reference condition and \(C=1\) the comparison condition. A simple model is:

$$ \log(\mu_{igt})=\log(s_{it})+\beta_{0g}+\beta_{Cg}C+\beta_{Tg}t+\beta_{CTg}(C t)+b_i. $$

Here \(s_{it}\) is the library-size normalization factor, and \(t\) is measured in a stated unit, such as months. Conditional on the donor effect, the model assumes the specified negative-binomial mean-variance relationship and a suitable dependence structure for repeated measurements. The random intercept accounts for donor-specific baseline shifts; it does not automatically guarantee that every aspect of within-donor correlation is adequately modeled.

In the reference condition, the log mean changes by \(\beta_{Tg}\) per time unit. In the comparison condition, the change is \(\beta_{Tg}+\beta_{CTg}\). Thus \(\beta_{CTg}\) is the difference in time slopes on the log-mean scale. Exponentiating it gives the ratio of the conditions’ per-unit multiplicative changes in expected expression. The condition log-mean contrast at time \(t\) is \(\beta_{Cg}+\beta_{CTg}t\), so \(\beta_{Cg}\) alone is the condition contrast only when \(t=0\).

If time is categorical, choose a reference time and include indicators for the other time points, along with condition-by-indicator products. Each interaction coefficient describes how the condition contrast at that time differs from its contrast at the reference time. A joint test of all these coefficients asks whether the condition contrasts are constant across the observed times. This flexible approach does not impose a straight-line trajectory.

Key interpretation: A significant interaction is evidence that the modeled condition contrast changes over time. It does not, on its own, establish that either condition has a nonzero change from baseline, identify which time points differ, or imply a causal effect.

Testing the Interaction

For numeric time with a single interaction coefficient, the gene-level null hypothesis is \(H_0:\beta_{CTg}=0\), against an alternative that it is nonzero. A Wald test divides the estimated coefficient by its standard error. A likelihood-ratio test instead compares the fitted model containing the interaction with a nested model that omits it. For categorical time, the null sets the full vector of interaction coefficients to zero; the joint test is essential because separate coefficient tests do not answer that overall question.

The comparison is meaningful only when the models use the same observations and the reduced model is nested in the full model. Both models should retain the same main effects and the same non-interaction terms. In particular, omitting the main effects while retaining an interaction changes the model’s interpretation and is generally not the intended test.

$$ z=\frac{\widehat{\beta}_{CTg}}{SE(\widehat{\beta}_{CTg})}, \qquad W=(R\widehat{\boldsymbol{\beta}})^\mathsf{T} \left[R\widehat{\mathrm{Var}}(\widehat{\boldsymbol{\beta}})R^\mathsf{T}\right]^{-1} (R\widehat{\boldsymbol{\beta}}). $$

The first statistic is used for a one-coefficient Wald test. The second gives a joint Wald statistic for linear restrictions \(R\boldsymbol{\beta}=0\); under the null and suitable large-sample conditions, it is compared with a chi-square distribution with degrees of freedom equal to the number of independent restrictions. In a likelihood-ratio test, twice the log-likelihood difference is compared with a chi-square distribution whose degrees of freedom equal the number of added parameters, subject to the model’s regularity and fitting assumptions.

Worked Example: Testing a Linear Slope Difference

Suppose a hypothetical repeated-measures pseudobulk analysis measures time in months and fits the numeric-time model above for one gene. The estimated interaction is \(\widehat{\beta}_{CTg}=-0.180\) log-mean units per month, with \(SE=0.070\). Assume the model converged, the design includes both conditions over the time range, and the stated covariance model is appropriate for the repeated donor observations.

The Wald statistic is:

$$ z=\frac{-0.180}{0.070}=-2.571. $$

The two-sided normal-reference p-value is approximately \(0.0101\). A 95% Wald confidence interval for the interaction is:

$$ -0.180\pm1.96(0.070) =-0.180\pm0.1372 =(-0.3172,-0.0428). $$

Because the interval excludes zero and the p-value is below \(0.05\), this hypothetical analysis provides evidence that the log-mean time slopes differ between conditions. The estimated ratio of per-month multiplicative changes is \(\exp(-0.180)=0.8353\), with interval \((\exp(-0.3172),\exp(-0.0428))=(0.7282,0.9581)\). In context, the estimated per-month multiplicative change in expected expression for the comparison condition is about 0.835 times that for the reference condition. This describes a relative slope difference, not a 16.5% absolute decrease in expression each month.

What the Interaction Means at Specific Times

An interaction is often most useful when translated into condition contrasts at interpretable time points. With numeric time, the contrast at time \(t\) is \(\beta_{Cg}+\beta_{CTg}t\). The change in that contrast from time 0 to time \(t\) is \(\beta_{CTg}t\). This is a difference-in-differences on the log-mean scale. For categorical time, the corresponding interaction coefficient is the difference between a time-specific condition contrast and the reference-time condition contrast.

Worked Example: Contrasts at Baseline and Follow-Up

Consider a second hypothetical model with time measured in visits, coded 0, 1, and 2. Suppose the fitted coefficients are \(\widehat{\beta}_{Cg}=-0.300\) and \(\widehat{\beta}_{CTg}=0.200\). Assume the model uses the same offset and repeated-measures structure at all visits, and that a linear time trend is a defensible summary over these three visits.

At baseline, the comparison-versus-reference condition log-mean contrast is \(-0.300+0.200(0)=-0.300\). At visit 2, it is \(-0.300+0.200(2)=0.100\). The contrast therefore changes by \(0.100-(-0.300)=0.400\) log-mean units, also obtained directly as \(2(0.200)=0.400\).

On the mean scale, the condition ratio is \(\exp(-0.300)=0.7408\) at baseline and \(\exp(0.100)=1.105\) at visit 2. The ratio of these two condition ratios is \(\exp(0.400)=1.492\). Thus, in this hypothetical model, the relative expression contrast at visit 2 is about 1.49 times its baseline value. These fitted contrasts describe the model; their uncertainty should be calculated from the covariance of the coefficient estimates before making inferential claims about either time-specific contrast.

Choosing Numeric or Categorical Time

The time coding determines what kinds of trajectory differences the interaction can detect. Numeric time with a linear interaction is parsimonious and answers whether average linear slopes differ over the modeled interval. It may be inappropriate if trajectories are curved or if time intervals are irregular but treated as equal increments. Categorical time allows each observed time to have its own mean pattern and tests whether condition contrasts vary across those times, but it uses more parameters and can be less precise.

A polynomial or other prespecified time basis can represent selected nonlinear patterns, but every added interaction term changes the hypothesis being tested. The scientific question and sampling schedule should guide the coding. Avoid selecting a time scale after inspecting which one produces the smallest p-value. When interpreting an interaction, report the time coding, reference levels, interaction test, and useful time-specific contrasts.

Worked Example: Joint Test With Categorical Time

In a third hypothetical analysis, time has three categories, with visit 0 as the reference. A model includes two condition-by-time coefficients, \(\widehat{\gamma}_1=0.400\) for visit 1 and \(\widehat{\gamma}_2=-0.100\) for visit 2. Suppose their estimated covariance matrix is diagonal, with standard errors 0.200 and 0.250. Assume this is the fitted covariance for the two coefficients, the model converged, and the donor covariance structure was prespecified and appropriate.

The joint null is \(H_0:\gamma_1=\gamma_2=0\). With the stated diagonal covariance, the Wald statistic is:

$$ W=\left(\frac{0.400}{0.200}\right)^2+ \left(\frac{-0.100}{0.250}\right)^2 =4.000+0.160=4.160. $$

There are two restrictions, so the reference distribution is chi-square with \(df=2\). The upper-tail p-value is approximately \(0.125\). At significance level \(\alpha=0.05\), this test does not provide sufficient evidence that the condition contrasts vary across the three visits. It does not show that the contrasts are identical; the estimates and their uncertainty remain relevant, and limited precision could obscure differences.

Derivation: Why the Interaction Is a Difference in Slopes

Ignoring the offset and donor term for the purpose of comparing the two fitted condition means, the model gives the reference-condition log mean as \(\beta_{0g}+\beta_{Tg}t\). For the comparison condition, substitute \(C=1\), giving \(\beta_{0g}+\beta_{Cg}+\beta_{Tg}t+\beta_{CTg}t\). Subtracting the reference-condition expression from the comparison-condition expression yields:

$$ (\beta_{0g}+\beta_{Cg}+\beta_{Tg}t+\beta_{CTg}t) -(\beta_{0g}+\beta_{Tg}t) =\beta_{Cg}+\beta_{CTg}t. $$

The condition contrast changes by \(\beta_{CTg}\) for each one-unit increase in time, since the contrast at \(t+1\) minus the contrast at \(t\) is \(\beta_{CTg}\). This algebra explains both the slope interpretation and why the main condition coefficient is tied to the chosen time origin. Centering numeric time changes the interpretation of the main condition coefficient, but not the fitted slope difference when the model is otherwise equivalent.

Diagnostics, Failure Modes, and Limitations

Before interpreting a test, check that the design supports it. Both conditions should be observed across the relevant times; if a condition-time combination is absent, some contrasts may not be estimable. Inspect the design matrix for rank deficiency and verify that the model converged. Sparse counts, few independent donors, or weak replication can make standard errors and asymptotic test references unreliable. A nominally significant interaction is not a substitute for checking uncertainty and biological plausibility.

Also assess whether the chosen time form is plausible. A linear interaction can miss a crossing or curved trajectory; a categorical interaction can be imprecise when many time indicators are fitted. Check fitted means and residual behavior using suitable diagnostics, and examine whether a small number of donors dominates the estimated interaction. These checks do not prove the model is correct, but can reveal major incompatibilities between the data and its assumptions.

  • Testing the condition main effect instead. That tests a contrast at the reference time, not whether trajectories differ.
  • Testing each categorical interaction separately. Separate tests do not replace the joint test of whether all interaction coefficients are zero.
  • Calling a nonsignificant test evidence of no interaction. A large p-value indicates insufficient evidence against the null under the model; it does not establish equality.
  • Ignoring repeated donors. Treating correlated measurements as independent can misstate uncertainty and invalidate inference.
  • Overinterpreting a linear slope. A linear interaction summarizes the specified scale and interval; it does not establish that the biological trajectory is truly linear.

The main alternatives are a numeric-time interaction for a compact slope comparison and a categorical-time interaction for a less restrictive comparison across observed visits. A likelihood-ratio test and a Wald test are alternative inferential procedures for the same specified nested-model question, though their finite-sample behavior can differ. Neither approach resolves design limitations, confounding, model misspecification, or low replication. For gene-wide inference, interaction p-values also require appropriate multiplicity control, as with other gene-level tests.

1
Define the target contrast.
Specify the conditions, time unit or categories, and whether the scientific question concerns a linear slope difference or any difference across visits.
2
Fit the repeated-measures count model.
Include the offset, condition and time main effects, the interaction terms, and a suitable donor-level covariance structure.
3
Test the appropriate null.
Use a one-coefficient test for a numeric linear interaction or a joint test for all categorical-time interaction coefficients.
4
Interpret and check.
Report the interaction estimate and uncertainty, translate it into contrasts at meaningful times, and assess design support, fit, and limitations.

Check Your Understanding

Answer each question, then compare with the solutions.

  1. With numeric time, what does \(\beta_{CTg}\) represent, and what does \(\exp(\beta_{CTg})\) represent?
  2. If time has four categories and one is the reference, how many condition-by-time coefficients are tested in the joint interaction test?
  3. A fitted numeric interaction is \(0.12\) log-mean units per week. What is the estimated ratio of the conditions’ per-week multiplicative changes?
  4. Why does a test of the condition main effect not generally test whether trajectories differ?

Solutions. (1) It is the difference between the conditions’ log-mean slopes per time unit; its exponential is the ratio of their per-unit multiplicative changes. (2) Three, one for each non-reference category. (3) \(\exp(0.12)=1.127\), approximately; the comparison condition’s per-week multiplicative change is estimated to be 1.127 times the reference condition’s. (4) The main condition coefficient is the condition contrast at the reference time (or time zero); the interaction tests whether that contrast changes over time.

Clinical Biostats