Ragged data matrix and some questions using "segment" function

Dear stan users:

I essentially I want to model a histogram-valued data (counts per bin per histogram) using multinomial likelihood. In this case I have 24 histograms, each with varying number of bins. I observed counts per bin per histogram.

This is my model ( right now just using some fake data to ensure that I get the expected results, especially dealing with ragged data structures) , I can import the data fine using the idea from the sparse and ragged matrix section in user manual

data {
  int<lower=0> N_counts;// #total number of observed counts
  int<lower=0> N_bounds;// #total number of observed histogram breaks
  int<lower=0> tt;//# of groups, in this case there are 24 hourly-blocks per day
  //vector[N_counts] counts;// observed counts 
  int<lower=0> counts[N_counts];
  int s_counts[tt];//group size for counts
  vector[N_bounds] int_bounds;//observed bounds
  int s_bounds[tt];//group size for bounds
  vector[2] mu;
  vector[2] sigma;
  vector[tt] lambda;

  real<lower=0> DUMMY;


generated quantities{
  int pos;
  int pos1;
  vector[6] result;
  pos = 1;
  pos1 = 1;
  for(t in 1:tt){ 
 for(b in 1:s_counts[t]){
     print("aa=",  log( lambda[t] * ( normal_cdf( segment (int_bounds, pos1, s_bounds[t]) [b+1], mu[1], sigma[1])
                               - normal_cdf(segment (int_bounds, pos1, s_bounds[t]) [b],  mu[1], sigma[1]) ) +
                      (1-lambda[t]) * (normal_cdf(segment (int_bounds, pos1, s_bounds[t]) [b+1], mu[2], sigma[2])- normal_cdf(segment (int_bounds, pos1, s_bounds[t]) [b],  mu[2], sigma[2])) ) ); 

    pos = pos +s_counts[t];
    pos1 = pos1 +s_bounds[t];

However, my bin probabilities, bb, for a particular bin at a given time period also have different size, using the print statement, I can see that bb has a size of N_counts.

I am stuck in getting using statements like segment (counts, pos, s_counts[t])) ~ multinomial ( bin probability) which I am able to obtain in bb, but it is printing out, I do not know how can I extract the bb elements.

Any help is greatly appreciated.


I think this didn’t get answered because it’s not clear what you’re asking.

If you have 24 histograms with different number of bins, then you’ll need a different simplex parameter for each one.

You might find the slicing syntax a[from:to] easier than segment().

For debugging, try to get one thing working at a time. Just do the extraction and make sure you’re getting that right by printing in a small (but not too small) case.