---
title: 'Quadratic Shrinkage: Theory & Applications'
url: https://www.emergentmind.com/topics/quadratic-shrinkage
type: topic
---

# Quadratic Shrinkage: Theory & Applications

Quadratic shrinkage denotes a family of estimators that modify a baseline estimator by shrinking it toward zero, a fixed target, a subspace, or a structured prior in order to reduce risk under a quadratic criterion. In the canonical multivariate normal mean model \(X\sim N_p(\mu,\sigma^2I_p)\) with loss \(L(\mu,\delta)=\|\delta-\mu\|^2\), the usual estimator \(X\) is minimax but, for \(p\ge3\), inadmissible; this observation motivates James–Stein rules, Bayes formulations based on superharmonic marginals, and a large body of extensions to predictive densities, matrix-valued estimators, quadratic-penalty methods, and high-dimensional models under general quadratic loss [1203.5617][2205.13840].

## 1. Classical normal-mean formulation

The classical setup observes
\[
X\sim N_p(\mu,\sigma^2I_p),
\]
and seeks an estimator \(\delta(X)\) of \(\mu\) under quadratic loss
\[
L(\mu,\delta)=\|\delta-\mu\|^2,
\qquad
R(\mu,\delta)=E_\mu\|\delta(X)-\mu\|^2.
\]
The usual estimator \(\delta_0(X)=X\) has constant risk
\[
R(\mu,X)=p\sigma^2,
\]
so its maximum risk over \(\mu\) equals \(p\sigma^2\), and it is minimax [1203.5617].

For \(p\ge3\), Stein (1956) showed that \(X\) is inadmissible. James and Stein (1961) proposed
\[
\delta_{JS}(X)=\left(1-\frac{(p-2)\sigma^2}{\|X\|^2}\right)X,
\]
with risk
\[
R(\mu,\delta_{JS})
=
p\sigma^2
-
E_\mu\!\left[\frac{(p-2)^2\sigma^4}{\|X\|^2}\right]
<
p\sigma^2
\qquad
\text{for all }\mu,\; p\ge3.
\]
Hence \(\delta_{JS}\) is minimax and strictly better than \(X\). If \(\sigma^2\) is unknown, one may replace \(\sigma^2\) by an estimate.

The essential feature is a data-dependent multiplicative factor smaller than \(1\). Risk reduction is obtained not by unbiasedness but by a bias–variance trade-off that lowers quadratic risk uniformly over \(\mu\).

## 2. Bayes representation and superharmonicity

The Bayes formulation makes shrinkage explicit. Under the conjugate Gaussian prior
\[
\mu\sim N_p(0,\nu I_p),
\]
the posterior of \(\mu\) given \(X\) is
\[
N_p\!\left(\frac{\nu}{\nu+\sigma^2}X,\frac{\sigma^2\nu}{\nu+\sigma^2}I_p\right),
\]
and the Bayes rule is the posterior mean
\[
\delta_B(X)
=
\left(1-\frac{\sigma^2}{\nu+\sigma^2}\right)X.
\]
As \(\nu\to\infty\), this tends to \(X\); as \(\nu\to0\), it tends to \(0\). One may choose \(\nu\) to optimize Bayes risk or replace the shrinkage factor by its unbiased estimate in an empirical-Bayes construction [1203.5617].

For a general prior \(\pi(\mu)\), possibly improper, let
\[
m_\pi(x)=\int e^{-\|x-\mu\|^2/(2\sigma^2)}\,\pi(\mu)\,d\mu.
\]
Brown (1971) and Stein (1974) showed that the Bayes-posterior mean can be written as
\[
\delta_\pi(x)=E_\pi[\mu\mid x]=x+\sigma^2\nabla_x\log m_\pi(x).
\]
Substituting this representation into the risk and using integration by parts yields
\[
R(\mu,X)-R(\mu,\delta_\pi)
=
4\sigma^4E_\mu\!\left[-\frac{\Delta\sqrt{m_\pi(X)}}{\sqrt{m_\pi(X)}}\right].
\]
A sufficient condition for minimaxity is either
\[
\Delta m_\pi(x)\le0
\quad\text{for all }x,
\]
so that \(m_\pi\) is superharmonic, or the weaker condition
\[
\Delta\sqrt{m_\pi(x)}\le0
\quad\text{for all }x.
\]

