11.3 Inference in logistic regression

The exact distribution of \(\hat {\boldsymbol \beta }\) is intractable, and inference rests on large-sample approximations. The matrix \(\mathcal {I}(\boldsymbol \beta ) = \mathbf {X}^\top \mathbf {WX}\), minus the second-derivative matrix, is the information matrix.

Theorem 11.3.1 (Large-sample distribution). For large samples, approximately, \[ \hat {\boldsymbol \beta } \sim N\bigl (\boldsymbol \beta , (\mathbf {X}^\top \mathbf {WX})^{-1}\bigr ) , \] and the variance is estimated by \((\mathbf {X}^\top \hat {\mathbf {W}}\mathbf {X})^{-1}\), with \(\hat {\mathbf {W}}\) evaluated at \(\hat {\boldsymbol \beta }\).

Derivation. Expand the score about the true \(\boldsymbol \beta \): \(\mathbf {0} = \mathbf {U}(\hat {\boldsymbol \beta }) \approx \mathbf {U}(\boldsymbol \beta ) - \mathcal {I}(\hat {\boldsymbol \beta } - \boldsymbol \beta )\), so \(\hat {\boldsymbol \beta } - \boldsymbol \beta \approx \mathcal {I}^{-1}\mathbf {U}(\boldsymbol \beta )\). Now \(\mathbf {U}(\boldsymbol \beta ) = \mathbf {X}^\top (\mathbf {y} - \boldsymbol \mu )\) has mean \(\mathbf {0}\) and, since the \(y_i\) are independent with variances \(w_i\), variance \(\mathbf {X}^\top \mathbf {WX} = \mathcal {I}\). It is a sum of many independent terms, so by the central limit theorem it is approximately normal. Then \(\mathcal {I}^{-1}\mathbf {U}\) is approximately normal with mean \(\mathbf {0}\) and variance \(\mathcal {I}^{-1}\mathcal {I}\mathcal {I}^{-1} = \mathcal {I}^{-1}\). □

So the Wald test of \(H_0\colon \beta _j = 0\) uses \(Z_0 = \hat \beta _j/\operatorname {se}(\hat \beta _j)\), compared with the standard normal, and an approximate \(100(1 - \alpha )\%\) confidence interval for \(\beta _j\) is \(\hat \beta _j \pm z_{\alpha /2}\operatorname {se}(\hat \beta _j)\). Exponentiating its end points gives an interval for the odds ratio \(e^{\beta _j}\).

Deviance and the likelihood ratio test

The largest log-likelihood any model can reach is that of the saturated model, which has a separate parameter for every group and so fits \(\hat \pi _i = y_i/m_i\) exactly. The deviance of a model measures how far it falls short: \[ D = 2\bigl [\ell (\text {saturated}) - \ell (\hat {\boldsymbol \beta })\bigr ] = 2\sum _{i=1}^n\Bigl [y_i\log \frac {y_i}{\hat \mu _i} + (m_i - y_i)\log \frac {m_i - y_i}{m_i - \hat \mu _i}\Bigr ] , \] the second form by substituting \(\pi _i = y_i/m_i\) and \(\pi _i = \hat \mu _i/m_i\) in the first expression for \(\ell \); a term with \(y_i = 0\) or \(y_i = m_i\) is taken as \(0\). The deviance plays the part of the residual sum of squares: it is never negative, and it is zero only for a perfect fit.

Theorem 11.3.2 (Likelihood ratio test). Let model \(M_0\) be obtained from model \(M_1\) by setting \(q\) of its coefficients to zero. If \(M_0\) is correct, then for large samples, approximately, \[ D(M_0) - D(M_1) = 2\bigl [\ell (\hat {\boldsymbol \beta }_1) - \ell (\hat {\boldsymbol \beta }_0)\bigr ] \sim \chi ^2_q . \]

