On the problem of inflation in transcriptome-wide association studies

Hundreds of thousands of loci have been associated with complex traits via genome-wide association studies (GWAS), but an understanding of the mechanistic connection between GWAS loci and disease remains elusive. Genetic predictors of molecular traits are useful for identifying the mediating roles of molecular traits and prioritizing actionable targets for intervention, as demonstrated in transcriptome-wide association studies (TWAS) and related studies. Given the widespread polygenicity of complex traits, it is imperative to understand the effect of polygenicity on the validity of these mediator-trait association tests. We found that for highly polygenic target traits, the standard test based on linear regression is inflated Eχtwas2>1. This inflation has implications for all TWAS and related methods where the complex trait can be highly polygenic—even if the mediating trait is sparse. We derive an asymptotic expression of the inflation, estimate the inflation for gene expression, metabolites, and brain image derived features, and propose a solution to correct the inflation.


Introduction
To explain the mechanism behind the hundreds of thousands of loci discovered via genome-wide association studies (GWAS), researchers have studied the role of molecular traits as mediators, leading to development of transcriptome-wide association studies (TWAS) and related methods for other molecular traits (Gamazon 1 et al., 2015;Gusev et al., 2016;Zhang et al., 2022).These mediator-trait association studies have been increasingly important in genetic studies.However, a recent paper (Leeuw et al., 2023) has reported statistical inflation in TWAS.This inflation could be falsely identifying loci as putatively causal; therefore, it is important that we identify the source of the inflation error and correct for it.Leeuw et al. (2023) attributed the observed inflation to the prediction error in the mediator.However, established error-in-variable literature (Fuller, 1987) indicates error in the independent variable does not cause inflation.The Leeuw et al. (2023) derivation conflates the true parameter (an unknown but fixed number) and the estimated parameter (a random variable).The proposed null hypothesis sets the estimated regression coefficient to be zero, which is an event of probability zero.This is naturally not a reasonable null hypothesis.
We show in this study that the inflation in TWAS and related methods is actually due to the polygenicity of the target trait.Most complex traits are highly polygenic (Yang et al., 2010;O'Connor et al., 2019;Boyle et al., 2017), and while the effect of polygenicity has been explored and leveraged in the context of GWAS with methods such as linkage disequilibrium (LD) score regression and related approaches (Bulik-Sullivan et al., 2015), the effect of polygenicity on TWAS has not been explored rigorously.
Here, we demonstrate analytically that when the target trait has a polygenic component, the distribution of the χ 2 twas statistics of the associations between the target trait and the gene expression (or other mediators) under the null has a positive noncentrality parameter, rather than the standard χ 2 with a mean = 1, causing the observed inflation.Further, we characterize the properties of this inflation, highlight that prediction error does not cause inflation, show the effect of prediction error in the power of the test, demonstrate the practical relevance of inflation in real TWAS, and, finally, propose a solution to correct it.

Results
TWAS and related methods nominate potential causal mediators (gene expression, protein levels, etc.) by testing the effect of the mediating trait T on a target trait Y assuming the following model: where β is the (fixed) effect of the mediator on the target trait, ϵ twas is the component unexplained by the mediator and is independent of it.γ k are the genetic effects on the mediator and X k are the genotype dosages.
γ k and ϵ twas are normally distributed and independent of each other.To allow for sparse architecture of the mediator, γ k can have a positive probability of being 0. The normality assumption can be replaced by finite variance assumption for all variables.However, traditional TWAS ignores the effect polygenicity of the target trait and prediction error may have on the association statistics.
To examine the effect of the polygenicity of the target trait Y on the association statistics, we explicitly modeled the direct genetic effects δ k on the target trait using equation (1).To account for prediction error, we took into account that the TWAS regression is performed against a noisy version of the mediator T .
where δ k , γ k , e γ,k , and ϵ are all normally distributed and independent of each other.As in traditional TWAS, γ k can have a positive probability of being 0 to include sparse architecture for the mediator.For simplicity and without loss of generality, we assumed that Y , T , and X k have mean = 0 and variance = 1.
To quantify the effect of polygenicity and prediction error, we examined the asymptotic distribution of the χ 2 statistic of the regression of the phenotype Y on the noisy predicted expression T .We call this statistic χ 2 twas to emphasize that this is the statistic that the standard TWAS method would calculate when the polygenicity and prediction error are present but ignored.As shown in the Methods section, the χ 2 twas statistic has a noncentral χ 2 distribution with the mean given by where N is the sample size, β is the effect of the mediator T on the target trait Y , and h 2 δ is the polygenic portion of Y , i.e., the heritability of the target trait explained by the genetic effects δ k .Φ is defined as: where M is the number of causal SNPs for the target trait, Σ is the limit of R for N → ∞, and γ = γ + e γ is the M-dimensional vector of noisy prediction weights.τ 2 is the precision of the prediction of the mediator, i.e., the signal to noise ratio of T : The precision τ 2 is also known as the reliability ratio in the error-in-variables literature (Fuller, 1987).

Polygenicity of target trait Y causes inflation under the null
Under the null, β 2 = 0 so that the expected χ 2 in equation ( 4) reduces to:   One immediate implication of this result is that the standard TWAS test-which assumes that Eχ 2 = 1 when in fact Eχ 2 = 1+N h 2 δ Φ-will yield inflated p-values.The degree of inflation is illustrated in Figure 1 using an example.As discussed below, we will correct this inflation by estimating Φ and using the correct null distribution (i.e., the noncentral χ 2 with noncentrality parameter N h 2 δ Φ) to calculate the p-values of the association.
Characterization of the inflation and definition of the mediator specific factor Φ The inflation in TWAS is determined by the noncentrality parameter, which under the null is N h 2 δ Φ.This term is positive when the target trait is polygenic, i.e. it has nonzero h 2 δ .It is linear in h 2 δ and N .The dependence on the number and LD of the causal SNPs is encapsulated in the factor Φ.
For better understanding of how Φ behaves, we investigated its properties.From the definition in equation ( 5), we can see that Φ is only a function of the number of causal SNPs of the target trait M , the prediction weights γ, and the LD matrix Σ. which is the large N limit of X ′ • X/N and thus is no longer dependent on N .Therefore, it does not depend on the sample size, heritability, or other properties of the target trait.
This indicates that Φ is a property of the mediator, which can be pre-calculated and applied to any target trait Y .
We also show that Φ is strictly positive, bounded below and above as follows: The lower bound is attained when Σ is the identity matrix, and the upper bound is achieved when the SNPs are perfectly correlated (Methods).When the SNPs are independent, we have: and hence the Eχ 2 twas ≈ 1 + N M h 2 δ under the null.
The factor Φ is independent of the precision of the predictor τ 2 : τ 2 is a scalar that affects both the numerator and denominator of Φ through the noisy weights γ and therefore has no net effect.This result is more evident when we assume polygenicity of the mediator trait, in which case, Φ is given by: Prediction error has no effect on the inflation under the null Next, we examined the effect of the precision of the prediction of the mediator on the inflation.Under the null, both T and estimated T are independent of the target trait and have exactly the same dependence on genotype data; thus, the expected χ 2 should not change.Indeed, we corroborated this by verifying that, under the null, the simulated expected χ 2 is constant across all values of the precision as shown in Figure 2a.
Under the alternative, we find that Eχ 2 twas increases monotonically with the precision of the prediction as shown in Figure 2b, indicating that prediction error reduces the power of the test.Details of the simulation are described in the following section.
Consistent with established literature (Fuller, 1987), our results show that prediction error in the mediator causes a loss of power, but it does not affect the inflation under the null-contrary to the conclusion in Leeuw et al. (2023).

Validation of inflation equation under the null
To assess how well this approximation works under the null, we simulated both the target trait Y and the mediating trait T from infinitesimal models with independent effect sizes according to the equations (1) to (3).We simulated genotype data using independent binomial random variables with probability of 0.4corresponding to the minor allele frequency of the SNPs-and therefore assumed no LD between SNPs.We used a range of values for the heritability of Y (h 2 δ : 0.1 -0.9), sample sizes (N : 100 -10, 000), and number of causal SNPs (M : 99 -6000).
For each combination of h 2 δ , N , and M , we simulated 1000 target traits Y and 99 mediating traits T s and T s unrelated to Y .Each of the 99 predicted mediating traits were simulated with different levels of precision (τ 2 : 0.10 -0.9).We then regressed each target trait on each mediating trait separately and averaged the square of the Z-scores across the 1,000 simulations, thereby obtaining an estimated Eχ 2 twas for each combination of h 2 δ , N , M , T k , and τ 2 .
This estimated value was well approximated by our theoretical expression under the independent SNP assumption as shown in Figure 3a, where data points fall in the vicinity of the identity line.Panels b-d of the figure corroborate the linear relationship between the expected χ 2 and 1/M , h 2 δ , and N , respectively, as predicted by our equation ( 7).

Validation of inflation equation under the alternative
To examine the effect of the prediction error on the test under the alternative, we used the same simulation setup we used for the null hypothesis with β 2 > 0. Hence, we simulated the target trait Y = β T + k δ k X k +ϵ and the mediating trait T = k γ k X k .We performed the association using the noisy version of the mediator, Similarly to the null case, the estimated Eχ 2 was well approximated by our theoretical epxression under the independent SNP assumption as shown in Figure 4a, where the data points fall in the vicinity of the identity  against 1/M, h 2 δ , and N. Dash-dotted lines in the figure show the predicted χ 2 statistics from equation ( 4) assuming independent SNPs-for most conditions, the match is good enough to obscure the dash-dotted line in the figure.We used M − 1 for M, which further improved the match.

Observed inflation in TWAS with actual genotype data
To assess the practical relevance of this inflation in a traditional TWAS, we computed the association statistic between actual predicted expression levels and null polygenic target traits.We used genotype data from unrelated White British individuals in the UK Biobank with sample sizes of 1,000, 5,000 and 10,000.We predicted expression levels of 7,131 genes in whole blood using the GTEx v8 prediction models (Barbeira et al., 2021).To generate null polygenic target traits, we simulated Y with heritability ranging from 0.01 to 0.99 using the same UK Biobank genotype data.We sampled the effect sizes δ k and independent error term ϵ from independent standard normal distributions.For each combination of gene, heritability value, and sample size, we generated 1,000 independent simulated traits.We regressed out the first five genetic principal components from the simulated trait to avoid capturing associations due to population structure, which we found to be sufficient to account for population structure in our simulations.Finally, we regressed the residuals of the simulated traits against predicted expression levels and estimated the expected χ 2 twas statistics, averaging the results across the 1000 simulations.
Figure 5 shows the resulting average χ 2 twas , which show a linear dependence on the heritability of the target trait and the sample size of the association consistent with equation (4).Panels (a) and (c) show the average χ 2 twas for all genes, whereas panels (b) and (d) show the average χ 2 twas for the top 10 genes with the highest inflation.Similar results were obtained for other mediators such as metabolites and brain features (Figure S7 and S8), showing the robustness of equation ( 4) to mediators with differing genetic architecture.
Next, we used the average χ 2 to estimate of the factor Φ gene for 7,131 genes expression in whole blood.
We regressed the average χ 2 twas on the product of sample size and heritability (Y ∼ N h 2 δ ).To borrow information across genes, we used a mixed effects model with a common intercept for all genes and random slopes for each gene.For a realistic scenario with biobank-level sample size of 500,000 individuals and heritability of 0.50, the noncentrality parameter (N h 2 δ Φ) would be 5, illustrated as an example in Figure 1.We also applied the same procedure to 1,156 plasma metabolites and 308 diffusion MRI brain features.

Proposed solution to the inflation problem
To account for the inflation in the PrediXcan software, the p-values will be calculated using the noncentral χ 2 distribution with noncentrality parameter N h 2 δ Φ.For example, in R the command would be pchisq(chi2, ncp=N*h2Y*phi, df=1, lower.tail=FALSE),where chi2 is the χ 2 twas statistic = Z 2 , ncp is the noncentrality parameter, df is the degrees of freedom, h2Y is the heritability of the target trait, and phi is the estimated noncentrality parameter for the gene (or the mediator more generally).
To help users of the PrediXcan framework implement the correction, we will make the corrected values the default output from the software, and we will share the estimated Φ in the same database as the prediction weights.When the user performs the the actual TWAS analysis with GWAS data, the software will estimate the heritability of the GWAS trait.The sample size will have to be provided.Then, the noncentrality parameter will be calculated using 1 + h 2 N Φ gene , where h 2 is the heritability of the GWAS trait and N is the sample size of the GWAS.

Discussion
We quantified the effect of polygenicity and prediction error of mediating traits in TWAS and related methods.To accomplish this, we derived a general closed-form asymptotic expression for the inflation, the validity of which was confirmed by simulations and TWAS associations based on real genetically predicted mediators of varying genetic architecture.
We found that the inflation is driven by the polygenicity of the target trait and not the polygenicity of the mediator trait.This is important to note because cis gene expression levels are mostly sparse ( Wheeler et al. (2016)).We also demonstrated that the inflation is not affected by the precision (or, equivalently, the uncertainty) of the prediction of the mediator trait under the null.
The closed-form solution for independent SNPs provided us with useful insights of the cause of the inflation: as the number of independent causal SNPs increases, the factor Φ trends to zero.However, the sample size increase offsets that trend.Polygenicity in the target trait is necessary for inflation to occur, since inflation is zero when the polygenic component of Y is zero.
We also show how to compute the correct p-values by estimating the parameters of the noncentral χ 2 distribution.We will share the estimated parameters in the same database as the prediction weights and provide updated PrediXcan software to calculate the correct p-values.
The fact that the uncertainty in prediction does not increase type I error means that when we leave out mediators based on their prediction performance-a common practice in TWAS-we can be less stringent with the threshold.Adding more noisy predictors reduces the power due to multiple testing but will not contribute to inflation.The decision on which parameters to prioritize will depend on the goal of the study.
For example, mediators with suspected functional roles but low heritability, and therefore predicted with larger uncertainty, need not be excluded from the TWAS analysis for fear of increasing the false positive rate.This could increase the discovery yield of TWAS and related methods without increasing the false positives.In summary, this study demonstrates the following key conclusions: 1) Polygenicity of the target trait induces inflation in the test statistics regardless of the genetic architecture of the mediating trait.2) Uncertainty in the prediction of the mediator does not cause inflation.3) Uncertainty in the prediction of the mediator reduces the power of the test.4) LD is not necessary for the inflation to occur (our simulations were done using independent SNPs).5) The inflation can be corrected by using the noncentral χ 2 distribution with noncentrality parameter N h 2 δ Φ, where the factor Φ can be pre-calculated independent of the GWAS.
We have attempted to use realistic assumptions and real genotype data and prediction models for mediators to illustrate the validity of our derivations.However, our study has several limitations: We made several assumptions to streamline the theoretical derivation.These are commonly used assumptions, but they may not hold in practice (e.g., the independence of the effect sizes and the noise terms).However, the results based on actual genotype data and prediction models for gene expression, metabolites, and brain features suggest that prediction results are likely to be robust to these assumptions.We assume an additive infinitesimal model for the target trait to simplify the theoretical derivation, but traits may deviate from this model in practice.Improvements on the estimation of the inflation factor Φ could be made for different genetic architectures.Even so, we expect that the infinitesimal model will be conservative.Our model does not include gene-gene and gene-environment interactions.Their effect on inflation is not clear and will need to be addressed in future work.Finally, our derivations used a linear regression framework, but many GWAS are performed using logistic regression.Linear regression results are a good approximation for logistic regression when the case control ratio is balanced.Hence, we expect our results to be broadly applicable for balanced designs but will need adjustments for unbalanced designs.

