Logistic Fluctuation for Binary-Outcome TMLE

Why Fluctuate a Binary-Outcome Prediction?

In TMLE, an initial outcome model supplies predicted outcome means, and a fluctuation model updates those predictions in a specified direction. For a binary outcome, the updated predictions must remain between zero and one. A logistic fluctuation achieves this by updating the log odds rather than adding a constant directly to a probability.

This tutorial focuses on fitting that one-dimensional logistic model. It does not by itself define the full treatment-effect estimator: the quantity of interest here is the fitted fluctuation parameter and the corresponding updated predictions. The earlier material on Clever Covariates in Binary-Treatment TMLE defines the direction of the update; below, that direction is placed inside a Bernoulli likelihood and estimated.

Definition: Let \(A_i\in\{0,1\}\) be treatment, \(W_i\) baseline covariates, and \(Y_i\in\{0,1\}\) the observed outcome for observation \(i\). Let \(\hat Q(a,W_i)\) be an initial estimate of the conditional mean outcome under treatment level \(a\), and let \(\hat g(W_i)\) estimate the conditional probability of treatment. For the average treatment effect, define the clever covariate \[ \hat H(A_i,W_i)=\frac{A_i}{\hat g(W_i)}-\frac{1-A_i}{1-\hat g(W_i)}. \] The logistic fluctuation submodel is \[ \operatorname{logit}\{\hat Q_\epsilon(A_i,W_i)\} =\operatorname{logit}\{\hat Q(A_i,W_i)\} +\epsilon\hat H(A_i,W_i), \] where \(\operatorname{logit}(q)=\log\{q/(1-q)\}\), \(\epsilon\) is a scalar parameter, and \(\hat Q_\epsilon\) is the updated prediction.

The initial prediction on the right-hand side is evaluated at the treatment actually observed for that individual. Thus a treated observation uses \(\hat Q(1,W_i)\) and a control observation uses \(\hat Q(0,W_i)\). For a treated observation, \(\hat H=1/\hat g(W_i)\); for a control observation, \(\hat H=-1/\{1-\hat g(W_i)\}\). The same scalar \(\epsilon\) updates all observations, while the direction and size of each update depend on their clever covariate values.

The offset is \(\operatorname{logit}\{\hat Q(A_i,W_i)\}\): it is included in the linear predictor with its coefficient fixed at one. The fitted covariate is \(\hat H(A_i,W_i)\), with no separately estimated intercept. Estimating an intercept or a coefficient for the offset would define a different model, not this fluctuation submodel.

Likelihood and Score Equation

Assume the observed outcomes are modeled as Bernoulli conditional on the observed treatment and covariates, with fitted means \(\hat Q_\epsilon(A_i,W_i)\). The log-likelihood for the scalar fluctuation parameter, up to terms not depending on \(\epsilon\), is

\ell(\epsilon)=\sum_{i=1}^{n}\left[ Y_i\log\{\hat Q_\epsilon(A_i,W_i)\} +(1-Y_i)\log\{1-\hat Q_\epsilon(A_i,W_i)\} \right].

Write \(\hat H_i=\hat H(A_i,W_i)\) and \(\hat Q_{\epsilon,i}=\hat Q_\epsilon(A_i,W_i)\). Differentiating the Bernoulli log-likelihood gives the score

U(\epsilon)=\frac{d\ell(\epsilon)}{d\epsilon} =\sum_{i=1}^{n}\hat H_i\{Y_i-\hat Q_{\epsilon,i}\}.

A finite maximum-likelihood estimate \(\hat\epsilon\) satisfies \(U(\hat\epsilon)=0\). Its second derivative is

\frac{d^2\ell(\epsilon)}{d\epsilon^2} =-\sum_{i=1}^{n}\hat H_i^2\hat Q_{\epsilon,i}(1-\hat Q_{\epsilon,i})\leq 0.

This establishes concavity of the log-likelihood. If at least one \(\hat H_i\ne 0\) and the fitted probabilities are strictly between zero and one, the second derivative is strictly negative at every finite \(\epsilon\). In that case, the log-likelihood is strictly concave, so it has at most one finite maximizer. Concavity does not guarantee that a finite maximizer exists: in a separated data configuration, the likelihood can approach its supremum only as the magnitude of \(\epsilon\) diverges.

A numerical optimizer can find \(\hat\epsilon\) by maximizing \(\ell(\epsilon)\), or by solving \(U(\epsilon)=0\). A Newton update from a current value \(\epsilon^{(t)}\) is

