Quadrature for Logistic-Normal Selection Models

Why Quadrature Appears in Logistic-Normal Selection Models

When dropout depends on the current outcome, the probability of dropout depends on a value that is not observed for people who drop out. If the outcome is modeled as normal and the dropout probability as logistic in that outcome, the dropout probability averaged over the outcome distribution is a logistic-normal integral. In general, this integral has no elementary closed form. Gaussian quadrature provides a deterministic numerical approximation.

The task here is to calculate the relevant integrals accurately and transparently. Quadrature evaluates the integrand at a small set of carefully chosen outcome values and combines those evaluations with weights. It does not replace the assumed outcome or dropout model, and it cannot make an unidentified or poorly specified model reliable.

Definition: Let $Y$ be a current outcome, with $Y\mid x\sim N(\mu,\sigma^2)$, where $x$ denotes the observed history and covariates, and $\sigma>0$. Let $D=1$ indicate dropout before the current outcome is observed, and $D=0$ indicate that it is observed. A logistic current-outcome selection equation is $$ P(D=1\mid Y=y,x)=\operatorname{expit}(\eta(y,x)),\qquad \eta(y,x)=\gamma_0+\gamma_y y+\gamma^\mathsf{T}x, $$ where $\operatorname{expit}(a)=1/(1+\exp(-a))$. The model specifies the conditional outcome density and the conditional dropout probability; their product gives the joint contribution when $y$ is observed.

This specification assumes the conditional normal distribution and logistic dropout equation are appropriate for the current occasion and the given history. It also assumes the recorded dropout indicator is defined consistently with the timing of that occasion. For a person with $D=1$, $Y$ is unobserved, so the observed-data contribution must average over its modeled distribution. For a person with $D=0$, the observed value can be substituted directly into the outcome density and retention probability.

From the Normal Outcome to a Standard-Normal Integral

Write the standardized outcome as $Z=(Y-\mu)/\sigma$. Under the normal outcome model, $Z\mid x\sim N(0,1)$ and $Y=\mu+\sigma Z$. The marginal dropout probability conditional on $x$ is therefore

q(x)=P(D=1\mid x) =\int_{-\infty}^{\infty} \operatorname{expit}(\eta(y,x))f(y\mid x)\,dy =E\!\left[\operatorname{expit}\{\eta(\mu+\sigma Z,x)\}\mid x\right].

For a dropout whose current outcome is missing, $q(x)$ is the observed-data contribution for the dropout indicator, conditional on the modeled history. More generally, the same change of variable applies to an integral of any function $h(y)$ against the normal outcome density.

Definition: The normalized $m$-node Gauss–Hermite approximation for $E\{h(\mu+\sigma Z)\}$, where $Z\sim N(0,1)$, is $$ E\{h(\mu+\sigma Z)\}\ \approx\ \sum_{j=1}^{m} w_j h(\mu+\sigma z_j), \qquad \sum_{j=1}^{m}w_j=1. $$ The $z_j$ and $w_j$ below are nodes and weights for a standard normal variable. They are already normalized for an expectation; do not multiply by another normalizing constant.

The standard-normal nodes are not the nodes for a normal variable with mean $\mu$ and standard deviation $\sigma$. Each standard-normal node $z_j$ must be transformed to the outcome scale as $y_j=\mu+\sigma z_j$. The weights remain unchanged. For the logistic-normal dropout probability, set $h(y)=\operatorname{expit}\{\eta(y,x)\}$.

Rule	Standard-normal nodes	Corresponding weights
Three nodes	$-\sqrt{3},\ 0,\ \sqrt{3}$	$1/6,\ 2/3,\ 1/6$
Five nodes	$-\sqrt{5+\sqrt{10}},\ -\sqrt{5-\sqrt{10}},\ 0,\ \sqrt{5-\sqrt{10}},\ \sqrt{5+\sqrt{10}}$	$(7-2\sqrt{10})/60,\ (7+2\sqrt{10})/60,\ 8/15,\ (7+2\sqrt{10})/60,\ (7-2\sqrt{10})/60$

The five-node values are approximately $-2.85697,\ -1.35563,\ 0,\ 1.35563,\ 2.85697$, with corresponding weights approximately $0.0112574,\ 0.222076,\ 0.533333,\ 0.222076,\ 0.0112574$. Gauss–Hermite rules integrate polynomials of degree at most $2m-1$ exactly, but a logistic function is not a polynomial. Thus, even though the rule is well defined, its approximation to this integral must be checked when accuracy matters.

Worked Calculations

Worked Example: A Three-Node Dropout Probability

Suppose $Z\sim N(0,1)$ and the dropout probability at standardized outcome $z$ is $\operatorname{expit}(-0.5+1.2z)$. Approximate the marginal dropout probability using the three-node rule.