Methods
Derivation of the distribution of the χ 2 statistic under the null Here we derive the formula (4) restated below for completeness.
We also restate the definitions of the variables used in the formula.To simplify notation and derivation, we use the vector form of the model.
where δ, γ, e γ , and ϵ all normally distributed and independent of each.γ is defined as γ + e γ .γ, δ are , and e γ is M -dimensional vector.Let M be the number of causal SNPs for target trait Y , N be the sample size, and h 2 β = β 2 be the heritability of the target trait explained by the mediator T .
We define the sample (M × M ) LD matrix as and its large N limit Σ.
Without loss of generality (WLOG), we assume Y , T , and X k have mean = 0 and variance = 1.For simplicity we ignore the difference between N and N − 1, since we are interested in the case where the sample size N is large.
We list here additional assumptions and facts used in the derivation.
• δ, ϵ, γ, e γ are all normally distributed, although we can relax this assumption and only need to require that they are well-behaved so that central limit theorems apply.γ k 's can have a positive probability of being 0 to allow sparse architecture for the mediator.
• We assume that we know X, γ, γ, which means that we are conditioning on these variables unless stated otherwise.This also implies that we are conditioning on T and T .
• We considered the case where the sample size is larger than the number of independent causal SNPs for Y , i.e., N > M .With GWAS studies increasing in sample size, this is the most likely scenario.
We defined τ 2 as the ratio of the variances of T and T , but we related it to the sample variances The estimated β in a TWAS is (Weisberg, 2005) The expected value of β is since the expected value of the last two terms in (18) equal to zero, i.e.E T ′ Next, we show that β tends to τ 2 β for large N .
We want to calculate the χ 2 twas statistic.For that we need the TWAS-estimated variance of β.From standard regression results (Weisberg, 2005), the variance of β is estimated as where the estimate of the variance of the error term v ar(ϵ twas ) is calculated as the residual some of squares (RSS) divided by N − 1 (Weisberg, 2005).
The variance of β used in TWAS, which is unaware of the polygenic term This variance, v ar twas ( β), is not the actual variance of β when the target trait is polygenic, which is the source of the inflation.
Hence, χ 2 twas is given by The χ 2 statistics in TWAS is given by Now, we take expectation of the χ 2 stat to get the mean of the distribution.
When we expand the square of the terms between the brackets, all the cross terms have expectation 0 as shown below.
So the expectation is given by the expected value of the squared terms Assuming a polygenic mediator, Multiplying both sides by τ 2 /N where Φ polygenic is defined as Derivation of Φ and Eχ 2 when SNPs are independent We calculate tr(R 2 ) using the asymptotic distribution of the eigenvalues of R 2 .When SNPs are independent and M and N are large but comparable, i.e.M/N converges to a number that is finite and greater than 0, the eigenvalues of R follow a Marchenko-Pastur distribution (Bai and Silverstein, 2010).The trace of R 2 can be approximated by the second moment of the Marchenko-Pastur distribution multiplied by M (see 28).
We use Lemma 3.1 in (Bai and Silverstein, 2010), which states that the kth moment is In the case of k = 2, eq.( 26) becomes Therefore, (4) becomes to the following when SNPs are independent.
Since these derivations are based on asymptotic approximations with N ≥ 1 and M ≫ 1, they do not distinguish between M or M − 1.

