# Modelling a jump component

Dear All! I would like to ask a question. How to implement the following model in Stan?

I have a time-series:

y_t = \mu_t + \xi_t J_t + \epsilon_t, \qquad \epsilon_t \sim \mathcal{N}(0,\exp(h_t/2)).

I have no problem with modelling the SV part and \mu_t, but I do not understand how to incorporate a jump component. Ideally, I would like to have the following specification:

\xi_t \sim \mathcal{N}(0, \sigma_{\xi}^2) - measures the intensity of the jump.
J_t = \Phi(\lambda_0 + \lambda_x X) - determines the probability of the jump, where X is observable data that gives me information about the probability of the jump. If I use

 bernoulli(Phi_approx(lam_1  + lam_x * x));


Ill-typed arguments to ‘~’ statement. No distribution ‘bernoulli’ was found with the correct signature.

Because Stan does not support discrete parameter. Could you please suggest how to translate the model in Stan model or propose a similar model that could be implemented in Stan? Below, I provide a code that I wrote so far, but as I wrote before it does not match the structure exactly

data {
int<lower=0> T_len;           // Number of observations
vector[T_len] y;              // Time series data
vector[T_len] x;              // Time series data
}

parameters {
vector[T_len] xi;                             // jump magnitude at time t
real<lower=0> sig_y;                          // Standard deviation of the y process
real<lower=0> sig_xi;                         // Standard deviation of the jump magnitude
real gamma_y;                                 // Intercept of the y process
real lam_1;                                   // Parameters for jump
real lam_x;                                   // Parameters for jump
}

transformed parameters {
vector[T_len] jump = Phi_approx(lam_1  + lam_x * x);          // Jump index at time t
}

model {
// Prior
sig_y    ~ cauchy(0, 5);               // Prior on the standard deviation of the y process
sig_xi   ~ cauchy(0, 5);               // Prior on the standard deviation of the jump magnitude
gamma_y  ~ normal(0, 10);              // Prior on the intercept of the y process
lam_1    ~ normal(0, 10);              // Prior on the intercept of the y process
lam_x    ~ normal(0, 10);              // Prior on the intercept of the y process

xi       ~ normal(0, sig_xi * jump);                       // Jump magnitude at time t

// Likelihood (speed up)
y ~ normal(gamma_y + xi, sig_y);            // Likelihood of the data given
}


Thank you very much!

Best,
Borys

Discrete parameters must be integrated out of Stan models.

For simplicity, let’s start with a likelihood that uses bernoulli() (incorrectly). I gather you essentially want

for (t in 1:T_len) {
xi[t] ~ normal(0, sig_xi);
if (bernoulli(jump[t])) {
y[t] ~ normal(mu[t] + xi[t], sig_y);
} else {
y[t] ~ normal(mu[t], sig_y);
}
}


(your code had gamma_y but I substituted mu[t] as that’s what’s shown in the equations.)

First, integrate out the jump size xi. (It’s continuous so integrating it out is not strictly required but (I expect) will improve sampling efficiency significantly). So, the likelihood is now

for (t in 1:T_len) {
if (bernoulli(jump[t])) {
y[t] ~ normal(mu[t], sqrt(sig_xi^2 + sig_y^2));
} else {
y[t] ~ normal(mu[t], sig_y);
}
}


and the model has no parameter xi.

Integrating out discrete choices requires us to use target += syntax instead of the sampling statement. Here’s how it looks. The following code is equivalent to the previous block, just different syntax that explicitly spells out log probability density calculations.

for (t in 1:T_len) {
if (bernoulli(jump[t])) {
target += normal_lpdf(y[t] | mu[t], sqrt(sig_xi^2 + sig_y^2));
} else {
target += normal_lpdf(y[t] | mu[t], sig_y);
}
}


Now we can handle the discrete Bernoulli choice.
Stan has a function specifically for that, log_mix.

for (t in 1:T_len) {
target += log_mix(jump[t],
normal_lpdf(y[t] | mu[t], sqrt(sig_xi^2 + sig_y^2)),
normal_lpdf(y[t] | mu[t], sig_y));
}

2 Likes

Dear nhuurre, thank you very much for your response. It’s exactly what I looked for. I’m sorry for the late reply, since I was extremely busy and did not work on this problem