Firth Logistic Models for Treatment Probabilities

Why Use Firth Logistic Regression for Treatment Probabilities?

In a treatment model, logistic regression estimates the probability of treatment given measured covariates or covariate history. With sparse treatment data, a covariate pattern may occur only among treated individuals, or only among untreated individuals. Ordinary maximum-likelihood logistic regression can then send one or more coefficients toward infinity: the likelihood improves as predicted probabilities approach zero or one, but no finite coefficient achieves its supremum. This is separation, a problem introduced in the earlier material on sparse treatment data.

Firth logistic regression changes the estimation criterion rather than imposing an arbitrary upper or lower bound on fitted probabilities. Its penalty is based on the information matrix and produces an adjusted score that counteracts first-order small-sample bias. In inverse probability weighting, this can yield usable finite fitted probabilities in separated samples. It does not, by itself, establish positivity in the population or guarantee stable weights.

Definition and Adjusted Score

Let \(A_i\in\{0,1\}\) be the treatment received by observation \(i\), and let \(x_i\) be a \(p\)-dimensional column vector of model covariates, including an intercept if one is fitted. Assume independent Bernoulli observations conditional on the fixed design vectors, with \(\operatorname{logit}\{P(A_i=1\mid x_i)\}=x_i^\top\beta\). Write \(p_i=\operatorname{expit}(x_i^\top\beta)\), \(w_i=p_i(1-p_i)\), and let \(X\) be the matrix with rows \(x_i^\top\). For finite \(\beta\), the Fisher information is \(I(\beta)=X^\top W X\), where \(W\) is diagonal with entries \(w_i\).

The full-column-rank condition on \(X\) ensures that \(I(\beta)\) is positive definite for every finite \(\beta\), since every \(w_i>0\). If the columns are linearly dependent, the determinant is zero and this penalty is not defined as written; the design must first be reparameterized or otherwise made identifiable.

To derive this, the ordinary score is \(U(\beta)=X^\top(A-p)\). For any parameter component \(\beta_j\), differentiation of the penalty uses \(\partial\log|I|/\partial\beta_j=\operatorname{tr}\{I^{-1}(\partial I/\partial\beta_j)\}\). Because \(\partial w_i/\partial\beta_j=w_i(1-2p_i)x_{ij}\), we have \(\partial I/\partial\beta_j=\sum_i w_i(1-2p_i)x_{ij}x_ix_i^\top\). Consequently,

where \(h_i=w_i x_i^\top I^{-1}x_i\), exactly the \(i\)th diagonal element of \(H\). Adding this derivative to the ordinary score gives the stated result. Thus Firth estimation can be computed by solving an ordinary logistic score plus a leverage-weighted correction; the correction is largest for observations with high leverage. \(\square\)

Worked Calculations

Proof. In the group-logit parameterization, the information matrix is diagonal with entries \(n_jp_j(1-p_j)\). The ordinary log-likelihood is a sum of independent group contributions, and the penalty is \(\tfrac12\sum_j\log\{n_jp_j(1-p_j)\}\). Terms involving \(n_j\) are constant in the probabilities. Thus the criterion for group \(j\) is \(s_j\log(p_j)+(n_j-s_j)\log(1-p_j)+\tfrac12\log(p_j(1-p_j))\). Its derivative is \((s_j+1/2)/p_j-(n_j-s_j+1/2)/(1-p_j)\). Setting it to zero gives \((s_j+1/2)(1-p_j)=(n_j-s_j+1/2)p_j\), hence \(p_j=(s_j+1/2)/(n_j+1)\). The second derivative is \(-(s_j+1/2)/p_j^2-(n_j-s_j+1/2)/(1-p_j)^2<0\); the solution is the unique maximum. Since \(0<s_j+1/2<n_j+1\), it is strictly between zero and one. \(\square\)

Using the Estimates in Inverse Probability Weights

For a single treatment decision, let \(\hat p_i\) be the fitted probability of treatment given the covariates in the denominator model. The unstabilized weight is \(1/\hat p_i\) for a treated observation and \(1/(1-\hat p_i)\) for an untreated observation. In longitudinal settings, fit the appropriate sequential treatment model at each time and substitute its Firth-fitted probabilities into the corresponding factors of the longitudinal weight. For example, the denominator contribution at time \(t\) is

A stabilized weight uses its specified numerator probability as well; changing how the denominator is estimated does not change that definition. Record which model produced each fitted probability, and compute the weight from the probability corresponding to the observed treatment, not automatically from \(\hat p_{it}\).

Firth estimates are not obtained by simply adding one-half to each observed treatment count in every arbitrary logistic regression. The simple \((s+1/2)/(n+1)\) calculation is exact for a saturated group model. For a general covariate design, the adjustment depends on the fitted probabilities and leverage values through the adjusted score.

Diagnostics, Alternatives, and Limitations

A finite coefficient is a useful computational outcome, not a complete diagnostic. Before using fitted probabilities for weights, check that the design matrix has full column rank, the optimizer converged, and the fitted probabilities correspond to the intended treatment and covariate history. Inspect predicted probabilities by treatment group and the resulting weight distribution. Large weights can remain even when all estimated probabilities are strictly between zero and one. A finite estimate does not create empirical support for treatment assignments absent from important covariate regions.

Ordinary maximum likelihood is the direct alternative when the model is estimable and sparse-data bias or separation is not a concern. Under separation it may diverge or fail to converge. Firth’s penalty specifically targets first-order bias reduction and is not equivalent to clipping probabilities at a chosen threshold: clipping leaves the fitted model unchanged and alters the weights afterward. Exact logistic methods are another possible small-sample approach, but they answer a different computational problem and can be costly. No estimation method repairs an incorrect treatment model, unmeasured confounding, or a violation of the causal positivity assumption.

• Rank deficiency: remove redundant columns or otherwise specify an identifiable model before applying the determinant penalty.
• Extreme weights: investigate the fitted probabilities and covariate support; do not assume finite Firth estimates imply harmless weights.
• Numerical warnings: check scaling, coding, convergence tolerances, and whether some fitted probabilities are so close to zero or one that calculations lose precision.
• Interpretation: Firth changes the finite-sample fitting criterion. It does not redefine the target causal estimand, and it does not guarantee improved causal inference if the model or assumptions are inadequate.