The nodes are $-\sqrt{3},0,\sqrt{3}$, and the weights are $1/6,2/3,1/6$. The two outer linear predictors are $-0.5-1.2\sqrt{3}\approx-2.57846$ and $-0.5+1.2\sqrt{3}\approx1.57846$. Their dropout probabilities are approximately $0.0705376$ and $0.8289865$; at the center the probability is $\operatorname{expit}(-0.5)\approx0.3775407$.

Substitution into the weighted sum gives

\begin{aligned} q &\approx \frac{1}{6}(0.0705376)+\frac{2}{3}(0.3775407)+\frac{1}{6}(0.8289865)\\ &\approx 0.0117563+0.2516938+0.1381644\\ &\approx 0.401614. \end{aligned}

The estimated dropout probability is about $0.4016$, or 40.16%. This is a quadrature approximation under the stated model, not a probability calculated from observed dropout frequencies.

Worked Example: The Contribution When the Outcome Is Observed

Suppose $Y\mid x\sim N(10,2^2)$, $D=0$, and the observed current outcome is $y=12$. Let the dropout equation be $\operatorname{logit}P(D=1\mid Y=y,x)=-1+0.2y$, with no additional covariate term. Find the joint contribution for the observed outcome and retention indicator.

Because $y$ is observed, no integral over the current outcome is needed. The outcome density is

f(12\mid x)=\frac{1}{2\sqrt{2\pi}} \exp\!\left\{-\frac{(12-10)^2}{2(2^2)}\right\} \approx 0.120985.

At $y=12$, the dropout predictor is $-1+0.2(12)=1.4$, so the dropout probability is $\operatorname{expit}(1.4)\approx0.802184$. The retention probability is $1-0.802184=0.197816$. Therefore the joint density contribution is

f(12\mid x)P(D=0\mid Y=12,x) \approx 0.120985(0.197816) \approx 0.023933.

This is a density contribution per unit of $Y$, not a standalone probability for the exact continuous value $Y=12$. The example also shows why quadrature is needed for a missing current outcome but not for this observed one.

Worked Example: Checking the Five-Node Approximation

Return to the dropout probability $\operatorname{expit}(-0.5+1.2z)$ for $Z\sim N(0,1)$. Apply the five-node rule and compare it with the three-node estimate.

At the center, the weighted contribution is $(8/15)\operatorname{expit}(-0.5)\approx0.2013550$. At the two inner nodes $z=\pm1.355626$, the predictors are approximately $-2.12675$ and $1.12675$, giving dropout probabilities approximately $0.1065$ and $0.7550$. Together, their weighted contribution is approximately $0.1913767$. At the outer nodes $z=\pm2.856970$, the predictors are approximately $-3.92836$ and $2.92836$, giving probabilities approximately $0.0193$ and $0.9493$; their combined weighted contribution is approximately $0.0109031$.

Adding the contributions using the five-node weights gives

q\approx 0.2013550+0.1913767+0.0109031 =0.4036348\approx0.403635.

The five-node estimate is about $0.4036$, compared with $0.4016$ from the three-node rule, a difference of approximately $0.00202$. The five-node result is a more refined approximation, but the difference warns that the three-node result may not be adequate for a calculation requiring tighter precision. Agreement between higher orders should be checked before deciding that numerical error is negligible.

Observed-Data Contributions and a Short Justification

The same model gives distinct contribution forms depending on whether the current outcome is observed. Ignoring any factors from earlier observed history that do not involve the current outcome or dropout indicator, the contribution for $D=0$ and observed $Y=y$ is the normal density times the conditional retention probability. For $D=1$ and missing $Y$, integrate out $Y$ from the joint model.

\begin{aligned} L_{\mathrm{obs}}(y,D=0\mid x) &=f(y\mid x)\{1-\operatorname{expit}(\eta(y,x))\},\\ L_{\mathrm{obs}}(D=1\mid x) &=\int_{-\infty}^{\infty} f(y\mid x)\operatorname{expit}(\eta(y,x))\,dy\\ &=E\!\left[\operatorname{expit}\{\eta(\mu+\sigma Z,x)\}\right]\\ &\approx \sum_{j=1}^{m}w_j\operatorname{expit}\{\eta(\mu+\sigma z_j,x)\}. \end{aligned}

The second equality follows directly from the substitution $y=\mu+\sigma z$: the normal density and differential transform together to the standard-normal density, so integrating the dropout probability against $f(y\mid x)$ is the expectation of that probability at $Y=\mu+\sigma Z$. Replacing this expectation with the weighted sum gives the quadrature approximation. This derivation also explains why the node transformation uses $\mu+\sigma z_j$, rather than $\mu+z_j$ or $\mu+\sigma^2 z_j$.

