Past exam of the mathematics course of the University of Cambridge / 2022 / iii / Paper 218 (Statistical Learning in Practice)

Put $\varphi =1/\gamma$ and $\theta =-1/\alpha$. The density can be writtenwhereThus this is an exponential dispersion family. Its derivative identities giveThe variance function is $V(\mu )={\mu }^2$. In the usual Gamma-GLM convention the canonical inverse link isit differs by a minus sign from the natural parameter $\theta =-1/\mu$.

Let $x_i^T$ be row $i$ of the design matrix, ${\eta }_i=x_i^T\beta =1/{\mu }_i$, and $\gamma =1/\varphi$. The full log-likelihood isDifferentiation gives the score functionandThe Hessian does not depend on $Y$, so the Fisher information matrix is

The score equation $X^T(\mu -Y)=0$ is independent of $\gamma$, because $\gamma >0$ is only a common factor. Fixing $\varphi =1$ or estimating it therefore gives the same maximum-likelihood estimator $\hat{\beta }$.

Under the usual full-rank and regularity conditions,with $W$ evaluated consistently at the fitted means. Consequently every standard error from "mod2" is $\sqrt{\hat{\varphi }}$ times the corresponding standard error computed with dispersion one in "mod1". In particular,

The null hypothesis isagainst the alternative that at least one coefficient is nonzero. The code uses the scaled reduction in exponential-family deviance,and compares it with ${\chi }_3^2$.

That chi-squared calibration is appropriate when the dispersion is known, and is asymptotically valid after consistent dispersion estimation. With unknown Gamma dispersion, the standard finite-sample GLM comparison instead usesagainst an $F_{3,57}$ distribution. Its p-value is approximately $0.034$, so the correctly calibrated test still rejects at the 5% level and gives evidence that component type or position contributes to mean failure time.

First recode the labels as $y_i=2Y_i-1\in \{-1,1\}$; with labels $0$ and $1$, the displayed constraints for class zero can never hold. Add an unpenalized intercept $b$ so that the separating hyperplane need not pass through the origin, and solveCentering and scaling predictor columns is also advisable because the Euclidean penalty depends on their units.

The support vector machine classifier isUnder the constraints, the two classes lie beyond the parallel support-vector-machine margin boundaries $b+x^T\beta =\pm 1$. Their distance is $2/\Vert \beta {\Vert }_2$, so minimizing the norm maximizes the geometric margin. The observations touching the margin are the support vectors.

The hard-margin constraints are feasible exactly when the classes are separated by a separating hyperplane. Failure after the corrections in part a therefore means the classes overlap.

Introduce slack variables of a support vector machine and solve the soft-margin problemsubject toLarge $C$ strongly penalizes violations and approaches the hard-margin solution when separation is possible. Small $C$ tolerates more violations in exchange for a wider, more strongly regularized margin.

The kernel trick replaces inner products $\varphi (x_i{)}^T\varphi (x_j)$ by evaluations of a positive-semidefinite kernel $k(x_i,x_j)$ without explicitly constructing the feature vectors. In the corresponding reproducing-kernel Hilbert space, the hard-margin problem isEquivalently, its dual isThere is no upper bound ${\alpha }_i\le C$ because this is the hard-margin, rather than soft-margin, problem.

The successful hard-margin fit proves that the transformed observations are separable. This creates complete separation in unpenalized logistic regression: scaling a separating coefficient vector continually raises the likelihood, so no finite maximum-likelihood estimate exists. The nearly singular observed information produces the enormous reported standard errors.

A ridge-penalized logistic regression, or equivalently a Bayesian logistic model with a proper Gaussian prior, gives finite, stable coefficients. Firth bias reduction is another standard remedy.

For word-count vector $x$, the fitted logistic model isHolding other counts fixed, one additional occurrence of "dollar" multiplies the fitted spam odds by $e^{1.859318}$.

The logistic classifier isThe coefficients maximize the Bernoulli logistic-regression model likelihood over the training emails.

The decision boundary is the hyperplaneAt threshold $q$ the classifier predicts spam whenso $q=1/2$ gives the original classifier. Varying $q$ translates the boundary parallel to itself: raising $q$ shrinks the region classified as spam and generally trades fewer false positives for more false negatives.

Writing ${\eta }_i={\beta }_0+x_i^T\beta$, the log-likelihood isWith 500 observations and 5000 word predictors, the augmented design matrix cannot have full column rank. There are nonzero coefficient directions that leave every ${\eta }_i$ unchanged, so any maximizer belongs to an affine family and is not unique. In addition, the high-dimensional features may completely separate the classes, in which case no finite maximizer exists at all.

Principal component analysis centers the word-count vectors, diagonalizes their sample covariance matrix, and projects onto the eigenvectors with the largest eigenvalues. Fitting logistic regression to $d<n$ principal-component scores removes exact collinearity and yields a lower-dimensional design.

Each principal component is an unsupervised linear combination of many words chosen to explain predictor variance; it is not a word selected for association with spam. The dimension can be chosen from a scree plot or cumulative explained variance, but cross-validating the downstream classification loss better targets prediction.

