When calculating the required sample size for a proportion ci, if n=50.13, the sample size we use is

Mendelian randomization (MR) study has become a powerful approach to assess the potential causal effect of a risk exposure on an outcome. Most current MR studies are conducted under the two-sample setting by combining summary data from two separate genome-wide association studies (GWAS), with one providing measures on associations between genetic markers and the risk exposure, and the other on associations between genetic markers and the outcome. We develop a power calculation procedure for the general two-sample MR study, allowing for the use of multiple genetic markers, and shared participants between the two GWAS. This procedure requires a few easy-to-interpret parameters and is validated through extensive simulation studies.

Keywords: genome-wide association studies, Mendelian randomization, power calculation, two-stage least squares estimator

Mendelian randomization (MR) analysis uses genetic variants as instrumental variables (IVs) to estimate the causal effect of a risk factor on an outcome based on observational studies (Davies, Holmes, & Davey Smith, 2018; Lawlor, Harbord, Sterne, Timpson, & Davey Smith, 2008; Smith & Ebrahim, 2003). By taking advantage of many robust findings from large-scale genome-wide association studies (GWAS) within the past decade, MR analysis is becoming a powerful procedure to elucidate causal relationships between risk exposures and disease outcomes (Hemani et al., 2018; MacArthur et al., 2017; Staley et al., 2016; Z. Zhu et al., 2016).

MR studies can be conducted under the one-sample setting, where measurements on IVs, risk exposure, and the outcome are available on all study participants. Most recent MR studies are often conducted under the two-sample setting by combining data from two separate GWAS, with one study providing measures on IVs and the risk factor, and the other on IVs and the outcome (Holmes et al., 2014; Lyall et al., 2017; Trajanoska et al., 2018). In some cases, due to the nature of GWAS consortia, the two GWAS studies used in the MR study can have nonnegligible overlaps, with a noticeable proportion of subjects participating in both studies (Burgess, Davies, & Thompson, 2016). This nonindependent two-sample setting can also occur in the subsample design in which the genetic study of the risk exposure is part of a larger GWAS of the outcome (Pierce & Burgess, 2013).

Several power calculation methods based on the two-stage least squares estimator (2SLS) have been developed for the design of MR studies, but they all focused on the one-sample MR study with a single IV (Brion, Shakhbazov, & Visscher, 2013; Burgess, 2014; Freeman, Cowling, & Schooling, 2013; Pierce, Ahsan, & Vanderweele, 2011). For MR study with multiple IVs, a common strategy is to convert them into a single IV called genetic risk score, which is a weighted sum of high-risk allele counts, and then apply the single IV power calculation procedure assuming weights are known without uncertainty (Burgess, 2014; Freeman et al.,2013).

Given that many MR studies are conducted under the general two-sample setting (with or without overlapping participants), with multiple IVs, we aim at developing a power calculation procedure for MR studies conducted under such more general circumstances. We consider both the outcome and risk exposure as continuous variables.

Let Y represent the continuous outcome, X represent the continuous exposure measure on the risk factor, and Z be the vector of random variables representing the set of k IVs. In MR studies, IVs are chosen genetic markers, such as single nucleotide polymorphisms (SNPs) with each element of Z taking value of 0, 1, or 2. Following standard conditions for the IVs analysis, we assume that they are connected with the following two regression models,

where β is the coefficient for the causal effect of interest, Π is the vector of regression coefficients for the k IVs, and u and v are two correlated error terms. We call (1) and (2) the first, and second stage regression model, respectively. To simplify the notation, we assume there is no intercept term or adjusted covariates in models (1) and (2), with E(X) = 0, E(Y) = 0, and E(Z) = 0. For more general model, we can replace (X, Y, Z) with their corresponding residuals after projecting them onto the space spanned by vectors of covariates (including the vector of 1s for the intercept term). We denote the variance-covariance matrix of Z as Qz, and let Var(u)=σu2, Var(v)=σv2, and Cov(u,v)=σuv. Based on (1) and (2), we have

with ω≜vβ+u. We define the variance of ω as σω2, which can be expressed as

