Past exam of the mathematics course of the University of Cambridge / 2001 / iii / Paper 32 (Statistical Theory)

Work in a regular statistical model, with fixed support, differentiability under the integral, and nonsingular Fisher information. Let $U=(U_{\psi },U_{\lambda })$ be the score function. The blocks are orthogonal statistical parameters at $\theta$ whenThis is global parameter orthogonality if it holds throughout the parameter domain; otherwise it is an assertion at a specified point.

For an independent size-$n$ sample, the regular maximum-likelihood estimator satisfieswhere $I$ is the information in one observation. Thus the fitted blocks are independent to first order in their limiting normal distribution. The efficient information for $\psi$ is generally the Schur complement $I_{\psi \psi }-I_{\psi \lambda }{I}_{\lambda \lambda }^{-1}{I}_{\lambda \psi }$; under orthogonality it is just $I_{\psi \psi }$. Estimating the nuisance block therefore has no first-order variance cost for the interest block, compared with knowing the nuisance parameter.

This also reduces local movement of the constrained nuisance fit. If $j=-{\ell }^{\prime \prime }$ is the observed information, differentiating its nuisance score equation givesUnder regular independent sampling and orthogonality, $j_{\lambda \lambda }=O_p(n)$ and the centered mixed block is $O_p(\sqrt{n})$. Thus this derivative is $O_p(n^{-1/2})$ locally, and an $O_p(n^{-1/2})$ change in $\psi$ changes the nuisance fit by only $O_p(n^{-1})$. Neither exact finite-sample independence nor vanishing observed cross-information follows from expected-information orthogonality.

For the displayed density family, $\lambda$ is scalar andThe mean-zero score identity gives ${\mathbb{E}}_{\theta }{\nabla }_{\psi }t=0$ for $\lambda \ne 0$. ButsoThis is score factorization implies parameter orthogonality. If zero belongs to the parameter domain, the conclusion extends there by continuity of the regular expected derivatives. At zero the density itself has no dependence on $\psi$, so identifiability and nonsingular information for that block must be treated separately.

Partition the statistical parameter as $\theta =(\psi ,\lambda )$, with $\psi$ of interest. For each fixed $\psi$, maximize over the nuisance parameter to obtain ${\hat{\lambda }}_{\psi }$. The profile likelihood and profile log-likelihood areTheir maximizer is the interest component of the joint maximum-likelihood estimator. At a regular nuisance maximum, ${\ell }_p^{\prime }(\psi )={\ell }_{\psi }(\psi ,{\hat{\lambda }}_{\psi })$ and its negative Hessian matrix is the Schur complementThe profile likelihood ratio is unchanged by a one-to-one nuisance reparametrization, and transforms naturally under a one-to-one interest reparametrization. Under the regular conditions of Wilks theorem, $2\{{\ell }_p(\hat{\psi })-{\ell }_p({\psi }_0)\}$ has a limiting chi-squared distribution with dimension of $\psi$ degrees of freedom. Profiling is maximization, rather than integration or conditioning; $L_p$ need not be a normalized data density, and nuisance fitting can cause appreciable higher-order distortion.

The modified profile likelihood corrects for nuisance curvature and sample-coordinate distortion. In suitable coordinates $(\hat{\psi },\hat{\lambda },a)$, with $a$ an ancillary statistic, setholding $\hat{\psi },a$ fixed when differentiating the data coordinates. Up to an interest-independent factor, the Barndorff-Nielsen modified profile likelihood isEquivalently, implicit differentiation of the nuisance score equation gives $J_{\psi }=|j_{\lambda \lambda }{|}^{-1}|{\ell }_{\lambda ;\hat{\lambda }}|$, andThe semicolon denotes sample-coordinate differentiation, not a second parameter derivative. In orthogonal coordinates, the Cox-Reid adjusted profile likelihood uses the information adjustment alone; it agrees with this modification when $J_{\psi }=1$, as happens below.

