Tutorials › Generalized Linear Models › BIC And Laplace Approximation

Selection and Predictive Comparison · Tutorial 803 of 981

BIC And Laplace Approximation

Learn how a Laplace approximation to marginal likelihood produces BIC, and why each additional parameter incurs a penalty that grows with the logarithm of the sample size.

Advanced 11 min read

What You'll Learn

  • Express a model’s marginal likelihood as an integral of its likelihood against a parameter prior
  • Derive the leading Laplace approximation around an interior maximum-likelihood estimate
  • Explain the regularity conditions needed for the BIC approximation
  • Obtain the \(k\log n\) penalty and identify what \(k\) counts
  • Calculate and interpret BIC differences in concrete GLM examples
  • Recognize why BIC is an approximation rather than an exact marginal likelihood

From Likelihood to Marginal Likelihood

A fitted model’s maximized likelihood measures how well its best-fitting parameter values describe the observed data. Bayesian model comparison instead uses the marginal likelihood: it averages the likelihood over parameter values, weighting them by a prior distribution. This averaging accounts for the fact that a model’s fit is not achieved at every possible parameter value.

For a candidate GLM \(M\) with parameter vector \(\theta\in\mathbb{R}^k\), let \(L_M(\theta)\) be its likelihood and let \(\pi_M(\theta)\) be a proper prior density. The marginal likelihood is

$$ m_M(y)=\int L_M(\theta)\pi_M(\theta)\,d\theta. $$

The integral can be difficult to compute. When the likelihood is concentrated around a well-behaved interior maximum, the Laplace approximation replaces the log-likelihood near that maximum by a quadratic function. This tutorial derives how the resulting approximation leads to the Bayesian information criterion (BIC).

Definition: For a model with \(k\) free parameters fitted to \(n\) independent observations, the Bayesian information criterion is $$ \operatorname{BIC}=-2\ell(\hat{\theta})+k\log n, $$ where \(\ell(\hat{\theta})\) is the maximized log-likelihood. When comparing models, smaller BIC is preferred. The criterion is an approximation to minus twice the log marginal likelihood, up to terms that are bounded as \(n\) grows under regularity conditions.

The Laplace Approximation

Assume the observations are independent, so the log-likelihood is a sum of \(n\) contributions. For a fixed-dimensional parameter \(\theta\), write \(\hat{\theta}\) for an interior maximum-likelihood estimate and define the observed information matrix

$$ H_n=-\left.\frac{\partial^2\ell(\theta)}{\partial\theta\,\partial\theta^\mathsf{T}}\right|_{\theta=\hat{\theta}}. $$

The derivation needs more than the existence of a maximum. A standard sufficient set of conditions is: the parameter space contains an open neighborhood of \(\hat{\theta}\); the log-likelihood has three continuous derivatives there; \(H_n\) is positive definite with eigenvalues of order \(n\); the prior density is continuous and positive at \(\hat{\theta}\); and the contribution to the integral outside a fixed neighborhood of \(\hat{\theta}\) is negligible compared with \(e^{\ell(\hat{\theta})}n^{-k/2}\). In addition, the local Taylor remainder must be small on the \(n^{-1/2}\) scale. These conditions hold in many regular, fixed-dimensional GLMs, but not automatically in every model.

Result (Laplace approximation): Under the conditions above, $$ m_M(y) = e^{\ell(\hat{\theta})}\pi_M(\hat{\theta})(2\pi)^{k/2}|H_n|^{-1/2}\{1+o(1)\}. $$ Here \(|H_n|\) is the determinant of the observed information matrix.

To see why, expand the log-likelihood around its maximum. Its gradient is zero at \(\hat{\theta}\), so for a small displacement \(u\),

$$ \ell(\hat{\theta}+u) = \ell(\hat{\theta})-\tfrac12u^\mathsf{T}H_nu+r_n(u), $$

where the Taylor remainder \(r_n(u)\) is negligible for displacements of order \(n^{-1/2}\). In the integral for \(m_M(y)\), make the change of variables \(z=H_n^{1/2}u\). The determinant of this change of variables gives \(du=|H_n|^{-1/2}dz\). The local part of the integral becomes