Specify the conditional models.
Record the normal outcome mean $\mu$, standard deviation $\sigma$, and dropout linear predictor $\eta(y,x)$ for the person’s observed history and covariates.

Identify what is observed.
For an observed current outcome with $D=0$, evaluate the normal density and retention probability at that outcome. For a dropout with missing $Y$, calculate the integral over $Y$.

Transform and evaluate.
For the missing-outcome integral, calculate each $y_j=\mu+\sigma z_j$, evaluate the logistic probability there, multiply by $w_j$, and sum.

Check numerical stability and accuracy.
Repeat with more nodes when needed, use stable evaluation of the logistic function for extreme predictors, and retain sufficient precision until the final reported result.

Diagnostics, Alternatives, and Limitations

A quadrature result should be treated as a numerical approximation, not as automatically exact because the outcome distribution is normal. Useful checks include increasing the node count and comparing the resulting integral at the precision needed for the analysis. Inspect the predictors $\eta(\mu+\sigma z_j,x)$ as well: a large $|\gamma_y|\sigma$ can make the logistic probability change sharply across the outcome distribution, so a low-order rule may miss important curvature. Calculations can also go wrong through incorrect node scaling, mismatched weights, omitted weights, or premature rounding.

For extreme predictors, evaluating $1/(1+\exp(-a))$ naively can overflow numerically. Use a stable logistic evaluation: for $a\geq0$, compute $1/(1+\exp(-a))$; for $a<0$, compute $\exp(a)/(1+\exp(a))$. If individual likelihood factors are combined computationally, log-scale calculations can further reduce underflow. These are numerical safeguards; they do not change the model.

Quadrature is a deterministic alternative to approximating the expectation by random draws. It is often efficient for a one-dimensional normal integral. Other numerical approaches may be preferable for difficult or higher-dimensional integrals; Monte Carlo integration is one such alternative, treated separately. A closed-form calculation is available for some different link-and-distribution combinations, but changing the link to obtain a convenient formula changes the selection model and is not a numerical shortcut for the logistic specification.

Finally, accurate integration does not establish that the normal outcome model or logistic selection equation is appropriate. Nor does it resolve whether the outcome-dependent dropout effect is sufficiently informed by observed data. Quadrature addresses the computation conditional on specified model parameters; it cannot remove uncertainty about assumptions or identify information absent from the observed data.

Common pitfall: Using standard-normal nodes directly as outcome values, or using $y_j=\mu+\sigma^2z_j$, gives the wrong integral. The correct transformation is $y_j=\mu+\sigma z_j$, because $\sigma$ is the standard deviation. Also, a close match between two low-order approximations is a useful diagnostic, not proof that the model is correct.

Check Your Understanding

Answer each question, then compare with the solutions below.

For $Y\mid x\sim N(8,3^2)$, what outcome-scale values correspond to standard-normal quadrature nodes $z=-\sqrt{3},0,\sqrt{3}$?
Why does a dropout with missing current outcome require an integral, while an observed outcome with $D=0$ does not?
In a three-node standard-normal rule, what are the weights, and why should they not be divided by $\sqrt{\pi}$?
A three-node estimate is $0.41$ and a five-node estimate is $0.39$. What is an appropriate next numerical check, and what does the discrepancy not tell you?

Solutions. (1) Use $y_j=\mu+\sigma z_j$: the values are $8-3\sqrt{3}\approx2.804$, $8$, and $8+3\sqrt{3}\approx13.196$. (2) For dropout, $Y$ is unobserved, so its modeled density must be integrated against the dropout probability. With $D=0$ and observed $y$, evaluate the density and retention probability at that value. (3) The weights are $1/6,2/3,1/6$, which already sum to 1 and are normalized for a standard-normal expectation. Dividing by $\sqrt{\pi}$ would apply an additional normalization incorrectly. (4) Increase the quadrature order and compare results until numerical differences are small enough for the required precision. The discrepancy indicates possible quadrature sensitivity; it does not establish whether the outcome or dropout model is correct.

Rule	Standard-normal nodes	Corresponding weights
Three nodes	\(-\sqrt{3},\ 0,\ \sqrt{3}\)	\(1/6,\ 2/3,\ 1/6\)
Five nodes	\(-\sqrt{5+\sqrt{10}},\ -\sqrt{5-\sqrt{10}},\ 0,\ \sqrt{5-\sqrt{10}},\ \sqrt{5+\sqrt{10}}\)	\((7-2\sqrt{10})/60,\ (7+2\sqrt{10})/60,\ 8/15,\ (7+2\sqrt{10})/60,\ (7-2\sqrt{10})/60\)