Solution (source code)

= Solution

Use the unscaled penalty convention consistent with the displayed bound:
$$
\boxed{\widehat f_\lambda=\arg\min_{f\in\mathcal H}\left\{\sum_{i=1}^n(Y_i-f(x_i))^2+\lambda\|f\|_{\mathcal H}^2\right\}.}
$$
For an objective with averaged squared loss, the corresponding penalty parameter would be $\lambda/n$. We retain the unscaled $\lambda$ throughout.

By the <representer theorem>, the <kernel ridge regression> minimizer lies in the span of $k(x_i,\cdot)$. A valid coefficient vector is $\widehat a=(K+\lambda I_n)^{-1}Y$, giving
$$
\widehat f_\lambda(\cdot)=\sum_i\widehat a_i k(x_i,\cdot),\qquad\widehat m=HY,\qquad H=K(K+\lambda I_n)^{-1}.
$$
The <kernel matrix> is a <positive semidefinite matrix>, so the inverse exists even if $K$ is singular. In the coefficient objective $\|Y-Ka\|_2^2+\lambda a^TKa$, this choice satisfies $K\{(K+\lambda I_n)a-Y\}=0$; there is no need to cancel $K$. The function itself is unique because of the strictly convex squared Hilbert norm. The matrix $H$ is the <kernel-ridge hat matrix> in the present normalization.

Let $m^0=(f^0(x_i))_{i=1}^n$. Since $Y=m^0+\varepsilon$, the zero mean and covariance assumption give the exact <bias-variance decomposition for linear prediction>
$$
\frac1n\mathbb E\|\widehat m-m^0\|_2^2=\frac1n\|(I_n-H)m^0\|_2^2+\frac{\sigma^2}{n}\operatorname{tr}(H^2).
$$
No normality of the noise is required for this identity.

Choose an <orthonormal eigenbasis> $v_i$ of $K$, with eigenvalues $d_i\geq0$. Its shrinkage factors are $d_i/(d_i+\lambda)$, so
$$
\operatorname{tr}(H^2)=\sum_i\frac{d_i^2}{(d_i+\lambda)^2}\leq\sum_i\min\left\{\frac{d_i}{4\lambda},1\right\}=\frac1\lambda\sum_i\min(d_i/4,\lambda).
$$
The first bound follows from $(d_i+\lambda)^2\geq4d_i\lambda$ and from each shrinkage factor being at most one.

To control the bias by the <Reproducing kernel Hilbert space> norm, for $d_i>0$ define
$$
g_i=\frac1{\sqrt{d_i}}\sum_j(v_i)_j k(x_j,\cdot)\in\mathcal H.
$$
The <reproducing property> gives $\langle g_i,g_\ell\rangle_{\mathcal H}=\delta_{i\ell}$ and
$$
v_i^Tm^0=\sqrt{d_i}\langle g_i,f^0\rangle_{\mathcal H}.
$$
If $d_i=0$, the corresponding unscaled kernel combination has Hilbert norm squared $v_i^TKv_i=0$, so $v_i^Tm^0=0$. There is therefore no bias component in a zero-eigenvalue direction. For the remaining directions,
$$
\begin{aligned}
\|(I_n-H)m^0\|_2^2
&=\sum_{d_i>0}\frac{\lambda^2d_i}{(d_i+\lambda)^2}|\langle g_i,f^0\rangle_{\mathcal H}|^2\\
&\leq\frac\lambda4\sum_{d_i>0}|\langle g_i,f^0\rangle_{\mathcal H}|^2\\
&\leq\frac\lambda4\|f^0\|_{\mathcal H}^2.
\end{aligned}
$$
The scalar bound is again $(d_i+\lambda)^2\geq4d_i\lambda$, and the final step is <Bessel inequality>. Combining the variance and bias bounds proves the <prediction-risk bound for kernel ridge regression>:
$$
\boxed{\frac1n\mathbb E\sum_i\{f^0(x_i)-\widehat f_\lambda(x_i)\}^2\leq\frac{\sigma^2}{n\lambda}\sum_i\min(d_i/4,\lambda)+\frac\lambda{4n}\|f^0\|_{\mathcal H}^2.}
$$