This criterion organizes a broad class of shrinkage priors. The improper harmonic prior
\[
\pi_H(\mu)\propto \|\mu\|^{-(p-2)}
\]
has marginal
\[
m_H(x)\propto \|x\|^{-(p-2)},
\]
which is superharmonic for \(p\ge3\), so the corresponding Bayes rule \(\delta_H\) dominates \(X\). Strawderman’s proper priors \(\pi_a\), obtained as scale-mixtures of Gaussians with
\[
\mu\mid s\sim N(0,sI),
\qquad
s\propto(1+s)^{a-2},
\qquad
0\le a<2,
\]
provide examples for which \(\sqrt{m_a(x)}\) is superharmonic even though \(m_a\) is not.

The same device yields multiple shrinkage estimators. Given targets or subspaces \(B_1,\dots,B_N\), one forms recentered superharmonic marginals \(m_i(x)\) and then a mixture marginal
\[
m_*(x)=\sum_{i=1}^N w_i\,m_i(x),
\]
which remains superharmonic. The resulting Bayes rule
\[
\delta_*(x)=x+\sigma^2\nabla\log m_*(x)
\]
is minimax and adapts to whichever target is most supported by the data.

An empirical-Bayes analogue of the conjugate-normal rule uses the unbiased estimate \((p-2)\sigma^2/\|X\|^2\) of the shrinkage factor and truncates it at zero, producing the positive-part estimator
\[
\delta_{JS}^+(X)=\left(1-\frac{(p-2)\sigma^2}{\|X\|^2}\right)_+X,
\]
described as minimax and admissible.

## 3. Second-order and higher-order polynomial shrinkage

A different use of the adjective “quadratic” appears in polynomial shrinkage estimators of a multivariate normal mean. With
\[
X\sim N_p(\theta,I_p),
\qquad
U=\|X\|^2,
\]
and balanced squared-error loss
\[
L_\omega(\delta,\theta)=\omega\|\delta-X\|^2+(1-\omega)\|\delta-\theta\|^2,
\]
one considers estimators of the form
\[
\delta(X)=(1-\phi(U))X.
\]
The second-degree, or quadratic, shrinkage function is
\[
\phi(U)=aU^{-1}-bU^{-2},
\]
so that
\[
\delta_b^{(2)}(X)=\left(1-\frac{a}{U}+\frac{b}{U^2}\right)X.
\]
In the balanced-loss setting the optimal first- and second-order coefficients are
\[
a=(1-\omega)(p-2),
\qquad
b=2(1-\omega)(p-6).
\]
The estimator is minimax relative to the MLE \(X\) if and only if
\[
E[\phi(U)^2U]\le 2(1-\omega)(p-2)E[\phi(U)].
\]
When \(\phi(U)=a/U\), this reduces to
\[
0\le a\le 2(1-\omega)(p-2)
\]
[2107.14021].

The explicit risk formula for the quadratic estimator is
\[
R_\omega(\delta_b^{(2)},\theta)
=
R_\omega(\delta_{JS},\theta)
-4b(1-\omega)E[U^{-2}]
+b^2E[U^{-3}].
\]
Choosing
\[
\hat b=2(1-\omega)(p-6)
\]
gives
\[
R_\omega(\delta_{\hat b}^{(2)},\theta)
=
(1-\omega)p
-(p-2)^2(1-\omega)^2E[U^{-1}]
-4(1-\omega)^2(p-6)E[U^{-2}]
\le
R_\omega(\delta_{JS},\theta)
\le
R_\omega(X,\theta).
\]
For \(p>6\), the quadratic rule therefore strictly dominates both the usual James–Stein estimator and the MLE.

