StataStan and indexes

Interfaces Other
1

Operating System: osx
StataStan version: 1.2.2
CmdStan version: 2.14.0

I have a question on how to pass data to Stan from Stata. The only written example I can find is the Bernoulli and 8 school examples, which only use one index (value passed as a global). Specifically, I’m trying to implement a risk aversion model already written up in stan:
ra_prospect.stan where stan is expecting two dimensional objects.

I get an error “dims declared=(5,140); dims found=(700)” for gamble, using a similar call from the Bernoulli example:

stan gamble gain cert loss , inline thisfile(“raprospect.do”) modelfile(“raprospect.stan”) cmd(“$cmdstandir”) globals(“all”) load mode

having set

global N=5
global T = 140
global numPars = 3

I assume this is simple to fix, but I couldn’t see an example e.g. the actual stata command for the IRT stata paper doesn’t appear to be included in the paper.

Many Thanks, Paul

Stan file below:

data {
  int<lower=1> N;
  int<lower=1> T;
  real<lower=0> gain[N,T];
  real cert[N,T];  
  real<lower=0> loss[N,T];  # absolute loss amount
  int<lower=0, upper=1> gamble[N,T];
}
transformed data {
}
parameters {
  vector[3] mu_p;  
  vector<lower=0>[3] sigma;  
  vector[N] rho_p;  
  vector[N] lambda_p;  
  vector[N] tau_p;  
}
transformed parameters {
  vector<lower=0,upper=2>[N] rho;
  vector<lower=0,upper=5>[N] lambda;
  vector<lower=0>[N] tau;
     
  for (i in 1:N) {
    rho[i]    = Phi_approx( mu_p[1] + sigma[1] * rho_p[i] ) * 2;
    lambda[i] = Phi_approx( mu_p[2] + sigma[2] * lambda_p[i] ) * 5; 
  }
  tau = exp( mu_p[3] + sigma[3] * tau_p );
}
model {
  # ra_prospect: Original model in Soko-Hessner et al 2009 PNAS
  # hyper parameters
  mu_p  ~ normal(0, 1.0); 
  sigma ~ cauchy(0, 5.0);
  
  # individual parameters w/ Matt trick
  rho_p    ~ normal(0, 1.0);   
  lambda_p ~ normal(0, 1.0);   
  tau_p    ~ normal(0, 1.0);
  
  for (i in 1:N) {
      for (t in 1:T) {
      real evSafe;    # evSafe, evGamble, pGamble can be a scalar to save memory and increase speed. 
      real evGamble;  # they are left as arrays as an example for RL models.
      real pGamble;
      
      # loss[i,t]=absolute amount of loss (pre-converted in R)
      evSafe   = pow(cert[i,t], rho[i]);    
      evGamble = 0.5 * (pow(gain[i,t], rho[i]) - lambda[i] * pow(loss[i,t], rho[i]) ); 
      pGamble  = inv_logit( tau[i] * (evGamble - evSafe) );
      gamble[i,t] ~ bernoulli( pGamble );
    }
  }
}
2

Sorry this has gone unanswered. Not many of us use StataStan—I don’t even have Stata installed on my machine. @robertgrant maintains it, and you could always open an issue on the GitHub issue tracker for StataStan, which might get recognied.

is saying that you declared a 2D, 5 x 140 data structure in the Stan model, but when it was doing IO, it found a 1D, 700-element data structure. I’m not sure how you’d normally pass two-dimensional arrays in Stata, but that’s what it should be accepting.

My last suggestion is to look on GitHub for tests.

3

Yes, sorry but I don’t get much time for monitoring forums and such nowadays. There is a matrix option in the -stan- command. Say you define a matrix in Stata, maybe with -mkmat-, and call it gain, then you can send its contents to Stan by invoking that option. It would look like this:
stan gamble cert loss , inline thisfile(“raprospect.do”) modelfile(“raprospect.stan”) cmd("$cmdstandir") globals(“all”) load mode matrix(“gain”)

but I don’t know your data structure. It sounds like it’s long, and gain is a variable. Maybe then you need to reshape gamble and cert and loss as well? Or maybe you don’t have to pass it as a matrix.

4

Many thanks - I’ll have a play with this.

The data are long, but Stan will need some variables to a matrix, so this is useful as no examples seemed to have dealt with such a situation.