$$ e^{\ell(\hat{\theta})}|H_n|^{-1/2} \int e^{-\frac12z^\mathsf{T}z+r_n(H_n^{-1/2}z)} \pi_M(\hat{\theta}+H_n^{-1/2}z)\,dz. $$

For each fixed \(z\), the remainder tends to zero and the prior density tends to \(\pi_M(\hat{\theta})\). The local curvature and tail conditions ensure that the transformed integrand is dominated by an integrable Gaussian bound locally and that the omitted tails are negligible. Thus dominated convergence gives the limiting integral \(\pi_M(\hat{\theta})\int e^{-z^\mathsf{T}z/2}dz\). The \(k\)-dimensional Gaussian integral is \((2\pi)^{k/2}\), yielding the stated approximation. \(\square\)

Why the Penalty Is \(k\log n\)

For independent observations, the curvature of the log-likelihood is typically proportional to sample size. More precisely, suppose \(H_n=nJ_n\), where the eigenvalues of \(J_n\) stay bounded away from zero and infinity. Then \(|H_n|=n^k|J_n|\). Taking logarithms of the Laplace approximation gives

$$ -2\log m_M(y) = -2\ell(\hat{\theta})+k\log n +\log|J_n|-2\log\pi_M(\hat{\theta})-k\log(2\pi)+o(1). $$

If the last three displayed terms remain bounded as \(n\) grows, the leading terms are exactly BIC. The penalty is \(k\log n\) because the likelihood’s concentration region has width of order \(n^{-1/2}\) in each of \(k\) parameter directions. Its volume is therefore of order \(n^{-k/2}\); on the \(-2\log\) scale this contributes \(k\log n\).

Result (BIC from Laplace approximation): For a regular, fixed-dimensional model with a proper prior positive near the maximum-likelihood estimate and information growing proportionally to \(n\), $$ -2\log m_M(y)=-2\ell(\hat{\theta})+k\log n+O_p(1). $$ Consequently, BIC retains the terms that grow with sample size and omits bounded terms that can depend on the prior and the model.

The order of the remainder matters. BIC is not generally equal to minus twice the log marginal likelihood: prior density, local curvature, and other terms contribute to the omitted part. The approximation is especially useful for comparing models when these terms do not dominate the growing likelihood and dimension terms.

Worked Example: Laplace Approximation for a Binomial Intercept

Suppose \(n=100\) independent binary outcomes follow a logistic regression with an intercept only. There are 30 successes. Parameterize the success probability as \(p=\{1+\exp(-\eta)\}^{-1}\), and put a uniform prior on \(p\) over \((0,1)\). Calculate the Laplace approximation to the marginal likelihood, using the likelihood for the ordered binary observations.

The maximum-likelihood estimate is \(\hat{p}=30/100=0.3\), so \(\hat{\eta}=\log(0.3/0.7)\). The maximized log-likelihood is

$$ \ell(\hat{\eta}) =30\log(0.3)+70\log(0.7) =-61.0864. $$

For an intercept-only logistic model, the observed information for \(\eta\) is \(np(1-p)\). At the estimate it is \(100(0.3)(0.7)=21\). A uniform prior density on \(p\) becomes a density on \(\eta\) by multiplying by \(dp/d\eta=p(1-p)\), so \(\pi_\eta(\hat{\eta})=0.21\). The one-dimensional Laplace approximation is therefore

$$ m(y)\approx e^{-61.0864}(0.21)\sqrt{\frac{2\pi}{21}} =e^{-61.0864}(0.1149). $$

Thus \(-2\log m(y)\) is approximately \(122.1728-2\log(0.1149)=126.50\). BIC gives \(-2\ell(\hat{\eta})+\log(100)=122.1728+4.6052=126.78\). Their difference is small here but not zero: BIC omits bounded terms from the prior and curvature.

Interpreting the Sample-Size Penalty

When comparing two models fitted to the same \(n\) observations, the log-likelihood improvement must compensate for the extra parameter penalty. If model 1 has \(k_1\) parameters and model 0 has \(k_0\), then

$$ \operatorname{BIC}_1-\operatorname{BIC}_0 = -2\{\ell_1(\hat{\theta}_1)-\ell_0(\hat{\theta}_0)\} +(k_1-k_0)\log n. $$

