Intermediate quantities

Modeling
1

I’m having a lot of difficulty setting up to define and use an intermediate quantity. In the program below, at the point *** I want to define two such: xb_obs, the matrix-vector product of x_obs and beta (ie x_obs * beta), and similarly for xb_cens. (I know that it’s possible to avoid them here, but I want to know how to use them for a case where they will be needed).

I’ve tried to say (at ***):

real xb_obs [N_obs] ;
xb_obs <- x_obs * beta ;

and

real[N_obs] xb_obs ;
xb_obs<-x_obs * beta ;

but in both cases rstan complains that variable “real” does not exist.

I’m sure that I’m missing something extremely obvious: can anyone tell me what it is?

Also, unrelated: let ind be a vector of indices. What’s the stan
equivalent of R’s x_obs[ind,] - ie rows provided by the indices, all columns?

Thanks very much!

============================ stan program
data {
int<lower=0> N_obs ;
int<lower=0> N_cens ;
int<lower=0> K ; // number of predictors
real y_obs[N_obs] ; //
real<lower=max(y_obs)> U;
matrix[N_obs,K] x_obs; // uncensored predictors data mx, with 1s in col 1
matrix[N_cens,K] x_cens; // censored predictors data mx
}

parameters {
vector[K] beta ;
real<lower=0> sigma ;
}

model {
beta[1] ~ normal(10,5) ;
beta[2] ~ normal(5,5) ;
beta[3] ~ normal(5,5) ;
beta[4] ~ normal(-5,5);
sigma ~ cauchy(0,2) ; # see eg Rethinking p. 249. also stan manual p.54

y_obs ~ normal(xb_obs,sigma);
increment_log_prob( normal_ccdf_log( U , xb_cens , sigma)) ;
}

2

Add it at the top of your model block. Variable definitions have to go at the top. Also, matrix * vector is still vector.

model {
  vector[N_obs] xb_obs = x_obs * beta;
  ... rest of stuff ...
}

If you want it to be an array you can convert it like: to_array_1d(x_obs * beta).

Otherwise, increment_log_prob(v) has been deprecated in favor of target += v, and I think normal_lccdf is the new version of normal_ccdf_log. If you’re using an older version of Stan, probly a good idea to go ahead and update.

Hope that helps!

3

Also, unrelated: let ind be a vector of indices. What’s the stan equivalent of R’s x_obs[ind,] - ie rows provided by the indices, all columns?

Should be just x_obs[ind], but with multiple indexing I find it’s always good to double check just to make sure I’m thinking right (since different people do it differently sometimes). Check out 26.1. Multiple Indexing in the manual.

4

Thanks very much for the extremely fast response — much appreciated!