There are a bunch of new features out now in the mrgsolve 1.0.0
release. Most of these are related to model syntax (what / how you write
your model in the .mod or .cpp file) rather
than model functionality.
All models will be able to use THETA(n) to refer to
THETAn. For example:
[ nmxml ]
run = 100
root = "cppfile"
[ main ]
double CL = THETA(1) * pow(WT/70, 0.75);You’ll notice that we referred to a nonmem run with
[ nmxml ] so this syntax is meant be one small step toward
making it easier to translate your model from the nonmem control stream
to mrgsolve format. This is pretty simple, but will mean less time
removing parentheses when copying over model code.
autodecPlugins are extra functionality that you can invoke when coding your model. For example, plugins exist for calculating time after dose, adding Rcpp functionality to your model etc and more.
One new plugin is called autodec. This plugin will
automatically declare C++ variables for you. For example, I can code the
snipped from the previous block using autodec plugin and
mrgsolve will discover that CL needs to be declared:
[ plugin ] autodec
[ nmxml ]
run = 100
root = "cppfile"
[ main ]
CL = THETA(1) * pow(WT/70, 0.75);The motivation behind this feature is to make it easier to code models with lots of variables when all your variables are simple double precision floating point numbers.
nm-varsThis plugin provided a nonmem-like syntax for certain model elements, including
Fn, Dn, Rn, and
ALAGnA(n), A_0(n) and DADT(n)T (used for SOLVERTIME in
$ODE)So your mrgsolve model translated from nonmem might look like this
see(modlib("nm-like"))$PROB Model written with some nonmem-like syntax features
$PLUGIN nm-vars autodec
$PARAM
THETA1 = 1, THETA2 = 21, THETA3 = 1.3, WT = 70, F1I = 0.5, D2I = 2
KIN = 100, KOUT = 0.1, IC50 = 10, IMAX = 0.9
$CMT @number 3
$PK
CL = THETA(1) * pow(WT/70, 0.75);
V = THETA(2);
KA = THETA(3);
F1 = F1I;
D2 = D2I;
A_0(3) = KIN / KOUT;
$DES
CP = A(2)/V;
INH = IMAX*CP/(IC50 + CP);
DADT(1) = -KA*A(1);
DADT(2) = KA*A(1) - (CL/V)*A(2);
DADT(3) = KIN * (1-INH) - KOUT * A(3);
$ERROR
CP = A(2)/V;Notice that not all nonmem syntax is supported; just
select data structures that use A and DADT to
refer to compartments and differential equations. Also notice that we
still require ; at the end of each line and we still need
to use pow(base, exponent) to calculate exponent of some
number.
The $CMT block has two new options: @number
and @prefix that let you quickly specify a numbered series
of compartments in the model similar to nonmem. The default
@prefix is A so that this code
$CMT @number 3will put compartments A1, A2,
A3 into your model.
modlib: nm-likeYou can see a model coded with these new syntax features in the
internal model library; the model is called nm-like
mod <- modlib("nm-like")
see(mod).
. Model file: nm-like.cpp
. $PROB Model written with some nonmem-like syntax features
.
. $PLUGIN nm-vars autodec
.
. $PARAM
. THETA1 = 1, THETA2 = 21, THETA3 = 1.3, WT = 70, F1I = 0.5, D2I = 2
. KIN = 100, KOUT = 0.1, IC50 = 10, IMAX = 0.9
.
. $CMT @number 3
.
. $PK
. CL = THETA(1) * pow(WT/70, 0.75);
. V = THETA(2);
. KA = THETA(3);
.
. F1 = F1I;
. D2 = D2I;
. A_0(3) = KIN / KOUT;
.
. $DES
. CP = A(2)/V;
. INH = IMAX*CP/(IC50 + CP);
.
. DADT(1) = -KA*A(1);
. DADT(2) = KA*A(1) - (CL/V)*A(2);
. DADT(3) = KIN * (1-INH) - KOUT * A(3);
.
. $ERROR
. CP = A(2)/V;Here, I’ve coded the THETAs in a parameter block; these
would ordinarily come into the model via $NMXML or
$NMEXT, which automatically import these parameters and
estimates. But note that we can still refer to THETA(1) and
others … this resolves to THETA1 regardless.
When mrgsolve loads your model, it checks the $ODE block
to make sure you have differential equations for every compartment in
your model. A long time ago, I started calling this check
audit and it was something that was on by default but could
be turned off via the call to mread().
Starting with version 1.0.0, audit can be controlled at
the $ODE block level like this
[ ode ] @audit
dxdt_A1 = ...
dxdt_A2 = ...to explicitly tell mrgsolve to audit the equations or
[ ode ] @!audit
dxdt_A1 = ...
dxdt_A2 = ...to disable the audit. The audit is on by default and should only be
turned off when you have written some special construct in
[ ode ] to code the differential equations (e.g. odes are
written programmatically).
The new feature is an ii spacer for sequences of event
objects.
Using event objects, we can create complicated dosing regimens from simpler ones. For example, we might have a week of 50 mg dosing followed by a week of 100 mg dosing and then up to 200 mg dosing. We can accomplish that with an event sequence
First create the individual pieces
a <- ev(amt = 50, ii = 24, addl = 6)
b <- ev(amt = 100, ii = 24, addl = 6)
c <- ev(amt = 200, ii = 24, addl = 6)Then put them in a sequence
data <- ev_seq(a, b, c)
data. Events:
. time amt ii addl cmt evid
. 1 0 50 24 6 1 1
. 2 168 100 24 6 1 1
. 3 336 200 24 6 1 1When they are sequenced, you’ll see that the second piece
(b) starts one dosing interval after the last dose in the
first piece (a).
We can put a 24 hour spacer between a and
b
seq(a, wait = 24, b). Events:
. time amt ii addl cmt evid
. 1 0 50 24 6 1 1
. 2 192 100 24 6 1 1Here, the last dose in a is given, we wait one dosing
interval (24 hours) then wait another 24 hours (via wait)
and then start b.
Rather than using wait, we can use ii to
specify the amount of time from the last dose in a to the
first dose in b. So if we want to wait 3 days between the
last dose in a and the first dose in b
ev_seq(a, ii = 3*24, b). Events:
. time amt ii addl cmt evid
. 1 0 50 24 6 1 1
. 2 216 100 24 6 1 1Notice that this same behavior can be achieved with wait
but using ii might be easier to interpret in some
cases.
Sometimes in your model, you have multiple OMEGA or SIGMA blocks like this
code <- '
$OMEGA @name first
1 2 3
$OMEGA @name second
4
'
mod <- mcode("collapse", code, compile = FALSE)
omat(mod). $first
. [,1] [,2] [,3]
. 1: 1 0 0
. 2: 0 2 0
. 3: 0 0 3
.
. $second
. [,1]
. 4: 4If I want to update these values but I only have a 4x4 matrix
mat. [,1] [,2] [,3] [,4]
. [1,] 0.1 0.0 0.0 0.0
. [2,] 0.0 0.2 0.0 0.0
. [3,] 0.0 0.0 0.3 0.0
. [4,] 0.0 0.0 0.0 0.4I can’t do it; mrgsolve wants a 3x3 matrix for the first slot and 1x1 matrix for the second.
A new function will collapse the model matrix into larger matrices
mod <- collapse_omega(mod)
omat(mod). $...
. [,1] [,2] [,3] [,4]
. 1: 1 0 0 0
. 2: 0 2 0 0
. 3: 0 0 3 0
. 4: 0 0 0 4mod <- update(mod, omat = mat)
omat(mod). $...
. [,1] [,2] [,3] [,4]
. 1: 1 0 0 0
. 2: 0 2 0 0
. 3: 0 0 3 0
. 4: 0 0 0 4You can’t split matrices up into smaller chunks, but if there are more than two blocks, you can select which blocks to join
code <- '
$OMEGA @name first
1 2 3
$OMEGA @name second
4
$OMEGA @name third
5 6 7 8
'
mod <- mcode("collapse2", code, compile = FALSE)
omat(mod). $first
. [,1] [,2] [,3]
. 1: 1 0 0
. 2: 0 2 0
. 3: 0 0 3
.
. $second
. [,1]
. 4: 4
.
. $third
. [,1] [,2] [,3] [,4]
. 5: 5 0 0 0
. 6: 0 6 0 0
. 7: 0 0 7 0
. 8: 0 0 0 8mod <- collapse_omega(mod, range = c(2, NA), name = "remainder")
omat(mod). $first
. [,1] [,2] [,3]
. 1: 1 0 0
. 2: 0 2 0
. 3: 0 0 3
.
. $remainder
. [,1] [,2] [,3] [,4] [,5]
. 4: 4 0 0 0 0
. 5: 0 5 0 0 0
. 6: 0 0 6 0 0
. 7: 0 0 0 7 0
. 8: 0 0 0 0 8We recently rolled out a feature where the user could call
simeta() and pass an integer that would indicate a single
ETA() to update, leaving the others alone. This proved to
be more than a little dangerous when $OMEGA had
off-diagonal elements. We will be deprecating this feature starting with
version 1.0.0. When you load a model that contains
simeta(n) or simeps(n), there will be a
warning about the syntax. The warning includes some code you can write
to silence the warning, but this will be temporary. In a future release,
we will always warn when this syntax us used and eventually generate an
error. It is recommended to always use simeta() or
simeps() when re-simulating random effects.
mrgsolve: mrgsolve.github.io | metrum research group: metrumrg.com