# How to specify ODE output when there are replicates

**URL:** <https://discourse.mc-stan.org/t/how-to-specify-ode-output-when-there-are-replicates/28626>\
**Category:** Modeling\
**Created:** [August 23, 2022, 4:59pm UTC](https://discourse.mc-stan.org/t/how-to-specify-ode-output-when-there-are-replicates/28626 "2022-08-23T16:59:53Z")\
**Posts on this page:** 4\
**Page:** 1

<div class="post-metadata">

**Author:** ![Xiang\_Zhang](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/xiang_zhang/32/13172_2.png) [@Xiang\_Zhang](https://discourse.mc-stan.org/u/Xiang_Zhang)\
**Post date:** [August 23, 2022, 4:59pm UTC](https://discourse.mc-stan.org/t/how-to-specify-ode-output-when-there-are-replicates/28626/1 "2022-08-23T16:59:53Z")

</div>

Hi,  
I have been confused for a while about the example 13.4 in the user guide. In the generated quantities section, the ODE solution is defined as below

```nohighlight
vector[2] y_sim[T] = ode_rk45(sho, y0, t0, ts, theta);

```

I understand that the solution of ODE is a vector of 2 and the ODE is solved T times. My first question is that this kind of specification of y\_sim is not considered as the optimal way (there will be a warning when the model is compiled). But it doesn’t work if it is array[T] vector[2] y\_sim.

In this example, there is one solution based on one initial condition. What if I have a number of different initial conditions, then y\_sim should be a array with [n\_conditions, T, 2]. But the compiling won’t work. I wonder if anyone could give me a hint.

---

<div class="post-metadata">

**Author:** ![Bob\_Carpenter](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/bob_carpenter/32/9230_2.png) [@Bob\_Carpenter](https://discourse.mc-stan.org/u/Bob_Carpenter)\
**Post date:** [August 23, 2022, 8:40pm UTC](https://discourse.mc-stan.org/t/how-to-specify-ode-output-when-there-are-replicates/28626/2 "2022-08-23T20:40:01Z")

</div>

> [@Xiang\_Zhang](#):
>
> My first question is that this kind of specification of y\_sim is not considered as the optimal way (there will be a warning when the model is compiled). But it doesn’t work if it is array[T] vector[2] y\_sim.

Can you please include the full Stan code generating that example and say which version of which Stan interface you’re using? My guess is that you’re using a version of RStan that’s using an old version of Stan that doesn’t support the new array syntax. The old array syntax has been deprecated and is scheduled for removal soon. I would suggest moving to CmdStanR or CmdStanPy, which are up to date with Stan and will support the `array[T]` syntax.

---

<div class="post-metadata">

**Author:** ![Xiang\_Zhang](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/xiang_zhang/32/13172_2.png) [@Xiang\_Zhang](https://discourse.mc-stan.org/u/Xiang_Zhang)\
**Post date:** [August 24, 2022, 8:04am UTC](https://discourse.mc-stan.org/t/how-to-specify-ode-output-when-there-are-replicates/28626/3 "2022-08-24T08:04:32Z")

</div>

Hi Bob,

I am using cmdstanR 0.5.3 (CmdStan 2.30.1).

The full model is

functions {  
vector sho(real t,  
vector y,  
real theta) {  
vector[2] dydt;  
dydt[1] = y[2];  
dydt[2] = -y[1] - theta \* y[2];  
return dydt;  
}  
}  
data {  
int\<lower=1\> T;  
array[T] vector[2] y;  
real t0;  
array[T] real ts;  
}  
parameters {  
vector[2] y0;  
vector\<lower=0\>[2] sigma;  
real theta;  
}  
model {  
vector[2] mu[T] = ode\_rk45(sho, y0, t0, ts, theta);  
sigma ~ normal(0, 2.5);  
theta ~ std\_normal();  
y0 ~ std\_normal();  
for (t in 1:T) {  
y[t] ~ normal(mu[t], sigma);  
}  
}

When I compile the example there is a warning:  
Compiling Stan program…  
Warning in ‘/var/folders/9x/86f88pr155s4lsxpdt623gx57n0ss1/T/Rtmpim4LXT/model-61687d701e5d.stan’, line 23, column 2: Declaration  
of arrays by placing brackets after a variable name is deprecated and  
will be removed in Stan 2.32.0. Instead use the array keyword before the  
type. This can be changed automatically using the auto-format flag to  
stanc

This is caused by the vector[2] mu[T] = ode\_rk45(sho, y0, t0, ts, theta);

When I replace this with  
array[T] vector [2] mu = ode\_rk45(sho, y0, t0, ts, theta);  
The warning is gone. This is good.

Now, if I have a number of replicates (n\_rep), and I want to fit the ODE model to each of them. I attempted to do the model below:

functions {  
vector sho(real t,  
vector y,  
real theta) {  
vector[2] dydt;  
dydt[1] = y[2];  
dydt[2] = -y[1] - theta \* y[2];  
return dydt;  
}  
}  
data {  
int\<lower=1\> T;  
array[T] vector[2] y;  
real t0;  
array[T] real ts;  
int\<lower=0\> n\_rep;  
}  
parameters {  
vector[2] y0;  
vector\<lower=0\>[2] sigma;  
real theta;  
}  
model {  
array[T,n\_rep] vector[2] mu = ode\_rk45(sho, y0, t0, ts, theta);  
sigma ~ normal(0, 2.5);  
theta ~ std\_normal();  
y0 ~ std\_normal();  
for (t in 1:T) {  
y[t] ~ normal(mu[t], sigma);  
}  
}

## But I got an error message Semantic error in ‘/var/folders/9x/86f88pr155s4lsxpdt623gx57n0ss1/T/Rtmpim4LXT/model-61682954e7c1.stan’, line 24, column 2 to column 65:

```
22: }
23: model {
24: array[T,n_rep] vector[2] mu = ode_rk45(sho, y0, t0, ts, theta);
       ^
25: sigma ~ normal(0, 2.5);
26: theta ~ std_normal();

```

* * *

Ill-typed arguments supplied to assignment operator =: lhs has type array[,] vector and rhs has type array vector

I couldn’t figure out myself and I hope you could provide me some hint.

---

<div class="post-metadata">

**Author:** ![Xiang\_Zhang](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/xiang_zhang/32/13172_2.png) [@Xiang\_Zhang](https://discourse.mc-stan.org/u/Xiang_Zhang)\
**Post date:** [August 24, 2022, 9:39am UTC](https://discourse.mc-stan.org/t/how-to-specify-ode-output-when-there-are-replicates/28626/4 "2022-08-24T09:39:41Z")

</div>

I think I have found the solution  
The key is to define mu as  
array[n\_rep, T] vector[2] mu;

Then  
for (i in 1:n\_rep) {  
mu[i] = ode\_rk45(sho, y0, t0, ts, theta);  
}

mu[i] is an array (dimension = T) of vectors (length = 2).

Hopefully it is helpful for others who are interested