\epsilon^{(t+1)} =\epsilon^{(t)} +\frac{U(\epsilon^{(t)})} {\sum_{i=1}^{n}\hat H_i^2\hat Q_{\epsilon^{(t)},i} (1-\hat Q_{\epsilon^{(t)},i})}.

The denominator is the negative of the second derivative, so the plus sign in this expression is essential. For numerical stability, software typically uses a likelihood optimizer rather than relying on a manually chosen number of Newton steps.

Worked Example: A Finite Fluctuation Estimate

Consider four hypothetical observations. Both treatment groups have \(\hat g(W_i)=0.5\), and every initial prediction at the observed treatment is \(\hat Q(A_i,W_i)=0.5\). The two treated outcomes are 1 and 0; the two control outcomes are 0 and 0. Find the fitted fluctuation parameter.

The offset is \(\operatorname{logit}(0.5)=0\). For treated observations, \(\hat H=1/0.5=2\); for controls, \(\hat H=-1/(1-0.5)=-2\). Therefore the updated means are \(q_T=\operatorname{expit}(2\epsilon)\) for treated observations and \(q_C=\operatorname{expit}(-2\epsilon)\) for controls, where \(\operatorname{expit}(x)=1/(1+\exp(-x))\).

The two treated observations contribute \(2(1-q_T)+2(0-q_T)=2-4q_T\) to the score. The two controls contribute \((-2)(0-q_C)+(-2)(0-q_C)=4q_C\). Setting the total score to zero gives

2-4q_T+4q_C=0 \quad\Longleftrightarrow\quad q_T-q_C=0.5.

Because \(\operatorname{expit}(x)+\operatorname{expit}(-x)=1\), here \(q_T+q_C=1\). Solving the two equations gives \(q_T=0.75\) and \(q_C=0.25\). Thus \(2\hat\epsilon=\operatorname{logit}(0.75)=\log(3)\), and

\hat\epsilon=\frac{\log(3)}{2}=0.5493.

The finite estimate is consistent with the score calculation: the treated predictions increase from 0.5 to 0.75, while the control predictions decrease from 0.5 to 0.25. The observations’ outcomes determine the likelihood fit; the updated values are model predictions, not replacements for the observed outcomes.

Numerical Fitting and a Score Check

With different offsets or clever covariates, the score equation will generally not have a simple closed-form solution. The next example shows how to check a numerical fit using the score and curvature. All values are hypothetical.

Worked Example: One Newton Update and the Fitted Value

Suppose there are three treated and three control observations, with \(\hat g(W_i)=0.5\) and initial predictions of 0.5 throughout. The treated outcomes are 1, 1, 0, and the control outcomes are 0, 0, 1. Starting at \(\epsilon^{(0)}=0\), calculate a Newton update and the fitted value.

As before, the clever covariates are \(2\) for treated observations and \(-2\) for controls. At \(\epsilon=0\), all updated means equal 0.5. The score is \[ U(0)=2\{(1-0.5)+(1-0.5)+(0-0.5)\} -2\{(0-0.5)+(0-0.5)+(1-0.5)\}=2. \] The denominator in the Newton update is \(6(2^2)(0.5)(0.5)=6\), so \[ \epsilon^{(1)}=0+\frac{2}{6}=0.3333. \]

The root can also be found exactly for this symmetric setup. At any \(\epsilon\), let \(q_T=\operatorname{expit}(2\epsilon)\) and \(q_C=\operatorname{expit}(-2\epsilon)=1-q_T\). The score reduces to \(U(\epsilon)=8-12q_T\). Setting it to zero gives \(q_T=2/3\), and hence \[ \hat\epsilon=\frac{\operatorname{logit}(2/3)}{2} =\frac{\log(2)}{2}=0.3466. \]

At the first Newton update, \(q_T=\operatorname{expit}(2/3)\approx0.6608\), so \(U(0.3333)=8-12(0.6608)\approx0.0709\), still slightly positive. The curvature denominator there is \(24q_T(1-q_T)\approx5.380\), giving a further update of about \(0.0709/5.380=0.0132\). The resulting value, about \(0.3465\), is close to the root \(0.3466\). A final fit should be checked by evaluating the score at the optimizer’s reported estimate, not by assuming that one update is sufficient.

When a Finite Estimate Does Not Exist

Strict concavity means there cannot be two distinct finite maximizers, but it does not rule out separation. In a separated configuration, increasing or decreasing \(\epsilon\) can keep improving the likelihood, while the fitted means approach zero or one. A numerical routine may then report a very large estimate, a convergence warning, or a score that is small only at an extreme value.

Worked Example: Separation in the Negative Direction