Derivation of useful results
Calculating By "fully polygenic T ", we mean that γ i ∼ iid N (0, σ 2 γ ) for each variant i across the genome.Here we want to show that γ R is symmetric so it can be eigenvalue-decomposed where C is the matrix of eigenvectors of R as columns and is orthonormal C ′ C = I.Λ is a diagonal matrix with the eigenvalues of R, λ k .
Since the eigenvalues of 2 are the squares of the eigenvalues of R, we have Similarly, we get Proof that σ 2 γ = 1/M From eq (30) and ( 31), we have

Generalization of inflation formula to sparse mediators
Next, we show that the formula for the TWAS χ 2 statistics in equation ( 24) holds for both polygenic and sparse mediators.However, when the mediator is a sparse mediator, the approximate identities -such as γ • R 2 • γ ≈ tr(R 2 )var(γ) and var(γ) ≈ 1 M γ 2 -will not work since they rely on the fact that there are a large number of variants to sum over and M representing the number of variants in the mediator is equal to the number of variants in the polygenic Y .
Instead, to include a sparse mediator into the scope, we can keep γ as is in the derivation and such result will apply to both polygenic and sparse mediators.Recall that the equation 24 is as follow and γ can take any value, i.e. for a sparse mediator, most entries in γ equals to zero.
Similar to the derivation for the polygenic case, The last approximation, (1 − 1−τ 2 N τ 2 ) ≈ 1, is obtained as suggested by the derivation of the polygenic case.
So, for any mediators, our main result still holds sparse and polygenic mediator whereas to accommodate sparse T , we define Φ R as to infinity.So we can approximate the ratio with the ratio of the expectations, using the corollary of the Slutsky's theorem (Ferguson, 1996).We use From equations ( 38) and ( 39), we have So, for any predictor γ, We note that the result on the scenario when variants are independent (equation 27) can be seen as a special case of equation 41.
To find an extreme of Φ we will find an extreme of γ′ • Σ 2 • γ with the condition that 1/τ 2 = γ′ • Σ • γ, using the Lagrange multiplier approach, i.e. we will find the extreme of the function where L is the Lagrange multiplier.We will also write the numerator in a more convenient form by using the eigenvalue decomposition of Σ.Since Σ is symmetric, it can be eigenvalue-decomposed where C is the matrix of eigenvectors of Σ as columns and is orthonormal Since the eigenvalues of Σ 2 are the squares of the eigenvalues of Σ, we have γ′ where L is the Lagrange multiplier.The solution is obtained by setting the derivatives with respect to λ k and L to 0.
We can show that this is indeed a minimum of the function f by calculating the Hessian of the function and verifying that all the eigenvalues are positive definite.f ′′ = 2γ 2 C,k > 0 for all k and cross derivatives are 0, so that the Hessian is diagonal so that its eigenvalues are 2γ 2 C,k > 0. Therefore, this solution corresponds to a minimum.
The minimum of the function corresponds to the case where all the eigenvalues are equal, therefore Σ is the identity matrix.
If we plug in the identity matrix in the equation for Φ, we get Using the eigenvalue decomposition of Φ, we can reexpress as: This upper bound is attained when the first eigenvalue equals M and all others are 0, i.e., all SNPs are perfectly correlated as seen below Since variance of T has to be non zero, γ 2 C,1 is > 0.

