Common Computational Errors in Profile Intervals

Why Profile-Interval Computation Needs Several Checks

A profile interval is not obtained by holding every fitted parameter fixed and varying only the parameter of interest. At each proposed target value, the nuisance parameters must be optimized again. Even with the correct objective, an optimizer can return a misleading answer if it starts poorly, searches outside the feasible parameter space, or stops before an endpoint has converged.

These are computational errors, but they change the inferential result: they can move an endpoint, exclude a valid part of the parameter space, or make an interval appear narrower than it should be. This tutorial develops checks for those failures. It treats the likelihood-ratio cutoff as given; choosing or calibrating that cutoff is a separate question.

Definition: Let $\psi$ be a scalar parameter of interest and let $\lambda$ denote nuisance parameters. At a fixed feasible value $\psi$, the profile log-likelihood is $$ \ell_p(\psi)=\sup_{\lambda\in\Lambda(\psi)}\ell(\psi,\lambda), $$ where $\Lambda(\psi)$ is the set of nuisance values allowed at that $\psi$. Write $\hat{\ell}=\sup_{\psi,\lambda}\ell(\psi,\lambda)$. For a specified cutoff $c>0$, the likelihood-ratio confidence set is $\{\psi: R(\psi)\le c\}$, where $R(\psi)=2\{\hat{\ell}-\ell_p(\psi)\}$. When $R$ is continuous at an interior endpoint of the confidence set, that endpoint satisfies $R(\psi)=c$.

The supremum in the definition matters. At every fixed $\psi$, the nuisance parameters are allowed to adapt to that constraint. In practice, an optimizer may only approximate this supremum, so the analyst needs to check both the optimization and the endpoint search. If the feasible parameter space has a boundary, the optimization and endpoint search must also respect it.

Reestimate Nuisance Parameters at Every Target Value

For each candidate $\psi$, maximize $\ell(\psi,\lambda)$ over the allowed nuisance values. Reusing the nuisance estimate from the unrestricted fit generally evaluates a slice through the likelihood, not the profile likelihood. The distinction can be made precise.

Result: Suppose $\lambda_0\in\Lambda(\psi)$ is a feasible nuisance value at the fixed target $\psi$. Then $$ \ell_p(\psi)\ge \ell(\psi,\lambda_0), \qquad 2\{\hat{\ell}-\ell_p(\psi)\}\le 2\{\hat{\ell}-\ell(\psi,\lambda_0)\}. $$ Thus fixing nuisance parameters instead of profiling them can only increase the computed likelihood-ratio statistic, provided the fixed value is feasible.

Proof. By definition, $\ell_p(\psi)$ is the supremum of $\ell(\psi,\lambda)$ over all $\lambda\in\Lambda(\psi)$. Since $\lambda_0$ belongs to that set, its log-likelihood cannot exceed the supremum: $\ell(\psi,\lambda_0)\le\ell_p(\psi)$. Subtracting both sides from the same value $\hat{\ell}$ reverses the inequality, and multiplying by two preserves it. $\square$

Worked Example: Reprofiling a Poisson Rate Ratio

Suppose independent counts $Y_1=9$ and $Y_0=4$ have Poisson means $e^{\alpha+\psi}$ and $e^\alpha$, respectively. The target is the log rate ratio $\psi$, and $\alpha$ is a nuisance parameter. Constants not depending on the parameters can be dropped, giving $$ \ell(\alpha,\psi)=9(\alpha+\psi)-e^{\alpha+\psi}+4\alpha-e^\alpha. $$

At fixed $\psi$, differentiate with respect to $\alpha$. The score equation is $13-e^\alpha(1+e^\psi)=0$, and the second derivative is $-e^\alpha(1+e^\psi)<0$, so this gives the unique maximum: $$ e^{\hat{\alpha}(\psi)}=\frac{13}{1+e^\psi}. $$

The total fitted mean is therefore 13 at every fixed $\psi$. Substitution gives the profile log-likelihood, up to an additive constant, as $9\psi-13\log(1+e^\psi)$. Its maximum occurs where $9-13e^\psi/(1+e^\psi)=0$, hence $\hat{\psi}=\log(9/4)$. At $\psi=0$, profiling gives the common fitted mean $e^{\hat{\alpha}(0)}=13/2=6.5$. The likelihood-ratio statistic is $$ R(0)=2\left\{9\log\left(\frac{9}{6.5}\right)+4\log\left(\frac{4}{6.5}\right)\right\}\approx 1.974. $$