For the Inverse Gaussian distribution, write $s={\sum }_i{y}_i$, $r={\sum }_i{y}_i^{-1}$ and $\bar{y}=s/n$. Terms independent of both parameters can be absorbed into $c(y)$, givingThe nuisance score function isIt is positive below $\bar{y}$ and negative above it. Thus its global maximum isindependent of $\psi$. Put $D=r-n/\bar{y}$. ThenThe nuisance observed information at this fit isSince the constrained nuisance estimator is literally the full nuisance estimator as a function of the data, $J_{\psi }=1$ in the specified sample coordinates. Hence the modified profile likelihood for inverse Gaussian shape isFor $n\ge 2$ and a nonconstant positive sample, the Cauchy-Schwarz inequality gives $D>0$. The corresponding fits are ${\hat{\psi }}_p=n/D$ and ${\hat{\psi }}_m=(n-1)/D$. A constant sample has $D=0$ and no finite profile maximum; under continuous independent sampling this event has probability zero when $n\ge 2$. In particular, the modified coefficient accounts for the fitted mean nuisance parameter.

Let $K(t)=\log \mathbb{E}{e}^{tY}$ be the one-observation cumulant-generating function, finite on an open interval containing zero. At a value $y$, find the real tilt $\hat{t}$ solvingIf the distribution is nondegenerate then $K^{\prime \prime }(t)>0$ in the interior, so the solution is unique whenever $y$ is in the range of $K^{\prime }$. The leading saddlepoint density approximation to the sample mean isTo derive it, use exponential tilting:Under this tilted distribution the mean is $K^{\prime }(t)$ and the variance is $K^{\prime \prime }(t)$, by differentiating its cumulant function $K(t+u)-K(t)$. For $S_n={\sum }_i{Y}_i$, changing measure gives the exact density identityAt $t=\hat{t}$, $ny$ is the tilted mean. A local normal approximation there gives $f_{S_n,\hat{t}}(ny)\simeq (2\pi n{K}^{\prime \prime }(\hat{t}){)}^{-1/2}$. The Jacobian determinant for $S_n=n\bar{Y}$ contributes the factor $n$, yielding the boxed formula.

This approximation retains the exponential large-deviation factor, while the curvature is evaluated at the tilt appropriate to each $y$. At the population mean, $\hat{t}=0$, so it reduces to the usual central normal approximation. Under suitable smooth nonlattice transform-inversion conditions, its relative error at fixed interior $y$ is $O(n^{-1})$; existence of a few moments alone does not establish this density accuracy. Boundary values and absent exponential moments need separate treatment. The leading density need not integrate exactly to one, so it can be divided by its integral when a normalized approximation is desired.

For the unit-scale Laplace distribution, integrating separately over the two half-lines gives the moment-generating functionfor $|t|<1$; the appropriate half-line integral diverges outside that interval. ThereforeIts derivatives needed for the saddlepoint density approximation areSet $D=\sqrt{1+y^2}$. Solving $K^{\prime }(\hat{t})=y$ and selecting the root in $(-1,1)$ gives the stable formulaincluding $\hat{t}=0$ at $y=0$. The other quadratic root is outside the domain. At the correct saddlepoint,Substitution gives the saddlepoint density of a Laplace sample mean:It is symmetric in $y$, as the underlying Laplace distribution requires. At zero its value is $\sqrt{n/(4\pi )}$, consistent with the local normal approximation having variance $2/n$. This is an approximation, not an exact finite-$n$ density or an automatically normalized one.