We consider the general two-sample setting, where we have measures on IVs and the outcome from a genetic association study with n1 + n2 subjects, and measures on IVs and the exposure from another association study with n2 + n3 subjects. The two studies have n2 overlapping subjects. We call the data from the first study the outcome data, the one from the second study the risk exposure data. We can break all data into three exclusive subsets. For the first subset, we have n1 subjects on which we have measures on (Y, Z), and denote them as a n1 × 1 vector Y1, and a n1 × k matrix Z1. For the second subset, we have n2 subjects on which we have measures on (Y, X, Z), and denote them as n2 × 1 vectors Y2 and X2, and a n2 × k matrix Z2. For the third subset, we have n3 subjects on which we have measures on (X, Z), and denote them as a n3 × 1 vector X3, and a n3 × k matrix Z3. Let N=n1+n2+n3 be the total sample size. When the outcome data and the risk exposure data are completely overlapped (i.e., n1 = n3 = 0), we call the MR study a one-sample study. When the two data sets have no overlapping samples (i.e., n2 = 0), we call the MR study an independent two-sample study. When the risk exposure data set is a subset of the outcome data set (i.e., n3 = 0), we call the MR study a subsample study.

For the power calculation, we focus on the standard 2SLS estimator, which can be expressed as,

β^={Π^′(Z1′Z1+Z2′Z2)Π^}−1Π^′(Z1′Y1+Z2′Y2),

where

Π^=(Z2′Z2+Z3′Z3)−1(Z2′X2+Z3′X3).

To derive the asymptotic variance of the 2SLS estimator, we assume ni/N→ρi,i=1,2,3, with

N=n1+n2+n3 being the total sample size. In the Appendix I we show that the asymptotic variance of β^ can be represented as

var(β^)={N(ρ1+ρ2)Π′QzΠ}−1{σω2+ρ1+ρ2ρ2+ρ3β02σv2−2ρ2ρ2+ρ3(β02σv2+β0σuv)}={N(ρ1+ρ2)r21−r2}−1{σω2σv2+ρ1+ρ2ρ2+ρ3β02−2ρ2ρ2+ρ3(β02+β0σuvσv2)}≡N−1σ2(β0),

(5)

where r2=Π′QzΠ/var(X), with Qz being the variance-covariance matrix for Z. Notice that r2 can be interpreted as the proportion of the total variation of X explained by Z. According to (4), we can also represent σ2(β0) as

σ2(β0)={(ρ1+ρ2)r21−r2}−1(σu2σv2+ρ1+ρ3ρ2+ρ3β02+2ρ3β0ρ2+ρ3σuvσv2).

(6)

For an one-sample MR study, (5) reduces to the one given by Burgess (2014). For an independent two-sample MR study, (5) becomes the same as the one given by Inoue and Solon (2010).

Denote the null hypothesis as H0:β=0, and the alternative hypothesis as H1:β=β0. We can reject H0 at the level of significance of α if

where zα is the α quantile of the standard normal distribution, and N−1σ˜2(β^) is the estimated variance, with consistent estimates being plugged in for ρi, σω2, σv2, and σuv.

For a given sample size N, with ni/N=ρi,i=1,2,3, the power function of the above test can be written as

power(β0)=Pr(|Nσ˜−1(β^)β^|>z1−α/2H1)≈Φ−(z1−α/2+Nσ−1(β0)β0)+Φ(−z1−α/2−Nσ−1(β0)β0).

(7)

To estimate the power, or equivalently to calculate σ2(β0), we need to assume values for r2, σω2, σv2, σuv. Here we provide practical guidance on choosing those population level parameters.

We assume the variance of the outcome (denoted as σY2), the variance of the exposure (denoted as σX2) are known from existing literature. Also, we assume we know r2, the proportion of variance of X attributed to IVs. If there are overlapping samples between the outcome data and the risk exposure data, we also have to assume a value for σuv, which might be estimated from an existing study with measures on (Y, X, Z). With these values, we can estimate remaining parameters needed for σ2(β0) as following,

We do not need to know σuv to calculate the power for an independent two-sample MR study as var (β^) defined by (5) does not involve with σuv when ρ2 = 0.

When we have access to summary statistics on each IV from a GWAS of Y (the outcome) and a GWAS of X (the exposure), we can use those data to estimate σX2, σY2, r2, and σuv. Here is a brief description of the procedure. Suppose the sample size is nY for the outcome GWAS, and nX for the exposure GWAS. Let the summary statistics (the coefficient estimate for the SNP-outcome association, and its standard error) on the lth IV from the outcome GWAS of be β^Yl and σYl, l = 1, …, k, with k being the number of considered IVs. Similarly, let β^Xl and σXl be the summary statistics for the lth IV from the risk exposure GWAS. Following Deng, Zhang, Song, and Yu (2019) we can estimate σY2 as

σY2≈1k∑l=1k{(nY−2)σYl2+β^Yl2}σZl2,