For one additional parameter, the larger model has smaller BIC precisely when its increase in maximized log-likelihood exceeds \(\tfrac12\log n\). This threshold grows with \(n\), though only logarithmically.

Worked Example: Comparing Models at Two Sample Sizes

Consider two candidate GLMs: a base model with \(k_0=2\) parameters and an expanded model with \(k_1=3\). Suppose the expanded model improves the maximized log-likelihood by 3, so \(\ell_1-\ell_0=3\). First take \(n=20\). The BIC difference is

$$ \operatorname{BIC}_1-\operatorname{BIC}_0 =-2(3)+\log(20) =-6+2.996=-3.004. $$

The expanded model has the smaller BIC. Now take \(n=1000\), with the same log-likelihood improvement:

$$ \operatorname{BIC}_1-\operatorname{BIC}_0 =-6+\log(1000) =-6+6.908=0.908. $$

Now the base model has the smaller BIC. This isolates the effect of sample size: for one extra parameter the penalty rises from \(\log(20)=2.996\) to \(\log(1000)=6.908\). Holding the total fit improvement fixed across sample sizes is only an illustration of the penalty, not a claim about how such improvements usually behave in repeated sampling.

Worked Example: Checking the Approximation in a Gaussian GLM

Let \(X_1,\ldots,X_4\) be independent \(N(\theta,1)\) observations, a Gaussian GLM with known variance and one unknown intercept. Suppose their sample mean is \(\bar{x}=0.5\), and use a prior \(\theta\sim N(0,1)\). The maximum-likelihood estimate is \(\hat{\theta}=0.5\), and the observed information is \(H_4=4\). The Laplace approximation, divided by the maximized likelihood, is

$$ \pi(0.5)\sqrt{\frac{2\pi}{4}} = \frac{e^{-0.5^2/2}}{\sqrt{2\pi}}\sqrt{\frac{2\pi}{4}} =0.4413. $$

In this example the marginal likelihood can also be integrated exactly. Relative to the maximized likelihood, the integrand is the prior density times \(\exp\{-4(\theta-0.5)^2/2\}\). Completing the square gives

$$ \frac{m(y)}{L(\hat{\theta})} = \frac{1}{\sqrt{5}}\exp\left\{-\frac{4(0.5)^2}{2(5)}\right\} = \frac{e^{-0.1}}{\sqrt{5}} =0.4047. $$

The Laplace and exact values are not identical at \(n=4\). BIC keeps the leading penalty \(\log(4)=1.386\) for the one estimated parameter; it does not retain every finite-sample or prior-dependent adjustment.

What to Count and a Common Pitfall

Here \(k\) is the number of free parameters in the likelihood under the model being fitted. For a GLM, this includes the regression coefficients. If a dispersion or other likelihood parameter is estimated, it must also be counted unless the comparison is explicitly using a convention under which that parameter is fixed or cancels. Use the same likelihood definition and parameter-counting convention across candidate models.

A common mistake is to treat BIC’s penalty as a fixed surcharge per parameter, independent of sample size. That is not what its derivation gives: each additional parameter contributes \(\log n\) to BIC, or \(\tfrac12\log n\) on the log-marginal-likelihood scale. Another mistake is to call the Laplace formula exact. It is an asymptotic approximation and relies on regularity: boundary estimates, nonidentifiability, or a singular information matrix can invalidate the usual derivation.

Key takeaway: BIC follows by approximating the marginal-likelihood integral near a regular likelihood maximum. The volume of the concentrated likelihood contributes a factor of order \(n^{-k/2}\), which becomes the \(k\log n\) penalty on the BIC scale.

Check Your Understanding

Use the definitions and derivation above to answer these questions.

  1. For a model with \(k=4\) parameters and \(n=200\) observations, write its BIC in terms of its maximized log-likelihood.
  2. Two models differ by one parameter. What increase in maximized log-likelihood is needed for the larger model to have lower BIC when \(n=50\)?
  3. In the Laplace approximation, why does the determinant of the information matrix produce a term involving \(k\log n\)?
  4. Name two regularity conditions that can fail and make the usual BIC derivation unreliable.
  5. Why can \(-2\log m_M(y)\) and BIC differ by a nonzero amount even when the Laplace approximation is appropriate?
Clinical Biostats