Hi All,
I am trying to fit a graded response model in rstan. Unfortunately, R keeps crashing with the error “R for Windows terminal front-end has stopped working.”
I get the following output before the crash happens:
> stanout <- stan('latent-factor-model.stan'
+ , data=standat
+ , iter=4000
+ , chains=1)
hash mismatch so recompiling; make sure Stan code ends with a blank line
In file included from C:/Users/CURTS025/programs/R/library/BH/include/boost/config.hpp:39:0,
from C:/Users/CURTS025/programs/R/library/BH/include/boost/math/tools/config.hpp:13,
from C:/Users/CURTS025/programs/R/library/StanHeaders/include/stan/math/rev/core/var.hpp:7,
from C:/Users/CURTS025/programs/R/library/StanHeaders/include/stan/math/rev/core/gevv_vvv_vari.hpp:5,
from C:/Users/CURTS025/programs/R/library/StanHeaders/include/stan/math/rev/core.hpp:12,
from C:/Users/CURTS025/programs/R/library/StanHeaders/include/stan/math/rev/mat.hpp:4,
from C:/Users/CURTS025/programs/R/library/StanHeaders/include/stan/math.hpp:4,
from C:/Users/CURTS025/programs/R/library/StanHeaders/include/src/stan/model/model_header.hpp:4,
from file46d031c2faa.cpp:8:
C:/Users/CURTS025/programs/R/library/BH/include/boost/config/compiler/gcc.hpp:186:0: warning: "BOOST_NO_CXX11_RVALUE_REFERENCES" redefined
# define BOOST_NO_CXX11_RVALUE_REFERENCES
^
<command-line>:0:0: note: this is the location of the previous definition
SAMPLING FOR MODEL 'latent-factor-model' NOW (CHAIN 1).
Gradient evaluation took 0 seconds
1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
Adjust your expectations accordingly!
Iteration: 1 / 4000 [ 0%] (Warmup)
My R code and Stan code are below.
##
## Simulated data
##
library(tidyverse)
library(rstan)
simgrm <- function(n_subj, n_items, n_cat){
item_prob <- function(t, a, k){
c(0, plogis(k - a*t), 1) %>%
diff()
}
item_sample <- function(t, a, k){
item_prob(t, a, k) %>% sample(1:length(.), 1, prob=.)
}
theta <- rnorm(n_subj)
alpha <- rnorm(n_items)
kappa <- t(sapply(1:n_items
, function(j){
sort(c(rnorm(n_cat - 1)))
}))
out <- data_frame(Respondent_Idx=1:n_subj) %>%
crossing(Item_Idx=1:n_items) %>%
mutate(Y=mapply(item_sample
, t=theta[Respondent_Idx]
, a=alpha[Item_Idx]
, k=as.list(data.frame(t(kappa[Item_Idx, ])))))
return(out)
}
set.seed(1984)
n_cat <- 3L
n_respondents <- 500L
n_items <- 5L
tmpdat <- simgrm(n_respondents, n_items, n_cat)
standat <- list(
n_data=nrow(tmpdat)
, n_respondents=length(unique(tmpdat$Respondent_Idx))
, n_items=length(unique(tmpdat$Item_Idx))
, n_cat=n_cat
, Respondent_Idx=tmpdat$Respondent_Idx
, Item_Idx=tmpdat$Item_Idx
, Y=tmpdat$Y
)
stanout <- stan('latent-factor-model.stan'
, data=standat
, iter=4000
, chains=1)
data {
int <lower=1> n_data;
int <lower=1> n_items;
int <lower=1> n_respondents;
int <lower=1> n_cat;
int <lower=1, upper=n_items> Item_Idx[n_data];
int <lower=1, upper=n_respondents> Respondent_Idx[n_data];
int <lower=1, upper=n_cat> Y[n_data];
}
parameters {
vector[n_respondents] theta;
real<lower=0> alpha[n_items];
ordered[n_cat - 1] kappa[n_items];
real mu_kappa;
real<lower=0> sigma_kappa;
}
model {
theta ~ normal(0.0, 1.0);
alpha ~ cauchy(0.0, 5.0);
mu_kappa ~ normal(0.0, 5.0);
sigma_kappa ~ cauchy(0.0, 5.0);
for (j in 1:n_items){
for (k in 1:(n_cat - 1)){
kappa[j, k] ~ normal(mu_kappa, sigma_kappa);
}
}
for (n in 1:n_data) {
Y[n] ~ ordered_logistic(theta[Respondent_Idx[n]] * alpha[Item_Idx[n]], kappa[Item_Idx[n]]);
}
}
Any help would be appreciated!
Thank you,
McKay
[edited to set off code]