Modeling Time Trends in Repeated Samples

Representing Time in a Repeated-Sample Model

When pseudobulk samples are collected repeatedly from the same donors, a gene’s expected expression can change over time while measurements from one donor remain correlated. A time term describes the pattern of the gene-level mean across visits; a random effect or other specified covariance structure accounts for dependence among observations from the same donor. These are distinct parts of the model: changing the time term does not, by itself, account for repeated-measure correlation.

The main choice here is how to represent time. Treating visits as categories estimates a separate mean at each visit without imposing a particular curve shape. Treating time as a continuous variable uses a numerical time scale and, in the simplest model, assumes a straight-line trend on the log-mean scale. The choice determines what the coefficients mean and what changes the analysis can detect.

Definition: For gene $g$, donor $d$, and visit $j$, let $Y_{gdj}$ be the pseudobulk count, $s_{dj}$ the sample-specific size factor, and $b_{gd}$ a donor random intercept. A negative-binomial mixed model with a log link can be written $$ Y_{gdj}\mid b_{gd}\sim \operatorname{NB}(\mu_{gdj},\alpha_g),\qquad \operatorname{Var}(Y_{gdj}\mid b_{gd})=\mu_{gdj}+\alpha_g\mu_{gdj}^{2}, $$ $$ \log(\mu_{gdj})=\log(s_{dj})+\beta_{0g}+f_g(t_j)+b_{gd},\qquad b_{gd}\sim N(0,\tau_g^2). $$ Here $\log(s_{dj})$ is an offset, $\alpha_g$ is the gene-specific dispersion, $f_g(t_j)$ is the chosen time effect, and the random intercept represents donor-level variation. Conditional on the random effects, observations are assumed to follow the stated count distribution; the chosen random-effects or covariance structure must also reasonably describe within-donor dependence.

This specification is an example, not the only possible covariance model. The covariance structure should be selected for the repeated design and kept comparable when assessing alternative time terms. The fixed time coefficients describe a shared trend, while donor-specific departures are represented by the random part of the model. If donors have substantially different time slopes, a common fixed slope plus random intercepts may not describe the data well.

Two Ways to Encode Time

For categorical time, choose a reference visit $t_0$. If the other visits are $t_1,\ldots,t_K$, define indicator variables for those visits and set $f_g(t_j)=\sum_{k=1}^{K}\beta_{kg}I(t_j=t_k)$. The intercept is the log mean for the reference visit when the random effect is zero, apart from the offset. Each $\beta_{kg}$ is a log-rate contrast against that reference. Exponentiating it gives a ratio of expected means, conditional on the same donor effect and size factor. This encoding does not require equally spaced visits or a monotone trend.

For continuous time, select a meaningful unit and often center time at a reference visit, defining $x_j=t_j-t_0$. A linear time term uses $f_g(t_j)=\beta_{1g}x_j$. Its coefficient is a change in log expected expression per time unit. Thus $\exp(\beta_{1g})$ is the expected-mean ratio for a one-unit increase in time, comparing observations with the same donor effect and size factor. Centering changes the intercept’s interpretation, not the slope.

A continuous linear term is parsimonious, but it constrains changes between visits to follow a straight line on the log scale. Categorizing time is more flexible but uses more coefficients and usually has less precision for individual comparisons. Neither representation is automatically preferable: the scientific question, visit schedule, sample size, and plausible trend shape matter.

Worked Example: Interpreting Categorical-Time Coefficients

Suppose a hypothetical study measures pseudobulk expression for one gene at days 0, 7, and 14. A fitted model uses day 0 as the reference and reports $\hat\beta_{7}=0.30$ with standard error $0.12$, and $\hat\beta_{14}=0.10$. Assume the model’s count, offset, and repeated-measure assumptions are adequate for this illustration.

The day-7 versus day-0 estimated mean ratio is $\exp(0.30)=1.350$. At the same size factor and donor effect, the estimated expected expression at day 7 is therefore 1.35 times the day-0 mean, or about 35.0% higher. Similarly, the estimated day-14 versus day-0 ratio is $\exp(0.10)=1.105$, about 10.5% higher.