An M-estimator based on an estimating equation chooses a consistent root ofIts population target $T(F)={\theta }_F$ solves ${\mathbb{E}}_F\psi (X,{\theta }_F)=0$. For a scalar parameter, suppose the root is locally unique andFor the contaminated distribution $F_{\epsilon }=(1-\epsilon )F+\epsilon {\delta }_x$, differentiate its population equation at zero contamination. The influence function satisfiesConsequentlyThe constant of proportionality is independent of $x$. A local expansion of the sample estimating equation similarly givesThe central limit theorem yields asymptotic normality with the sandwich variance of an M-estimatorThus the first-order variance of the unscaled estimator is $V/n$. For a vector parameter the same differentiation gives $\mathrm{IF}=A^{-1}\psi$ and $V=A^{-1}B{A}^{-\mathsf{T}}$, with $A=-\mathbb{E}{\partial }_{\theta }\psi$ and $B=\mathbb{E}\psi {\psi }^{\mathsf{T}}$.

For the location problem, take the target location to be zero and let $Z$ have the standard normal distribution, with density $\varphi$ and distribution function $\Phi$. An odd score has $\mathbb{E}\psi (Z)=0$. For a differentiable score, integration by parts givesThe last expression also gives the appropriate derivative of the population equation for bounded monotone scores with jumps, by differentiating the shifted normal density. Assume $A>0$ and put $h=\psi /A$. ThenMultiplying a score by a positive constant changes neither its estimator nor this normalized influence function.

First establish feasibility. The bound impliesso necessarily $C\ge \sqrt{\pi /2}$. For $C>\sqrt{\pi /2}$, let $K>0$ solveThis solution is unique: the ratio tends to $\sqrt{\pi /2}$ as $K\downarrow 0$, tends to infinity as $K\to \infty$, and has positive derivative becauseTake the Huber score ${\psi }_K(x)=\max (-K,\min (x,K))$. Its derivative is one on $(-K,K)$ and zero outside, so $A=p_K$. HenceIt is odd, nondecreasing, meets the bound, and satisfies $\mathbb{E}[Z{h}_K(Z)]=1$.

To prove optimality, put $a=1/p_K$. For each fixed $x$, the minimum of $u^2-2axu$ over $|u|\le C$ is attained by projecting $ax$ onto $[-C,C]$, which is exactly $h_K(x)$. Therefore every competing normalized influence function $h$ satisfiesTaking expectations, the right-hand side vanishes because both functions satisfy $\mathbb{E}[Zh(Z)]=1$. Thus $\mathbb{E}{h}^2\ge \mathbb{E}{h}_K^2$. This proves the optimal bounded influence function for normal location, even over the larger class of all bounded normalized functions, and hence over the required odd nondecreasing class. Its minimum asymptotic variance isThe numerator is $\mathbb{E}[Z^2{1}_{\{|Z|\le K\}}]+K^2\mathrm{Pr}(|Z|>K)$.

At the included boundary $C=\sqrt{\pi /2}$, equality in the feasibility bound forces $h(x)=C\mathrm{sgn}(x)$ almost everywhere. This is the influence function of the sample median, with $V=\pi /2$. It is obtained from the rescaled scores ${\psi }_K/K$ as $K\downarrow 0$. The boundary solution is the sign score, or equivalently the median estimator. Literal substitution of $K=0$ in the clipped score gives the identically zero function and no identifiable estimating equation; the printed clipped formula therefore requires this limiting interpretation at equality.

Let a group $G$ act on the sample space. A statistic $T$ is a maximal invariant when it is invariant and separates group orbits:Together these conditions say that $T(x)=T(y)$ exactly when $x$ and $y$ belong to the same orbit.

The basic result is that invariant tests factor through a maximal invariant. Specifically, a test function $\phi (x)\in [0,1]$ is invariant exactly when it can be written $b(T(x))$. To prove the nontrivial direction, define $b(t)=\phi (x)$ for any representative with $T(x)=t$. Two representatives are in the same orbit by maximality, and invariance gives them the same test value, so $b$ is well-defined. Conversely $b(T(gx))=b(T(x))$ proves invariance. Measurability is taken on the quotient statistic space; for the finite rank spaces used below it is automatic.