The paper’s figures plot the risk ratios
\[
R_\omega(\delta_{JS})/R_\omega(X),
\quad
R_\omega(\delta_b^{(2)})/R_\omega(X),
\quad
R_\omega(\delta_c^{(3)})/R_\omega(X)
\]
against \(\lambda=\|\theta\|^2\) for \(p=8,12,14,18\) and \(\omega=0.1,0.4\); all curves lie below \(1\), and each higher-order polynomial pushes the ratio down further, especially when \(\lambda\) and \(\omega\) are small. The same Stein-lemma argument extends to third- and higher-order rules with
\[
\phi(U)=a/U-b/U^2+c/U^3-\cdots,
\]
but one needs \(p>2m+2\) to control \(E[U^{-j}]\) and guarantee domination. The practical gain diminishes because \(E[U^{-j}]\to0\) for large \(\lambda\) and high \(j\), and the complexity and variance of estimating higher-order inverse moments become unattractive.

## 4. Predictive density, matrix, and singular-value generalizations

Quadratic shrinkage extends beyond point estimation of a single normal mean. In the normal linear model,
\[
X\mid\beta\sim N_m(A\beta,\sigma^2I_m),
\qquad
Y\mid\beta\sim N_n(B\beta,\sigma^2I_n),
\]
independent, the Bayes predictive density under \(\pi(\beta)\) can be written in the shrinkage form
\[
p_\pi(y\mid x)=\frac{m_\pi(w;v_w)}{m_\pi(x;v_x)}\,p_U(y\mid x),
\]
where \(p_U\) is the uniform-prior predictor. Diagonalizing the covariances reduces the problem to the canonical mean setting, and the same superharmonic-marginal conditions ensure minimaxity. In nonparametric regression, the Gaussian sequence model
\[
y_i=\theta_i+e_i,
\qquad
e_i\sim N(0,1/n),
\qquad
i=1,\dots,n,
\]
with \(\theta\) constrained in an ellipsoid, has minimax KL-risk asymptotically attained by Bayes predictive densities under Gaussian priors, in direct analogy with Pinsker’s theorem for \(L_2\)-risk [1203.5617].

A matrix-valued version arises when estimating an \(n\times p\) mean matrix \(\Theta\) from an \(n\times p\) normal data matrix \(X\) under the matrix quadratic loss
\[
L_{\mathrm{matrix}}(\widehat\Theta,\Theta)
=
(\widehat\Theta-\Theta)^T(\widehat\Theta-\Theta).
\]
The usual MLE is \(\widehat\Theta_{\mathrm{MLE}}=X\), with risk
\[
R_{\mathrm{matrix}}(X,\Theta)=n\sigma^2I_p.
\]
Applying vector James–Stein shrinkage separately to each column gives
\[
\widehat\theta_{a,(j)}
=
x_{(j)}-a\,g_{(j)}(x_{(j)}),
\qquad
g_{(j)}(x_{(j)})
=
\sigma^2\frac{n-2}{\|x_{(j)}\|^2}x_{(j)},
\]
or equivalently
\[
\widehat\Theta_a=XD,
\qquad
d_j=1-a\,\frac{\sigma^2(n-2)}{\|x_{(j)}\|^2}.
\]
Abu-Shanab, Kent and Strawderman showed that if the columnwise shrinkers satisfy the vector cross-product inequality, then \(\widehat\Theta_a\) strictly dominates \(X\) under the matrix loss for
\[
0<a<\frac{2}{p}.
\]
The same argument extends to positive-part James–Stein, empirical Bayes, subspace shrinkage, and superharmonic-prior rules, and also to unknown \(\sigma^2\) and correlated columns after whitening [1101.3412].

General quadratic loss can also be encoded by a positive-definite matrix \(Q\succ0\). In the multivariate Gaussian sequence and matrix-denoising settings,
\[
L_Q(\hat\theta,\theta)
=
\sum_{i=1}^\infty
(\hat\theta_i-\theta_i)^\top Q(\hat\theta_i-\theta_i),
\]
and
\[
\ell_Q(\hat M,M)
=
\mathrm{tr}\bigl((\hat M-M)^\top Q(\hat M-M)\bigr).
\]
The Efron–Morris singular-value shrinker
\[
\hat M_{\rm EM}
=
X\Bigl(I_p-(n-p-1)(X^\top X)^{-1}\Bigr)
\]
shrinks each singular value and dominates the MLE under Frobenius loss. Under arbitrary \(Q\)-loss it satisfies the oracle inequality
\[
R_Q(M;\hat M_{\rm EM})
\le
\min_C R_Q(XC,M)+2(p+1)\,\mathrm{tr}(Q).
\]
This leads to a multivariate Pinsker linear shrinker on Sobolev ellipsoids and to a blockwise Efron–Morris construction that is stated to be exactly adaptive minimax over multivariate Sobolev ellipsoids under the corresponding quadratic loss [2205.13840].

