Create parameter constraint index vector inside the transformed data block

Modeling
1

I would like to implement a zero-based constraints of a Q dimensional parameter vector such that a subset (i.e. D parameters) is set to zero while the remaining elements (P= Q-D parameters) are assigned random samples. Currently, I create the index vectors in R and pass them on to rstan as data

G <- 5 
Q <- G*G
gamma_zeros <- c(1:G, (1:(G-1))* G +1)
gamma_nonzeros <- setdiff(1:Q, gamma_zeros)
D <- length(gamma_zeros)
P <- length(gamma_nonzeros)

These would be used in the transformed parameter block as such

data { 
  int<lower=1> Q;  // dimension of parameter vector
  int<lower=1> D; // total number of zeros constrains 
  int<lower=1> P; // total number of non zeros constraints
  int<lower=1> gamma_zeros[D];  // zero constraint index
  int<lower=1> gamma_nonzeros[P];  // non zero index
}

parameters { 
  vector[P] gamma_tilde;
}

transformed parameters {
  vector[Q] gamma;

 for(d in 1:D) { 
    gamma[gamma_zeros[d]] = 0;
  }

 for(p in 1:P) {
   gamma[gamma_nonzeros[p]] = gamma_tilde[p];
  }
} 


model {
gamma_tilde ~ student_t(4, 0, 4);
// ..
}

I tired to create the index vectors in stan’s model, but I was not able to do so. I was expecting that stan would be able to append arrays together and convert real to integer types but that does not seem to be the case. Is there an easy way to create integer index arrays in stan?

Thanks

2

Either use a loop to fill in the elements or use the R-like syntax start:end.

3

Let’s say Q=16 and I want to create these two integer arrays

gamma_zeros = {1,  2,  3,  4,  5,  9, 13}
gamma_nonzeros = { 6,  7,  8, 10, 11, 12, 14, 15, 16}

The for loop can help me create the sequence of integers

  for(q in 1:Q) { 
    s[q] = q;
  }

but I am not sure the R-like syntax start:end would work since I would need to append integer arrays (i.e. c(1:G, (1:(G-1))* G +1)at one point and such an operation is not implemented in stan. For example, say

append_array(s[q][1:G], s[q][1:(G-1)]* G +1)
4

There is an append_array function in 2.18 but it only takes two arguments. You are probably best to just create your array via a loop in the transformed data block.

5

Yes, that is what I ended up doing.

  int<lower=1> D= 2*G-1; // total number of zeros constrains 
  int<lower=1> P= G*G-D; // total number of non zeros 
  int<lower=1> gamma_zeros[D];  // zero constraint index array
  int<lower=1> gamma_nonzeros[P];  // nonzero index array
  int<lower=1> dd;
  int<lower=1> pp;
  dd=1;
  pp=1;
  for(i in 1:G) { 
    for(j in 1:G) {
   	if (i==1 || j==1) {
    	  gamma_zeros[dd] = (i-1)*(G)+j;
    	  dd += 1;}
   	else {
          gamma_nonzeros[pp]= (i-1)*(G)+j;
          pp += 1;}
    }
  }
6

You can look at my implementation of the R %in% function in Stan here – Stan equivalent of R's %in% function

Look down in the comments for an even faster version from a real programmer.

This can help with accessing ragged arrays.