Probit Marginalization With Random Slopes

Why Random Slopes Change the Marginal Mean

In a probit GLMM, a random intercept shifts a cluster’s conditional index, while a random slope lets the predictor’s conditional association vary across clusters. Averaging over both effects does more than shrink the fixed effects by a common factor: the variance of the random contribution changes with the predictor. Consequently, the population-averaged probit index is generally nonlinear in the predictor, even though every cluster’s conditional index is linear.

The exact probit marginalization result for Gaussian random intercepts gives the key identity. With random slopes, the same identity still applies, but the variance being marginalized is a function of the predictor. The result below is specific to a probit link and jointly Gaussian random effects; it should not be assumed for other links.

Model and assumptions: For cluster $i$ and observation $j$, let $Y_{ij}$ be binary and let $x_{ij}$ be a scalar predictor. The conditional model is $$ P(Y_{ij}=1\mid b_{0i},b_{1i},x_{ij})=\Phi(\beta_0+\beta_1x_{ij}+b_{0i}+b_{1i}x_{ij}), $$ where $\Phi$ is the standard normal cumulative distribution function. Assume the cluster-specific random effects have a bivariate normal distribution, $$ \begin{pmatrix}b_{0i}\\b_{1i}\end{pmatrix} \sim N\left[ \begin{pmatrix}0\\0\end{pmatrix}, \begin{pmatrix}\tau_{00}&\tau_{01}\\\tau_{01}&\tau_{11}\end{pmatrix} \right]. $$ The covariance matrix must be positive semidefinite: $\tau_{00}\geq 0$, $\tau_{11}\geq 0$, and $\tau_{00}\tau_{11}-\tau_{01}^{2}\geq 0$. Clusters are independent, and conditional on the random effects and covariates, the binary outcomes are independent Bernoulli observations. The marginal mean averages over the random-effects distribution while holding $x_{ij}$ fixed.

Deriving the Marginal Probit Mean

At a fixed predictor value $x$, the random contribution is $Z=b_{0i}+xb_{1i}$. A linear combination of jointly Gaussian variables is Gaussian, with mean zero and variance $$ \operatorname{Var}(Z)=\tau_{00}+2x\tau_{01}+x^2\tau_{11}. $$ Write this variance as $q(x)$. To marginalize the conditional probability, introduce an independent standard normal variable $\epsilon$. Since $\Phi(t)=P(\epsilon\leq t)$, averaging the conditional probability is equivalent to evaluating $P(\epsilon\leq \beta_0+\beta_1x+Z)$. The inequality can be rewritten as $\epsilon-Z\leq \beta_0+\beta_1x$. The variable $\epsilon-Z$ is normal with mean zero and variance $1+q(x)$, because $\epsilon$ and $Z$ are independent. Standardizing this normal variable gives the marginal probability.

\begin{aligned} P(Y_{ij}=1\mid x_{ij}=x) &=E\{\Phi(\beta_0+\beta_1x+Z)\}\\ &=P(\epsilon-Z\leq \beta_0+\beta_1x)\\ &=\Phi\left(\frac{\beta_0+\beta_1x}{\sqrt{1+q(x)}}\right)\\ &=\Phi\left(\frac{\beta_0+\beta_1x} {\sqrt{1+\tau_{00}+2x\tau_{01}+x^2\tau_{11}}}\right). \end{aligned}

The added $1$ is the variance of the standard normal variable implicit in the probit link. The denominator is always positive: $q(x)$ is a variance and therefore nonnegative. Define the marginal probit index $m(x)$ as the argument of $\Phi$:

m(x)=\frac{\beta_0+\beta_1x}{\sqrt{V(x)}}, \qquad V(x)=1+\tau_{00}+2x\tau_{01}+x^2\tau_{11}.

Thus the marginal mean remains a probit probability, but $m(x)$ is generally not a straight line. If $\tau_{11}=0$, positive semidefiniteness requires $\tau_{01}=0$, and the denominator is constant; this is the random-intercept-only case. If $\tau_{11}>0$, the predictor-dependent variance generally makes the marginal index curve.

Worked Example: Computing Marginal Probabilities

