find the autocovariance function of the process $Y_t$ Consider the processes $X_t = \phi X_{t-1} + v_t$ and $Y_t = \phi Y_{t-1} + X_t + e_t$, in which $|\phi| < 1$ and $v_1$ and $e_t$ are non-correlated random errors with zero mean and variances equal to $\sigma^2$. Based on these informations, find the autocovariance function of the process $Y_t$.

I already know that if $|\phi| < 1$ than $Y_t$ should follow a stationary process and it’s autocovariance is:
$\begin{equation} 
E[Y_{t+h}Y_t] - \mu^2
\end{equation}$
but I do not know how to go further
 A: Note that all the processes have zero mean, i.e. $E[Y_t]=E[X_t]=E[v_t]=E[\epsilon_t]=0$.
Let $h>0$:
$$\begin{align}E[Y_{t+h}Y_t]&=E[(\phi Y_{t+h-1}+X_{t+h}+\epsilon_{t+h})Y_t]=\phi E[Y_{t+h-1}Y_t]+E[X_{t+h}Y_t]+E[\epsilon_{t+h}Y_t]\\&=\phi r_y(h-1)+\phi^h r_{xy}(0)+E[\epsilon_{t+h}]E[Y_t]=\phi r_y(h-1)+\phi^h r_{xy}(0)\end{align}$$
Note that: $r_{xy}(h)=E[X_{t+h}Y_t]=E[\phi^hX_tY_t+(\sum_{j=0}^{h-1} v_{t+h-j}\phi^j)Y_t]=\phi^h r_{xy}(0)$
For $h=0$:
$$\begin{align}r_y(0)&=E[Y_t^2]=E[(\phi Y_{t-1}+X_t+\epsilon_t)^2]\\&=\phi ^2r_y(0)+E[X_t^2]+E[\epsilon_t^2]+2\phi \underbrace{E[Y_{t-1}X_t]}_{r_{xy}(1)}+2E[X_t\epsilon_t]+2E[Y_{t-1}\epsilon_t]\\&=\phi^2r_y(0)+r_x(0)+r_e(0)+2\phi^2 r_{xy}(0)\rightarrow r_y(0)=\frac{r_x(0)+\sigma^2+2\phi^2 r_{xy}(0)}{1-\phi^2}\end{align}$$
Similarly, we can find that $$r_x(0)=\frac{\sigma^2}{1-\phi^2}, \ \ r_x(h)=\phi r_x(h-1)$$
It just remains to find $r_{xy}(0)$:
$$\begin{align}r_{xy}(0)&=E[X_tY_t]=[(\phi X_{t-1}+v_t)(\phi Y_{t-1}+X_t+\epsilon_t)]\\&=\phi^2E[X_{t-1}Y_{t-1}]+\phi r_x(1)+E[X_tv_t]\\&=\phi^2r_{xy}(0)+\phi r_x(1)+\sigma^2\end{align}$$
