Sequential Mediation Across Molecular Layers

Why Order Matters in Multi-Layer Mediation

A genetic variant may be associated with several molecular measurements, but an association among genotype, methylation, expression, and protein abundance does not by itself establish a sequence of effects. A sequential mediation model makes the proposed order explicit: genotype affects methylation, methylation affects expression, and expression affects protein. It also allows shorter routes, such as a direct genotype-to-protein effect, so that the contribution of the complete chain can be distinguished from other routes.

The goal is not simply to multiply regression coefficients. A product describes a causal pathway only when the model’s causal and ordering assumptions are defensible. Without them, the same calculations describe a decomposition of fitted associations, not evidence that one molecular layer caused a change in the next.

A Recursive Model for Ordered Molecular Layers

Let \(G\) denote genotype, coded as allele dosage \(0\), \(1\), or \(2\). Let \(M\), \(E\), and \(P\) denote methylation, expression, and protein abundance. For illustration, consider these linear structural equations, conditional on measured covariates \(C\):

\begin{aligned} M &= a_M G + \boldsymbol{\gamma}_M^{\mathsf T}C+\varepsilon_M,\\ E &= a_E G+d_{EM}M+\boldsymbol{\gamma}_E^{\mathsf T}C+\varepsilon_E,\\ P &= a_P G+d_{MP}M+d_{EP}E+\boldsymbol{\gamma}_P^{\mathsf T}C+\varepsilon_P. \end{aligned}

The arrows represented by this model are \(G\to M\), \(G\to E\), \(G\to P\), \(M\to E\), \(M\to P\), and \(E\to P\). It is recursive because each molecular layer is specified using only genotype, covariates, and earlier layers. Coefficients are slopes in the units of their outcome per unit of their predictor, holding the other terms in that equation fixed. For example, \(d_{EM}\) is the modeled change in expression per unit increase in methylation, conditional on \(G\) and \(C\).

Definition: The sequential indirect effect along \(G\to M\to E\to P\) in this additive linear model is the contribution transmitted through all three arrows, \(a_Md_{EM}d_{EP}\). The coefficients \(a_M\), \(d_{EM}\), and \(d_{EP}\) must be defined for the same model and scale conventions.

This full-chain contribution is not the only indirect route in the model. The route \(G\to E\to P\) contributes \(a_Ed_{EP}\), and the route \(G\to M\to P\) contributes \(a_Md_{MP}\). The direct route \(G\to P\) contributes \(a_P\). The words “direct” and “indirect” refer to routes in the specified model; they do not establish that the corresponding arrows are causal.

A causal interpretation requires, among other conditions, correct temporal ordering; well-defined measurements and interventions; adequate overlap in relevant exposures; and no unmeasured confounding of the modeled relations after conditioning on \(C\). For molecular mediation, this includes confounding of the genotype–mediator and mediator–outcome relations, and no mediator–outcome confounder caused by genotype that is left unaccounted for. Genotype is not automatically free of confounding: ancestry, relatedness, and selection can induce associations that need to be addressed. The linear product interpretation also relies on the stated additive model, including no unmodeled interactions that change path effects across genotype or mediator levels.

Deriving the Path Contributions

The path decomposition follows by substituting the earlier equations into the protein equation. To focus on the effect of changing \(G\), hold \(C\) fixed and consider the part of each equation that changes with \(G\). The change in methylation per unit change in genotype is \(a_M\). Substituting this into the expression equation gives a change in expression of \(a_E+d_{EM}a_M\). Substitution into the protein equation then gives:

\begin{aligned} \frac{\partial P}{\partial G} &=a_P+d_{MP}a_M+d_{EP}(a_E+d_{EM}a_M)\\ &=a_P+a_Md_{MP}+a_Ed_{EP}+a_Md_{EM}d_{EP}. \end{aligned}

Each term corresponds to a distinct directed route: \(a_P\) for \(G\to P\), \(a_Md_{MP}\) for \(G\to M\to P\), \(a_Ed_{EP}\) for \(G\to E\to P\), and \(a_Md_{EM}d_{EP}\) for the complete sequence \(G\to M\to E\to P\). This algebra is a model-based decomposition. It is a causal total effect only if the structural coefficients support causal interpretation and the model’s assumptions hold.

In a strict chain with only \(G\to M\), \(M\to E\), and \(E\to P\), the other routes are excluded, and the modeled total effect is \(a_Md_{EM}d_{EP}\). That restriction should come from scientific justification, not from deleting paths merely because their coefficients are inconvenient or imprecisely estimated.

Worked Example: Calculating the Complete Sequential Effect

Suppose a hypothetical study fits compatible linear models, with \(M\), \(E\), and \(P\) each measured in standard-deviation units. The estimated slopes are \(a_M=-0.20\) SD of methylation per allele, \(d_{EM}=-0.50\) SD of expression per SD of methylation, and \(d_{EP}=0.40\) SD of protein per SD of expression. Calculate the contribution along \(G\to M\to E\to P\).

The product is:

a_Md_{EM}d_{EP} =(-0.20)(-0.50)(0.40) =0.100(0.40) =0.0400.

The modeled sequential contribution is \(0.0400\) SD of protein abundance per additional allele. The sign is positive because the genotype-to-methylation and methylation-to-expression slopes are both negative, while the expression-to-protein slope is positive. This calculation isolates the complete route; it does not include any direct or shorter routes. It is a causal indirect effect only under the identification assumptions stated above.