Suppose the conditional model has $\beta_0=-0.4$ and $\beta_1=0.8$, with $x$ measured in predictor units. Let $\tau_{00}=0.36$, $\tau_{01}=0.06$, and $\tau_{11}=0.16$. The covariance matrix is valid: its determinant is $0.36(0.16)-0.06^2=0.054>0$, and both variances are positive.

At $x=0$, $V(0)=1+0.36=1.36$, so $\sqrt{V(0)}=1.1662$. The marginal index is $-0.4/1.1662=-0.3430$, giving a marginal probability of approximately $\Phi(-0.3430)=0.3658$.

At $x=1$, calculate every variance term: $$ V(1)=1+0.36+2(1)(0.06)+1^2(0.16)=1.64. $$ Therefore $\sqrt{V(1)}=1.2806$, and $m(1)=0.4/1.2806=0.3123$. The marginal probability is approximately $\Phi(0.3123)=0.6226$.

At $x=2$, $V(2)=1+0.36+2(2)(0.06)+2^2(0.16)=2.24$. Hence $\sqrt{V(2)}=1.4967$, $m(2)=1.2/1.4967=0.8018$, and the marginal probability is approximately $\Phi(0.8018)=0.7887$. These are population-averaged predictions at the specified predictor values, not predictions for a cluster whose random effects are known.

How the Marginal Slope Varies

Differentiate the marginal index to see explicitly why its slope changes with $x$. Since $V'(x)=2(\tau_{01}+x\tau_{11})$, differentiating $(\beta_0+\beta_1x)V(x)^{-1/2}$ gives

\begin{aligned} m'(x) &=\frac{\beta_1V(x)-(\beta_0+\beta_1x)(\tau_{01}+x\tau_{11})} {V(x)^{3/2}}\\ &=\frac{\beta_1(1+\tau_{00})-\beta_0\tau_{01} +x(\beta_1\tau_{01}-\beta_0\tau_{11})} {V(x)^{3/2}}. \end{aligned}

The second equality follows by expanding the numerator: the terms involving $x^2\beta_1\tau_{11}$ cancel, leaving a numerator that is at most linear in $x$. Even so, the denominator varies with $x$, so $m'(x)$ is not generally constant. The marginal probability has derivative $\phi(m(x))m'(x)$, where $\phi$ is the standard normal density. This probability-scale change also depends on the current index, so it is not a constant effect per predictor unit.

Worked Example: Comparing Marginal Index Slopes

Use the parameters from the preceding example. The numerator of $m'(x)$ is $1.112+0.112x$, since $\beta_1(1+\tau_{00})-\beta_0\tau_{01}=0.8(1.36)-(-0.4)(0.06)=1.112$ and $\beta_1\tau_{01}-\beta_0\tau_{11}=0.8(0.06)-(-0.4)(0.16)=0.112$.

At $x=0$, $V(0)^{3/2}=1.36^{3/2}\approx1.586$, so $m'(0)=1.112/(1.36^{3/2})\approx0.7011$ index units per predictor unit. At $x=1$, $V(1)^{3/2}=1.64(1.2806)=2.100$, so $m'(1)=(1.112+0.112)/2.100=0.5828$ index units per predictor unit. The marginal index is increasing at both points, but its local slope is smaller at $x=1$. The conditional fixed-part slope, $\beta_1=0.8$, is constant; it is not the marginal slope.

When $\tau_{11}>0$, the variance grows quadratically in the tails. In particular, as $x$ tends to positive infinity, $m(x)$ approaches $\beta_1/\sqrt{\tau_{11}}$; as $x$ tends to negative infinity, it approaches $-\beta_1/\sqrt{\tau_{11}}$. These are limits of the marginal index, not claims that every cluster’s conditional index has the same direction of change. A cluster’s conditional index has slope $\beta_1+b_{1i}$, which can be negative even when $\beta_1>0$.

Worked Example: Interpreting Variation in Cluster Slopes

In the numerical model above, the random slope has standard deviation $\sqrt{0.16}=0.4$, so $b_{1i}\sim N(0,0.16)$. The conditional slope in cluster $i$ is $0.8+b_{1i}$. It is negative when $b_{1i}<-0.8$, which has probability $\Phi(-0.8/0.4)=\Phi(-2)=0.0228$, or about 2.28%.