## 5. Quadratic penalties and structured regression shrinkage

Quadratic shrinkage may be implemented as a penalty term rather than as a multiplicative factor. In sparse Laplacian shrinkage for high-dimensional regression, one minimizes
\[
\mathcal L(\beta)
=
\frac{1}{2n}\|Y-X\beta\|_2^2
+
\sum_{j=1}^p \rho(|\beta_j|;\lambda_1,\gamma)
+
\frac{\lambda_2}{2}\beta^TL\beta,
\]
where \(\rho\) is the minimax concave penalty and \(L=D-A\) is the graph Laplacian derived from a signed adjacency matrix. The quadratic form is
\[
\beta^TL\beta
=
\sum_{j<k}|a_{jk}|
\bigl(\beta_j-\mathrm{sgn}(a_{jk})\beta_k\bigr)^2.
\]
This term promotes local smoothness among coefficients associated with connected predictors; unlike a global ridge penalty, the shrinkage is local on the graph. Under the stated sub-Gaussian, restricted-eigenvalue, penalty-size, and minimum-signal conditions, the estimator is selection-consistent and equal to the oracle Laplacian shrinkage estimator with probability at least \(1-3\delta\). Proposition 5.2 further gives a generalized grouping property within graph cliques [1112.3450].

A related quadratic penalization appears in shrinkage for categorical regressors and group means. If \(\bar Y_g\) denotes the sample mean in group \(g\), shrinkage toward first-stage targets \(\delta_g^0\) is obtained by minimizing
\[
L(\delta;\lambda)
=
\sum_{g=1}^G n_g(\bar Y_g-\delta_g)^2
+
\lambda\sum_{g=1}^G w_g(\delta_g-\delta_g^0)^2.
\]
Writing \(\Lambda=\mathrm{diag}(n_1,\dots,n_G)\) and \(W=\mathrm{diag}(w_1,\dots,w_G)\), the unique minimizer is
\[
\hat\delta(\lambda)
=
(\Lambda+\lambda W)^{-1}
\bigl[\Lambda\bar Y+\lambda W\delta^0\bigr].
\]
Under the local-to-zero framework
\[
\sqrt n(\mu-\delta^0)\to d\in\mathbb R^G,
\]
shrinkage introduces asymptotic bias
\[
-\lambda_0(P+\lambda_0W)^{-1}Wd
\]
but reduces variance. For every bounded \(d\in\mathbb R^G\), provided \(G\ge4\), the plug-in shrinkage estimator has asymptotic risk strictly below OLS uniformly over \(\|d\|\le M\). In Monte Carlo experiments with \(G=4\) and \(n=400\), the plug-in quadratic-shrinkage estimator never exceeds OLS risk and often improves by up to \(30\%\) [1901.01898].

## 6. High-dimensional, spectral, and heteroscedastic extensions

