Next, causal estimate was conducted by applying weighted generalized linear regression method using correlated IVs [27] to distinguish causality model from pleiotropy model. First, we assume the estimate of association for IV k = 1, …, K with the gene’s expression (X) is βkX with standard error σkX. The estimate of association for IV k with BMD can be expressed as βkY with standard error σkY. The correlation between IVs k1 and k2 can be defined as ρk1k2. Then we can perform a weighted generalized linear regression of the βkY parameters on the βkX parameters using the σkY−2 parameters as inverse-variance weights and taking into account the correlation between the IVs. If we define Ωk1k2 = σk1Yσk2Yρk1k2, bXY can be estimated from a weighted generalized linear regression as bXY = (βkXTΩβkX)−1βkXT Ω−1βkY [39], withβkXT being the transposed vector of βkX. The standard error of the estimate is se(bXY)=(βkXTΩ−1βkX)−1 [27]. After pruning for SNPs with linkage disequilibrium r2 > 0.8, all variants with a PeQTL < 5 × 10−8 in the cis-eQTL region of each probe available in both BMD and