-
Notifications
You must be signed in to change notification settings - Fork 1
Expand file tree
/
Copy pathcomponents.qmd
More file actions
347 lines (252 loc) · 14.3 KB
/
components.qmd
File metadata and controls
347 lines (252 loc) · 14.3 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
# Model components {#sec-model-components}
```{r, include = FALSE}
source("setup.R")
```
This chapter details the different components of a model in `mrgsolve`. Each
component listed here is maintained within a "model object". This is an
updatable S4 object in `R` that contains all of the basic information required
to properly configure and simulate from the model.
## Parameter list {#sec-component-param}
The parameter list is an updatable set of name-value pairs. Referencing the name
of an item in the parameter list will substitute the current value associated
with that name. While the name "parameter" may have a certain connotation in
the modeling world, in `mrgsolve` a "parameter" could be any category of numeric
data: covariates (e.g. `WT`, `AGE`, `SEX`), flags, other numeric data that we
commonly call "parameter" (e.g. `CL` or `VC`).
The parameter list is declared in the code block `$PARAM`. While there may be
multiple `$PARAM` blocks in a model, these are condensed to a single parameter
list stored in the model object. The names and numbers of all parameters in the
model must be declared at the time that the model is compiled. Also, a default
value for each parameter must be declared at model compile time, but the value
of each parameter may be updated in one of several ways.
The parameters in a model object can be queried or updated with the
`param()` function.
See also: @sec-block-param, `?param` in the `R` help system after
loading `mrgsolve`.
### Central role of parameters in planning simulations
The data items in the parameter list are more than just values associated with a
name. When an name is added to the parameter list, that name becomes a key word
that `mrgsolve` will start to recognize in input data sets or when manipulating
the model object.
For example, when you want to include a covariate in the model, say weight
(`WT`), you'll include a column in the data set called `WT` that will indicate
the weight of this or that patient. It is crucial that you also list `WT` in
`$PARAM` with some default value. It helps if that value is sensible too. When
`mrgsolve` receives the data set prior to simulating, the `WT` column is matched
up with the `WT` parameter name. As `mrgsolve` works its way through the input
data set (from person to person or from time to time), the value of `WT` is
updated so that the symbol `WT` in `$MAIN` or `$ODE` or `$TABLE` always points
to the value of `WT`. If the `WT` name is not in the parameter list, it won't
matter if it is in the data set or not. Only listing a name in `$PARAM` gets it
"into the game".
Understanding the parameter update mechanism is very important for planning
complicated simulations with `mrgsolve`. Please see the information in
@sec-datasets and in @sec-topic-parameter-update.
## Compartment list {#sec-component-init}
Like the parameter list, the compartment list is a series of name-value pairs.
The compartment list defines the number, names, and initial values of each
compartment in the model. The names, numbers, and order of the compartment
in a model is established at the time of model compile and changes to the
compartment list require re-compilation of the model.
Compartments are declared in one of two code blocks: `$INIT` and `$CMT`. Nominal
initial values must be supplied for each compartment. The main difference
between `$INIT` and `$CMT` is that `$CMT` assumes a default initial value of 0
for each compartment; thus only compartment names are entered. When using
`$INIT`, both names and values must be explicitly stated for each compartment.
The initial values for each compartment can be queried with the `init()`
function. There are several different ways to set the initial conditions in a
model; @sec-topic-init illustrates several of these.
See also: @sec-topic-init and `?init` in the `R` help system after
loading `mrgsolve`.
## Simulation time grid {#sec-component-stime}
The `mrgsolve` model object stores the parameters for the series of time points
to be output for a a simulation. This is the default output time grid that will
be used if not over-ridden by another mechanism.
The elements of the simulation time grid are: `start`, `end`, `delta` and `add`.
`start`, `end`, `delta` are passed to `seq()` as `from`, `to`, and `by`,
respectively. `add` is any arbitrary vector of additional times to simulate.
The simulation time grid in a model object may be queried with the `stime()`
function or by printing the model object to the `R` console.
See also @sec-data-set for discussion of the simulation time grid
and input data sets and @sec-component-tgrid and @sec-topic-tgrid for
using time grid objects .
### `tgrid` objects {#sec-component-tgrid}
A `tgrid` object has `start`, `end`, `delta` and `add` attributes. This object
is independent of the model object. `tgrid` objects may be created and combined
to create complex sampling designs.
See @sec-topic-tgrid for examples and usage.
## Solver settings
`mrgsolve` uses the `DLSODA` solver like the one from `ODEPACK`. Several of the
settings for that solver are stored in the model object and passed to the solver
when the problem is started. Settings include: `atol`, `rtol`, `maxsteps`,
`hmax`, `hmin`, `ixpr`, `mxhnil`.
Solver settings can be changed through the `update()` method for your mrgsolve
model object. For example, to change `rtol` to `1e-5`, write
```{r, eval = FALSE}
mod <- update(mod, rtol = 1e-5)
```
where `mod` is your mrgsolve model object.
### `atol`
Absolute tolerance parameter. Adjust this value lower when you see state variables
(compartments) that are becoming very small and possibly turning negative. For
example:
```{r, message = FALSE}
mod <- modlib("viral1", end = 144)
out <- mrgsim_e(mod, ev(amt = 1000)) %>% filter(V < 0)
out
```
Adjusting `atol` to `1E-20` or `1E-30` will prevent this.
```{r}
mrgsim_e(mod, ev(amt = 1000), atol = 1E-20) %>% filter(time %in% out$time)
```
Setting `atol` in this way will apply to _every_ compartment in the model.
### `rtol`
Relative tolerance parameter. Adjust this value lower when you want more
precision around the calculation of state variables as the system advances.
Setting `rtol` in this way will apply to _every_ compartment in the model.
### Custom `rtol` and `atol`
Starting with mrgsolve 1.6.1, each compartment can have a customized `rtol` or
`atol` value. These customized `rtol` and `atol` values are _only_ used when
advancing the ODE system.
See the following functions for an overview of the interface for setting these
values:
- `custom_tol()`, `custom_rtol()`, `custom_atol()`
- `reset_tol()` `reset_rtol()`, `reset_atol()`
- `use_custom_tol()`, `use_scalar_tol()`
- `get_tol()`, `get_tol_list()`
Release notes are here: [https://github.com/metrumresearchgroup/mrgsolve/releases/tag/1.6.0](https://github.com/metrumresearchgroup/mrgsolve/releases/tag/1.6.0)
### `maxsteps`
This is the maximum number of steps the solver will take when advancing from one
time to the next. If the solver can't make it in `maxsteps` it will stop and
give an error message like this:
```{c,eval=FALSE}
DLSODA- At current T (=R1), MXSTEP (=I1) steps
taken on this call before reaching TOUT
In above message, I =
[1] 2000
In above message, R =
[1] 0.0004049985
DLSODA- ISTATE (=I1) illegal.
In above message, I =
[1] -1
DLSODA- Run aborted.. apparent infinite loop.
Error in (function (x, data, idata = null_idata, carry.out = character(0), :
error from XERRWD
```
You might see this when you have to integrate along time between records in a
data set. There isn't necessarily a problem, but the solver might have to
advance over many doses to get to the next record and it only has a limited
number of steps it can take between those records before it stops with this
error.
When you see this, increase `maxsteps` to 50000 or larger.
But keep in mind that sometimes the solver can't make it to the next record
because there are issues with the model. It might take thousands of steps to
make it 24 hours down the road. In that case, go back to the model code and
look for problems in how it is coded.
### `hmax`
The __maximum__ step size. By default, the solver will take steps of different
sizes based on what is happening in the simulation. Setting `hmax` tells the
solver not to take a step larger than that value. So in a model where `time` is
in hours, reducing `hmax` to `0.1` will prevent the solver from taking a step
larger than `0.1` hours as it tries to advance to the next time. The will slow
down the simulation a bit. But sometimes helpful when the solver starts taking
large steps. We don't recommend using this routinely; for most applications, it
should be reserved for troubleshooting situations. If your model doesn't give
the results that you want without setting `hmax`, we'd recommend a new setup
where this isn't needed.
### `hmin`
The __minimum__ step size. Only set this if you know what you're doing.
### `ixpr`
A flag to enable printing messages to the R console when the solver
switches between non-stiff and stiff solving modes. Rarely used.
### `mxhnil`
The maximum number of messages printed when the model is solving. If you
have a lot of messages, keep working on your model code.
## Functions \ {#sec-model-functions}
There are four `C++` functions that `mrgsolve` creates and manages: `PREAMBLE`,
`MAIN`, `ODE`, `TABLE`. Each function is created from an entire code block in
the model specification file. The user is responsible for writing correct `C++`
code in each of these blocks. mrgsolve will parse these blocks and augment this
code with the necessary elements to create the `C++` function.
These functions may be specified in any order in the model specification file,
but there is a __specific calling order__ for these functions. Recognizing
and understanding this calling order will help understand how the different
pieces of the model specification fit together.
Just prior to starting the problem, mrgsolve calls `$PREAMBLE`. Then, during
advance from time `T1` to `T2`, first `$MAIN` is called, then `$ODE` is called
repeatedly as the solver finds the values of state variables at `T2`, and, once
the solution is found, `$TABLE` is called to calculate derived quantities at
`T2` and to specify variables that should be included in the model output. So,
it is helpful to write model specification files in the order:
1. `$PREAMBLE` - called __only once__ just prior to processing the first record
of the data set
1. `$MAIN` - __before__ advancing the system
1. `$ODE` - the system __advances__ to `T2`
1. `$TABLE` - __after__ advancing the system
But the order in which they are coded will not affect model compilation or the
simulation result.
### The `$PREAMBLE` function
The `PREAMBLE` function gets called only once, just prior to processing the
first record of the data set. This function is composed of `C++` code and is
used to initialize variables and get them set up prior to starting on
the problem.
See @sec-block-preamble for details.
### The `$MAIN` function
The `MAIN` function gets called at least once before the the solver advances
from the current time (`T1`) to the next time (`T2`). In the `MAIN` function,
the user may:
* Set initial conditions for any compartment
* Derive new variables to be used in the model
* Write covariate models
* Add between-subject variability to quantities to structural model
parameters (e.g. `CL` or `VC`).
In addition to getting called once per record, the `MAIN` function may be
called several times prior to starting the simulation run. The `MAIN`
function is also called whenever the user queries the compartment list.
mrgsolve allows you access compartment amounts in the `MAIN` function. But
it is important to remember the calling order of these model functions
(@sec-model-functions): because `MAIN` is called _before_ the system
advances, the compartment amounts inside this function will reflect the
pre-advance values. This is in contrast to accessing compartment amounts
inside the `TABLE` function, which will reflect the values _after_ the
system advances. While there are some use cases where it is useful
to check the pre-advance compartment amounts in `MAIN`, most applications
should interact with compartment amounts after the system advances in
`TABLE`.
The `MAIN` block contains code that is equivalent to the code you'd write in
NONMEM's `$PK` block. You can use either `$MAIN` or `$PK` in your
mrgsolve model.
See @sec-block-main and @sec-block-pk for details.
### The `$ODE` function
The `ODE` function is where the user writes the model differential equations.
Any derived quantity that depends on a state variable and is used to
advanced the system must be calculated inside `$ODE`. But, this
function is called repeatedly during the simulation run, so any calculation
that __can__ be moved out of `$ODE` (for example: to `$MAIN`) should be.
The `ODE` block contains code that is equivalent to the code you'd write in
NONMEM's `$DES` block. You can use either `$ODE` or `$DES` in your mrgsolve
model.
See @sec-block-ode and @sec-block-des for details.
### The `$TABLE` function
The `TABLE` function is called __after__ the solver advances in time. The
purpose of `TABLE` is to allow the user to interact with the values of the
state variables after advancing, potentially derive new variables, and to
insert different outputs into the table of simulated results.
The `TABLE` block contains code that is equivalent to the code you'd write in
NONMEM's `$ERROR` block. You can use either `$TABLE` or `$ERROR` in your
mrgsolve model.
See @sec-block-table and @sec-block-error for details.
## Random effect variances
The `mrgsolve` model object keeps track of a arbitrary number of block matrices
that are used to simulate variates from multivariate normal distributions.
Users can specify `OMEGA` matrices for simulating between-subject random effects
(one draw per individual) or `SIGMA` matrices for simulating within-subject
random effects (one draw per observation).
See @sec-matrix-chapter for complete details on how to work with `OMEGA` and
`SIGMA` in your mrgsolve model.
### `OMEGA`
The matrices are specified in `$OMEGA` blocks in the model specification file.
`OMEGA` may be queried or updated with the `omat()` function.
### `SIGMA`
The matrices are specified in `$SIGMA` blocks in the model specification file.
`SIGMA` may be queried or updated by the `smat()` function.