with σZl2 being the variance of the lth IV, which can be estimated from an existing set of reference genomes, such as the 1000 Genomes Project (Genomes Project Consortium et al., 2015). σX2 can be estimated similarly. r2 can be estimated as

where β^X=(β^X1,…,β^Xk)′ and Σij=ρijσzl−1σzj−1, with ρij is the correlation coefficients between the ith IV and the jth IV. Again, we can obtain an estimate of Σ using a set of reference genomes. More details are given in Deng et al. (2019).

If the two GWAS share ns samples, we can estimate σuv using the idea of Liu and Lin (2018). Similar idea has been suggested by others, for example, Province and Borecki (2013), Bulik-Sullivan et al. (2015), and X. Zhu et al. (2015). For any given null SNP that is not associated with both Y and X, denote its summary statistics as (β^Yl∗, σ^Yl∗), and (β^Xl∗, σ^Xl∗). According to Liu and Lin (2018) we have

Cov(β^Yl∗σYl∗,β^Xl∗σXl∗)≈ns2nXnYCor(X,Y)=ns2nXnYσX2β02+σuvσXσY.

We can estimate Cov(β^Yl∗σYl∗,β^Xl∗σXl∗) using summary statistics on a set of independent null SNPs. With this estimated covariance, we are able to obtain an estimate of σuv.

In two-sample MR analysis, the variance for the causal effect estimate sometimes is obtained under the “NO Measurement Error” (NOME) assumption (Bowden et al., 2016), that is, the uncertainty in estimating Π is ignored. The same assumption was made in several power calculation procedures for the one-sample MR study (Burgess, 2014; Freeman et al., 2013). The NOME assumption in the general two-sample setting implicitly assumes that n3 is much larger than n1 and n2, or

Therefore, under NOME assumption, σ2(β0) can be approximated as

σNOME2(β0)={(ρ1+ρ2)r21−r2}−1σω2σv2.

(8)

Note that σNOME2(β0) is still a function of β0 as σω2 is influenced by β0 according to (4).

The corresponding power function can be approximated as

power*(β0)=Φ{−z1−α/2+NσNOME−1(β0)β0}+Φ{−z1−α/2−NσNOME−1(β0)β0}.

(9)

According to (8), it is clear that we only need to specify σX2, σY2, and r2 for the power calculation, as the value of σuv is not required in (8). By comparing (5) with (8), we can see that power*(β0) is also appropriate when β0 is close to 0.

Although power*(β0) is easy to use, but its validity relies on the NOME assumption, which might not be true in practice. It asserts that the power is only determined by the sample size of the outcome data, assuming the exposure data is sufficiently large. When n2 = 0, it is clear that σNOME2(β0)<σ2(β0) for any β0 ≠ 0. Therefore, power*(β0) overestimates the true power for an independent two-sample MR study. For a two-sample study with shared samples between the outcome data and exposure data, power*(β0) can either underestimate, or overestimate the true power depending on those modeling parameters influencing σ2(β0). Later we will demonstrate the difference between power(β0) and power*(β0) under various settings.

We first conducted simulation studies to verify the accuracy of the power function power(β0) given by (7). We considered a set of k = 5 independent IVs, each with the minor allele frequency of 0.3. We denoted (centralized) IVs as Z=(z1,…,zk)′. The two models for X and Y were X=∑j=1kπjzj+v, and Y=β0X+u, respectively. The two error terms u and v were both N(0, 2), with a correlation coefficient of 0.5. We chose sample sizes such that n=n1+n2=n2+n3, with n ∈ {1,000, 2,000, …, 10,000} and varied the overlapping percentage, which was defined as ψ=n2/(n1+n2), by changing n2. The value for Π=(π1,…,πk)′ was chosen to ensure r2=Π′QzΠ/var(X)=0.20, by assuming all πj,j=1,…,k, were the same. We let β0∈{−0.08,−0.05,0.0,0.05,0.08}.

Under any given (n,ψ,Π,β0) configuration, we simulated 10,000 data sets. By summarizing results over those simulated data sets under each condition, we obtained empirical size or power of the MR analysis at the significance level α = .05. Analytic power was calculated by (7) using true model parameters. Results are reported in Tables 1–3, with ψ = 1.0, 0.4, and 0.0 corresponding to the one-sample, two-sample with overlapping participants, and independent two-sample scenario, respectively. It is clear from Tables 1–3 that the theoretic power calculated by power(β) matches quite well with the empirical one in all considered cases, especially for studies with relatively large sample sizes.