N
(Note that we are using E R to indicate that R is being integrated out, not conditioned on as in the rest of the derivation.Here we are not taking expectation with respect to R but using the fact that R and its quadratic forms are converging to their expected values.) From Theorem 3.1 in (Haff, 1979), for a random variable S ∼ W M (Σ, N ), 43) and (44) (45) Estimation of Φ for genes, metabolites, and brain features.
We downloaded prediction weights for gene expression in whole blood, metabolites in plasma, and brain features from the PredictDB website (Barbeira et al., 2021;Liang et al., 2021).We used the weights for the 7,131 genes, 1,156 metabolites, and 308 brain features to calculate predicted values using genotype data from unrelated White British UK Biobank individuals.We used HapMap3 subset of SNPs.We removed SNPs and individuals with missingness greater than 1%.Batches of 1,000, 5,000, 10,000 individuals were selected for the analysis.After filtering, the sample sizes were 999, 4,992, and 9,979.We computed the genetically predicted gene expression and metabolite levels using the PrediXcan software (Barbeira et al., 2018).For brain features, we used the BrainXcan tool (Liang et al., 2022).We simulated the target trait Y using an infinitesimal model Y = X • δ + ϵ, where X was the genotype matrix from the UK Biobank and the M-dimensional vector δ and N-dimensional vector ϵ were sampled from independent standard normal distributions.All values were scaled to achieve the desired heritability.We calculated Eχ 2 for each gene, metabolite, and brain feature for a range of values of h 2 δ and N .Heritability ranged from 0.01 to 0.99.We calculated Φ for each gene, metabolite, and brain feature as the slope of the linear regression of Eχ 2 on the product of the sample size and heritability N h 2 δ .To borrow information across genes, we ran a mixed effects model where we allowed for random slope but constant intercept.We used the slope of the linear regression as an estimate of Φ for each gene, metabolite, and brain feature.