If a null family is transitive under the group action, the null law of $T$ is parameter-free. Indeed, if $P_{{\theta }^{\prime }}$ is the law of $gX$ when $X\sim P_{\theta }$, then for every measurable set $B$,Thus invariant procedures use exactly the information in $T$, and its common null law calibrates tests without estimating the transformed nuisance distribution. This explains the role of maximal invariants in nonparametric statistics.

For simultaneous strictly increasing bijections of the real line, the labeled ranks of observations are a maximal invariant on samples without ties. Increasing maps preserve every comparison. Conversely, two such samples with the same ranks can be matched by a piecewise linear increasing bijection through their corresponding ordered observations, extended with positive-slope tails. This proves ranks as a maximal invariant under increasing transformations. Under independent sampling from any common continuous distribution, every labeled ordering is equally likely by exchangeability. Hence rank tests have distribution-free null calibration even without a global transitivity assertion for the entire class of continuous distributions.

Use the usual two mutually independent samples, and rank their $N=n+m$ pooled observations increasingly. Continuity gives no ties with probability one. Let $R_i$ be the pooled rank of $X_i$ and define the Wilcoxon rank sum test statisticThe last identity counts the ranks contributed by the other sample after subtracting the within-$X$ ranks.

When $F_X\ge F_Y$, the $X$ observations are stochastically smaller, so reject for small $W_X$, equivalently small $U_X$. For instance,The inequality follows from $F_X\ge F_Y$ and $\int F_{Y}d{F}_Y=1/2$ for continuous $F_Y$. A common-uniform quantile coupling also makes this direction transparent: $F_X^{-1}(u)\le F_Y^{-1}(u)$, and $U_X$ is nondecreasing in each $X_i$.

Under the common-distribution null hypothesis, the pooled observations are exchangeable. Every selection of $n$ ranks for the $X$ labels has probability ${\left(\genfrac{}{}{0ex}{}{N}{n}\right)}^{-1}$. Thus the exact null law isChoose a lower-tail critical value with probability at most the desired size, optionally randomizing at its boundary to attain the size exactly. This law is independent of the unknown common distribution. Common increasing transformations leave the pooled ranks of observations unchanged, so the statistic is a function of the maximal invariant from part (i). Although a rank sum is not itself maximal, it is an invariant test statistic. These facts justify the test within nonparametric statistics.

Sampling ranks without replacement givesso ${\mathbb{E}}_0{U}_X=nm/2$ with the same variance. For $n,m\to \infty$ with $n/N\to \rho \in (0,1)$, the normal approximation for a rank sum statesEquivalently replace $W_X-n(N+1)/2$ by $U_X-nm/2$. The large-sample one-sided test rejects when this standardized statistic is below the lower $\alpha$ normal distribution quantile. This is the two-sample rank-sum test, rather than the paired Wilcoxon signed-rank test.

An Edgeworth expansion refines a central normal approximation using higher cumulants. Let $Y_i$ be independent and identically distributed, with mean $\mu$, positive variance ${\sigma }^2$, and cumulants ${\kappa }_r$. DefineThe logarithm of its characteristic function has the expansionExpanding the exponential adds the squared cubic term ${\gamma }_3^2(it{)}^6/(72n)$. Fourier inversion turns $(it{)}^r{e}^{-t^2/2}$ into $\varphi (z)H_r(z)$, where the probabilists' Hermite polynomials satisfyThus the density Edgeworth expansion through order $n^{-1}$ isHere $H_3=z^3-3z$, $H_4=z^4-6z^2+3$, and $H_6=z^6-15z^4+45z^2-15$. Since $(\varphi H_{r-1}{)}^{\prime }=-\varphi H_r$, the corresponding cumulative distribution function expansion isThese remainder statements require adequate moments and smooth nonlattice conditions allowing transform inversion and error control, not merely a formal cumulant expansion. The density version needs the corresponding local smoothness. A symmetric distribution has ${\gamma }_3=0$, removing the first correction. Truncation can produce negative densities or nonmonotone distribution functions far from the center; it is not automatically a uniformly reliable relative tail approximation. A lattice distribution needs its lattice correction.

