Why Nonlinear Outcomes Break Ordinary Two-Stage Substitution

Why Substitution Is Tempting—and Why It Fails

In linear instrumental-variables analysis, a familiar procedure uses the first stage to predict treatment and then uses the prediction in a second-stage calculation. It is tempting to carry the same idea over to a probit or logit outcome model: predict treatment, insert that prediction into the nonlinear outcome equation, and interpret the resulting coefficient as an IV effect. The difficulty is not merely that the outcome is binary. The key issue is that a nonlinear function of an average is generally not the average of the nonlinear function.

For a binary treatment, the first-stage fitted value is often a predicted treatment probability. A person with fitted value \(0.2\), for example, did not receive “0.2 units” of a binary treatment. Rather, \(0.2\) represents a conditional average or probability under the first-stage model. Substituting that value into a nonlinear outcome response function does not usually reproduce the average outcome probability across people with treatment values zero and one.

To isolate the mathematical problem, let \(Z\) be a binary instrument and \(D\) a binary treatment. Suppose, for the algebra in the next example, that the conditional response restriction is \(P(Y=1\mid D,Z)=g(\alpha+\beta D)\), where \(g\) is a nonlinear response function such as the standard normal cumulative distribution function \(\Phi\) or the logistic function. This restriction explicitly says that, conditional on treatment, the response probability does not additionally depend on \(Z\). It is a simplifying condition for demonstrating the substitution problem; it is not a complete set of assumptions for causal identification.

The Mixture Versus the Plug-In Probability

Write \(p_z=P(D=1\mid Z=z)\). Because \(D\) is binary, conditioning on its two possible values gives the actual outcome probability at instrument value \(z\):

Ordinary predictor substitution instead evaluates the response function at the fitted treatment probability: \(g(\alpha+\beta p_z)\). These are different expressions. The first averages two response probabilities, weighting them by the probabilities of the two treatment states. The second treats the average treatment value as if it were itself a treatment level to be passed through the nonlinear response function.

This is an instance of the general expectation identity that fails for nonlinear functions: \(E[g(V)\mid Z]\ne g(E[V\mid Z])\) in general. Jensen’s inequality can sometimes describe the direction of the difference if the function is convex or concave over the entire relevant range. But probit and logit response functions are not globally convex or globally concave as functions of their index, so there is no universal direction for the substitution error.

Worked Binary-Outcome Examples

A Short Derivation for Continuous Treatment

The same issue is not confined to binary treatment. A simple probit calculation makes the averaging operation explicit. Suppose, as a hypothetical distributional setup, that \(D\mid Z\sim N(\mu,\sigma^2)\) and \(Y=1\{\alpha+\beta D+\varepsilon>0\}\), where \(\varepsilon\sim N(0,1)\) is independent of \(D\) conditional on \(Z\). Then \(P(Y=1\mid D,Z)=\Phi(\alpha+\beta D)\). To obtain the probability conditional only on \(Z\), write \(D=\mu+\sigma U\), where \(U\sim N(0,1)\) is independent of \(\varepsilon\). The latent index becomes \(\alpha+\beta\mu+\beta\sigma U+\varepsilon\), whose random part is normal with variance \(1+\beta^2\sigma^2\). Thus:

Plugging the conditional mean \(\mu\) into the response function gives \(\Phi(\alpha+\beta\mu)\), which omits the variance term in the denominator. Unless that omission has no effect in a special case, the two probabilities differ. The derivation depends on the stated normality and independence conditions; it is an illustration of averaging, not a general correction formula for IV estimation.

Why Linear 2SLS Is Different

In a linear conditional-mean model, the response function is affine: \(g(v)=v\), or more generally a linear expression in the predictors. For any random variable \(V\), \(E[V\mid Z]\) is exactly the conditional mean. Consequently, taking a conditional expectation of a linear expression replaces the random treatment by its conditional mean without the nonlinear averaging discrepancy shown above.

This algebraic fact does not mean that every informal “run a regression twice” procedure has correct standard errors or causal interpretation. Linear 2SLS has a specific estimator and inference theory, and it uses instruments to construct variation in treatment. In a nonlinear second stage, simply inserting \(\hat D\) into a probit or logit likelihood does not inherit those linear properties. It may yield a numerical coefficient, but that coefficient generally does not estimate the structural treatment coefficient in the model for \(Y\mid D\).

The distinction is also important for interpretation. A coefficient in a probit or logit model is a change on the latent index or log-odds scale, respectively; it is not itself a change in outcome probability. The nonlinear probability scale depends on the covariates and the treatment value. A plug-in procedure introduces the additional problem that its fitted probabilities need not equal the relevant averages over treatment values.

Assumptions, Diagnostics, and Failure Modes

The mixture identity requires a response model that specifies \(P(Y=1\mid D,Z)\), not merely \(P(Y=1\mid D)\). If the instrument has a direct association with the outcome conditional on treatment, the two treatment-specific probabilities may depend on \(Z\), and the simple expression using \(g(\alpha)\) and \(g(\alpha+\beta)\) does not follow. In a causal IV analysis, relevance, independence, and exclusion are distinct design conditions; the algebraic substitution problem neither establishes nor replaces them. A correctly specified outcome response function and an identified causal estimand are additional concerns.

Several practical checks can reveal that a nonlinear substitution analysis is on shaky ground:

• Inspect the treatment scale. With binary \(D\), fitted values between zero and one are probabilities or conditional means, not observed treatment states.
• Compare mixture and plug-in predictions. If the response model and treatment distribution are specified, calculate both quantities at representative covariate and instrument values. A difference directly displays the nonlinear averaging gap.
• Check conditional response assumptions. A formula that averages treatment-specific probabilities must use probabilities conditional on the relevant instrument and covariates, or justify why those probabilities are invariant to them.
• Separate model fit from identification. Good prediction of \(Y\) does not by itself show that the treatment coefficient has a causal IV interpretation.
• Do not treat standard errors as the only issue. Correcting a variance formula cannot make an estimator consistent for a target it does not generally estimate.

A common mistake is to reason, “The first stage estimates the expected treatment, so putting that expectation into the outcome model gives the expected outcome.” That step silently assumes \(E[g(V)\mid Z]=g(E[V\mid Z])\). The correction is to average the response probabilities over the conditional distribution of treatment, or to use an estimation strategy justified for the nonlinear structural model and the causal estimand of interest.

How to Choose an Approach

The choice should begin with the estimand and assumptions, not with whichever model is easiest to fit. Nonlinear IV methods differ in what effect they identify and in what functional-form or distributional restrictions they require. The central lesson here is narrower but essential: an ordinary nonlinear regression on first-stage fitted treatment is not justified simply by analogy with linear 2SLS.