Forecast Prediction distribution

Modeling
,
#1

I have the following forecast parameterization.

Forecast_h \sim normal(\bar{y}, \hat{\sigma}_{h} )

From a frequentist perspective \hat{\sigma}_{h} is:

\hat{\sigma}_{h} =\sqrt{\frac{1}{T-1}\sum_{t=1}^T e_t^2}\sqrt{1 + \frac{1}{T}}

I am struggling trying use this information to map to a Bayesian formula. In the frequentist specification all of the parameters are determined. How do I modify by stan program to be baysian?

data {
  int<lower=0> N;
  int<lower=0> N_horizon;
  vector[N] y;
}


parameters {
  real<lower=0> sigma;
}

model {
  y ~ normal(mean(y), sigma);
}

generated quantities {

  vector[N_horizon] forecast;
  for (T in 1:N_horizon){
    forecast[T] = normal_rng(mean(y), sigma * sqrt(1+(1/T)));
  }

}
#2

This is one option which does not use any estimated parameters. All this code chunk does is generate a quantity. This is not exactly what I am looking for, but this is the frequentist idea. I am still open to any Bayesian parameterization. Also, is there a way to pull the generated quantities with there no model is fit. I am using CmdStanPy

data {
  int<lower=0> N;
  int<lower=0> N_horizon;
  vector[N] y;
}

/*
parameters {
  real<lower=0> sigma;
}

model {
  y ~ normal(mean(y), sigma);
}
*/

generated quantities {
  vector[N_horizon] forecast;
  for (T in 1:N_horizon){
    forecast[T] = normal_rng(
    mean(y),sqrt((1/(T - 1))*sum(pow(mean(y) - y,2)))*(sqrt(1 + (1/T))));
  }

}
#3

It turns out that my confusion came from a misunderstanding of the subscripts. It makes more sense now.

This is my solution

data {
  int<lower=0> T;
  int<lower=0> horizon;
  vector[T] y;
}

parameters {
  real<lower=0> sigma;
}

model {
  y ~ normal(mean(y), sigma);
}

generated quantities {
  vector[horizon] forecast;
  for (h in 1:horizon){
    forecast[h] = normal_rng(mean(y),sigma * sqrt(1 + (1/T)));
  }

}
3 Likes