Empirical and theoretic power comparison for one-sample Mendelian randomization study under the significance level α = .05

β0 = 0β0 = .05β0 = −.05β0 = .08β0 = −.08 n Empi.Analy.Empi.Analy.Empi.Analy.Empi.Analy.Empi.Analy.

1,000	4.99	5.00	16.24	12.40	8.06	12.40	29.49	24.41	17.88	24.41
2,000	5.31	5.00	24.55	20.09	15.73	20.09	47.10	43.21	38.82	43.21
3,000	4.98	5.00	30.83	27.78	23.87	27.78	61.00	59.13	57.18	59.13
4,000	4.93	5.00	38.38	35.26	31.76	35.26	72.63	71.56	70.67	71.56
5,000	5.01	5.00	45.07	42.38	39.19	42.38	80.58	80.74	80.55	80.74
6,000	5.34	5.00	51.32	49.06	46.39	49.06	86.82	87.25	87.79	87.25
7,000	5.22	5.00	55.93	55.24	54.47	55.24	90.90	91.72	92.55	91.72
8,000	4.82	5.00	61.73	60.87	60.02	60.87	94.27	94.71	95.34	94.71
9,000	4.92	5.00	66.28	65.97	65.45	65.97	96.10	96.67	97.21	96.67
10,000	4.86	5.00	72.28	70.54	68.94	70.54	97.56	97.93	98.29	97.93

Empirical and theoretic power comparison for the independent two-sample Mendelian randomization study under the significance level α = .05

β0 = 0β0 = .05β0 = −.05β0 = .08β0 = −.08 n Empi.Analy.Empi.Analy.Empi.Analy.Empi.Analy.Empi.Analy.

1,000	4.63	5.00	10.42	12.01	11.56	12.76	20.83	22.74	24.67	25.82
2,000	4.80	5.00	18.77	19.29	20.39	20.82	39.64	40.18	45.00	45.71
3,000	4.75	5.00	25.35	26.58	28.72	28.85	54.63	55.40	62.53	62.11
4,000	4.70	5.00	33.21	33.72	36.64	36.63	66.80	67.72	74.24	74.51
5,000	5.11	5.00	40.22	40.57	43.65	44.00	77.11	77.20	83.16	83.36
6,000	4.66	5.00	47.90	47.03	49.38	50.86	84.92	84.22	89.26	89.40
7,000	4.97	5.00	52.64	53.04	57.26	57.16	88.88	89.27	93.51	93.38
8,000	4.82	5.00	58.19	58.59	62.87	62.86	92.79	92.81	95.83	95.94
9,000	4.91	5.00	63.42	63.65	68.26	67.97	95.52	95.25	97.60	97.55
10,000	5.17	5.00	68.30	68.22	72.30	72.51	97.05	96.90	98.46	98.54

For studies with relatively small sample size (e.g., for the one-sample MR study with n < 5,000 in Table 1), their empirical powers are noticeably different from corresponding analytic ones. This is caused by the bias in the TS2SLS estimate β^ when using relatively weak IVs. It is well known that β^ is not consistent when dealing with weak IVs (e.g., Stock, Wright, & Yogo, 2002). For example, the empirical mean of relative bias, defined as (β^−β0)/β0, observed in 10,000 simulated one-sample MR studies with n = 1,000 under β0 = 0.05 is around 11%. It drops to less than 2% when the sample size increases to 6,000. We can further explain the observed difference between empirical and analytic powers using the following bias formula given by Deng et al. (2019),

E(β^−β0)≈(−ρ1ρ1+ρ2β0+ρ2ρ1+ρ2σuvσv2)b(μ2,k),

where b(μ2,k)>0, with its definition given in Deng et al. (2019). Since σuvσv2=0.5, and |β0|≤0.08 in all settings considered in Tables 1–3, it is clear that we have E(β^−β0)>0 for both the one sample MR study (i.e., ρ1=ρ3=0,ρ2=1) considered in Table 1, and the two-sample study with overlapping subjects (i.e., ρ1=ρ3=3/8,ρ2=1/4) considered in Table 2. Therefore, we have |E(β^)|>|β0| when β0>0, and |E(β^)|<|β0| when β0<0. So, in Tables 1 and ​2 the empirical power for β0>0 appears to be larger than the analytic ones, and smaller for β0<0. For the independent two-sample MR study considered in Table 3 (i.e., ρ1=ρ3=0.5,ρ2=0), we always have |E(β^)|<|β0|. Thus, the empirical power tends to be smaller than the analytic one regardless of the sign of β0. As the sample size increases, the bias becomes negligible and the empirical power matches well with the analytic one in all cases.

