Hi folks,
I am trying to fit models using the Stan code produced from the blavaan package. The issue is that after I fit a first model and then do it for a second model, sampling is done not for the second model but for the first model. I do not understand why this is happening, so I would like to get some solutions from others. I first attach the code I used.
library(rstan)
options(mc.cores = 2)
rstan_options(auto_write = TRUE)
library(blavaan)
# 1. Data simulation
# True model
model.true <- '
# measurement equation
fac1 =~ x1 + 0.5*x2 + 0.5*x3
fac2 =~ x4 + 0.5*x5 + 0.5*x6
fac3 =~ x7 + 0.5*x8 + 0.5*x9
# structural equation
fac1 ~~ 0.5*fac2
fac1 ~~ 0.5*fac3
fac2 ~~ 0.5*fac3
'
set.seed(322)
sim.dat <- simulateData(model.true, sample.nobs = 500)
# 2. Models for comparison
model.1 <- '
# measurement equation
fac1 =~ x1 + x2 + x3
fac2 =~ x4 + x5 + x6
fac3 =~ x7 + x8 + x9
# structural equation
fac1 ~~ fac2
fac1 ~~ fac3
fac2 ~~ 0*fac3
'
model.2 <- '
# measurement equation
fac1 =~ x1 + x2 + x3
fac2 =~ x4 + x5 + x6
fac3 =~ x7 + x8 + x9
# structural equation
fac2 ~~ fac1
fac2 ~~ fac3
fac1 ~~ 0*fac3
'
With this simulated dataset and the two models, I fit models using Stan code produced from blavaan. I first get the Stan code from blavaan.
# 3. Fitting models using Stan code from blavaan
# Model 1
set.seed(322)
blav.fit.1 <- bsem(model.1, data=sim.dat, n.chains=2, mcmcfile="comparison", do.fit=FALSE)
Please note that, then, I renamed the two files (semstan.rda and sem.stan into semm1.rda and semm1.stan in the comparison folder). This procedure should be done to reproduce my issue!
load('comparison/semm1.rda')
set.seed(322)
stan.fit.1 <- stan('comparison/semm1.stan', data=stantrans$data, chains=2,
iter=1500, warmup=500, init=stantrans$inits)
The fitting is successfully completed. Then, for the second model,
# Model 2
set.seed(322)
blav.fit.2 <- bsem(model.2, data=sim.dat, n.chains=2, mcmcfile="comparison", do.fit=FALSE)
Please note that I again renamed the two files (semstan.rda and sem.stan into semm2.rda and semm2.stan in the comparison folder).
As soon as I run the following code…
load('comparison/semm2.rda')
set.seed(322)
stan.fit.2 <- stan('comparison/semm2.stan', data=stantrans$data, chains=2,
iter=1500, warmup=500, init=stantrans$inits)
In the ‘Viewer’, the message appears such that sampling happens for model ‘semm1’ not ‘semm2’!!!:
Click the Refresh button to see progress of the chains
starting worker pid=22092 on localhost:11616 at 21:53:37.407
starting worker pid=7916 on localhost:11616 at 21:53:37.785
SAMPLING FOR MODEL 'semm1' NOW (CHAIN 1).
Chain 1:
Chain 1: Gradient evaluation took 0.017 seconds
Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 170 seconds.
Chain 1: Adjust your expectations accordingly!
Chain 1:
Chain 1:
Chain 1: Iteration: 1 / 1500 [ 0%] (Warmup)
SAMPLING FOR MODEL 'semm1' NOW (CHAIN 2).
Chain 2:
Chain 2: Gradient evaluation took 0 seconds
Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
Chain 2: Adjust your expectations accordingly!
Chain 2:
Chain 2:
Chain 2: Iteration: 1 / 1500 [ 0%] (Warmup)
Chain 2: Iteration: 150 / 1500 [ 10%] (Warmup)
Chain 1: Iteration: 150 / 1500 [ 10%] (Warmup)
Chain 2: Iteration: 300 / 1500 [ 20%] (Warmup)
Chain 1: Iteration: 300 / 1500 [ 20%] (Warmup)
Chain 2: Iteration: 450 / 1500 [ 30%] (Warmup)
Chain 1: Iteration: 450 / 1500 [ 30%] (Warmup)
Chain 2: Iteration: 501 / 1500 [ 33%] (Sampling)
Chain 1: Iteration: 501 / 1500 [ 33%] (Sampling)
Chain 2: Iteration: 650 / 1500 [ 43%] (Sampling)
Chain 1: Iteration: 650 / 1500 [ 43%] (Sampling)
Chain 2: Iteration: 800 / 1500 [ 53%] (Sampling)
Chain 1: Iteration: 800 / 1500 [ 53%] (Sampling)
Chain 2: Iteration: 950 / 1500 [ 63%] (Sampling)
Chain 1: Iteration: 950 / 1500 [ 63%] (Sampling)
Chain 2: Iteration: 1100 / 1500 [ 73%] (Sampling)
Chain 1: Iteration: 1100 / 1500 [ 73%] (Sampling)
Chain 2: Iteration: 1250 / 1500 [ 83%] (Sampling)
Chain 1: Iteration: 1250 / 1500 [ 83%] (Sampling)
Chain 2: Iteration: 1400 / 1500 [ 93%] (Sampling)
Chain 1: Iteration: 1400 / 1500 [ 93%] (Sampling)
Chain 2: Iteration: 1500 / 1500 [100%] (Sampling)
Chain 2:
Chain 2: Elapsed Time: 30.131 seconds (Warm-up)
Chain 2: 37.617 seconds (Sampling)
Chain 2: 67.748 seconds (Total)
Chain 2:
Chain 1: Iteration: 1500 / 1500 [100%] (Sampling)
Chain 1:
Chain 1: Elapsed Time: 33.04 seconds (Warm-up)
Chain 1: 37.227 seconds (Sampling)
Chain 1: 70.267 seconds (Total)
Chain 1:
I do not know why this is happening. Is this okay and a naturally-happening process where the semm2 is being sampled actually? Or, is this really a problem? My guess is that the functionrstan_options(auto_write = TRUE) should have something to do with it, but I do not know.
So, could anyone help me in figuring out this problem?
Thanks!