MR-RAPS for Weak Instruments and Balanced Pleiotropy

Why MR-RAPS is useful

In summary-data Mendelian randomization, a genetic variant may be only weakly associated with the exposure. Treating its estimated exposure association as if it were known exactly can distort the causal-effect estimate. At the same time, variants can have direct effects on the outcome that are not mediated by the exposure: horizontal pleiotropy. MR-RAPS (robust adjusted profile score) addresses these issues by accounting for uncertainty in exposure associations and allowing balanced, heterogeneous direct effects. Its robust score can reduce the influence of variants with unusually large residuals.

These protections depend on assumptions. In particular, MR-RAPS does not make arbitrary pleiotropy harmless: its central model requires direct effects to average to zero and to be independent of instrument strength. The method also needs an appropriate model for the uncertainty and dependence in the association estimates.

Model and assumptions

For variant \(j\), let \(\gamma_j\) be its true association with the exposure, \(\hat{\gamma}_j\) its estimated exposure association, and \(\hat{\Gamma}_j\) its estimated outcome association. Let \(\beta\) denote the causal effect per unit increase in the exposure, and let \(\alpha_j\) denote the variant's direct effect on the outcome. A common two-sample summary-data model is

Here \(e_{Xj}\) and \(e_{Yj}\) are association-estimation errors, with variances \(s_{Xj}^2\) and \(s_{Yj}^2\). Under the usual independent two-sample approximation, these errors are independent. The variants are assumed independent of one another; correlated variants require a model that accounts for their covariance.

For the normal random-effects version, the residual at a proposed causal effect is

The term \(\beta^2s_{Xj}^2\) is essential: uncertainty in the exposure association contributes to uncertainty in the outcome residual. The term \(\tau^2\) allows residual variation beyond the two association standard errors, as expected when direct effects vary across variants. Under the stated model, \(r_j\) has mean zero and variance \(v_j\). A weak instrument has relatively large exposure-association uncertainty compared with its true exposure association; MR-RAPS does not simply discard that uncertainty by treating \(\hat{\gamma}_j\) as exact.

The adjusted profile score

With a squared-error score, the adjusted estimating equation for \(\beta\) is

The first term resembles a residual-weighted association score. The second term is the adjustment: it corrects for the fact that the estimated exposure association also appears in the residual. Without it, the score is generally not centered at zero when exposure-association estimates are noisy. Solving an estimating equation means finding the value of \(\beta\) for which the sum of contributions is zero; it is not a calculation that can generally be replaced by a single ratio.

MR-RAPS generalizes this equation by replacing the linear residual score with a robust score. Define \(Z_j=r_j/\sqrt{v_j}\), and let \(\psi\) be the derivative of a chosen robust loss function. The normal-adjusted score takes the form

Here \(\psi'\) is the derivative of \(\psi\), where it exists, and the expectation in the adjustment is taken under a standard normal variable. This particular adjustment is justified by a normal residual model. If direct effects are assumed only to have mean zero and variance \(\tau^2\), but are not assumed normal, the residual need not be normal and this expectation is not automatically justified. In that case, the adjustment must be appropriate to the residual distribution rather than borrowed unchanged from the normal model.

Why the adjustment is needed

The key point can be seen from a short derivation. Conditional on the true \(\gamma_j\), suppose the association errors are independent normal variables and \(\alpha_j\) is also independent and normal with variance \(\tau^2\). Then \(\hat{\gamma}_j\) and \(Z_j\) are jointly normal, with

For jointly normal variables, the normal integration-by-parts identity gives \(E\{X\psi(Z)\}=E(X)E\{\psi(Z)\}+\operatorname{Cov}(X,Z)E\{\psi'(Z)\}\). Here \(Z_j\) is standard normal and a symmetric score such as Huber's has \(E\{\psi(Z_j)\}=0\). Consequently,

Adding the adjustment in \(U_\psi\) cancels this expectation, centering the score at the true \(\beta\) under the model. This justification uses the stated joint normality; mean zero and variance \(\tau^2\) alone do not establish the identity. In the squared-error case, \(\psi(z)=z\) and \(E\{\psi'(Z)\}=1\), recovering the displayed adjusted score.

Robust loss choices

A robust loss grows less rapidly than squared error for large standardized residuals. Its score \(\psi(z)\) therefore limits how much a variant with an extreme residual can influence the estimate. Robustness is not automatic protection against invalid instruments: it changes the influence of residual outliers, while the balanced-pleiotropy and other identification assumptions remain important.

For Huber's loss, the usual tuning constant \(c=1.345\) gives about 95% asymptotic efficiency relative to squared-error estimation under a normal location model. That efficiency statement is not the proportion of observations inside the linear region: \(P(|Z|\le1.345)\) for \(Z\sim N(0,1)\) is about \(0.821\), or 82.1%. For Tukey's bisquare, the tuning constant is often set near \(4.685\) for a similar normal-model efficiency target. Such tuning choices describe behavior under specified conditions; they do not guarantee valid causal inference in a particular analysis.

Estimating overdispersion and interpreting results

The overdispersion parameter \(\tau^2\) is estimated rather than assumed to be known. In practice, MR-RAPS procedures estimate it jointly with the causal effect using a profile-likelihood or robust scale procedure appropriate to the chosen implementation. A larger fitted \(\tau^2\) increases \(v_j\), so residuals are standardized against more variation than association-sampling error alone. The adjustment term also depends on \(v_j\). Thus changing \(\tau^2\) can change both residual influence and score weighting; its effect on the final estimate need not be in one direction.

Inference commonly uses an estimated standard error for the robust estimate and a confidence interval such as \(\hat{\beta}\pm1.96SE\) when the relevant large-sample approximation is adequate. A 95% confidence procedure has 95% coverage across repeated samples under its assumptions; it does not mean there is a 95% probability that a fixed causal effect lies in the particular computed interval. The causal interpretation also depends on instrument validity, population structure control, and the balanced-pleiotropy model. Robust fitting does not prove those assumptions.

Common pitfalls

• Dropping exposure uncertainty. Using only \(s_{Yj}^2+\tau^2\) in the residual variance treats \(\hat{\gamma}_j\) as exact. Include \(\beta^2s_{Xj}^2\).
• Calling balanced pleiotropy no pleiotropy. The model permits nonzero \(\alpha_j\); it requires their distribution to be centered and suitably independent of instrument strength.
• Confusing robust with assumption-free. A bounded score limits the influence of large residuals, but cannot repair directional pleiotropy or confounding of genetic associations.
• Using the normal adjustment without its justification. The identity involving \(E\{\psi'(Z)\}\) requires a suitable joint normal model (or another valid derivation for the actual residual distribution).
• Misreading Huber's tuning constant. At \(c=1.345\), about 82.1% of standard normal values fall in the linear region. The commonly cited 95% is an efficiency property under a normal location model.
• Treating a score at one candidate as the estimate. The estimate solves the estimating equation; an individual score evaluation only shows whether that candidate satisfies it.