I’ve been going back and forth on this one. Any thoughts or advice are welcomed.
Right now, Torsten returns a matrix. This matrix stores quantities of interest at different times, in the same way the returned object of integrate_ode_*
does. Except the latter function returns a 2D array.
In most cases, Torsten computes the elements of the matrix one at a time. When it solves linear ODEs, it calls matrix_exp()
and computes the matrix one row at a time. Still, no real gain using a matrix (versus an array of vectors; I save a conversion step I would need to do for a 2D array from eigen::vector to std::vector for linear ODEs. That said with the matrix, I need to convert std::vectorstd::vector to eigen::vector for nonlinear ODEs).
My initial thought was that a matrix would be easier to work with in later parts of the model. Usually, I’m interested in one column and I want to compare it to data in my likelihood function. I understand sampling statements are more efficient with vectors. So it pays off to have cHat
(the column of interest) be a vector and it easier to extract it from a matrix.
transformed parameters {
...
x = PKModelTwoCpt(theta, time, amt, rate, ii, evid, cmt, addl, ss);
cHat = col(x, 2) ./ V1;
...
}
}
However, I only want certain elements of cHat
: namely predictions at times where a measurement occurs and I have data. So the object that appears in the likelihood function is cHatObs
, not cHat
:
transformed parameters {
...
cHat = col(x, 2) ./ V1;
for(i in 1:nObs){
cHatObs[i] = cHat[iObs[i]]; ## predictions for observed data records
}
}
model {
logCObs ~ normal(log(cHatObs), sigma);
}
While cHatObs
should be a vector, the type of cHat
doesn’t really matter, since I’ll be going through it one element at a time. It might be more difficult to construct from a 2D array or an array of (row-major) vectors.
All that said, if Torsten returns a 2D array, that would be more consistent with what the ODE integrators do. Was there any reason when designing the ODE solvers for returning a 2D array instead of a matrix?
thank you!