Now make the common computational error of fixing the baseline mean at its unrestricted estimate, $e^{\hat{\alpha}}=4$, rather than reestimating it at $\psi=0$. The fixed-$\alpha$ fit has both means equal to 4. The profile fit improves its log-likelihood over that fixed fit by $$ 13\log(6.5/4)-(13-8)=13\log(6.5/4)-5\approx 1.312. $$ Consequently the fixed-nuisance likelihood-ratio statistic exceeds the correctly profiled statistic by $2[13\log(6.5/4)-5]\approx 2.623$. This is not a harmless change in computation: it can alter whether a target value lies inside the interval.

This example also suggests a practical diagnostic. When a nuisance estimate is available in closed form, compare the software profile likelihood against that value at several target points. Otherwise, rerun the nuisance optimization from multiple starting values and check that the attained log-likelihoods agree to a stated tolerance.

Starting Values and Numerical Optimization

A poor starting value can lead to very large steps, overflow, failure of line searches, or premature termination. A convergence message alone does not establish that the relevant constrained maximum has been found. In canonical GLMs, many likelihoods have useful concavity properties, but a numerical routine still has to represent parameters and steps stably, especially near boundaries or when the likelihood is nearly flat.

Worked Example: A Poor Start in a Constrained Poisson Fit

Return to the two-count model and fix $\psi=0$. Ignoring constants, the nuisance log-likelihood is $\ell_0(\alpha)=13\alpha-2e^\alpha$. Its score is $13-2e^\alpha$, so the exact maximizing value is $\hat{\alpha}(0)=\log(6.5)\approx1.872$.

A raw Newton update is $\alpha_{\mathrm{new}}=\alpha-\ell_0'(\alpha)/\ell_0''(\alpha)=\alpha+13/(2e^\alpha)-1$. If it starts at $\alpha=\log(0.01)\approx-4.605$, then $13/(2e^\alpha)=13/0.02=650$, giving $\alpha_{\mathrm{new}}\approx644.395$. Exponentiating such a value will overflow in ordinary floating-point arithmetic. The update has failed numerically even though the likelihood has a simple, finite maximum.

A robust implementation can use a step-halving or line-search rule, reject steps that decrease the objective, and initialize from a sensible estimate such as $\log(13/2)$. It should still verify the score or compare likelihood values; a stable starting value is a safeguard, not a proof of convergence.

Warm starts can help when evaluating a sequence of nearby target values: use the nuisance estimate at one value as a starting point at the next. However, they should not replace checks. A poor solution at one point can be propagated along the sequence. Repeating selected points from independent starts is a useful way to detect that problem.

Respect Parameter Boundaries

Optimization must be over the model's actual parameter space. A numerical transformation can be convenient, but it does not remove the boundary from the statistical problem. For instance, a logit parameterization represents only probabilities strictly between zero and one; it cannot attain $p=0$ or $p=1$ at any finite logit value. A boundary maximum may therefore appear as an optimizer that keeps moving toward an extreme rather than converging to a finite transformed estimate.

Worked Example: A Binomial Maximum at Zero

Let $X\sim\operatorname{Binomial}(10,p)$ and observe $X=0$, with $0\le p\le1$. The likelihood is $L(p)=(1-p)^{10}$, so its constrained maximum occurs at $\hat p=0$ and $L(\hat p)=1$. For $0\le p<1$, the likelihood-ratio statistic is $$ R(p)=2\log\left(\frac{L(0)}{L(p)}\right)=-20\log(1-p). $$

For illustration, use the supplied cutoff $c=3.841$. Solving $R(p)=c$ gives $1-p=\exp(-3.841/20)$, and hence the upper endpoint is $p=1-\exp(-3.841/20)\approx0.175$. The lower end of the accepted set is the feasible boundary $p=0$; it is not a root of $R(p)=c$.

The calculation shows how to handle the likelihood and endpoint at the boundary, not how to calibrate a confidence level for this nonregular case. The cutoff must be justified for the inferential problem at hand. Computationally, the essential checks are that the optimizer admits $p=0$, that no invalid value such as $p<0$ is searched, and that the one-sided endpoint equation is solved inside the feasible interval.

A boundary check should be explicit: identify the feasible parameter set before optimization, determine whether its endpoints are included, and compare the likelihood at a boundary with the best interior value. Do not report a large finite logit as if it were the exact estimate $p=0$, or treat a missing finite root as proof that the endpoint is absent. The likelihood may attain its maximum at the boundary, or the confidence set may meet the boundary without a likelihood-ratio equality there.

