FInding ARMA model/coefficients with correlation function I want to find a stationary ARMA difference equation (i.e the ARMA coefficients) such that $\rho(1)=\frac23$, $\rho(2)=\frac16$ and $\rho(k)=0$ for all $k=3,4,...$, where $\rho$ is the auto correlation function.
I would like to use Yule Walker equations, but I think they only apply to AR models. Any ideas?
 A: Since the autocorrelation function cuts off at lags higher than 2 you have a MA(2) process $y_t = w_t - \theta_1 w_{t-1} - \theta_2 w_{t-2}$.  In general, to find the moving average coefficients, you equate the autocorrelations expressed in terms of the $\theta_i$'s to each $\rho_k$.  For an MA(1) process you end up with a quadratic equation with up to two solutions.  For an MA(2) process you end up with a non-linear system of two equations with up to four solutions.  You can solve these equations numerically, e.g. using the rootSolve R-package.
Suppose that 
$$
\theta(z) = (1-r_1 z)(1-r_2 z) \tag{1}
$$
is the MA-polynomial corresponding to a solution you've found and suppose $|r_1|>1$ such that the model is non-invertible ($\theta(z)$ has one root $1/r_1$ inside the unit circle).  Then a reparameterized invertible model 
$$
\theta'(z) = (1-r_1' z)(1-r_2 z)
$$ 
where $r_1'=1/r_1$ is also a solution since this MA(2)-process has the same autocorrelation function as (1) (you may want to verify this for yourself).  
In general you can always find an invertible MA$(q)$ model through this form of reparameterization unless $\theta(z)$ has one or more unit roots.  This can be easily proved using autocovariance generating functions. 