Figure 1 :
Figure 1: Example of standard and noncentral χ 2 distributions.The probability densities of a standard and noncentral χ 2 with noncentrality parameter (NCP) = 5 are shown.We would get a noncentrality parameter of 5 under realistic parameters such as Φ = 2 • 10 −5 , N = 500, 000, and h 2 δ = 0.50.For an observed χ 2 = 9 a standard χ 2 would yield p = 2.6 • 10 −3 (shaded area under the blue curve to the right of the threshold), whereas with the noncentral χ 2 would yield p = 0.22 (shaded area under the pink curve), illustrating the miscalibration induced by the polygenicity of the target trait.

Figure 2 :
Figure 2: Dependence of χ 2 statistics on the precision τ 2 of the prediction under the null and alternative.We calculated expected association χ 2 statistics from simulations with different heritability values h 2 δ of the target trait, number of causal SNPs M, and sample sizes N. We took the average over 1,000 simulations, as well as over h 2 δ and M. Dash-dotted lines indicate the predicted χ 2 statistics from equation (7) under the null and equation (4) under the alternative.

Figure 3 :
Figure 3: Expected χ 2 vs. heritability, sample size, and inverse of the number of causal SNPs.We calculated the expected association χ 2 statistics from simulations with different number of causal SNPs M, heritability values h 2 δ , and sample sizes N. We took the average over 1,000 simulations and 99 independent mediators T .Panel (a) shows the average χ 2 vs 1 + h 2 δ N/M. .Panels (b-d) show the average χ 2 against 1/M, h 2 δ , and N. The horizontal dashed line at 1 indicates where a calibrated χ 2 statistics should be.The error bars represent the 95% confidence intervals from the simulations.Dash-dotted lines in the figure show the predicted χ 2 statistics from equation (7)under the null-the expected χ 2 mostly obscures the dash-dotted line, indicating the linear relationship is consistent with the theoretical approximation; We used M − 1 for M, which further improved the match.

