Time series model

Hi there,
I’m new to coding time series models with Stan. I’ve got below R code to generate the data and a Stan model that fits the data. I think it’s all ok, but it would be great if someone who’s worked with time series models before could have a look.
I’m not quite sure whether i’ve done the priors correctly.

R code to generate the data:
Each outcome is drawn from a normal distribution of mean my.mus for that trial (with sd=my.sd).
These means also drift over time, following a normal distribution (with noise my.k).

my.outcomes = array(NA,dim=ntr)
my.mus      = array(NA,dim=ntr)
my.mus[1]   = my.mu1
for (itr in 1:ntr){
  my.outcomes[itr] = round(rnorm(1,my.mus[itr],my.sd));
  if (itr<ntr){
    my.mus[itr+1] = rnorm(1,my.mus[itr],my.k);

From looking at the Stan manual, this is not quite an autoregressive model I think as the outcome on the next trial does not depend on the outcome of the last trial, but instead there is a relationship over time between the (unobserved) means

Here is the Stan code:

data {
  int<lower=0> ntr;
  int outcomes[ntr];

parameters {
  real mag_mean[ntr];
  real<lower=0> mag_sd;
  real<lower=0> k;

model {
  // Priors
  mag_mean[1] ~ normal(50,20); // priors for other trials are given in loop (?)
  mag_sd      ~ normal(0,5);
  k           ~ normal(0,5);
  // Model
  for (itr in 1:(ntr)){
    outcomes[itr]     ~ normal(mag_mean[itr],mag_sd);
    if (itr<ntr){
      mag_mean[itr+1]   ~ normal(mag_mean[itr],k);

Many thanks

Are there multiple observations at each time point? If not, then your priors for the variance parameters are going to be doing a lot of work to render the model identifiable. Have you considered a Gaussian Process model?

Actually there was a previous post on a similar topic, sorry I had missed that before:

which seems to tell me that I should vectorise and use the ‘matt trick’ … trying this now

Thanks for the quick reply.
I’m looking at the Gaussian Process model chapter section 10.1 in current manual. It seems that these models have y and x? But my model only has observations y, nothing else is observed. Am I misunderstanding these models? (Well possible as I did not understand the text well at all given my lack of knowledge about this)

Why dont you write it as Kalman Filter? This Handook could help you alot. And with the gaussian_dlm_obs_lpdf function is really ease to write your model. I think ctsem package could do it for you. But the right person to ask is @Charles_Driver.

In the SUG, x is analogous to time. There’s a worked-through example of a GP regression (with inference on the latent means, which I suspect is what you’re after) here.

You should also take a look at the lgpr for a possibly more user-friendly interface to GPs.

1 Like