Thanks @bbbales2

I do have same number of observations for each person. I am still confused about the *N_obs* in the *stan_data* list. Would *N_obs = 5* for this example or would it need to be changed for each individual depending on if the observation is more than or equal to censored value?

And if so would coding be changed as

```
data {
int<lower=1> J; // number of individuals
int<lower=0> N_obs[J]; // number of uncensored observations for individual j = 1, ..., J
int<lower=0> N_cens[J]; //number of censored points for individual j = 1, ..., J
real y_obs[N_obs, J]; //observed measurements
real<upper=min(y_obs)> L; //censored point
}
parameters {
real<lower =0>sigma_raw;
real<lower=0> sigma[J];
}
transformed parameters {
real mu[N_obs, J];
for (j in 1:J){
mu[, j] = ode_solution;
}
}
model {
sigma_raw ~ exponential(1);
for (j in 1:J){
sigma[j] ~ exponential(sigma_raw);
y_obs[, j]~ normal(mu[,j], sigma[j]);
target += N_cens[j] * normal_lcdf(L | mu[,j], sigma[j]);
}
}
```

And data would be like:

```
stan_data <- list(J = 2,
N_obs = matrix(data = c(3,2),
N_cens = matrix(data = c(2, 3), 1, 2),
y_obs = matrix(data = c(y1_obs, y2_obs), 5, 2),
L = 1)
```

I wanted to make sure that I get this right as the observations corresponds to time series measurements and thus *mu[, j]* needs to be computed for the time points corresponding to jth individual.

Also, there are some columns in the data matrix with no censoring points (ie. all observations >= 1).

Otherwise you’ll need to do something like:

```
int<lower=0> N_obs_total; // Number of total observations
int<lower=0> which_individual[N]; // Indicator of which person owns which observation
real y_obs[N]; // All observed measurements pushed together in a flat array
```

Did you mean to say *N* instead of *N_obs_total* here?

I assume you’ll add ODE parameters and such?

Yes I do have ODE parameters, some of which are fixed and some of them are individual specific.

Thanks