Empirical and theoretic power comparison for the overlapping two-sample Mendelian randomization study under the significance level α = .05

β0 = 0β0 = .05β0 = −.05β0 = .08β0 = −.08 n Empi.Analy.Empi.Analy.Empi.Analy.Empi.Analy.Empi.Analy.

1,000	4.92	5.00	13.31	12.16	10.19	12.62	25.10	23.37	21.58	25.24
2,000	4.91	5.00	20.38	19.60	18.31	20.52	42.23	41.34	42.19	44.68
3,000	4.91	5.00	28.66	27.05	26.19	28.41	58.31	56.84	59.04	60.89
4,000	4.82	5.00	35.91	34.32	34.11	36.07	70.11	69.22	72.78	73.31
5,000	4.54	5.00	41.89	41.27	42.67	43.34	78.23	78.60	82.18	82.31
6,000	5.12	5.00	49.18	47.82	48.31	50.13	85.19	85.44	88.36	88.54
7,000	5.08	5.00	54.67	53.90	55.16	56.37	90.33	90.27	92.69	92.73
8,000	4.90	5.00	59.36	59.48	63.01	62.04	93.15	93.60	95.74	95.46
9,000	4.59	5.00	64.84	64.56	66.82	67.16	96.06	95.84	97.32	97.21
10,000	4.83	5.00	69.27	69.13	71.65	71.71	97.16	97.33	98.69	98.31

It can be noticed from Table 1 that the theoretic power for the one-sample MR study is a symmetric function of β0. This is expected as σ2(β0) is a symmetric function of β0 according to (6) with ρ1 = ρ3 = 0. Therefore, the power function given by (7) is also symmetric. However, this symmetric property does not hold if ρ3 > 0 (e.g., for an independent, or a partially overlapping two-sample MR study shown in Tables 2 and ​3).

For simulations shown in Tables 1–3, we always let the outcome data and the exposure data have the same sample size (i.e., =n1+n2=n2+n3), and σuv/σv2=0.5. Thus, Formula (6) can be simplified as

N−1σ2(β0)={nr21−r2}−1{σu2σv2+2(1−ψ)(β02+12β0)}.

It is clear in this case that N−1σ2(β0) increases as ψ decreases from 1 to 0 when β0>0. As a result, for β0>0, the power function decreases as the overlapping percentage decreases. This trend would reverse if β02+12β0<0, which is true as long as β0∈(−0.5,0). We can confirm those patterns by comparing theoretic powers among Tables 1–3 (with ψ changing from 100%, to 40%, and to 0%).

In Figure 1 we show the power curves generated by the power function power(β) under α = .05. The power was calculated for various r2 and sample size n, with ψ = 0.4, and β = .05. As expected, it is evident from Figure 1 that the power for detecting the causal effect increases as the sample size, or r2 increases.

Power curve under the significance level of 0.05 and true effect β0 = .05, for a two-sample Mendelian randomization study with its risk exposure and outcome data having the same sample size (n) and 40% overlap percentage. Curves are for instrumental variables with different strength (r2)

In Figure 2 we compared power(β) with power*(β) under the subsample design, in which we set r2 = 0.04, n3 = 0, and n1 + n2 = 5,000, with ψ=n2/(n1+n2) being 10%, 20%, 60%, or 100%. It appears that the two functions match well around β = 0, but do not agree with each other for detecting relatively large causal effect, especially for the subsample design with a small overlap percentage ψ. It can also be noticed that under the subsample design power*(β) approximates power(β) very well for β0∈[0, 0.4] when ψ = 20% (see Figure 2c), while there are noticeable differences between the two functions under other ψ values (see Figure 2). This happens due to parameters we chose for the power calculation in Figure 2. Since we let σuv/σv2=0.5, r2 = 0.04, and ρ3 = 0 for all cases in Figure 2, by (5) and (8) we have

N−1σ2(β0)−N−1σNOME2(β0)=24{(ψ−1−2)β02−β0}.

Comparison between the theoretic power function power(β) and its approximation power*(β) for a Mendelian randomization study under the subsample design. The sample size for the outcome data is fixed at n1 + n2 = 5,000. The sample size for the exposure data (n2) varies from 5,000 to 500

This difference remains to be very small for all β0∈[0,0.4] when ψ=20%, but it is not the cases for other ψ values.

