Information-Based Covariance Estimation

Why Information Gives a Covariance

A fitted GLM coefficient vector is uncertain because a different sample would usually produce a different maximum-likelihood estimate. Large-sample likelihood theory describes that uncertainty through the information in the data: more information means less variability, while dependence among coefficients appears as off-diagonal entries in the covariance matrix. The practical difficulty is that the information at the true coefficient is unknown. This tutorial derives the inverse-information approximation and explains why information evaluated at the fitted coefficient can be used instead.

We consider a coefficient vector $\beta\in\mathbb{R}^p$, where $p$ is fixed, and observations with log-likelihood $\ell_n(\beta)$. The design may be fixed or random, but the likelihood must be correctly specified and the regularity conditions for the usual GLM likelihood asymptotics must hold. In particular, the true coefficient $\beta_0$ is an interior, identifiable parameter, the maximum-likelihood estimator $\hat\beta$ is consistent, and the information grows at order $n$ without becoming singular in any coefficient direction.

Definition: The Fisher information for $\beta$ in the sample of size $n$, evaluated at $\beta$, is $$ I_n(\beta)=E_\beta\!\left[s_n(\beta)s_n(\beta)^\mathsf{T}\right], $$ where $s_n(\beta)=\partial\ell_n(\beta)/\partial\beta$ is the score vector. Under the usual differentiability and integrability conditions for a correctly specified regular likelihood, this also equals $E_\beta[-\partial^2\ell_n(\beta)/\partial\beta\partial\beta^\mathsf{T}]$.

Information is a matrix, not a single measure of precision. A large diagonal entry indicates substantial information about one coefficient; off-diagonal entries encode how information about coefficients is shared. The covariance matrix will therefore generally have nonzero off-diagonal entries as well.

For independent observations in a regular GLM with known dispersion $\phi$, a common expression for coefficient information is

I_n(\beta)=\sum_{i=1}^n x_i x_i^\mathsf{T} \frac{\left(d\mu_i/d\eta_i\right)^2}{\phi V(\mu_i)}, \qquad \eta_i=x_i^\mathsf{T}\beta,\quad \mu_i=g^{-1}(\eta_i).

Here $x_i$ is the predictor vector, $g$ is the link, and $V(\mu_i)$ is the variance function. This formula makes the dependence on the coefficient explicit: changing $\beta$ changes the fitted means and generally changes the information weights.

From the Normal Limit to Inverse Information

The earlier material on asymptotic normality gives, under regularity conditions, a limit of the form

\sqrt{n}(\hat\beta-\beta_0)\xrightarrow{d}N\!\left(0,J^{-1}\right), \qquad \frac{I_n(\beta_0)}{n}\longrightarrow J,

where $J$ is a finite, positive-definite $p\times p$ matrix. The scaling by $n$ matters. The total information $I_n(\beta_0)$ grows with sample size, whereas the information per observation approaches $J$.

Theorem: If $I_n(\beta_0)/n\to J$ for a positive-definite matrix $J$, and the maximum-likelihood estimator has the asymptotic normal limit above, then its large-sample normal approximation has covariance $I_n(\beta_0)^{-1}$. In particular, $n I_n(\beta_0)^{-1}\to J^{-1}$.

To see the scaling directly, the covariance matrix in the limiting distribution of $\sqrt n(\hat\beta-\beta_0)$ is $J^{-1}$. Dividing that random vector by $\sqrt n$ divides its covariance matrix by $n$. Thus the corresponding approximation for $\hat\beta$ is $J^{-1}/n$. Since $I_n(\beta_0)/n\to J$, its inverse is approximately $J^{-1}$, so $I_n(\beta_0)^{-1}$ is the same first-order covariance approximation.

\operatorname{Cov}_{\mathrm{approx}}(\hat\beta) =\frac{1}{n}J^{-1} \approx I_n(\beta_0)^{-1}.

The symbol $\operatorname{Cov}_{\mathrm{approx}}$ emphasizes that this is a covariance in the asymptotic normal approximation. Convergence in distribution alone does not guarantee convergence of the actual finite-sample covariance matrix; that stronger claim requires additional moment conditions. Nor is inverse information generally the exact finite-sample covariance of the estimator.

Replacing the Unknown Information