In high-dimensional mean estimation with unknown covariance, Wang et al. consider i.i.d. observations
\[
X_i=\Sigma_p^{1/2}\epsilon_i+\mu,
\qquad i=1,\dots,n,
\]
with \(E(\epsilon_i)=0\), \(\mathrm{Cov}(\epsilon_i)=I_p\), finite fourth moments, and no normality assumption. Accuracy is measured by
\[
L_A(\delta,\mu)=(\delta-\mu)'A(\delta-\mu),
\]
where \(A\) is known positive definite. Their estimator
\[
\delta_{\alpha,\beta}=\alpha\bar X+\beta e
\]
shrinks the sample mean toward \(e=(1,\dots,1)'\), with data-driven coefficients \(\hat\alpha^*,\hat\beta^*\) formed from four \(U\)-statistics \(Y_{1,n},\dots,Y_{4,n}\). Under Assumption 3.1 and the asymptotic condition \(n\,p^{-2}\,\|\mu\|^2\to0\), one has
\[
\hat\alpha^*-\alpha^*=o_p(1),
\qquad
\hat\beta^*-\beta^*=o_p(1).
\]
For the percentage relative improvement
\[
\mathrm{PRIAL}
=
\frac{E[(\bar X-\mu)'Q(\bar X-\mu)]-E[(\hat\mu^*-\mu)'Q(\hat\mu^*-\mu)]}
{E[(\bar X-\mu)'Q(\bar X-\mu)]},
\]
Corollary 3.1 states that if \(s_n=n\pi_2/p\to0\), then \(\mathrm{PRIAL}\to1\); if \(s_n\to C\in(0,\infty)\), then \(\mathrm{PRIAL}\to\) some \(C_1\in(0,1)\); and if \(s_n\to\infty\), then \(\mathrm{PRIAL}\to0\). Simulations and leukemia microarray analysis are reported to show practical improvement over the sample mean and competing shrinkage rules, especially when \(\Sigma_p\) is non-diagonal and training sizes are very small [1211.1456].

In spectral matrix estimation for partial coherencies, one observes
\[
\hat S\sim W_p^C(K,\Sigma)/K,
\]
and a standard shrinkage estimator is
\[
\hat S(\rho)=(1-\rho)\hat S+\rho T,
\qquad
T=(\mathrm{tr}\,\Sigma/p)I_p.
\]
Schneider-Luftman and Walden study the quadratic-loss criterion
\[
R_{QL}(\rho)=E\,\mathrm{tr}\{[\hat S(\rho)^{-1}-I_p]^2\},
\]
derive oracle rules QLa and QLb for spectral-matrix shrinkage, and also formulate precision-matrix shrinkage
\[
\hat P(\alpha,\beta)=\alpha\hat S^{-1}+\beta I_p.
\]
A central methodological point is that Hilbert–Schmidt shrinkage can “over-shrink” and wipe out true large partial coherencies. In EEG-derived simulations, HS-based shrinkage can increase error when \(K\) only slightly exceeds \(p\), whereas oracle QLb gives the best PRISE, \(68\%\) for \(K=12\). Among the feasible methods, QLP-est retains almost all of the oracle QLP gain, \(52\%\to52\%\) at \(K=12\) and \(34\%\to32\%\) at \(K=14\), with much less variance across frequency; it is recommended as the full-estimation method [1511.07030].

A further generalization treats families of distributions with quadratic variance function
\[
\mathrm{Var}(X_i)=A\mu_i^2+B\mu_i+C.
\]
Estimators
\[
\hat\mu_i=(1-\lambda_i)X_i+\lambda_i m
\]
then have oracle weight
\[
\lambda_i^*
=
\frac{\mathrm{Var}(X_i)}{\mathrm{Var}(X_i)+(m-\mu_i)^2},
\]
or, when a sample-size factor \(\tau_i\) appears,
\[
\lambda_i^*
=
\frac{(A\mu_i^2+B\mu_i+C)/\tau_i}{(A\mu_i^2+B\mu_i+C)/\tau_i+(m-\mu_i)^2}.
\]
Xie, Kou and Brown construct semiparametric URE-based shrinkage estimators
\[
\hat\mu_i(\mathbf b,m)=(1-b_i)X_i+b_im,
\]
with monotonicity \(b_i\le b_j\) whenever \(\mathrm{Var}(X_i)\le\mathrm{Var}(X_j)\), and show that the URE rule is asymptotically first-order risk-optimal in its class. Their simulations report \(20\)–\(50\%\) lower MSE than classical James–Stein and empirical-Bayes estimators when variances and means are correlated or heteroscedastic [1603.09144].

Across these formulations, quadratic shrinkage is unified less by a single estimator than by a common risk logic: a controlled introduction of bias to secure lower quadratic risk. The resulting theory is dimension-sensitive and model-sensitive. James–Stein dominance requires \(p\ge3\); second-order polynomial improvement requires \(p>6\); higher-order polynomial rules require \(p>2m+2\); and, in spectral estimation, methods based on inverse-power traces may be unstable even when their oracle versions are attractive. The main analytical tools are minimaxity, Bayes posterior means, superharmonicity, oracle inequalities, and explicit bias–variance decompositions.

Source: https://www.emergentmind.com/topics/quadratic-shrinkage