I am trying to fit a model to annual (time series) data with a time-dependent external regressor.

In my simple/simplistic scenario (no autoregressive terms, no moving-average terms, no seasonality), changes in response y are exclusively driven by changes in the temperature T, which in turn is expected to have some linear dependence on time t. My goal is to fit data to the model, and then predict y into the future based on future predictions for T.

In the context of ARIMA models, this would correspond to an ARIMA(0, 0, 0) model, with the single external regressor T.

I had expressed the model in the following form:

The likelihoods for y and T are coupled through \mu_t = \alpha + \mu_T.

Translating this into Stan and fitting the model results in errors of the form

Chain 1: Rejecting initial value:

Chain 1: Error evaluating the log probability at the initial value.

Chain 1: Exception: normal_lpdf: Location parameter[1] is nan, but must be finite! (in âmodel90841be2e428_8f83da75be386a1006ae4f72518f8aefâ at line 37)

I am not sure I understand whatâs going on. I had *hoped* to fit simultaneously T to the temperature data and y to the response data. In terms of forecasting, sampling from T for future time points and using these samples to draw samples for y would allow me to carry through uncertainties in the regressor T to samples of y.

Would somebody be able to point out where Iâm going wrong?

Reproducible R code example:

```
model_code <- "
data {
int<lower=1> N; // Number of time points
vector[N] t; // Time
vector[N] y; // Response
vector[N] T; // Temperature
}
parameters {
real<lower=0> sigma_y;
real<lower=0> sigma_T;
real alpha;
real beta0;
real beta1;
}
transformed parameters {
vector[N] mu_y;
vector[N] mu_T;
mu_y = alpha + mu_T;
mu_T = beta0 + beta1 * t;
}
model {
// Priors
sigma_y ~ cauchy(0, 2.5);
sigma_T ~ cauchy(0, 2.5);
alpha ~ normal(mean(y), 10);
beta0 ~ std_normal();
beta1 ~ normal(0, 10);
// Likelihood
T ~ normal(mu_T, sigma_T);
y ~ normal(mu_y, sigma_y);
}
"
data <- list(
N = 13L,
t = c(
2007, 2008, 2009, 2010, 2011, 2012, 2013,
2014, 2015, 2016, 2017, 2018, 2019),
y = c(
49.6417910447761, 46.2388059701492, 55.1940298507463, 54.4776119402985, 55.910447761194, 45.3432835820895, 52.5074626865672, 55.1940298507463, 44.6268656716418, 52.865671641791, 55.910447761194, 51.7910447761194, 59.1343283582089),
T = c(
39.3, 34.7, 40.1, 39.2, 37.55, 34.8, 40.85,
40.1, 36.8, 39, 41.35, 40, 41.2))
fit <- stan(model_code = model_code, data = data, chain = 1)
```