The expression $I_n(\beta_0)^{-1}$ cannot usually be calculated because $\beta_0$ is unknown. A practical solution is to evaluate the model’s information formula at $\hat\beta$. The justification is that a consistent coefficient estimate is close to the true coefficient, and information changes continuously with the coefficient under suitable regularity conditions.

Theorem: Suppose $p$ is fixed, $\hat\beta\xrightarrow{p}\beta_0$, and the normalized information satisfies $I_n(\hat\beta)/n\xrightarrow{p}J$ for a positive-definite matrix $J$. Then $$ n I_n(\hat\beta)^{-1}\xrightarrow{p}J^{-1}. $$ Consequently, $I_n(\hat\beta)^{-1}$ consistently estimates the first-order covariance scale $J^{-1}/n$.

A sufficient condition for the theorem’s information convergence is that the normalized information is continuous in the coefficient near $\beta_0$, with its variation controlled uniformly in $n$. For example, suppose

\frac{I_n(\beta_0)}{n}\longrightarrow J \quad\text{and}\quad \left\|\frac{I_n(\hat\beta)-I_n(\beta_0)}{n}\right\| \xrightarrow{p}0.

Then the triangle inequality gives

\left\|\frac{I_n(\hat\beta)}{n}-J\right\| \leq \left\|\frac{I_n(\hat\beta)-I_n(\beta_0)}{n}\right\| + \left\|\frac{I_n(\beta_0)}{n}-J\right\| \xrightarrow{p}0.

This proves convergence of the fitted, normalized information to $J$. Because $J$ is positive definite, its smallest eigenvalue is strictly positive. Matrices sufficiently close to $J$ are therefore nonsingular, and matrix inversion is continuous in a neighborhood of $J$. Applying this continuity yields

\left(\frac{I_n(\hat\beta)}{n}\right)^{-1} \xrightarrow{p}J^{-1}.

Since $\left(I_n(\hat\beta)/n\right)^{-1}=n I_n(\hat\beta)^{-1}$, the stated result follows. This is the plug-in justification: consistency of $\hat\beta$, together with regular information behavior, allows the unknown coefficient in the information formula to be replaced by its estimate. $\square$

In practice, the calculation is often performed by computing the model’s information matrix at the fitted coefficients and inverting that matrix. If the likelihood is expressed using information per observation rather than total information, its inverse must be divided by $n$ to obtain the covariance approximation for $\hat\beta$.

Worked Covariance Calculations

Worked Example: Bernoulli Intercept-Only Model

Suppose $Y_1,\ldots,Y_{100}$ are independent Bernoulli observations with success probability $p$. Use the logit parameter $\beta=\log\{p/(1-p)\}$. The information for $\beta$ is $I_{100}(\beta)=100p(1-p)$. If 30 of the 100 observations are successes, calculate the fitted-information covariance and standard error.

The maximum-likelihood estimate is $\hat p=30/100=0.3$. Its corresponding coefficient estimate is $\hat\beta=\log(0.3/0.7)=\log(3/7)$, approximately $-0.8473$. Substituting $\hat p$ into the information gives

I_{100}(\hat\beta)=100(0.3)(0.7)=21.

For a scalar parameter, the inverse-information covariance estimate is

\widehat{\operatorname{Cov}}_{\mathrm{approx}}(\hat\beta) =21^{-1} =\frac{1}{21} \approx 0.04762.

The corresponding standard error is $\sqrt{1/21}\approx 0.2182$. This is an asymptotic, model-based estimate of coefficient uncertainty; it is not an assertion that the finite-sample variance equals $1/21$ exactly.

Worked Example: Two-Coefficient Logistic Regression

Consider a logistic regression with predictor vector $x=(1,z)^\mathsf{T}$, with 50 observations at $z=0$ and 50 at $z=1$. Suppose the group data contain 10 successes at $z=0$ and 30 at $z=1$. The fitted probabilities are therefore $0.2$ and $0.6$. Find the fitted-information covariance matrix for the intercept and slope.

For logistic regression, each observation contributes $p_i(1-p_i)x_i x_i^\mathsf{T}$. The group at $z=0$ has weight $0.2(0.8)=0.16$ per observation and predictor vector $(1,0)^\mathsf{T}$, so its total contribution is

50(0.16) \begin{pmatrix} 1&0\\ 0&0 \end{pmatrix} = \begin{pmatrix} 8&0\\ 0&0 \end{pmatrix}.