Derivation. The saturated log-likelihood cancels in the difference, which gives the equality. For the distribution, consider first the case where \(M_0\) fixes all \(p\) coefficients at their true value \(\boldsymbol \beta \). Expanding \(\ell \) about \(\hat {\boldsymbol \beta }\), where the gradient is zero, \[ \ell (\boldsymbol \beta ) \approx \ell (\hat {\boldsymbol \beta }) - \tfrac 12(\hat {\boldsymbol \beta } - \boldsymbol \beta )^\top \mathcal {I}(\hat {\boldsymbol \beta } - \boldsymbol \beta ) . \] So \(2[\ell (\hat {\boldsymbol \beta }) - \ell (\boldsymbol \beta )] \approx \mathbf {v}^\top \mathbf {v}\) with \(\mathbf {v} = \mathcal {I}^{1/2}(\hat {\boldsymbol \beta } - \boldsymbol \beta )\), which by Theorem 11.3.1 is approximately \(N(\mathbf {0}, \mathbf {I}_p)\). A sum of \(p\) squared independent standard normals is \(\chi ^2_p\) (Lemma 3.7.1 with \(\mathbf {A} = \mathbf {I}_p\)). When \(M_0\) fixes only \(q\) coefficients, the same expansion applies to both fits, and the \(p - q\) coefficients estimated in both cancel, leaving \(q\) squared approximately independent standard normals. □

This is the analogue of the extra-sum-of-squares \(F\)-test. In particular, comparing a model with the intercept-only model tests whether the regressors matter at all. The likelihood ratio test is generally more reliable than the Wald test, which can behave badly when a coefficient is large.

Goodness of fit

When the \(m_i\) are reasonably large, the saturated model is a sensible alternative, and the deviance of the fitted model is approximately \(\chi ^2_{n-p}\) if the model is correct. A large deviance relative to \(n - p\) signals lack of fit. A similar statistic is Pearson’s chi-square, \[ X^2 = \sum _{i=1}^n \frac {(y_i - \hat \mu _i)^2}{m_i\hat \pi _i(1 - \hat \pi _i)} , \] the sum of the squared Pearson residuals \((y_i - \hat \mu _i)/\sqrt {m_i\hat \pi _i(1 - \hat \pi _i)}\). For ungrouped binary data (\(m_i = 1\)) neither statistic has a \(\chi ^2\) distribution, and goodness of fit is checked by grouping the observations by their fitted probabilities (the Hosmer–Lemeshow test) or by plotting residuals.

Example 11.3.3. For the insecticide data of Example 11.2.3:

(a): test whether dose affects mortality, by the Wald and likelihood ratio tests;
(b): find a \(95\%\) confidence interval for the odds ratio;
(c): check the fit;
(d): estimate the dose that kills half the insects, the \(\mathrm {ED}_{50}\), with its standard error.

Solution. It is easiest to work with the centred dose \(u = x - 3\), for which \(\operatorname {logit}\hat \pi = \hat \alpha _0 + \hat \beta _1u\) with \(\hat \alpha _0 = -3.036 + 3(1.012) = 0\). The weights \(w_i = 20\hat \pi _i(1 - \hat \pi _i)\) are \[ 2.062, \quad 3.910, \quad 5.000, \quad 3.910, \quad 2.062 , \] symmetric about \(u = 0\), so \(\sum w_iu_i = 0\) and \(\mathbf {X}^\top \hat {\mathbf {W}}\mathbf {X}\) is diagonal, with entries \(\sum w_i = 16.944\) and \(\sum w_iu_i^2 = 2(4 \times 2.062 + 3.910) = 24.314\).

