Hello team,

I have been trying to fit an ARIMAX model with order (0,1,1) in Stan with an exogenous variable. I have the Stan code for ARIMA, and I want to reuse this code to fit ARIMAX passing an exogenous variable to it.

```
data {
int<lower=0> T; //no. of observations
real y[T]; // observed outputs
}
parameters {
real phi; //AR coeff
real theta; //MA coeff
real mu; //intercept coeff
real<lower=0> sigma; // noise scale
}
model {
vector[T] nu; // prediction for time t
vector[T] err; // error for time t
nu[1] = mu + phi * mu; // assume err[0] == 0
err[1] = y[1] - nu[1];
for (t in 2:T) {
nu[t] = mu + phi * (y[t-1]) + theta * err[t-1];
err[t] = y[t] - nu[t];
}
mu ~ cauchy(-0.07, 0.1);
phi ~ normal(0.85, 0.1);
theta ~ normal(0, 1);
sigma ~ cauchy(0, 1);
err ~ normal(0, sigma);
}
```

Is there any way to pass the exogenous variable to this stan code and if yes, how can I reuse this code to fit the ARIMAX model? Is there some other code that I can use to fit the ARIMAX model with order (0,1,1) in Stan?

All suggestions are welcome.