Using fitting trivial model dX= \theta dW as example, what I’d like to do in Stan is like this
data {
int N; /* nb. of data points */
int M; /* nb. of auxiliary data between each pair of data points */
real<lower=0> T; /* end time */
vector[N] scal_w; /* observed data X */
}
transformed data {
real dt = T/(N*M); /* fixed time step */
}
parameters {
real<lower=0.0> sigma; /* data simulated from given sigma */
real<lower=0.0> eps; /* obs error */
real sn[M*N]; /* Wiener process' std normal steps */
}
model {
vector[M] x;
real x0;
sigma ~ normal(1.0, 0.5);
eps ~ cauchy(0, 1);
sn ~ normal(0, 1);
x0 = 0.0;
for (i in 1:N) {
/* nested loop equivalent to SDE integrator */
/* that outputs x, given sn, dt, and sigma */
for (j in 1:M) {
x[j] = (j == 1 ? x0 : x[j - 1]) + sigma * sqrt(dt) * sn[(i - 1) * M + j];
}
scal_w[i] ~ normal(x[M], eps);
x0 = x[M];
}
}
as proposed in http://apps.olin.wustl.edu/faculty/chib/papers/elerianchibshephard2001.pdf. The idea is to generate missing/auxiliary data between observations, with SDE integrator describing the recursive conditional likelihood p(X_{i+1}|X_i, W_{i+1}, \theta). This is a full-posterior approach compared with filter-based approaches.