Decomposing a Model With More Than One Route

A full sequential model often needs shorter or direct routes because molecular regulation need not proceed through every measured layer. The decomposition helps make clear what a reported “mediated effect” includes. In particular, a product involving only \(G\), \(M\), and \(P\) is not the same quantity as the contribution along \(G\to M\to E\to P\).

Worked Example: Separating the Sequential Route From Other Routes

Consider hypothetical coefficients in the same units as the preceding example: \(a_M=-0.20\), \(d_{EM}=-0.50\), \(d_{EP}=0.40\), \(a_E=0.10\), \(d_{MP}=0.15\), and \(a_P=0.02\). Calculate each route’s contribution and the model-implied total effect per allele.

The full sequential route is:

G\to M\to E\to P:\quad a_Md_{EM}d_{EP}=(-0.20)(-0.50)(0.40)=0.0400.

The shorter methylation route contributes \((-0.20)(0.15)=-0.0300\) SD of protein per allele. The expression route contributes \((0.10)(0.40)=0.0400\) SD per allele, and the direct route contributes \(a_P=0.0200\) SD per allele. Adding the four contributions gives:

$$ 0.0400-0.0300+0.0400+0.0200=0.0700. $$

Thus the model-implied total change is \(0.0700\) SD of protein per additional allele. The complete sequential contribution is \(0.0400\) SD per allele, not the total \(0.0700\). These are algebraic results for the specified model; causal labels require the assumptions to be credible.

When Observational Measurements Do Not Identify the Sequence

Identification asks whether the available data and assumptions determine a target quantity uniquely. A fitted regression coefficient can be calculated from observed data even when the causal effect represented by that coefficient is not identified. The distinction matters especially between mediators: an unmeasured common cause of methylation and expression can make their observed association compatible with multiple causal effects of methylation on expression.

For a simple demonstration, suppose genotype has already been conditioned on, and let \(m\) and \(e\) be the residual methylation and expression measurements. Assume the observed data have \(\operatorname{Var}(m)=1\), \(\operatorname{Var}(e)=1\), and \(\operatorname{Cov}(m,e)=0.60\). A regression of \(e\) on \(m\) has slope \(0.60\). But that slope does not, by itself, reveal the causal coefficient from \(M\) to \(E\).

Worked Example: Two Causal Models With the Same Observed Association

In one hypothetical model, write \(e=0.60m+v\), with \(\operatorname{Cov}(m,v)=0\) and \(\operatorname{Var}(v)=0.64\). This gives \(\operatorname{Var}(e)=0.60^2(1)+0.64=1\), and \(\operatorname{Cov}(m,e)=0.60\). The causal coefficient of \(m\) on \(e\) in this structural model is \(0.60\).

In another hypothetical model, write \(e=0.20m+v\), but allow an unmeasured common cause to produce \(\operatorname{Cov}(m,v)=0.40\). Choose \(\operatorname{Var}(v)=0.80\). Then \(\operatorname{Cov}(m,e)=0.20(1)+0.40=0.60\), and \(\operatorname{Var}(e)=0.20^2(1)+2(0.20)(0.40)+0.80=1\). The residual covariance structure is valid because \(0.80-(0.40)^2=0.64>0\). The causal coefficient is now \(0.20\), despite the same observed variances and covariance.

The observed regression slope is \(0.60\) in both cases, but the structural effect differs. If a downstream protein equation has a causal expression-to-protein coefficient of \(0.50\), the corresponding \(M\to E\to P\) contribution is \(0.60(0.50)=0.30\) in the first model and \(0.20(0.50)=0.10\) in the second. The observed association alone does not choose between these values. A sequential effect therefore cannot be identified from these measurements without additional assumptions or information that rules out the alternative confounding structure.

Practical Checks Before Reporting a Path Product

Before interpreting a sequential product, verify that the proposed sequence is supported by biological timing and study design, that all coefficients refer to the same defined population and compatible covariate adjustment, and that their units multiply to the units reported for the effect. In multi-omic data, measurements taken at one time point may not reveal which molecular layer changed first. Feedback, shared regulation, cell composition, technical artifacts, and unmeasured factors can also produce cross-layer associations that resemble a chain.

A product may be statistically estimable while its causal interpretation remains uncertain. Conversely, a biologically plausible pathway is not identified merely because its estimated product is nonzero. Report the route and its assumptions explicitly, and distinguish a model-implied association decomposition from a causal mediation claim.

Common pitfall: “The genotype is associated with methylation, methylation is associated with expression, and expression is associated with protein, so their coefficient product proves the sequential mechanism.” The correction is that coefficient multiplication gives the specified model’s path contribution; causal interpretation additionally requires the proposed order, adequate adjustment for confounding, and the other identification conditions to hold. Observational associations alone may not identify the direction or the causal effect at each mediator transition.

Check Your Understanding

Use the model definitions above, and distinguish algebraic path contributions from causal conclusions.

In the stated model, which product represents the complete route from genotype through methylation and expression to protein?
If \(a_M=0.30\), \(d_{EM}=0.20\), and \(d_{EP}=-0.50\), calculate the sequential contribution and state its units if all three molecular variables are measured in standard deviations.
Why is \(a_Md_{MP}\) a different path contribution from \(a_Md_{EM}d_{EP}\)?
What does it mean for a sequential causal effect to be unidentified, even when a regression slope can be estimated?
Name one reason why the observed ordering of molecular measurements may not establish the causal ordering of the pathway.