Hi,
we are two Italian statistical student of Sapienza Università di Roma. We are working about Testing Hypotheses via a mixture estimation model, based on the work of Kamary, Mergensen, Robert and Rousseau.
In particular we would like to implement this procedure for a mixture of two autoregressive processes AR(1) and AR(2). Our aim is to estimate the weights of the mixture and analyze if their behaviour is coherent with the posterior probability of the model computed through the Bayes Factor.
We tried to implement ourself an Metropolis within Gibbs algorithm for all the parameters of the mixture and it seems to work quite well when the true model is AR(1) but not when the true model is AR(2).
Furthermore we would like to implement the same procedure with the package Rstan. But we had a lot of problems in the code.
Could you help us with same advises on the implementation of our code.
We joint our code in pdf (r-markdown) with a short summary, an attempt of the implementation with Rstan package and the paper based on.
library(rstan)
y_2 = arima.sim(list(order=c(1,0,0),ar = c(0.5)) ,n=10000)
mix <- '
data {
int<lower=0> N;
real y[N];
}
parameters {
real<lower=0,upper=1> theta;
real alpha;
real mu;
real beta;
real gamma;
real<lower=0> sigma;
real<lower=0> tau;
}
model {
theta ~ beta(1,1);
mu ~ normal(0.5,6);
sigma ~ inv_gamma(2, 0.1);
for(n in 3:N)
target += log_mix(theta, normal_lpdf(y[n] | mu + beta * y[n-1], sigma),
normal_lpdf(y[n] |alpha + beta * y[n-1] + gamma * y[n-2], tau));
}
'
dat = list(y_2)
test = stan(model_code = mix, model_name = “example”,
data = dat, iter = 10000, chains = 2,
verbose = F)
Robert,Rousseau, Kameny,Mergensen “Testing hypotesis via mixture estimation model”.pdf (1.1 MB)
attempt1.pdf (1.5 MB)
[edit: escaped code]