Why Test Tissue-Specific Molecular Effects?
A genetic variant can be associated with a molecular trait in several tissues, but similar-looking associations do not establish that the genetic effect is the same in each tissue. A variant may change expression strongly in one tissue and weakly in another; the direction may even differ. Testing this question requires a direct comparison of tissue-specific effects. Finding an association in one tissue but not another is not, by itself, evidence that the effects differ.
Measurements from different tissues may also reflect different cell mixtures. For example, a bulk-tissue expression measurement combines expression across the cell types present in the sample. If genotype is associated with the proportion of a cell type, an apparent genotype-by-tissue difference could partly reflect this composition difference rather than a different effect within comparable cellular mixtures. The model and interpretation therefore need to make clear whether composition is adjusted for.
A Genotype-by-Tissue Interaction Model
Let \(Y_{it}\) be a molecular trait measured for donor \(i\) in tissue \(t\), and let \(G_i\) be the donor's genotype dosage, commonly coded as 0, 1, or 2 copies of an allele. With three tissues, take tissue 1 as the reference and define indicators \(I(t=2)\) and \(I(t=3)\). A linear model is
Here \(\mathbf{Z}_{it}\) contains other prespecified covariates, such as ancestry principal components or technical covariates. The vector \(\mathbf{Q}_{it}\) can contain measured cell-composition variables, with one component omitted when proportions sum to one. The donor term \(u_i\) accounts for a shared donor-level contribution when multiple tissues are measured on the same person. The residual \(\varepsilon_{it}\) represents remaining tissue-specific variation.
Under this coding, the genotype slope in tissue 1 is \(\beta_G\); in tissue 2 it is \(\beta_G+\delta_2\); and in tissue 3 it is \(\beta_G+\delta_3\). Each slope is the adjusted difference in expected trait value per additional allele copy, holding the included covariates fixed. The interaction coefficients are differences in slopes from tissue 1. A different reference tissue changes the individual coefficient values but not the fitted tissue-specific slopes or the hypothesis that all slopes are equal.
If donors contribute multiple tissue measurements, those observations are correlated and must not be treated as independent. A random-intercept model is one option when its covariance assumptions are reasonable. Generalized least squares or a model with an appropriate within-donor covariance structure can also be used. If each donor contributes only one tissue, a donor random intercept is not needed. The linear model assumes the conditional mean is adequately represented by the specified terms; standard model-based tests additionally rely on a correctly specified covariance model. Large-sample Wald inference may be used with a consistent covariance estimate, while small-sample settings require suitable finite-sample methods.
Worked Example: Interpreting an Adjusted Tissue Interaction
Suppose a hypothetical study measures a molecular trait in blood and liver from the same donors. The trait is standardized, so its unit is one standard deviation. In a donor-correlated model adjusted for ancestry covariates and measured cell proportions, blood is the reference tissue. The fitted genotype coefficient is \(\hat{\beta}_G=0.12\), and the genotype-by-liver interaction estimate is \(\hat{\delta}_{\text{liver}}=0.30\), with \(SE(\hat{\delta}_{\text{liver}})=0.11\).
The estimated blood slope is \(0.12\) standard deviations per additional allele copy. The estimated liver slope is
standard deviations per allele copy. The interaction test evaluates \(H_0:\delta_{\text{liver}}=0\). Its Wald statistic is \(z=0.30/0.11=2.73\), giving an approximate two-sided \(p\)-value of 0.0064 under a large-sample normal approximation. In this hypothetical analysis, the adjusted genetic slopes differ between liver and blood. This is evidence of a tissue difference in the conditional association, not proof that the variant causes the trait through a particular biological mechanism.
The interpretation depends on the model assumptions: the donor-level correlation has been accounted for, the trait is adequately modeled as linear in dosage and covariates, and the reported standard error is valid for that model. Because cell proportions were included, these slopes describe associations conditional on measured composition. They need not equal total effects if genotype also changes composition.
Testing Heterogeneity Across Several Tissues
With \(K\) measured tissues, let \(\boldsymbol{\beta}\) denote the vector of estimated tissue-specific genotype slopes and let \(\widehat{\operatorname{Var}}(\hat{\boldsymbol{\beta}})\) be their estimated covariance matrix. Testing whether the slopes are all equal is a joint test of \(K-1\) independent contrasts. Write the contrasts as \(C\boldsymbol{\beta}\), where each row of \(C\) subtracts one tissue slope from another. The Wald statistic is
If the model is correctly specified and regularity conditions hold, \(W\) is approximately chi-squared under the null hypothesis of equal slopes, with degrees of freedom equal to the rank of \(C\), usually \(K-1\). Using the covariance among tissue estimates matters: estimates from the same donors are generally not independent. A significant omnibus test indicates some heterogeneity, not which tissues differ; prespecified pairwise contrasts or follow-up estimates are needed to describe its location. If many traits or tissue contrasts are tested, multiplicity also needs to be considered.
Worked Example: A Joint Test of Three Tissue Slopes
Suppose a hypothetical fitted model gives three tissue-specific slopes, and two independent contrasts relative to tissue 3 are estimated as \(0.30\) and \(0.15\) trait units per allele copy. Their estimated covariance matrix is
The off-diagonal term is the covariance between the two contrasts. For example, if the contrasts are \(\hat{\beta}_1-\hat{\beta}_3\) and \(\hat{\beta}_2-\hat{\beta}_3\), it includes the covariance induced by their shared tissue-3 estimate. The determinant is
Thus the inverse covariance matrix is approximately
Substituting the contrast estimates gives
The arithmetic is \(W=0.000855/0.000193\approx4.430\), using the unrounded inverse covariance matrix. There are two independent contrasts, so the reference distribution is approximately \(\chi^2\) with \(df=2\). The upper-tail \(p\)-value is approximately \(0.109\). At a 0.05 significance level, these data do not provide sufficient evidence to reject equal slopes across the three tissues. This result is not evidence that the slopes are exactly equal; the estimates and their uncertainty should still be reported.
Cell Composition and the Meaning of Adjustment
A bulk-tissue molecular measurement can change because expression within cells changes, because the mixture of cell types changes, or both. Including measured composition variables can help distinguish these explanations, but it changes the estimand. The genotype slope conditional on composition asks how the expected molecular trait differs by genotype among samples with the same measured composition. If genotype influences composition, adjustment removes part of a pathway that may contribute to the overall genotype–trait association. Therefore, adjusted and unadjusted interaction estimates answer different questions.
Composition adjustment is not automatically a cure for confounding. It requires reasonably measured composition variables, an appropriate functional form, and no important unmodeled confounding of the conditional association. If a proportion is measured with substantial error, residual composition differences may remain. If composition is downstream of genotype, conditioning on it can also introduce bias under some causal structures. Researchers should describe which variables were adjusted for and avoid calling a conditional association a within-cell causal effect unless the design and assumptions justify that interpretation.
When the scientific target is an effect within a particular cell type, a bulk-tissue interaction model may not answer the question. It estimates a tissue-level conditional association, not necessarily a cell-intrinsic effect. The experimental unit and dependence structure also matter: repeated measurements from the same donor provide correlated information, not independent replicates.
Tissue-Specific Mediation Paths
A tissue-specific pathway may be represented by a mediator \(M_t\), such as an expression trait, and an outcome \(Y_t\). In a simple linear mediation specification, \(a_t\) is the association of genotype with the mediator in tissue \(t\), and \(b_t\) is the association of the mediator with the outcome conditional on genotype and the modeled covariates. Under the linear path model, the indirect effect in tissue \(t\) is \(D_t=a_tb_t\). Comparing \(D_t\) across tissues asks whether the estimated mediated path differs; it is not the same test as comparing genotype slopes for the molecular outcome.
The product \(a_tb_t\) is nonlinear in its estimated coefficients. A first-order delta-method variance uses the gradient of the product and the covariance matrix of the coefficient estimates. For comparing tissues, the gradient must include both paths in both tissues, and the full covariance should be used when estimates are correlated. The simplified sum of squared terms below assumes the four coefficient estimates are mutually independent.
Worked Example: Comparing Two Tissue-Specific Indirect Effects
Consider a hypothetical analysis of two tissues. In tissue A, the estimated genotype-to-mediator path is \(\hat a_A=0.50\) with standard error \(0.10\), and the mediator-to-outcome path is \(\hat b_A=0.40\) with standard error \(0.08\). In tissue B, \(\hat a_B=0.10\) with standard error \(0.10\), and \(\hat b_B=0.60\) with standard error \(0.08\). Units are such that the products are measured in outcome units per allele copy.
The estimated indirect effects are \(\hat D_A=(0.50)(0.40)=0.20\) and \(\hat D_B=(0.10)(0.60)=0.06\). Their difference is \(\hat D=\hat D_A-\hat D_B=0.14\). For the coefficient vector \((a_A,b_A,a_B,b_B)\), the gradient of \(D=a_Ab_A-a_Bb_B\) is \((b_A,a_A,-b_B,-a_B)\), which at the estimates is \((0.40,0.50,-0.60,-0.10)\).
Assuming the four coefficient estimates are mutually independent, the delta-method variance pairs each gradient component with the standard error of the coefficient with respect to which it is a derivative:
The standard error is \(\sqrt{0.006864}=0.08285\). A normal-approximation statistic is \(z=0.14/0.08285=1.69\), with a two-sided \(p\)-value of about 0.091. At \(\alpha=0.05\), this analysis does not provide sufficient evidence that the indirect effects differ between these tissues. The point estimates differ, but their estimated difference is uncertain. If path estimates share donors or otherwise covary, the independence assumption is inappropriate; the full coefficient covariance matrix must be used in the delta-method calculation.
Derivation: Why the Delta-Method Gradient Has These Terms
For a smooth function \(g(\boldsymbol{\theta})\) of an estimated coefficient vector \(\hat{\boldsymbol{\theta}}\), a first-order Taylor expansion around the true coefficient vector gives
Taking the variance of the linear term yields the delta-method approximation
For \(g=a_Ab_A-a_Bb_B\), differentiating with respect to \(a_A,b_A,a_B,b_B\) gives \(b_A,a_A,-b_B,-a_B\), respectively. Under independence, the covariance matrix is diagonal, so the quadratic form reduces to the sum of each squared gradient component times the variance of its corresponding coefficient. This explains why the pairing must follow coefficient order; using a coefficient's standard error with a different coefficient's derivative gives the wrong variance.
Common Mistakes and Practical Checks
- Comparing significance labels instead of testing a difference. One tissue can have \(p<0.05\) and another \(p>0.05\) even when their estimated slopes are not significantly different. Test the interaction or an explicit slope contrast.
- Ignoring covariance between tissues. When donors contribute several tissues, contrasts share information. Use the covariance from the fitted joint model rather than treating tissue estimates as independent.
- Calling a composition-adjusted coefficient a total effect. Adjustment conditions on measured composition. If genotype changes composition, that conditional estimate excludes a possible composition-mediated component.
- Mis-pairing delta-method terms. Each derivative is paired with the variance of its own coefficient. Write the parameter order and gradient explicitly before calculating.
- Overstating a mediation result. A product of fitted paths does not establish a causal mechanism without the required identification assumptions.
A useful report gives the tissue-specific slope estimates, their units and uncertainty, the exact heterogeneity contrasts, the covariance-aware test and degrees of freedom, and the covariates used, including any composition adjustment. This makes clear both what differed and what population-level or conditional question the model addressed.
Check Your Understanding
Answer each question using the model and interpretations developed above.
- With tissue 1 as reference, how do you express the genotype slope in tissue 3 using the interaction-model coefficients?
- Why is a significant genotype association in one tissue and a nonsignificant association in another not, on its own, evidence of heterogeneity?
- What does a genotype slope adjusted for cell proportions estimate if genotype may affect those proportions?
- For three tissues, how many independent contrasts are needed to test the joint null hypothesis that all slopes are equal?
- In a mediation comparison \(a_Ab_A-a_Bb_B\), what is the derivative with respect to \(b_B\), and which coefficient's standard error must be paired with it?