In Figure 3 we compared power(β) with power*(β) under an independent two-sample design, in which we set r2 = 0.02, n2 = 0, n1 = 5,000, and n3 = 0.25n1, 0.5n1, 1.0n1, or 4.0n1. It is clear from Figure 3 that power*(β) approaches power(β) as n3/n1 increases. This is expected as the NOME assumption becomes appropriate when n3 is much larger than n1.

Comparison between the theoretic power function power(β) and its approximation power*(β) for an independent two-sample Mendelian randomization study. The sample size for the outcome data (n1) is fixed at 5,000. The sample size for the exposure data (n3) varies from 1,250 to 20,000

We used the MR study conducted by Geng and Huang (2018) as an example to demonstrate the application of the power calculation procedure. Geng and Huang (2018) recently conducted a MR analysis to evaluate the causal effect of maternal obesity on birth size. They considered three risk factors representing maternal central obesity, and three outcome measures on the different birth size. For the purpose of illustration, we focused on the risk exposure defined by mother’s BMI adjusted hip circumference (HIPadjBMI), and chose birth length as the outcome. There was no overlapping sample in this MR study. A total of 41 independent SNPs were selected as IVs for HIPadjBMI. Summary statistics on the SNP-exposure and SNP-outcome associations were provided in tables given by Geng and Huang (2018). We first estimated σX2, σY2, and r2 using those summary statistics as described in Section 2.2. There was no need to estimate σuv as we were considering an independent two-sample MR study. With estimated values σX2=1.33, σY2=1.14, and r2 = 1.3 %, and by fixing the sample size for the risk exposure data at nX = 117,321, the same as that of Geng and Huang (2018), we conducted the power calculation under various sample sizes of the outcome data, and causal effect sizes (Figure 4).

Power calculation for the Mendelian randomization study of maternal obesity and birth length. The sample size of the exposure data is fixed at 117,321. Variances of the outcome and exposure, and the proportion of exposure variance explained by all selected instrumental variables are estimated from existing genome-wide association studies

We develop an easy to use power calculation procedure for general MR studies with multiple IVs. The procedure can be used for various MR study designs, including the one-sample, subsample, and two-sample with or without overlapping subjects between the outcome data and the exposure data. The proposed procedure estimates the power without making the NOME assumption. Unless the causal effect is small, or the sample size of the exposure data is relatively large compared to that of the outcome data, we show that calculating power under the NOME assumption is not appropriate. For an independent two-sample MR study, the power calculation procedure just requires the knowledge on the variance of the outcome, the variance of the exposure, as well the proportion of variance in the exposure that can be explained by selected IVs. For MR study under other designs, due to the complication introduced by shared samples between the outcome data and exposure data, the power calculation procedure requires additional knowledge of the unobserved confounding effect on the relationship between the outcome and the risk exposure, which can be estimated using summary statistics generated from GWAS of the outcome, and GWAS of the risk exposure if the two GWAS have overlapping participants. R code for the power calculation and instructions on how to use it are provided in Appendix II.

The proposed power calculation procedure is only applicable to MR studies with strong IVs, as it relies on the asymptotic normal distribution of the 2SLS estimate. For weak IVs, it is well known that the distribution of 2SLS estimate cannot be properly approximated by the normal distribution even with very large sample size. For MR studies with weak IVs, a different test statistic should be adopted, and therefore a different power calculation procedure needs to be developed.

The power calculation procedure requires the knowledge of several population level parameters, including the variance of the outcome, and risk exposure, as well as the proportion of variation in the risk exposure that can be explained by IVs (r2). With summary statistics from existing GWAS, we usually can reliably estimate variances for the outcome, and risk exposure. But we might overestimate r2 if summary data from the discovery GWAS is used for estimating r2, due to the selection bias (Yu et al., 2007; Zhong & Prentice, 2008). To be more conservative in power calculation, we recommend using the adjusted estimate of the regression coefficient (Zhong & Prentice, 2008) to obtain a less biased r2 estimate.

We only focus on MR studies with continuous outcomes as the property of 2SLS estimate under the linear regression model is well understood. For the study of a binary outcome under a logistic regression model, a similar two-stage regression model-based estimate can be derived. However, due to noncollapsibility of the logistic regression model, this type of estimate is not consistent. Further investigations are needed to develop power calculation procedures for MR studies with binary outcomes.

Click here to view.(26K, docx)

DATA AVAILABILITY STATEMENT

Data sharing is not applicable to this article as no new data were created or analyzed in this study.

SUPPORTING INFORMATION

Additional supporting information may be found online in the Supporting Information section.

CONFLICT OF INTERESTS