Take one hypothetical treated observation and one hypothetical control observation, each with initial prediction 0.5 and treatment probability estimate 0.5. Suppose the treated outcome is 0 and the control outcome is 1. Then \(\hat H_T=2\), \(\hat H_C=-2\), and the updated means are \(q_T=\operatorname{expit}(2\epsilon)\) and \(q_C=\operatorname{expit}(-2\epsilon)\).

The log-likelihood is \[ \ell(\epsilon)=\log(1-q_T)+\log(q_C). \] As \(\epsilon\to-\infty\), \(q_T\to0\) and \(q_C\to1\), so \(\ell(\epsilon)\to0\), its supremum. For any finite \(\epsilon\), both probabilities are strictly between zero and one, and each log term is negative; therefore the supremum is not attained at a finite parameter. In this example the likelihood improves as \(\epsilon\) tends to negative infinity. Reversing the two outcomes makes the likelihood improve as \(\epsilon\) tends to positive infinity. The direction depends on the outcome pattern and the signs of the clever covariates.

A practical fit should be inspected for convergence warnings, unusually large \(|\hat\epsilon|\), fitted probabilities extremely close to zero or one, and a score that remains appreciably different from zero. Also inspect the range of \(\hat g(W_i)\): probabilities close to zero or one create large clever covariates, which can make the likelihood highly sensitive to a small number of observations. These issues can signal weak overlap, separation, or numerical instability; they should not be hidden by reporting a converged-looking number without a diagnostic check.

What the Fluctuation Is—and Is Not

The logistic fluctuation is a logistic regression with a particular structure: the initial logit predictions are an offset, the clever covariate is the only fitted predictor, and there is no free intercept. In software, the specification is equivalent to fitting a binomial model with the offset \(\operatorname{logit}\{\hat Q(A_i,W_i)\}\) and predictor \(\hat H_i\), constrained to have no additional intercept. The estimated slope is \(\hat\epsilon\).

This differs from fitting a new, unrestricted logistic outcome model. An unrestricted regression would estimate its own intercept and potentially other coefficients; it would not generally produce the targeted update above. It also differs from leaving the initial predictions unchanged, which is equivalent to setting \(\epsilon=0\) whether or not that value maximizes the fluctuation likelihood. The fluctuation uses the observed outcomes to estimate a specific one-parameter correction.

The update is bounded because it applies the expit function to a real-valued updated logit. However, bounded predictions do not guarantee that the treatment-specific comparisons are stable: extreme clever covariates can still dominate the score. Nor does fitting this one-dimensional likelihood alone establish that the initial outcome model, treatment model, or causal assumptions are adequate. The update is a component of TMLE, not a substitute for checking the analysis assumptions.

Common pitfall: It is incorrect to treat the clever covariate as an ordinary outcome-model predictor with a freely estimated intercept, or to say that a large fitted \(\hat\epsilon\) must be positive under separation. The correction is to keep the offset fixed, fit only the coefficient of \(\hat H\), and inspect the likelihood and score: depending on the outcome pattern, the likelihood supremum may require \(\epsilon\to+\infty\) or \(\epsilon\to-\infty\).

For a reliable implementation, confirm that the initial means supplied to the logit are strictly between zero and one, that the clever covariate uses the intended treatment probability estimates, and that the fitted parameter approximately solves \(\sum_i\hat H_i\{Y_i-\hat Q_{\hat\epsilon,i}\}=0\). If predictions must be bounded away from zero and one for computation, any such adjustment should be explicit and consistent with the analysis; it changes the supplied offset values.

Check Your Understanding

Try each question, then compare your reasoning with the solutions.

For a control observation with \(\hat g(W_i)=0.8\), what is its clever covariate?
In the logistic fluctuation, which quantity is the offset, and which coefficient is estimated?
At a finite maximum-likelihood estimate, what equation should the fitted score satisfy?
Why does a large-magnitude fluctuation estimate require investigation rather than automatic acceptance as a stable fit?

Solutions. (1) The control clever covariate is \(-1/(1-0.8)=-5\). (2) The offset is \(\operatorname{logit}\{\hat Q(A_i,W_i)\}\), held fixed with coefficient one; the fitted coefficient is the scalar \(\epsilon\) multiplying \(\hat H(A_i,W_i)\), with no additional intercept. (3) It should satisfy \(\sum_i\hat H_i\{Y_i-\hat Q_{\hat\epsilon,i}\}=0\), up to numerical tolerance. (4) A large magnitude may indicate separation, extreme clever covariates caused by treatment probabilities near zero or one, or numerical instability. Check convergence, fitted probabilities, and the score rather than interpreting the estimate in isolation.