How $\hat{y}_0=x_0^T\beta+\sum^N_{i=1}l_i(x_0)\epsilon_i$? I am currently reading "Elements of Statistical Learning II" and I am not quite sure about one thing in section 2.5 on p.24 and p.26
So at the end of p.24 they write the following:
Suppose that we know that the relationship between $Y$ and $X$ is linear,
$$Y=X^T\beta + \epsilon$$
where $\epsilon \sim N(0,\sigma^2)$ and we fit the model by least squares to the training data. For an arbitrary test point $x_0$, we have $\hat{y}_0=x_0^T\hat{\beta}$, which can be written as $\hat{y}_0=x_0^T\beta+\sum^N_{i=1}l_i(x_0)\epsilon_i$, where $l_i(x_0)$ is the $i$th elements of $X(X^TX)^{-1}x_0$.
I'm struggling to derive how did they arrive to this equation. I know that $\hat{\beta}=(X^TX)^{-1}X^TY$ and so $x_0^T\hat{\beta}=x_0^T(X^TX)^{-1}X^TY$
Since $Y=X^T\beta+\epsilon$ then $x_0^T\hat{\beta}=x_0^T(X^TX)^{-1}X^TY=x_0^T(X^TX)^{-1}X(X^T\beta+\epsilon)$
Well, i'm not sure how exactly to arrive to the equation they do. Can someone show me how?

Let's start with $\underset{N\times 1}{y}=\underset{N\times d}{X}.\underset{d\times 1}{\beta}+\underset{N\times 1}{\epsilon}$, note the dimension of the matrices / vectors underneath.
Hence,
$\sum\limits^N_{i=1}l_i(x_0)\epsilon_i$
$=\langle\;\underset{N\times d}{X}\underset{d \times d \quad}{(X^TX)^{-1}}\underset{d\times 1}{x_0}\quad,\quad \underset{N\times 1}{\epsilon}\;\rangle$, an inner product of two length-$N$ vectors
$=\left(X(X^TX)^{-1}x_0\right)^T.\epsilon$, since $\langle a, b \rangle = a^Tb$
$=\left(x_0^T(X^TX)^{-T}X^T\right).(y-X\beta)$, since $y=X\beta+\epsilon$
$=x_0^T(X^TX)^{-1}X^Ty-x_0^T(X^TX)^{-1}.(X^TX)\beta$, since $(X^TX)^T=X^TX$
$=x_0^T\hat{\beta}-x_0^T\beta$, where $\hat{\beta}=(X^TX)^{-1}X^Ty$
$=\hat{y}_0-x_0^T\beta$