Endpoint Convergence Requires a Bracket

An endpoint is usually found by solving $g(\psi)=0$, where $g(\psi)=R(\psi)-c$. A small change between successive iterates is not enough to show that the desired root has been located: the iteration may have stalled, followed the wrong root, or stopped at a point where the profile likelihood itself was computed inaccurately. When $g$ is continuous and a sign-changing bracket is available, bisection gives a transparent convergence guarantee.

Result: If $g$ is continuous on $[a,b]$ and $g(a)g(b)<0$, bisection produces nested intervals containing a root. After $k$ bisection steps, the interval width is $(b-a)/2^k$; its midpoint is therefore within $(b-a)/2^{k+1}$ of a root in that interval.

Proof. Continuity and the opposite signs imply, by the intermediate value theorem, that $g$ has at least one root in $[a,b]$. At each step, evaluate $g$ at the midpoint $m=(a+b)/2$. If $g(m)=0$, the root has been found. Otherwise, $g(m)$ has the same sign as exactly one of the two endpoint values, so the half-interval with opposite endpoint signs still contains a root. Repeating this argument gives nested intervals, each half as wide as its predecessor. Thus after $k$ steps the width is $(b-a)/2^k$; the midpoint's distance from any root contained in that interval is at most half the width. $\square$

Worked Example: Bracketing a Poisson Profile Endpoint

For the Poisson data above, the profile likelihood is equivalent, up to constants, to the binomial likelihood for 9 successes out of 13 with $p=e^\psi/(1+e^\psi)$. Thus $$ R(\psi)=2\left\{9\log\left(\frac{9/13}{p}\right)+4\log\left(\frac{4/13}{1-p}\right)\right\}. $$ For a cutoff $c=3.841$, define $g(\psi)=R(\psi)-3.841$. At $\psi=2$, $p=e^2/(1+e^2)$, and direct substitution gives $R(2)\approx3.252$, so $g(2)\approx-0.589$. At $\psi=2.3$, substitution gives $R(2.3)\approx4.84$, so $g(2.3)>0$. The upper endpoint is therefore bracketed between 2 and 2.3.

The midpoint is 2.15. Evaluating the same expression gives $R(2.15)\approx4.017$, so the root is in $[2,2.15]$. The next midpoint is 2.075; its statistic is approximately 3.62, so the root is in $[2.075,2.15]$. The bracket width has fallen from 0.3 to 0.075 after two steps. Further bisection narrows the endpoint uncertainty by exactly a factor of two at each step.

In a real computation, retain the bracket and report the chosen endpoint tolerance. Recompute the profile likelihood accurately at both bracket ends; otherwise, rounding or optimization error can produce a false sign change. A root finder that reports a single number without a bracket or an endpoint tolerance has not provided the same convergence evidence.

Common Computational Mistakes

Holding nuisance estimates fixed. This substitutes a likelihood slice for the profile likelihood. Reoptimize nuisance parameters at every target value; the fixed-nuisance likelihood-ratio statistic can be too large.
Trusting one starting value or a success flag. A routine may stop after an unstable step or at a poor solution. Check the attained log-likelihood, gradients or score where available, and agreement across starts.
Ignoring the feasible boundary. A transformed parameterization can hide an attainable boundary maximum. Optimize on the constrained parameter space and compare boundary and interior likelihoods.
Stopping when iterates barely move. Small steps can indicate stalling, not endpoint convergence. Use a sign-changing bracket and a stated tolerance in the target parameter.
Confusing the endpoint equation with cutoff calibration. Solving $R(\psi)=c$ accurately does not establish that $c$ has the desired coverage, particularly in a nonregular boundary problem.

The central computational principle is to keep the inferential definitions visible in the numerical workflow: maximize over the correct nuisance set, evaluate only feasible parameters, and verify endpoint roots rather than relying on optimizer messages. These checks make it easier to distinguish a genuine feature of the likelihood from a numerical artifact.

Check Your Understanding

Use the definitions and checks in this tutorial to answer the following.

At a fixed target value, why can the profile likelihood-ratio statistic not exceed the statistic obtained by using any one feasible fixed nuisance value?
In the two-count Poisson model, what equation determines the nuisance maximum at a fixed $\psi$, and what fitted total mean does it imply?
Why can a logit parameterization fail to represent a binomial maximum at $p=0$ exactly?
What conditions make a sign-changing bisection bracket a reliable endpoint check, and how does its width change after each step?
Why does successful numerical solution of the endpoint equation not, by itself, justify the likelihood-ratio cutoff?