vignettes/ParameterTables.Rmd
ParameterTables.RmdParameter tables are a basic and essential look into a model estimate. Due to Nonmem’s flexible model definition approach, getting an informative parameter table can however be an involved process.
Nonmem parameters like THETA(1), OMEGA(2,2), SIGMA(1,1) etc. must be connected with physiologically interpretable parameters such as “Clearance” or “Between-subject variability on central volume” to be informative.
Physiologically interpretable parameters may be transformations of Nonmem parameters.
A manual or even semi-automated approach to establishing these relationships are manually involved to quality check.
NMdata offers two distinct approaches to connecting
Nonmem parameters to physiologially interpretable ones. Which one will
be most useful depends on your work style and how you organize your
control streams.
Flexible, user-defined interpretation of delimited comments in
the parameter sections ($THETA, $OMEGA, and $SIGMA)
(NMreadParText()). This means, if you added comments in
whatever consistent way in these sections, NMreadParText()
is likely to be able to get it into a formatted
data.frame.
Automatic translation of Nonmem code related to direct use of the
parameters (NMrelate()). This may be a pretty good approach
trying to report a Nonmem model without structured labeling in the
parameter sections. If you keep a coding standard, it may allow you to
skip manual parameter translation altogether. And lastly, it may be a
powerful QC tool in any case.
Notice, this vignette zooms in on how to read, label and format
“parameter estimates”. It does not deal with reading input and output
data (see the NMscanData()
vignette). Also, it does not go into how to get a nice print of the
table in a report or a slide deck. For that next step, I have personally
found R packages flextable and pmtables
useful. The scope is just to collect the information.
NMdata contains the following functions to read
estimates from Nonmem:
NMreadExt() A feature-rich processor of ext files.
Provides final parameter estimates with uncertainties and everything
else provided in the ext file, iterations, objective function value and
termination statusNMreadPhi() Reads individual posthoc estimates into a
data.frame
NMreadCov() Reads an estimated variance-covariance
matrix and formats as a matrixNMreadShk() to read shrinkage tablesNMreadExt() is central to the topic of this vignette,
and NMreadShk() is also used. Another interesting function
is NMrelate() which processes the control stream code to
automatically identify parameter names and relate them to parameters
(such as CL to $THETA(1) and
$OMEGA(1,1)). Read on, if interested.
Before getting into how to translate the Nonmem parameters, it is
worth mentioning how to even just read the parameter estimates with
NMdata. A lot of the relevant information is stored in the
.ext file and can be read using
NMreadExt().
file.mod <- system.file("examples/nonmem/xgxr133.mod",package="NMdata")
pars.ext <- NMreadExt(file.mod)This does provide a basic parameter table.
pars.ext[,.(par.name,FIX,value,se)] | par.name | FIX | value | se |
|---|---|---|---|
| THETA(1) | 0 | 0.8028730 | 0.0642182 |
| THETA(2) | 0 | 4.3165400 | 0.0586268 |
| THETA(3) | 0 | 2.5524100 | 0.0916355 |
| THETA(4) | 0 | 5.1554800 | 0.2547830 |
| THETA(5) | 0 | 2.1515900 | 0.0647584 |
| THETA(6) | 0 | 0.4427680 | 0.1405420 |
| OMEGA(1,1) | 1 | 0.0000000 | 0.0000000 |
| OMEGA(2,1) | 1 | 0.0000000 | 0.0000000 |
| OMEGA(2,2) | 0 | 0.1737760 | 0.0312273 |
| OMEGA(3,1) | 1 | 0.0000000 | 0.0000000 |
| OMEGA(3,2) | 0 | -0.0478580 | 0.0254373 |
| OMEGA(3,3) | 0 | 0.2542050 | 0.0433948 |
| OMEGA(4,1) | 1 | 0.0000000 | 0.0000000 |
| OMEGA(4,2) | 1 | 0.0000000 | 0.0000000 |
| OMEGA(4,3) | 1 | 0.0000000 | 0.0000000 |
| OMEGA(4,4) | 1 | 0.0000000 | 0.0000000 |
| OMEGA(5,1) | 1 | 0.0000000 | 0.0000000 |
| OMEGA(5,2) | 1 | 0.0000000 | 0.0000000 |
| OMEGA(5,3) | 1 | 0.0000000 | 0.0000000 |
| OMEGA(5,4) | 1 | 0.0000000 | 0.0000000 |
| OMEGA(5,5) | 1 | 0.0000000 | 0.0000000 |
| SIGMA(1,1) | 0 | 0.0812994 | 0.0054321 |
| SIGMA(2,1) | 1 | 0.0000000 | 0.0000000 |
| SIGMA(2,2) | 1 | 0.0000000 | 0.0000000 |
NMreadExt() does return more than that. These are all
the columns available:
colnames(pars.ext)
#> [1] "model" "TABLENO" "NMREP" "table.step" "par.type"
#> [6] "parameter" "par.name" "i" "j" "FIX"
#> [11] "value" "cond" "eigCor" "partLik" "se"
#> [16] "seStdDevCor" "stdDevCor" "termStat" "est" "iblock"
#> [21] "blocksize"where est is a deprecated copy of the value
column.
NMreadExt() can return other information than the
parameter estimates. This vignette is about creating parameter tables so
for now, let this be a teaser for what else NMreadExt()
provides (which is all the info in .ext files):
The above parameter table is not useful to many, maybe not even to the person who wrote the model. The parameter sections of the control stream contain consistently formatted information we can use.
NMreadSection(file.mod)[c("THETA","OMEGA","SIGMA")]#> $THETA (.1) ; 1 : LTVKA ; Absorption rate [1/h] ; log
#> $THETA (3) ; 2 : LTVV2 ; Central volume [L] ; log
#> $THETA (1) ; 3 : LTVCL ; Clearance [L/h] ; log
#> $THETA (4) ; 4 : LTVV3 ; Peripheral volume [L] ; log
#> $THETA (-1) ; 5 : LTVQ ; Intercomparmental clearance [L/h] ; log
#> $THETA .1 ; 6 : AGEEFF ; Age effect on clearance []; log
#>
#> $OMEGA 0 FIX ; 1 : BSV.KA ; KA Between-subject variability
#> $OMEGA BLOCK(2)
#> 0.1 ; 2 : BSV.V2 ; V2 Between-subject variability
#> 0.01 ; 2-3 : Cov.BSV.V2-CL ; V2-CL Covariance
#> 0.1 ; 3 : BSV.CL ; CL Between-subject variability
#> $OMEGA 0 FIX ; 4 : BSV.V3 ; V3 Between-subject variability
#> $OMEGA 0 FIX ; 5 : BSV.Q ; Q Between-subject variability
#>
#> $SIGMA 0.1 ; SigP - Prop err
#> $SIGMA 0 FIX ; SigA - Add errNMdata::NMreadParsText() provides a very flexible way to
read this information into R. One has to provide the format used in the
control stream in a simple way. Take a look at the annotations shown
above. Following an optional
$THETA/$OMEGA/$SIGMA, each line
contains the initial value specification used by Nonmem, and then
comments containing a counter, a variable name (symbol), and a label
string. $THETA also includes units and transformation
coding (log means the physiological parameter is the
exponentail of the $THETA). SIGMA is coded
simpler, only holding a symbol and a label separated by a dash. We
mimick that in the format arguments to
NMreadParsText().
partab <- NMreadParsText(file.mod,format="%init;%idx:%symbol;%label[%unit];%trans",
format.sigma="%init;%symbol-%label",
as.fun="data.table")
partab| init | idx | symbol | label | unit | trans | par.type | i | j | parameter | model |
|---|---|---|---|---|---|---|---|---|---|---|
| (.1) | 1 | LTVKA | Absorption rate | 1/h | log | THETA | 1 | NA | THETA1 | xgxr133 |
| (3) | 2 | LTVV2 | Central volume | L | log | THETA | 2 | NA | THETA2 | xgxr133 |
| (1) | 3 | LTVCL | Clearance | L/h | log | THETA | 3 | NA | THETA3 | xgxr133 |
| (4) | 4 | LTVV3 | Peripheral volume | L | log | THETA | 4 | NA | THETA4 | xgxr133 |
| (-1) | 5 | LTVQ | Intercomparmental clearance | L/h | log | THETA | 5 | NA | THETA5 | xgxr133 |
| .1 | 6 | AGEEFF | Age effect on clearance | log | THETA | 6 | NA | THETA6 | xgxr133 | |
| 0 FIX | 1 | BSV.KA | KA Between-subject variability | NA | NA | OMEGA | 1 | 1 | OMEGA(1,1) | xgxr133 |
| 0.1 | 2 | BSV.V2 | V2 Between-subject variability | NA | NA | OMEGA | 2 | 2 | OMEGA(2,2) | xgxr133 |
| 0.01 | 2-3 | Cov.BSV.V2-CL | V2-CL Covariance | NA | NA | OMEGA | 3 | 2 | OMEGA(3,2) | xgxr133 |
| 0.1 | 3 | BSV.CL | CL Between-subject variability | NA | NA | OMEGA | 3 | 3 | OMEGA(3,3) | xgxr133 |
| 0 FIX | 4 | BSV.V3 | V3 Between-subject variability | NA | NA | OMEGA | 4 | 4 | OMEGA(4,4) | xgxr133 |
| 0 FIX | 5 | BSV.Q | Q Between-subject variability | NA | NA | OMEGA | 5 | 5 | OMEGA(5,5) | xgxr133 |
| 0.1 | 1 | SigP | Prop err | NA | NA | SIGMA | 1 | 1 | SIGMA(1,1) | xgxr133 |
| 0 FIX | 2 | SigA | Add err | NA | NA | SIGMA | 2 | 2 | SIGMA(2,2) | xgxr133 |
Notice a few things here.
%
followed by alphanumeric characters.Non-alphanumeric characters can be used as delimiters. In this case
; is used as the delimiter most of the time, but not
concistently. After the counter a : was used. And in fact,
the brackets around the units ([1/h]) are delimiters too.
In this case [ is a delimiter on its own, and
] form a multi-charater delimitor together with the
following ;.
We did not specfiy any spaces, and spaces around delimitors are
ignored altogeter by default. You musst use the
spaces.split argument if you want to respect spaces as
delimitors.
In fact, we did not specify two distinct formats for
$THETA and $OMEGA even though
$OMEGA lines do no include units and transformations.
NA’s are filled in where information is missing. However,
you cannot skip delimitors in a line. You can skip contents, but the
structure must be consistent. As an example, see how the unit is skipped
on THETA(A) - AGEEFF.
$THETA,
$OMEGA and $SIGMA
format always applies to $THETA and is
inherited to use for $OMEGA if format.omega is
not supplied. The format used for $OMEGA is inherited for
$SIGMA if format.sigma is not provided.
idx was included, helping the
modeler keeping track of the parameter indexes in the control stream. As
of NMdata 0.1.9 is this by default not used and all i and
j values are automatically identified.In this example, the correlation between two $OMEGA
elements is estimated. NMreadParsText() needs the
off-diagonal elements in $OMEGA and/or $SIGMA
to be on their own lines to connect comments to them. If a diagonal
element and an off-diagonal element is on the same line,
NMreadParsText() will assign the comments to the diagonal
element. It may not be worth adding comments to off-diagonal elements
because they are simply covariances of already labeled diagonal
elements. A simple postprocessing step can therefore be used to add
them.
See more examples in ?NMreadParsText() to better
understand how you can tweak the formats to match your needs.
And we can format a parameter table with labels merging with the
.ext file:
partab2 <- mergeCheck(pars.ext,
partab,
by="parameter",
all.x=TRUE,
common.cols="drop.y",
quiet=TRUE
)
partab2[,value.trans:=value]
partab2[trans=="log",value.trans:=exp(value)]
partab2[,Estimate:=NMcalc::signif2(value.trans,3)][
,RSE:=scales::percent(se/value,accuracy=.1)]
partab.print <- partab2[,.(par.name,symbol,FIX,label,unit,Estimate,RSE=RSE)]
partab.print| par.name | symbol | FIX | label | unit | Estimate | RSE |
|---|---|---|---|---|---|---|
| THETA(1) | LTVKA | 0 | Absorption rate | 1/h | 2.23 | 8.0% |
| THETA(2) | LTVV2 | 0 | Central volume | L | 74.9 | 1.4% |
| THETA(3) | LTVCL | 0 | Clearance | L/h | 12.8 | 3.6% |
| THETA(4) | LTVV3 | 0 | Peripheral volume | L | 173 | 4.9% |
| THETA(5) | LTVQ | 0 | Intercomparmental clearance | L/h | 8.60 | 3.0% |
| THETA(6) | AGEEFF | 0 | Age effect on clearance | 1.56 | 31.7% | |
| OMEGA(1,1) | BSV.KA | 1 | KA Between-subject variability | NA | 0 | NA |
| OMEGA(2,1) | NA | 1 | NA | NA | 0 | NA |
| OMEGA(2,2) | BSV.V2 | 0 | V2 Between-subject variability | NA | 0.174 | 18.0% |
| OMEGA(3,1) | NA | 1 | NA | NA | 0 | NA |
| OMEGA(3,2) | Cov.BSV.V2-CL | 0 | V2-CL Covariance | NA | -0.0479 | -53.2% |
| OMEGA(3,3) | BSV.CL | 0 | CL Between-subject variability | NA | 0.254 | 17.1% |
| OMEGA(4,1) | NA | 1 | NA | NA | 0 | NA |
| OMEGA(4,2) | NA | 1 | NA | NA | 0 | NA |
| OMEGA(4,3) | NA | 1 | NA | NA | 0 | NA |
| OMEGA(4,4) | BSV.V3 | 1 | V3 Between-subject variability | NA | 0 | NA |
| OMEGA(5,1) | NA | 1 | NA | NA | 0 | NA |
| OMEGA(5,2) | NA | 1 | NA | NA | 0 | NA |
| OMEGA(5,3) | NA | 1 | NA | NA | 0 | NA |
| OMEGA(5,4) | NA | 1 | NA | NA | 0 | NA |
| OMEGA(5,5) | BSV.Q | 1 | Q Between-subject variability | NA | 0 | NA |
| SIGMA(1,1) | SigP | 0 | Prop err | NA | 0.0813 | 6.7% |
| SIGMA(2,1) | NA | 1 | NA | NA | 0 | NA |
| SIGMA(2,2) | SigA | 1 | Add err | NA | 0 | NA |
The left join we did with mergeCheck includes all
parameters from the .ext file. In this case we want to skip
all the off-diagonals in $OMEGA and $SIGMA
because they are not estimated. Let’s just look at the ones that were
estimated. In many cases some fixed parameters can be relevant to
include in a parameter table, and then you need to modify this step.
partab.print <- partab.print[FIX==0,.(par.name,symbol,label,unit,Estimate,RSE)] We need to explain that RSE’s were calculated on the log scale for the transformed parameters. We could also want to add confidence intervals based on the estimated standard error on the log scale. But we are getting close.
Shrinkage estimates can be taken from the .ext file
using NMreadShk(). Notice, depending on the Nonmem version
multiple different shrinkage estimates are provided (and fetched by
NMreadShk). Here, we just include the SD-based
\\eta-shrinkages.
shk <- NMreadShk(file.mod)
partab.print2 <- mergeCheck(partab.print,shk[,.(par.name,Pval=signif2(Pval,2),etabar,EtaShkSD=percent(EtaShkSD/100))],by="par.name",all.x=TRUE)
#> Column(s) added: Pval, etabar, EtaShkSD
partab.print2| par.name | symbol | label | unit | Estimate | RSE | Pval | etabar | EtaShkSD |
|---|---|---|---|---|---|---|---|---|
| THETA(1) | LTVKA | Absorption rate | 1/h | 2.23 | 8.0% | NA | NA | NA |
| THETA(2) | LTVV2 | Central volume | L | 74.9 | 1.4% | NA | NA | NA |
| THETA(3) | LTVCL | Clearance | L/h | 12.8 | 3.6% | NA | NA | NA |
| THETA(4) | LTVV3 | Peripheral volume | L | 173 | 4.9% | NA | NA | NA |
| THETA(5) | LTVQ | Intercomparmental clearance | L/h | 8.60 | 3.0% | NA | NA | NA |
| THETA(6) | AGEEFF | Age effect on clearance | 1.56 | 31.7% | NA | NA | NA | |
| OMEGA(2,2) | BSV.V2 | V2 Between-subject variability | NA | 0.174 | 18.0% | 0.97 | -0.0014793 | 7.3% |
| OMEGA(3,2) | Cov.BSV.V2-CL | V2-CL Covariance | NA | -0.0479 | -53.2% | NA | NA | NA |
| OMEGA(3,3) | BSV.CL | CL Between-subject variability | NA | 0.254 | 17.1% | 0.99 | -0.0006435 | 10.1% |
| SIGMA(1,1) | SigP | Prop err | NA | 0.0813 | 6.7% | NA | NA | NA |
A simple function addCor() is available to calculate
estimated correlations of random effects. addCor() simply
applies stats::cov2cor() to get the correlations from the
variance-covariance matrices $SIGMA and
$OMEGA.
Partab4 <- NMdata::addCor(partab2)
Partab4[FIX==0,.(par.name,symbol,FIX,label,unit,Estimate,RSE=RSE,corr)]| par.name | symbol | FIX | label | unit | Estimate | RSE | corr |
|---|---|---|---|---|---|---|---|
| THETA(1) | LTVKA | 0 | Absorption rate | 1/h | 2.23 | 8.0% | NA |
| THETA(2) | LTVV2 | 0 | Central volume | L | 74.9 | 1.4% | NA |
| THETA(3) | LTVCL | 0 | Clearance | L/h | 12.8 | 3.6% | NA |
| THETA(4) | LTVV3 | 0 | Peripheral volume | L | 173 | 4.9% | NA |
| THETA(5) | LTVQ | 0 | Intercomparmental clearance | L/h | 8.60 | 3.0% | NA |
| THETA(6) | AGEEFF | 0 | Age effect on clearance | 1.56 | 31.7% | NA | |
| OMEGA(2,2) | BSV.V2 | 0 | V2 Between-subject variability | NA | 0.174 | 18.0% | 1.0000000 |
| OMEGA(3,2) | Cov.BSV.V2-CL | 0 | V2-CL Covariance | NA | -0.0479 | -53.2% | -0.2277024 |
| OMEGA(3,3) | BSV.CL | 0 | CL Between-subject variability | NA | 0.254 | 17.1% | 1.0000000 |
| SIGMA(1,1) | SigP | 0 | Prop err | NA | 0.0813 | 6.7% | 1.0000000 |
Use the function NMreadInits() to break down the initial
values, bounds and “FIX” information as shown in the control stream.
This can easily be merged onto the parameter table as needed.
NMreadInits(file.mod)| par.type | parameter | par.name | i | j | iblock | blocksize | init | lower | upper | FIX |
|---|---|---|---|---|---|---|---|---|---|---|
| THETA | THETA1 | THETA(1) | 1 | NA | 1 | 1 | 0.10 | NA | NA | 0 |
| THETA | THETA2 | THETA(2) | 2 | NA | 2 | 1 | 3.00 | NA | NA | 0 |
| THETA | THETA3 | THETA(3) | 3 | NA | 3 | 1 | 1.00 | NA | NA | 0 |
| THETA | THETA4 | THETA(4) | 4 | NA | 4 | 1 | 4.00 | NA | NA | 0 |
| THETA | THETA5 | THETA(5) | 5 | NA | 5 | 1 | -1.00 | NA | NA | 0 |
| THETA | THETA6 | THETA(6) | 6 | NA | 6 | 1 | 0.10 | NA | NA | 0 |
| OMEGA | OMEGA(1,1) | OMEGA(1,1) | 1 | 1 | 1 | 1 | 0.00 | NA | NA | 1 |
| OMEGA | OMEGA(2,2) | OMEGA(2,2) | 2 | 2 | 2 | 2 | 0.10 | NA | NA | 0 |
| OMEGA | OMEGA(3,2) | OMEGA(3,2) | 3 | 2 | 2 | 2 | 0.01 | NA | NA | 0 |
| OMEGA | OMEGA(3,3) | OMEGA(3,3) | 3 | 3 | 2 | 2 | 0.10 | NA | NA | 0 |
| OMEGA | OMEGA(4,4) | OMEGA(4,4) | 4 | 4 | 4 | 1 | 0.00 | NA | NA | 1 |
| OMEGA | OMEGA(5,5) | OMEGA(5,5) | 5 | 5 | 5 | 1 | 0.00 | NA | NA | 1 |
| SIGMA | SIGMA(1,1) | SIGMA(1,1) | 1 | 1 | 1 | 1 | 0.10 | NA | NA | 0 |
| SIGMA | SIGMA(2,2) | SIGMA(2,2) | 2 | 2 | 2 | 1 | 0.00 | NA | NA | 1 |
The problem with the otherwise nicely looking table generated above is that we spotted a typo in one label. A “t” is missing in “Intercompartmental”. How are we going to fix that? The model was finalized so we don’t want to edit the control stream. For this specific problem we could fix it with a simple text replacement. But such patchwork can become messy with a large model. While the method of providing all the details in control stream may look like a great solution at first, it has problems
A practical issue is the tedious process keeping all the candidate models up to date as the information for the parameter tables improves. The parameter table work is often part of a model reporting step following model development. Since the control streams are created in the model development step, the parameter table fine tuning is difficult to integrate. Briefly, we need a way to update the parameter tables after finalizing the control stream.
tab.offline <- read.csv("derived/partab_offline.csv")
mergeCheck(partab2[FIX==0],tab.offline,by="symbol",common.cols="drop.x",all.x=T)
#> Column(s) added: label, unit| model | TABLENO | NMREP | table.step | par.type | parameter | par.name | i | j | FIX | value | cond | eigCor | partLik | se | seStdDevCor | stdDevCor | termStat | est | iblock | blocksize | init | idx | symbol | trans | value.trans | Estimate | RSE | label | unit |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| xgxr133 | 2 | 1 | IMP | THETA | THETA1 | THETA(1) | 1 | NA | 0 | 0.8028730 | 22.32860 | 0.129310 | 0.6574150 | 0.0642182 | 0.0000000 | 0.000000 | 0 | 0.8028730 | NA | NA | (.1) | 1 | LTVKA | log | 2.2319441 | 2.23 | 8.0% | Absorption rate | 1/h |
| xgxr133 | 2 | 1 | IMP | THETA | THETA2 | THETA(2) | 2 | NA | 0 | 4.3165400 | 0.12931 | 0.401499 | -0.8741820 | 0.0586268 | 0.0000000 | 0.000000 | 0 | 4.3165400 | NA | NA | (3) | 2 | LTVV2 | log | 74.9289252 | 74.9 | 1.4% | Central volume | L |
| xgxr133 | 2 | 1 | IMP | THETA | THETA3 | THETA(3) | 3 | NA | 0 | 2.5524100 | 2.88730 | 0.429168 | -0.3923980 | 0.0916355 | 0.0000000 | 0.000000 | 0 | 2.5524100 | NA | NA | (1) | 3 | LTVCL | log | 12.8380061 | 12.8 | 3.6% | Clearance | L/h |
| xgxr133 | 2 | 1 | IMP | THETA | THETA4 | THETA(4) | 4 | NA | 0 | 5.1554800 | 0.00000 | 0.481360 | -2.6760700 | 0.2547830 | 0.0000000 | 0.000000 | 0 | 5.1554800 | NA | NA | (4) | 4 | LTVV3 | log | 173.3790087 | 173 | 4.9% | Peripheral volume | L |
| xgxr133 | 2 | 1 | IMP | THETA | THETA5 | THETA(5) | 5 | NA | 0 | 2.1515900 | 0.00000 | 0.773260 | 3.3035300 | 0.0647584 | 0.0000000 | 0.000000 | 0 | 2.1515900 | NA | NA | (-1) | 5 | LTVQ | log | 8.5985192 | 8.60 | 3.0% | Intercompartmental clearance | L/h |
| xgxr133 | 2 | 1 | IMP | THETA | THETA6 | THETA(6) | 6 | NA | 0 | 0.4427680 | 0.00000 | 0.870593 | 0.0673492 | 0.1405420 | 0.0000000 | 0.000000 | 0 | 0.4427680 | NA | NA | .1 | 6 | AGEEFF | log | 1.5570111 | 1.56 | 31.7% | Age effect on clearance | |
| xgxr133 | 2 | 1 | IMP | OMEGA | OMEGA(2,2) | OMEGA(2,2) | 2 | 2 | 0 | 0.1737760 | 0.00000 | 0.000000 | 0.6164530 | 0.0312273 | 0.0374549 | 0.416864 | 0 | 0.1737760 | 2 | 2 | 0.1 | 2 | BSV.V2 | NA | 0.1737760 | 0.174 | 18.0% | Between-subject variability - Central volume | |
| xgxr133 | 2 | 1 | IMP | OMEGA | OMEGA(3,2) | OMEGA(3,2) | 3 | 2 | 0 | -0.0478580 | 0.00000 | 0.000000 | 1.2065900 | 0.0254373 | 0.1135310 | -0.227702 | 0 | -0.0478580 | 2 | 2 | 0.01 | 2-3 | Cov.BSV.V2-CL | NA | -0.0478580 | -0.0479 | -53.2% | NA | NA |
| xgxr133 | 2 | 1 | IMP | OMEGA | OMEGA(3,3) | OMEGA(3,3) | 3 | 3 | 0 | 0.2542050 | 0.00000 | 0.000000 | 0.7464250 | 0.0433948 | 0.0430343 | 0.504188 | 0 | 0.2542050 | 2 | 2 | 0.1 | 3 | BSV.CL | NA | 0.2542050 | 0.254 | 17.1% | Between-subject variability - Clearance | |
| xgxr133 | 2 | 1 | IMP | SIGMA | SIGMA(1,1) | SIGMA(1,1) | 1 | 1 | 0 | 0.0812994 | 0.00000 | 1.107170 | 16.4170000 | 0.0054321 | 0.0095256 | 0.285130 | 0 | 0.0812994 | 1 | 1 | 0.1 | 1 | SigP | NA | 0.0812994 | 0.0813 | 6.7% | NA | NA |
$PK or
$PRED code
In the last example, all we used from the control stream was the
symbol, i.e. an abbreviated parameter name. The rest was
merged in from a table maintained in a csv file.
labs.auto <- NMrelate(file.mod)
labs.auto[,symbol:=label]
labs.auto[par.type=="OMEGA",symbol:=paste0("BSV.",symbol)]
labs.auto| model | par.name | par.type | i | j | nrep.LHS | nrep.par | LHS | label | code | symbol |
|---|---|---|---|---|---|---|---|---|---|---|
| xgxr133 | THETA(1) | THETA | 1 | NA | 1 | 1 | LTVKA | LTVKA | LTVKA=THETA(1) | LTVKA |
| xgxr133 | THETA(2) | THETA | 2 | NA | 1 | 1 | LTVV2 | LTVV2 | LTVV2=THETA(2) | LTVV2 |
| xgxr133 | THETA(3) | THETA | 3 | NA | 1 | 1 | LTVCL | LTVCL | LTVCL=THETA(3) | LTVCL |
| xgxr133 | THETA(4) | THETA | 4 | NA | 1 | 1 | LTVV3 | LTVV3 | LTVV3=THETA(4) | LTVV3 |
| xgxr133 | THETA(5) | THETA | 5 | NA | 1 | 1 | LTVQ | LTVQ | LTVQ=THETA(5) | LTVQ |
| xgxr133 | THETA(6) | THETA | 6 | NA | 1 | 1 | AGEEFF | AGEEFF | AGEEFF=THETA(6) | AGEEFF |
| xgxr133 | OMEGA(1,1) | OMEGA | 1 | 1 | 1 | 1 | KA | KA | KA=EXP(MU_1+ETA(1)) | BSV.KA |
| xgxr133 | OMEGA(2,2) | OMEGA | 2 | 2 | 1 | 1 | V2 | V2 | V2=EXP(MU_2+ETA(2)) | BSV.V2 |
| xgxr133 | OMEGA(3,3) | OMEGA | 3 | 3 | 1 | 1 | CL | CL | CL=EXP(MU_3+ETA(3)) | BSV.CL |
| xgxr133 | OMEGA(4,4) | OMEGA | 4 | 4 | 1 | 1 | V3 | V3 | V3=EXP(MU_4+ETA(4)) | BSV.V3 |
| xgxr133 | OMEGA(5,5) | OMEGA | 5 | 5 | 1 | 1 | Q | Q | Q=EXP(MU_5+ETA(5)) | BSV.Q |
| xgxr133 | SIGMA(1,1) | SIGMA | 1 | 1 | 1 | 2 | SIGP | SIGP | SIGP=SIGMA(1,1) | SIGP |
| xgxr133 | SIGMA(1,1) | SIGMA | 1 | 1 | 2 | 2 | Y | Y - SIGMA(1,1) | Y=F+F*ERR(1)+ERR(2) | Y - SIGMA(1,1) |
| xgxr133 | SIGMA(2,2) | SIGMA | 2 | 2 | 1 | 2 | SIGA | SIGA | SIGA=SIGMA(2,2) | SIGA |
| xgxr133 | SIGMA(2,2) | SIGMA | 2 | 2 | 2 | 2 | Y | Y - SIGMA(2,2) | Y=F+F*ERR(1)+ERR(2) | Y - SIGMA(2,2) |
mergeCheck(labs.auto[,.(par.name,symbol)],tab.offline,by="symbol",all.x=T)
#> Column(s) added: label, unit| par.name | symbol | label | unit |
|---|---|---|---|
| THETA(1) | LTVKA | Absorption rate | 1/h |
| THETA(2) | LTVV2 | Central volume | L |
| THETA(3) | LTVCL | Clearance | L/h |
| THETA(4) | LTVV3 | Peripheral volume | L |
| THETA(5) | LTVQ | Intercompartmental clearance | L/h |
| THETA(6) | AGEEFF | Age effect on clearance | |
| OMEGA(1,1) | BSV.KA | Between-subject variability - Absorption rate | |
| OMEGA(2,2) | BSV.V2 | Between-subject variability - Central volume | |
| OMEGA(3,3) | BSV.CL | Between-subject variability - Clearance | |
| OMEGA(4,4) | BSV.V3 | Between-subject variability - Peripheral volume | |
| OMEGA(5,5) | BSV.Q | Between-subject variability - Intercompartmental clearance | |
| SIGMA(1,1) | SIGP | Proportional residual error | |
| SIGMA(1,1) | Y - SIGMA(1,1) | NA | NA |
| SIGMA(2,2) | SIGA | Additive residual error | ng/mL |
| SIGMA(2,2) | Y - SIGMA(2,2) | NA | NA |
file2.mod <- system.file("examples/nonmem/xgxr133.mod",package="NMdata")
pars.ext <- NMreadExt(file2.mod)
labs.auto <- NMrelate(file2.mod)
labs.auto[,symbol:=label]
labs.j <- labs.auto[par.type%in%c("OMEGA","SIGMA"),.(par.type.j=par.type,j,symbol.j=symbol)]
labs.ij <- labs.j[,labs.auto[par.type%in%c("OMEGA","SIGMA"),.(i,par.type.i=par.type,symbol.i=symbol)],by=labs.j]
labs.ij <- labs.ij[par.type.i==par.type.j & (i!=j)]
labs.ij[par.type.i==par.type.j,par.type:=par.type.i]
labs.ij[,symbol:=sprintf("Cov.%s-%s",symbol.i,symbol.j)]
labs.ij[,par.name:=sprintf("%s(%d,%d)",par.type.i,i,j)]
labs.ij[,label:=symbol]
cols.drop <- cc(par.type.j, symbol.j, par.type.i, symbol.i)
labs.ij[,(cols.drop):=NULL]
labs.auto[par.type=="OMEGA"&i==j,symbol:=paste0("BSV.",symbol)]
labs.auto2 <- rbind(labs.auto,labs.ij,fill=T)
### overwrite auto-gen labels with hard-coded ones
tab.labs.diag <- mergeCheck(
## labs.auto[par.type=="THETA",!("label")]
labs.auto2[,!("label")]
## labs.auto[par.type=="THETA"|(i==j),!("label")]
,
tab.offline
,by="symbol"
,all.x=T)
#> Column(s) added: label, unit
tab.labs.offdiag <- labs.auto2[par.type!="THETA"&i!=j]
tab.all.labs <- rbind(
tab.labs.diag
## ,tab.labs.offdiag
,fill=T)
## mege onto estimates (ext)
pars.aug <- mergeCheck(pars.ext[FIX==0],
tab.all.labs[!grepl("SIGMA",symbol)]
,
by="par.name",all.x=TRUE,common.cols="drop.y")
#> Column(s) added: nrep.LHS, nrep.par, LHS, code, symbol, label, unit
### back transformations
pars.aug[,value.trans:=value]
pars.aug[grepl("LTV",symbol),value.trans:=exp(value)]
## pars.aug
pars.aug[,.(par.name,symbol,FIX,label,unit,Estimate=value.trans,RSE=se/value)]| par.name | symbol | FIX | label | unit | Estimate | RSE |
|---|---|---|---|---|---|---|
| THETA(1) | LTVKA | 0 | Absorption rate | 1/h | 2.2319441 | 0.0799855 |
| THETA(2) | LTVV2 | 0 | Central volume | L | 74.9289252 | 0.0135819 |
| THETA(3) | LTVCL | 0 | Clearance | L/h | 12.8380061 | 0.0359016 |
| THETA(4) | LTVV3 | 0 | Peripheral volume | L | 173.3790087 | 0.0494198 |
| THETA(5) | LTVQ | 0 | Intercompartmental clearance | L/h | 8.5985192 | 0.0300979 |
| THETA(6) | AGEEFF | 0 | Age effect on clearance | 0.4427680 | 0.3174168 | |
| OMEGA(2,2) | BSV.V2 | 0 | Between-subject variability - Central volume | 0.1737760 | 0.1796986 | |
| OMEGA(3,2) | Cov.CL-V2 | 0 | Covariance - Clearance and central volume | -0.0478580 | -0.5315162 | |
| OMEGA(3,3) | BSV.CL | 0 | Between-subject variability - Clearance | 0.2542050 | 0.1707079 | |
| SIGMA(1,1) | SIGP | 0 | Proportional residual error | 0.0812994 | 0.0668157 |
Any of the code above can be run on multiple models at once.
NMdata functions generally return tables including a
model column for this purpose. Let’s get a quick overview
of a few models.
files <- system.file(file.path("examples/nonmem",c("xgxr132.mod","xgxr133.mod")),package="NMdata")
labs.auto <- NMrelate(files)
labs.auto[,symbol:=label]
labs.auto[par.type=="OMEGA",symbol:=paste0("BSV.",symbol)]
labs.auto <- labs.auto[!grepl("SIGMA",symbol)]
NMreadExt(files,return="pars")[FIX==0]|>
mergeCheck(labs.auto[,.(model,par.name,symbol)],by=c("model","par.name"),all.x=T) |>
mergeCheck(tab.offline,by="symbol",all.x=T) |>
dcast(par.name+symbol+label+unit~model,value.var="value")
#> Column(s) added: symbol
#> Column(s) added: label, unit| par.name | symbol | label | unit | xgxr132 | xgxr133 |
|---|---|---|---|---|---|
| OMEGA(2,2) | BSV.V2 | Between-subject variability - Central volume | 0.1754540 | 0.1737760 | |
| OMEGA(3,2) | NA | NA | NA | NA | -0.0478580 |
| OMEGA(3,3) | BSV.CL | Between-subject variability - Clearance | 0.3092220 | 0.2542050 | |
| SIGMA(1,1) | SIGP | Proportional residual error | 0.0809318 | 0.0812994 | |
| THETA(1) | LTVKA | Absorption rate | 1/h | 0.8134230 | 0.8028730 |
| THETA(2) | LTVV2 | Central volume | L | 4.3236000 | 4.3165400 |
| THETA(3) | LTVCL | Clearance | L/h | 2.4710300 | 2.5524100 |
| THETA(4) | LTVV3 | Peripheral volume | L | 5.4265500 | 5.1554800 |
| THETA(5) | LTVQ | Intercompartmental clearance | L/h | 2.2274200 | 2.1515900 |
| THETA(6) | AGEEFF | Age effect on clearance | 0.4197470 | 0.4427680 |