Thus the fixed-part index, obtained by setting both random effects to zero, increases with $x$, while a small fraction of clusters have a decreasing conditional index under this model. The calculation concerns the random-slope distribution, not an observed study result. Correlation between the random intercept and slope affects their joint distribution and the marginal variance across $x$, but does not change this marginal probability for $b_{1i}$ alone.

Assumptions, Diagnostics, and Alternatives

The formula is exact under the stated model, not an approximation. Its interpretation depends on the predictor scale, the Gaussian random-effects assumption, and the covariance parameters. Centering $x$ can make the intercept refer to a more meaningful predictor value and may reduce the numerical dependence between estimated intercept and slope. Centering changes the parameterization and the reported intercept-slope covariance; it does not remove the need to model the covariance structure appropriately.

Useful checks include whether there is enough within-cluster predictor variation to estimate a random slope, whether the fitted covariance matrix is positive semidefinite, and whether the model-implied predictions are plausible across the observed predictor range. A variance estimate near zero can indicate little support for the corresponding random-effect variation, though it is not by itself proof that the model is wrong. Also inspect sensitivity to the assumed random-effects distribution and whether a few clusters dominate the estimated slope variation. Extrapolation deserves particular caution: the limiting behavior in the tails is a mathematical property of the model, not a guarantee that the model describes data beyond the observed range.

A conditional prediction for a cluster with specified random effects uses $\Phi(\beta_0+\beta_1x+b_{0i}+b_{1i}x)$. Setting the random effects to zero instead gives a prediction for a cluster at the center of the assumed random-effects distribution, not the population average. The marginal prediction averages over the random effects and uses the derived denominator. A GEE specifies a marginal mean directly and targets population-averaged associations without requiring this particular random-effects distribution. These approaches answer different modeling questions; agreement should not be presumed.

Common pitfall: It is incorrect to call $\beta_1$ the marginal probit slope, or to say that every cluster’s conditional index increases whenever $\beta_1>0$. The marginal slope is $m'(x)$, which depends on $x$ and the covariance parameters; the cluster-specific slope is $\beta_1+b_{1i}$. State which scale and which prediction target an interpretation concerns.

When the Result Is Useful

This calculation makes explicit how a conditional GLMM translates into population-averaged predictions. It is particularly useful when the scientific target is a population probability at particular predictor values but the fitted model includes cluster-level heterogeneity in both baseline risk and predictor association. It also explains why reporting only the conditional coefficient can conceal meaningful changes in the marginal relationship.

The marginal probability curve is still monotone wherever $m'(x)$ has a fixed sign, because $\Phi$ is increasing. But it need not have a constant probit slope, and a positive conditional fixed slope alone does not establish that the marginal slope is positive at every predictor value. Check the derivative or the predictions over the range relevant to the analysis.

Check Your Understanding

Answer each question, then compare with the concise solutions below.

For random intercept and slope variances $\tau_{00}$, $\tau_{11}$ and covariance $\tau_{01}$, what is the variance of the random contribution at predictor value $x$?
Why does the marginal probit denominator contain a $1$ in addition to the random-effects variance?
In the numerical example, verify $V(1)$ and state the marginal probability at $x=1$.
What is the cluster-specific conditional slope, and why does $\beta_1>0$ not ensure that all such slopes are positive?
When $\tau_{11}>0$, what is the limit of the marginal probit index as $x$ tends to positive infinity?

Solutions. (1) It is $\tau_{00}+2x\tau_{01}+x^2\tau_{11}$, the variance of $b_{0i}+xb_{1i}$. (2) The probit link’s latent standard normal contributes variance $1$; adding the independent random contribution gives total variance $1+q(x)$. (3) $V(1)=1+0.36+0.12+0.16=1.64$; the marginal probability is $\Phi(0.4/\sqrt{1.64})\approx0.6226$. (4) It is $\beta_1+b_{1i}$; a random slope below $-\beta_1$ makes it negative. (5) The limit is $\beta_1/\sqrt{\tau_{11}}$, since the leading terms in numerator and denominator are $\beta_1x$ and $x\sqrt{\tau_{11}}$, respectively.