A Laplace approximation instead approximates an integral dominated by a sharp maximum. Suppose $h$ has a unique interior maximizer $\hat{\theta }\in {\mathbb{R}}^d$, with $A=-h^{\prime \prime }(\hat{\theta })$ positive definite, and the integral is localized near that point. For smooth amplitude $g$,Set $\theta =\hat{\theta }+z/\sqrt{n}$. A Taylor expansion givesExpanding $g$ and integrating the resulting Gaussian integral givesUnder sufficient smoothness and localization the relative correction is $O(n^{-1})$: odd terms of order $n^{-1/2}$ integrate to zero over the limiting symmetric Gaussian region. Higher terms use derivatives of $h,g$ and Gaussian moments.

For a likelihood function with ${\ell }_n=nh$ and a smooth positive prior density, this gives the integrated likelihood approximationwhere $j=-{\ell }_n^{\prime \prime }$ is the observed information. Ratios of such integrals give posterior expectation by Laplace approximation; the same method is useful when integrating nuisance variables. Boundary maxima, singular Hessian matrices, nonlocalized tails and several dominant modes require altered formulas, with contributions summed when appropriate.

Both methods start with a quadratic Gaussian approximation and then retain higher-order terms. The Edgeworth series corrects a sampling distribution through its cumulants and Fourier inversion; Laplace's method corrects an integral through local derivatives at its dominant point. Neither should be confused with simply replacing a likelihood by a normalized posterior density.

Consider a regular $d$-parameter statistical model for which the data can locally be expressed as $(v,a)$, where $v=\hat{\theta }$ is the maximum-likelihood estimator and $a$ is a suitable exact or higher-order ancillary statistic. Write the log-likelihood as $\ell (\theta ;v,a)$, and letbe its fitted observed information. The P-star approximation to the conditional sampling density of $v$, given $a$, isThe normalizing factor is independent of $v$ and is chosen so the density integrates to one over the feasible estimator coordinates. Its leading Gaussian value is $(2\pi {)}^{-d/2}$, but exact normalization can change it. Both the fitted likelihood and the fitted information vary with $v$; the observed dataset is not held fixed while integrating over this sampling coordinate.

To see the structure, expand the likelihood difference for $v$ close to $\theta$. Its leading term is minus one half of the information-weighted quadratic displacement, while $|j(v){|}^{1/2}$ supplies the local density scale. This recovers the first-order multivariate normal distribution of an efficient maximum-likelihood estimator. The formula is invariant under smooth one-to-one parameter changes: at a likelihood maximum the score vanishes, so the Hessian transforms as a quadratic form. Its square-root determinant supplies exactly the Jacobian determinant required for the estimator density to transform.

There is a direct derivation in a full regular exponential family. Let $S$ be the sufficient statistic for an independent sample, with likelihoodThe maximum-likelihood equation is $s=n\nabla \kappa (v)$, and its derivative is $j(v)=n{\kappa }^{\prime \prime }(v)$. Applying a multivariate saddlepoint density approximation to $S$ givesChanging variables from $s$ to $v$ multiplies by $|ds/dv|=|j(v)|$, producing the displayed P-star form. In more general models, ancillary conditioning supplies the relevant local sample coordinates; an ancillary is not an arbitrary extra statistic.

For an illustrative exact case, take an independent sample from an exponential distribution with mean $\theta$, so $v=\bar{Y}$. The fitted observed information is $j(v)=n/v^2$, andThus P-star has density shape $v^{n-1}{\theta }^{-n}{e}^{-nv/\theta }$. Normalizing giveswhich is exactly the gamma distribution sampling density of the mean. Here sample ratios provide ancillary coordinates and are independent of the mean.