Figure 4 :
Figure 4: Expected χ 2 under the alternative.We calculated expected χ 2 statistics by averaging 1000 simulations of target trait Y with different combinations of number of causal SNPs M, polygenic portion of the heritability of Y h 2 δ , and sample sizes N. Panel (a) shows the expected χ 2 against the theoretically predicted value in equation (4) assuming independent SNPs.Dashed gray line shows the identity line.Panel (b-d) show the expected χ 2

Figure 5 :
Figure 5: Expected χ 2 in TWAS.We calculated the expected χ 2 statistics in TWAS using simulated target traits with different heritability values h 2 δ of the target trait and sample sizes N. We averaged over 1000 simulated Y for a given gene, heritability, and sample size to estimate the Eχ 2 .Panel (a) shows the expected χ 2 averaged over genes as a function of the heritability h 2 δ , with different lines representing different sample sizes N. Panel (b) shows the expected χ 2 statistics for each gene against the heritability h 2 δ of the target trait averaged over different sample sizes.Colored lines in this panel represent selected genes; we highlight the 10 genes with the highest inflation.The lines in this panel are closely linear as predicted by the formula in equation (4).Panel (c) shows the expected χ 2 against the sample size N, with different colors showing different heritability values h 2 δ .Panel (d) shows the χ 2 statistics for each gene against the sample size N averaged over different heritability values h 2 δ .The lines in this panel are also closely linear, consistent with the formula in equation (4).