(a) \(\operatorname {se}(\hat \beta _1) = 1/\sqrt {24.314} = 0.2028\), and the Wald statistic is \(Z_0 = 1.012/0.2028 = 4.99\), far beyond \(1.96\). For the likelihood ratio test, the intercept-only model fits \(\hat \pi = 50/100 = 0.5\) everywhere, so \(\hat \mu _i = 10\). Dose \(3\) contributes nothing to its deviance, and doses \(5\) and \(4\) contribute the same as doses \(1\) and \(2\) (with \(y\) and \(20 - y\) exchanged), so \begin {align*} D(M_0) &= 2 \times 2\Bigl [\bigl (2\log \tfrac {2}{10} + 18\log \tfrac {18}{10}\bigr ) + \bigl (6\log \tfrac {6}{10} + 14\log \tfrac {14}{10}\bigr )\Bigr ] \\ &= 4(7.361 + 1.646) = 36.03 . \end {align*}

The fitted model has \[ D(M_1) = 0.0566 + 0.1114 + 0 + 0.1114 + 0.0566 = 0.336 , \] so \(D(M_0) - D(M_1) = 35.69\), against \(\chi ^2_{0.05,1} = 3.84\). Both tests reject \(H_0\colon \beta _1 = 0\) decisively.

(b) A \(95\%\) interval for \(\beta _1\) is \(1.012 \pm 1.96(0.2028) = (0.615, 1.410)\), and for the odds ratio \((e^{0.615}, e^{1.410}) = (1.85, 4.09)\).

(c) \(D(M_1) = 0.336\) on \(n - p = 5 - 2 = 3\) degrees of freedom, far below \(\chi ^2_{0.05,3} = 7.81\), and Pearson’s \(X^2 = 0.337\) agrees. The logistic model fits well.

(d) \(\pi = \frac 12\) where the logit is \(0\): \(\hat \alpha _0 + \hat \beta _1u = 0\), so \(\hat u_{50} = -\hat \alpha _0/\hat \beta _1 = 0\) and \(\widehat {\mathrm {ED}}_{50} = 3\). By the delta method, as in Theorem 10.1.5, with \(\operatorname {Cov}(\hat \alpha _0, \hat \beta _1) = 0\) and \(\hat \alpha _0 = 0\), \[ \widehat {\operatorname {Var}}(\hat u_{50}) = \frac {\widehat {\operatorname {Var}}(\hat \alpha _0)}{\hat \beta _1^2} + \frac {\hat \alpha _0^2\,\widehat {\operatorname {Var}}(\hat \beta _1)}{\hat \beta _1^4} = \frac {1/16.944}{1.012^2} + 0 = 0.0576 , \] so \(\operatorname {se}(\widehat {\mathrm {ED}}_{50}) = 0.240\), and an approximate \(95\%\) interval is \(3 \pm 0.47\), from \(2.53\) to \(3.47\).

Example 11.3.4. In a study of a vaccine, \(10\) of \(100\) vaccinated people became infected, and \(30\) of \(100\) unvaccinated people. Fit \(\operatorname {logit}\pi = \beta _0 + \beta _1z\), where \(z = 1\) for the vaccinated, and estimate the odds ratio.

Solution. By Theorem 11.2.4, with the unvaccinated as group \(0\), \[ \hat \beta _0 = \log \frac {30}{70} = -0.847, \qquad \hat \beta _1 = \log \frac {10 \times 70}{30 \times 90} = \log 0.259 = -1.350 , \] and \[ \operatorname {se}(\hat \beta _1) = \sqrt {\tfrac {1}{30} + \tfrac {1}{70} + \tfrac {1}{10} + \tfrac {1}{90}} = \sqrt {0.1587} = 0.398 . \] The Wald statistic is \(-1.350/0.398 = -3.39\). A \(95\%\) interval for \(\beta _1\) is \(-1.350 \pm 0.781\), so the odds ratio is \(0.26\), with interval \((0.12, 0.57)\). Vaccination cuts the odds of infection to about a quarter. (In terms of probabilities, from \(0.30\) to \(0.10\): the odds ratio and the ratio of probabilities differ when the probabilities are not small.)

[next] [prev] [prev-tail] [front] [up]

Questions on this section

Stuck on something here? Ask below and it stays attached to this topic.