The group at $z=1$ has weight $0.6(0.4)=0.24$ and predictor vector $(1,1)^\mathsf{T}$, giving

50(0.24) \begin{pmatrix} 1&1\\ 1&1 \end{pmatrix} = \begin{pmatrix} 12&12\\ 12&12 \end{pmatrix}.

Adding the group contributions gives fitted information

I_{100}(\hat\beta)= \begin{pmatrix} 20&12\\ 12&12 \end{pmatrix}.

Its determinant is $20(12)-12(12)=96$. The inverse is consequently

I_{100}(\hat\beta)^{-1} =\frac{1}{96} \begin{pmatrix} 12&-12\\ -12&20 \end{pmatrix} = \begin{pmatrix} 0.125&-0.125\\ -0.125&0.2083 \end{pmatrix}.

The diagonal entries approximate the intercept and slope variances, and the negative off-diagonal entries indicate negative covariance between the fitted coefficients. The off-diagonal terms should not be discarded when joint uncertainty is relevant.

Worked Example: Poisson Regression With Exposure

Suppose two independent Poisson counts have means $e_0\lambda_0$ and $e_1\lambda_1$, with exposures $e_0=20$ and $e_1=10$. Model the rates by $\log(\lambda_i)=\beta_0+\beta_1 z_i$, where $z_0=0$ and $z_1=1$. If the observed counts are $Y_0=10$ and $Y_1=20$, calculate the fitted-information covariance matrix.

The fitted rates are $\hat\lambda_0=10/20=0.5$ and $\hat\lambda_1=20/10=2$. Thus $\hat\beta_0=\log(0.5)$, and $\hat\beta_1=\log(2)-\log(0.5)=\log 4$. For a Poisson log-rate model, the information contribution from a group is $e_i\lambda_i x_i x_i^\mathsf{T}$. The fitted weights are $20(0.5)=10$ and $10(2)=20$. Hence

I_2(\hat\beta) =10\begin{pmatrix}1&0\\0&0\end{pmatrix} +20\begin{pmatrix}1&1\\1&1\end{pmatrix} =\begin{pmatrix}30&20\\20&20\end{pmatrix}.

The determinant is $30(20)-20(20)=200$, so

I_2(\hat\beta)^{-1} =\frac{1}{200} \begin{pmatrix}20&-20\\-20&30\end{pmatrix} =\begin{pmatrix}0.10&-0.10\\-0.10&0.15\end{pmatrix}.

This calculation illustrates how exposure and fitted rates affect the information weights. With only two aggregated counts, the numerical covariance approximation should not be treated as a precision guarantee; the large-sample interpretation concerns a regime in which information grows and the regularity conditions hold.

Common Mistake and Practical Interpretation

A frequent error is to use the inverse of information per observation as the covariance estimate for the coefficient vector. If $J_n(\beta)=I_n(\beta)/n$, then $J_n(\hat\beta)^{-1}$ estimates $J^{-1}$, the covariance scale for $\sqrt n(\hat\beta-\beta_0)$. The coefficient estimate itself has covariance on the smaller scale $J^{-1}/n$. Equivalently, use $I_n(\hat\beta)^{-1}$ when $I_n$ is total sample information.

A second error is to describe inverse information as the exact variance of a finite-sample MLE. The argument here uses an asymptotic normal approximation and consistency of the fitted information. Small samples, weakly identified coefficients, or a nearly singular information matrix can make that approximation unreliable even when the formula can be computed.

Key takeaway: Under regular GLM likelihood conditions, the large-sample covariance approximation for $\hat\beta$ is the inverse of total information. Because the true coefficient is unknown, evaluate the information at $\hat\beta$; consistency and continuity justify that substitution. Always check whether the information used is total or per observation.

Check Your Understanding

Use the definitions and scaling arguments above to answer these questions.

If $I_n(\beta_0)/n\to J$, what is the limiting covariance matrix for $\sqrt n(\hat\beta-\beta_0)$?
Why does positive definiteness of $J$ matter when replacing fitted information by its inverse?
A reported matrix is the inverse of information per observation. What scaling is needed to use it as a covariance approximation for $\hat\beta$?
In the logistic example, what does the negative off-diagonal entry of the covariance matrix describe?
Why does the inverse-information approximation not, by itself, prove that the actual finite-sample covariance equals the reported matrix?