The authors declare that there are no conflict of interests.

• Bowden J, Del Greco MF, Minelli C, Davey Smith G, Sheehan NA, & Thompson JR (2016). Assessing the suitability of summary data for two-sample Mendelian randomization analyses using MR-Egger regression: The role of the I2 statistic. International Journal of Epidemiology, 45(6), 1961–1974. 10.1093/ije/dyw220.//www.ncbi.nlm.nih.gov/pubmed/27616674 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Brion MJ, Shakhbazov K, & Visscher PM (2013). Calculating statistical power in Mendelian randomization studies. International Journal of Epidemiology, 42(5), 1497–1501. 10.1093/ije/dyt179. //www.ncbi.nlm.nih.gov/pubmed/24159078 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Bulik-Sullivan B, Finucane HK, Anttila V, Gusev A, Day FR, Loh PR, … Neale BM (2015). An atlas of genetic correlations across human diseases and traits. Nature Genetics, 47(11), 1236–1241. 10.1038/ng.3406. //www.ncbi.nlm.nih.gov/pubmed/26414676 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Burgess S. (2014). Sample size and power calculations in Mendelian randomization with a single instrumental variable and a binary outcome. International Journal of Epidemiology, 43(3), 922–929. 10.1093/ije/dyu005. //www.ncbi.nlm.nih.gov/pubmed/24608958 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Burgess S, Davies NM, & Thompson SG (2016). Bias due to participant overlap in two-sample Mendelian randomization. Genetic Epidemiology, 40(7), 597–608. 10.1002/gepi.21998. //www.ncbi.nlm.nih.gov/pubmed/27625185 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Davies NM, Holmes MV, & Davey Smith G. (2018). Reading Mendelian randomisation studies: A guide, glossary, and checklist for clinicians. BMJ, 362, k601. 10.1136/bmj.k601. //www.ncbi.nlm.nih.gov/pubmed/30002074 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Deng L, Zhang H, Song L, & Yu K. (2019). Approximation of bias and mean-squared error in two-sample Mendelian randomization analyses [published online ahead of print October 25, 2019]. Biometrics, 10.1111/biom.13169. //www.ncbi.nlm.nih.gov/pubmed/31651042 [PMC free article] [PubMed] [CrossRef]
• Freeman G, Cowling BJ, & Schooling CM (2013). Power and sample size calculations for Mendelian randomization studies using one genetic instrument. International Journal of Epidemiology, 42(4), 1157–1163. 10.1093/ije/dyt110. //www.ncbi.nlm.nih.gov/pubmed/23934314 [PubMed] [CrossRef] [Google Scholar]
• Geng TT, & Huang T. (2018). Maternal central obesity and birth size: A Mendelian randomization analysis. Lipids in Health and Disease, 17(1), 181. 10.1186/s12944-018-0831-4. //www.ncbi.nlm.nih.gov/pubmed/30064420 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Genomes Project Consortium, Auton A, Brooks LD, Durbin RM, Garrison EP, Kang HM, … Abecasis GR (2015). A global reference for human genetic variation. Nature, 526(7571), 68–74. 10.1038/nature15393. //www.ncbi.nlm.nih.gov/pubmed/26432245 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Hemani G, Zheng J, Elsworth B, Wade KH, Haberland V, Baird D, … Haycock PC (2018). The MR-Base platform supports systematic causal inference across the human phenome. eLife, 7, 10.7554/eLife.34408. //www.ncbi.nlm.nih.gov/pubmed/29846171 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Holmes MV, Lange LA, Palmer T, Lanktree MB, North KE, Almoguera B, … Keating BJ (2014). Causal effects of body mass index on cardiometabolic traits and events: A Mendelian randomization analysis. American Journal of Human Genetics, 94(2), 198–208. 10.1016/j.ajhg.2013.12.014. //www.ncbi.nlm.nih.gov/pubmed/24462370 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Inoue A, & Solon G. (2010). Two-sample instrumental variables estimators. The Review of Economics and Statistics, 92(3), 557–561. [Google Scholar]
• Lawlor DA, Harbord RM, Sterne JA, Timpson N, & Davey Smith G. (2008). Mendelian randomization: Using genes as instruments for making causal inferences in epidemiology. Statistics in Medicine, 27(8), 1133–1163. 10.1002/sim.3034. //www.ncbi.nlm.nih.gov/pubmed/17886233 [PubMed] [CrossRef] [Google Scholar]
• Liu Z, & Lin X. (2018). Multiple phenotype association tests using summary statistics in genome-wide association studies. Biometrics, 74(1), 165–175. 10.1111/biom.12735. //www.ncbi.nlm.nih.gov/pubmed/28653391 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Lyall DM, Celis-Morales C, Ward J, Iliodromiti S, Anderson JJ, Gill JMR, … Pell JP (2017). Association of body mass index with cardiometabolic disease in the UK Biobank: A Mendelian Randomization Study. JAMA Cardiol, 2(8), 882–889. 10.1001/jamacardio.2016.5804. //www.ncbi.nlm.nih.gov/pubmed/28678979 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• MacArthur J, Bowler E, Cerezo M, Gil L, Hall P, Hastings E, … Parkinson H. (2017). The new NHGRI-EBI Catalog of published genome-wide association studies (GWAS Catalog). Nucleic Acids Research, 45(D1), D896–D901. 10.1093/nar/gkw1133. //www.ncbi.nlm.nih.gov/pubmed/27899670 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Pierce BL, Ahsan H, & Vanderweele TJ (2011). Power and instrument strength requirements for Mendelian randomization studies using multiple genetic variants. International Journal of Epidemiology, 40(3), 740–752. 10.1093/ije/dyq151. //www.ncbi.nlm.nih.gov/pubmed/20813862 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Pierce BL, & Burgess S. (2013). Efficient design for Mendelian randomization studies: Subsample and 2-sample instrumental variable estimators. American Journal of Epidemiology, 178(7), 1177–1184. 10.1093/aje/kwt084. //www.ncbi.nlm.nih.gov/pubmed/23863760 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Province MA, & Borecki IB (2013). A correlated meta-analysis strategy for data mining “OMIC” scans. Pacific Symposium on Biocomputing, 236–246. //www.ncbi.nlm.nih.gov/pubmed/23424128 [PMC free article] [PubMed]
• Smith GD, & Ebrahim S. (2003). ‘Mendelian randomization’: Can genetic epidemiology contribute to understanding environmental determinants of disease? International Journal of Epidemiology, 32(1), 1–22. //www.ncbi.nlm.nih.gov/pubmed/12689998 [PubMed] [Google Scholar]
• Staley JR, Blackshaw J, Kamat MA, Ellis S, Surendran P, Sun BB, … Butterworth AS (2016). PhenoScanner: A database of human genotype-phenotype associations. Bioinformatics, 32(20), 3207–3209. 10.1093/bioinformatics/btw373. //www.ncbi.nlm.nih.gov/pubmed/27318201 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Stock J, Wright J, & Yogo M. (2002). A survey of weak instruments and weak identification in generalized method of moments. Journal of Business and Economic Statistics, 20(4), 518–529. [Google Scholar]
• Trajanoska K, Morris JA, Oei L, Zheng HF, Evans DM, Kiel DP, Ohlsson C, Richards JB, & Rivadeneira F. (2018). Assessment of the genetic and clinical determinants of fracture risk: Genome wide association and mendelian randomisation study. BMJ, 362, k3225. 10.1136/bmj.k3225. //www.ncbi.nlm.nih.gov/pubmed/30158200. [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Yu K, Chatterjee N, Wheeler W, Li Q, Wang S, Rothman N, & Wacholder S. (2007). Flexible design for following up positive findings. American Journal of Human Genetics, 81(3), 540–551. 10.1086/520678. //www.ncbi.nlm.nih.gov/pubmed/17701899 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Zhong H, & Prentice RL (2008). Bias-reduced estimators and confidence intervals for odds ratios in genome-wide association studies. Biostatistics, 9(4), 621–634. 10.1093/biostatistics/kxn001. //www.ncbi.nlm.nih.gov/pubmed/18310059 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Zhu X, Feng T, Tayo BO, Liang J, Young JH, Franceschini N, … Redline S. (2015). Meta-analysis of correlated traits via summary statistics from GWASs with an application in hypertension. American Journal of Human Genetics, 96(1), 21–36. 10.1016/j.ajhg.2014.11.011. //www.ncbi.nlm.nih.gov/pubmed/25500260 [PMC free article] [PubMed] [CrossRef] [Google Scholar]
• Zhu Z, Zhang F, Hu H, Bakshi A, Robinson MR, Powell JE, … Yang J. (2016). Integration of summary data from GWAS and eQTL studies predicts complex trait gene targets. Nature Genetics, 48(5), 481–487. 10.1038/ng.3538. //www.ncbi.nlm.nih.gov/pubmed/27019110 [PubMed] [CrossRef] [Google Scholar]