To compare day 14 directly with day 7, subtract the two log contrasts: $\hat\beta_{14}-\hat\beta_{7}=0.10-0.30=-0.20$. The estimated ratio is $\exp(-0.20)=0.819$, so the day-14 mean is estimated to be 18.1% lower than the day-7 mean. Suppose the fitted covariance matrix gives this contrast a standard error of $0.12$. A Wald statistic is $z=-0.20/0.12=-1.67$, with a two-sided normal-approximation $p$-value of about $0.096$. This is not strong evidence against equal day-14 and day-7 means at a 0.05 threshold. A confidence interval is also useful: $-0.20\pm1.96(0.12)=(-0.435,0.035)$ on the log scale, or approximately $(0.647,1.036)$ after exponentiating. It includes a ratio of 1.

The numerical pattern does not establish a smooth decline or any particular shape between visits. The categorical model estimates visit-specific means; it does not say what happened between the visits.

What a Continuous-Time Slope Means

The interpretation of a continuous slope follows directly from the log link. For two times $t_a$ and $t_b$, suppose the size factor is the same and the donor effect is held fixed. The model gives $\log(\mu(t_b))-\log(\mu(t_a))=\beta_{1g}(t_b-t_a)$. Exponentiating both sides yields the expected-mean ratio. This result explains why the slope applies multiplicatively to expression, rather than adding a fixed number of counts at every time unit.

\frac{\mu(t_b)}{\mu(t_a)} =\exp\{\beta_{1g}(t_b-t_a)\}, \qquad \text{percent change per unit} =100\{\exp(\beta_{1g})-1\}\%.

The percent change is not exactly $100\beta_{1g}\%$, although that is a close approximation when the slope is small. For a time interval of $m$ units, the fitted ratio is $\exp(m\beta_{1g})$, not $m\exp(\beta_{1g})$. This distinction matters when translating a daily slope into a two-week contrast.

Worked Example: Converting a Daily Slope to a Fitted Mean

Consider a hypothetical gene-level fit with $t_0=0$ days, size factor 1, donor random effect set to zero, fitted mean 20 counts at day 0, and estimated continuous-time slope $\hat\beta_1=0.028$ per day. Suppose its standard error is $0.010$. The model assumptions include a suitable negative-binomial mean-variance relationship, an adequate repeated-measure structure, and a linear trend on the log-mean scale.

The estimated mean ratio per day is $\exp(0.028)=1.028$, corresponding to an estimated increase of $100(1.028-1)=2.84\%$ per day. At day 14, the fitted mean is $20\exp(0.028\times14)=20\exp(0.392)\approx29.60$ counts.

If time were recorded in weeks instead, the same model would use a slope of $0.028\times7=0.196$ per week. At two weeks the fitted mean is $20\exp(0.196\times2)=20\exp(0.392)\approx29.60$ counts. Changing units changes the numerical slope but not the fitted value for a given elapsed time.

For the null hypothesis of no linear time trend, $H_0:\beta_1=0$, the Wald statistic is $z=0.028/0.010=2.80$. Its two-sided normal-approximation $p$-value is about $0.0051$. A 95% Wald interval for the slope is $0.028\pm1.96(0.010)=(0.0084,0.0476)$ per day. Exponentiating the endpoints gives daily ratios of about $1.0084$ to $1.0488$, or an increase of approximately 0.84% to 4.88% per day under the model. The interval procedure is designed to cover the fixed slope in 95% of repeated samples under its assumptions; it is not a 95% probability statement about this particular interval.

Comparing a Linear Trend With Visit-Specific Means

With three equally spaced visits, a linear continuous-time model is a restricted version of a categorical-time model. The categorical model allows two independent contrasts against a reference; the linear model constrains those contrasts to lie on a straight line in log mean across the equally spaced times. One can assess whether that restriction is adequate by comparing the models with a likelihood-ratio test, provided they use the same observations, offset, count distribution, and random-effects or covariance specification.

Worked Example: Testing Whether a Linear Shape Is Adequate

For a hypothetical gene measured at days 0, 7, and 14, fit a linear-time model and a categorical-time model with the same repeated-measure structure. Suppose the categorical model has a log likelihood of $-120.00$, while the linear model has a log likelihood of $-122.51$. The categorical model adds one time parameter, so the models differ by one degree of freedom.