In suitable regular models the unnormalized leading construction has relative error of order $n^{-1}$, and higher-order ancillary constructions with normalization can attain relative order $n^{-3/2}$. Such accuracy needs the stated smoothness and ancillary hypotheses; it is not a universal consequence of asymptotic normality. The exact scale-model example is likewise a special property, not a general identity. The formula is useful for refined conditional inference, score distributions, and nuisance-likelihood adjustments, including modified profile likelihood. It describes a sampling density of the estimator, rather than a posterior density of the parameter.

A finite-dimensional exponential family has common dominating measure and densityHere $\eta$ is the natural parameter of an exponential family, and normalization determines the cumulant function of an exponential family $\kappa$. A full regular family lets $\eta$ vary over an open natural parameter domain. A curved exponential family restricts $\eta$ to a lower-dimensional submanifold, while a minimal representation removes redundant affine combinations of the statistic.

Differentiating the normalization under regularity givesConsequently $\kappa$ is convex and its Hessian is the Fisher information matrix in natural coordinates. For an independent size-$n$ sample, $S={\sum }_{i}T(Y_i)$ is a sufficient statistic: the joint density factors into ${\prod }_{i}h(y_i)$ times $\exp \{{\eta }^{\mathsf{T}}S-n\kappa (\eta )\}$, so the conditional law given $S$ loses its parameter dependence. An interior maximum-likelihood estimator solvesmatching empirical and theoretical sufficient-statistic means. Under minimality the Hessian is positive definite, so the log-likelihood is strictly concave; existence of an interior maximizer still needs the observed statistic to be in the appropriate interior convex support.

For a full open natural family, the canonical sufficient statistic is a complete statistic: if ${\mathbb{E}}_{\eta }g(S)=0$ throughout the open domain, exponential tilting about an interior point makes the Laplace transform of the integrable signed measure $g(s)P_{{\eta }_0}(ds)$ zero on a neighborhood. Transform uniqueness makes that measure zero, so $g(S)=0$ almost surely. This requires integrability under the family. It is not an automatic property of curved restrictions. The sufficient-statistic reduction, explicit information and cumulant identities make these families especially useful for exact and higher-order inference. Normal distributions, Poisson distributions and gamma distributions supply familiar examples.

A transformation model is built differently: a group acts on the sample and parameter spaces, and its probability laws satisfyIn a transitive model all laws are transforms of a fixed reference law. The location family is generated by translations, with density $f_{\mu }(y)=f_0(y-\mu )$. The positive-scale family has density $f_{\sigma }(y)={\sigma }^{-1}{f}_0(y/\sigma )$, and the location-scale family combines them asThe Jacobian factor is essential. For independent observations its power is ${\sigma }^{-n}$. The base density need not have moments, and a transformation model need not have fixed support.

An equivariant estimator respects the parameter action. With a well-defined equivariant fit $T(y)$ valued in the transformation group, the normalized sampleis a maximal invariant. Equivariance gives $a(gy)=a(y)$; if normalized samples agree, $y^{\prime }$ is obtained from $y$ by $T(y^{\prime })T(y{)}^{-1}$. In a transitive transformation model its law is parameter-free, so it is an ancillary statistic. For location and positive scale, centered standardized residuals provide this construction on nondegenerate samples. If the parameter action has a nontrivial stabilizer, normalization must retain its quotient action rather than pretending that a unique group element is determined.

This geometry motivates conditioning on an ancillary and using the P-star approximation for the fitted parameter. Ancillarity means parameter-free distribution; it does not by itself assert independence from the estimator. In a normal location model with known scale, the sample mean and centered residuals are independent by their joint Gaussian law and zero covariance. In an exponential scale model, the sample sum and its normalized proportions are independent, again an additional model property.

The two classes overlap but are distinct. A normal location-scale family is both an exponential family and a transformation model; an exponential mean family is both a one-dimensional natural exponential family and a scale model. A general Cauchy distribution location family is a transformation model without a fixed finite-dimensional regular exponential-family sufficient-statistic reduction. Exponential-family structure emphasizes sufficient statistics and natural-parameter convexity; transformation structure emphasizes equivariance, group orbits and ancillary reduction.