Formula for an ARIMA(1,1,1) solving for y I have an Arima(1,1,1) model with predictors var1+var2+var3, but am struggling with how to write the equation. The problem is that on all of the sources I see a variation of the following is given.
$$\left( 1 - \sum_{i=1}^p \phi_i L^i\right)
(1-L)^d y_t = \delta + \left( 1 + \sum_{i=1}^q \theta_i L^i \right) \varepsilon_t . $$
Is there a way to write the equation solving for $y_t$? I find that the above equation is difficult to understand. Though I have the model saved in R the users of the forecast want the coefficients so they can plug them into excel and I don't have a good way to explain the formula solving for $y_t$
 A: Arima() fits a so-called regression with ARIMA errors. Note that this is different from an ARIMAX model. In your particular case, you regress your focal variable on three predictors, with an ARIMA(1,1,1) structure on the residuals:
$$ y_t = \beta_1x_{1t}+\beta_2x_{2t}+\beta_3x_{3t}+\epsilon_t $$
with $\epsilon_t\sim\text{ARIMA}(1,1,1)$.
To write down the formulas for $\epsilon_t$, we use the backshift operator. For ARIMA background, see here. A general ARIMA(1,1,1) model with AR parameter $\phi$ and MA parameter $\theta$ has the following form (note that some packages flip the sign of $\theta$, and that one sometimes fits a nonzero intercept $c$):
$$ (1-\phi B)(1-B)\epsilon_t = (1+\theta B)\eta_t, $$
where $\eta_t\sim N(0,\sigma^2)$ are IID innovations. Some calculations show that
$$ \begin{align*} & \epsilon_t-\phi\epsilon_{t-1}-\epsilon_{t-1}+\phi\epsilon_{t-2} \\
= & (1-\phi B)(\epsilon_t-\epsilon_{t-1}) \\
= & (1-\phi B)(1-B)\epsilon_t \\
= & (1+\theta B)\eta_t \\
= & \eta_t+\theta\eta_{t-1},
\end{align*}
$$
or
$$ \epsilon_t = (1+\phi)\epsilon_{t-1}-\phi\epsilon_{t-2}+\eta_t+\theta\eta_{t-1}.$$