Figure 6 :
Figure 6: Distribution of Estimated Φ.Inflation factors Φ for gene expression, metabolites, and brain features (diffusion MRI) are shown in the log10 scale.The factor Φ for each mediator is estimated using the average χ 2 statistics of the association between genetically predicted mediator and 1,000 simulated target traits for each combination of heritability of target trait h 2 δ and sample size N.The slope of the regression of Eχ 2 on Nh 2 δ is used to estimate Φ.Most values (78% of genes, 94% of metabolites, 100% of brain features) are in the range of 10 −5 .
γ and using this for the final version of Φ We can further simplify the expression of Φ by noting that the numerator and denominator of γ′ • R 2 • γ γ′ • R • γ converge to their expectations as N increases.The sample covariance R follows a Wishart distribution:

Figure S8 :
Figure S8: Inflation in BrainXcan.Expected χ 2 statistics in TWAS calculated with simulated target traits with different heritability values of the target trait h 2 δ and sample sizes N.Estimated χ 2 is calculated averaging over 1000 simulated Y for a given brain feature, heritability, and sample size.a shows the average Eχ 2 over 308 brain features χ 2 vs h 2 δ .b) shows the Eχ 2 for each of the top 20 most inflated brain features vs h2.c) shows the average Eχ 2 over 308 brain features vs N. d) shows the Eχ 2 for each brain feature vs N.