The likelihood-ratio statistic is $2\{-120.00-(-122.51)\}=5.02$. Under the null restriction that the log means follow a linear trend, this is compared approximately with a chi-square distribution with $df=1$. The resulting $p$-value is about $0.025$. At a 0.05 threshold, the categorical model fits significantly better by this test, suggesting that the linear restriction may not capture the visit pattern adequately.

This comparison does not prove which biological process produced the observed pattern. It indicates that the additional categorical flexibility improves fit more than expected from random variation under the linear model, subject to the model assumptions and the validity of the likelihood-ratio approximation.

Assumptions, Diagnostics, and Common Pitfalls

A time coefficient is only as useful as the model supporting it. Check that the time coding matches the actual visit schedule; that the offset represents sequencing depth or the specified library-size adjustment; and that the chosen random-effects or covariance structure captures repeated observations adequately. In the negative-binomial model, the mean and dispersion assumptions should be plausible for the gene. With small numbers of donors, large-sample Wald and likelihood-ratio approximations may be less reliable.

For a continuous term, compare observed or estimated visit-level patterns with the fitted log-linear trend. Systematic curvature, a reversal in direction, or a sharp change isolated to one visit can indicate that the linear form is too restrictive. For categorical time, inspect estimates and uncertainty at each visit and avoid reading a noisy sequence of estimates as a precise smooth curve. Extreme estimates, wide intervals, or strong sensitivity to one donor may signal sparse information or influential observations; these issues warrant investigation rather than automatic interpretation.

The models answer somewhat different questions. Categorical time asks whether and how expected expression differs among specified visits, without assuming a functional pattern. Continuous linear time estimates a compact average rate of change over the chosen interval and has fewer coefficients, but can conceal nonlinearity. Splines or other flexible continuous functions can represent curvature while retaining time ordering, but require additional choices and enough information to estimate the added flexibility. The appropriate comparison is not simply which model gives the smaller $p$-value; it is whether the trend representation matches the scientific question and observed data.

Pitfall: treating a slope as an additive count change. A coefficient of $0.028$ does not mean 0.028 more counts each day. It is a log-mean slope; exponentiating it gives the multiplicative change per day.
Pitfall: assigning a numerical scale to visits without justification. Coding visits as 0, 1, 2 assumes equal spacing for a linear trend. If actual gaps differ, use the elapsed time in meaningful units.
Pitfall: claiming a categorical model shows a continuous trajectory. It estimates means at observed visits only. No between-visit shape is implied.
Pitfall: interpreting a time association causally. A fitted gene-level trend describes modeled expression over time; it does not alone establish that time or an intervention caused the change.

Key takeaway: Categorical time makes few shape assumptions and estimates visit contrasts; continuous linear time summarizes change per time unit but imposes a straight line on the log-mean scale. Both require an adequate count model and a defensible account of within-donor dependence.

Check Your Understanding

Try each question, then compare your reasoning with the solutions below.

A model uses a continuous slope of $0.05$ per week. What is the fitted expected-mean ratio over four weeks, and what is the corresponding percent change?
Why can a categorical-time model be preferable when visits occur at days 0, 3, and 20 and no simple trend shape is justified?
For three equally spaced visits, a categorical model has one more time coefficient than a linear-time model. What does a likelihood-ratio test of these models assess?
A colleague says that a significant continuous-time coefficient proves expression rises by the same number of counts at every visit. Identify and correct the error.

Solutions. (1) The four-week ratio is $\exp(4\times0.05)=\exp(0.20)=1.221$, a fitted increase of $100(1.221-1)=22.1\%$ over four weeks. (2) Categorical time allows a separate mean at each visit and does not impose equal spacing or a straight-line pattern; it is suitable when a particular shape is not justified, though it uses more parameters. (3) The test assesses whether the categorical model’s extra flexibility improves fit beyond the linear restriction; a small $p$-value is evidence that the linear log-mean shape is inadequate, conditional on the models and assumptions. (4) The coefficient is a change in log expected expression per time unit. The model implies a multiplicative change in the mean, not a constant additive number of counts.