Just one final update on this—sorry it took me almost a year to post a proper closure on the topic.
@martinmodrak your insights here were much appreciated. It turns out I had some of the parameter units wrong, and the efficiency formula now does yield values bounded between 0 and 1, so there’s no need to constrain a transformed parameter, as you emphasised.
The units are:
K_a and K_e in 1 / minute
D in mg
We can calculate the mass of administered drug, m(t), that was absorbed into the body since time t = 0 as
m(t) = D (1 - e^{-K_at})
So, when t is large, m_a(t) \equiv D when a drug is 100% absorbable. So the efficiency at the peak time, t_p, becomes easily calculated, i.e. the mass of drug that is still in circulation, m_c(t_p), relative to m(t_p). Because
m_c(t_p) = \frac{DK_a\left(e^{-K_et_p} - e^{-K_at_p}\right)}{K_a - K_e}, then
\epsilon(t_p) = m_c(t_p)/m(t_p) becomes independent of D, and is bounded between 0 and 1 no matter what the combinations of K_a and K_e are.
So now \epsilon can either be calculated inside the generated quantities block or outside in R.