For fixed $A$, differentiating under the constraint ${\sum }_i{z}_i=0$ givesBecause the columns $u_j$ are orthogonal,The fitted value is therefore $\overline{X}+{\Pi }_A(X_i-\overline{X})$, where ${\Pi }_A$ is the orthogonal projection onto the column space of $A$.

The residual sum of squares is minimized by choosing this column space to be the span of the $d$ leading eigenvectors ofThus $\hat{A}$ may be taken to have those orthonormal eigenvectors as columns. Their nonzero scales are immaterial because inverse scaling of the scores leaves $A{z}_i$ unchanged.

For student $i$ in school $j$, "lme1" is the random-intercept linear mixed modelwhereindependently. The estimates areThe fixed intercept is the population-average expected post-study score for an untreated reference student with PTHK zero. The random intercept $b_j$ is school $j$'s deviation from that population intercept.

Treatment was assigned by school, so observations within one school are clustered data and plausibly correlated. A random intercept models shared unobserved school characteristics and prevents the effective amount of treatment-level information from being exaggerated by treating all students as independent.

The fitted school standard deviation is $0.209$, so the estimated between-school variation is nonzero. Accounting for it raises the TV and SC standard errors from about $0.065$ in "lm1" to about $0.105$ in "lme1", reflecting that only 28 independently randomized schools identify those effects.

Conditionally on the school effect, the response variance is ${\hat{\sigma }}^2=1.6531$. Marginally,and two distinct students in the same school have estimated covariance $0.0437$.

There is no single coefficient estimate for the random effect because the 28 school deviations are modeled as draws from a mean-zero distribution. The summary reports the estimated distributional variance; individual empirical Bayes predictions can be extracted separately.

Model "lme2" replaces $b_j$ by a random intercept and random PTHK slope:with a fitted bivariate normal covariance matrix for $(b_{0j},b_{1j})$. This adds a random-slope variance and an intercept-slope covariance.

The likelihood-ratio statistic isAn ordinary chi-squared reference is unreliable because the null random-slope variance is on the boundary and its correlation is unidentified there; the reported correlation of one also signals a nearly singular fit. A valid practical test is a parametric bootstrap: simulate many datasets from fitted "lme1", refit both models by maximum likelihood to each, recompute the likelihood-ratio statistic, and estimate the p-value by the fraction at least $8.4$. The AIC favors "lme2" slightly, whereas its BIC is larger, so the descriptive criteria do not agree.

For a model with $k$ fitted parameters and maximized likelihood $L(\hat{\theta })$, the Akaike information criterion isBackward selection starts with the full model, deletes the single variable giving the smallest AIC when that AIC is lower than the current value, and repeats until no deletion improves it. The first deletion is "indus", which gives AIC $1597.8$.

Adding variables can reduce training bias and increase the maximized likelihood, but it also increases estimation variance and optimism. The $2k$ penalty estimates this optimism, so minimizing AIC implements a bias-variance tradeoff aimed at expected out-of-sample Kullback–Leibler performance.

For $Y\sim N_n(X\beta ,{\sigma }^{2}I)$,

At the maximum-likelihood estimates,soThere are $p$ regression coefficients and one variance parameter. Adding twice this parameter count gives

The common terms cancel, andThis is negative exactly when

There are three regression coefficients and $97$ residual degrees of freedom, so $n=100$. The reported residual standard error uses the unbiased divisor:Hence the maximum-likelihood variance estimate is ${\hat{\sigma }}^2=\mathrm{RSS}/100$. Substituting this, $n=100$, and the four fitted parameters—three coefficients plus ${\sigma }^2$—into the formula in part c(ii) gives the AIC.

The prior isThe software samples from the Bayesian posterior; the displayed "Estimate" values are posterior means computed from the Markov-chain Monte Carlo draws.

With unit error variance, the Gaussian likelihood has log-likelihoodFull column rank gives $\hat{\beta }=(X^{T}X{)}^{-1}{X}^{T}Y$ and the least-squares identityMultiplying the likelihood by ${\prod }_j\varphi ({\beta }_j)$ and absorbing all terms independent of $\beta$ into $C$ gives

For a Laplace prior $\varphi (u)\propto e^{-\lambda |u|}$, maximizing the posterior is equivalent to minimizingso the posterior mode is the Lasso. For a Gaussian prior $\varphi (u)\propto e^{-\lambda u^2/2}$, the mode minimizesand is the ridge regression estimator $(X^{T}X+\lambda I{)}^{-1}{X}^{T}Y$. Gaussian conjugacy makes the posterior normal, so its mode and posterior mean coincide at this ridge estimate.

The nearly identical predictor columns create severe multicollinearity. Their sum is well identified, as shown by the excellent fitted response, but their difference is weakly identified, allowing ordinary least squares to choose huge opposite coefficients with huge standard errors. The independent $N(0,1)$ priors impose the ridge penalty from part c, shrinking that unstable difference toward zero and producing the stable estimates near $(1,1)$.

The reported 95% credible interval for ${\beta }_1$ isConditional on the model, prior, and observed data, its posterior probability is 0.95. A frequentist 95% confidence interval instead has 95% coverage under repeated sampling before the data are observed; it does not assign a sampling probability to the fixed parameter after observing this dataset.