Design evaluation
PK model from the library of models
PFIM provides a library of pre-implemented standard
PK/PD structural models. Using a library model is preferable to a
user-defined ODE when the structural form is standard, because it avoids
transcription errors and uses pre-computed closed-form sensitivity
equations, which are faster and more numerically stable than ODE-based
sensitivity computation.
The list modelFromLibrary is passed to
Evaluation() or Optimization() in place of
modelEquations (the two arguments are mutually exclusive).
The key "PKModel" is used for pharmacokinetic models.
The string "Linear1InfusionSingleDose_ClV" decodes
as:
| Token | Meaning |
|---|---|
Linear |
First-order (linear) elimination |
1 |
One-compartment disposition |
Infusion |
IV infusion route; requires Tinf in
Administration() |
ClV |
Parameterised by clearance \(Cl\) and volume \(V\) |
The predicted concentration during and after infusion follows superposition across repeated doses:
\[C_c(t) = \begin{cases} \dfrac{D}{T_{inf} \cdot Cl}\!\left(1 - e^{-\frac{Cl}{V}t}\right) & 0 \leq t \leq T_{inf} \\[6pt] C_c(T_{inf})\,e^{-\frac{Cl}{V}(t - T_{inf})} & t > T_{inf} \end{cases}\]
Model parameters
The one-compartment model has two structural parameters to estimate, both log-normally distributed:
| Parameter | Description | \(\mu\) | \(\omega\) |
|---|---|---|---|
| \(V\) | Volume of distribution (L) | 50 | \(\sqrt{0.26} \approx 0.51\) |
| \(Cl\) | Elimination clearance (L/h) | 5 | \(\sqrt{0.34} \approx 0.58\) |
omega is expressed as the standard deviation on the log
scale. Unlike Example 01, no parameters are fixed here (all have
inter-individual variability).
Residual error model
A Combined1 (proportional + additive) error model is
used for the PK response:
\[\mathrm{SD}(\varepsilon) = \sigma_{\mathrm{inter}} + \sigma_{\mathrm{slope}} \cdot f(\theta, \xi)\]
The additive term dominates at trough concentrations (between doses), while the proportional term prevails near \(C_{\max}\).
Administration parameters
The dosing schedule encodes the full multiple-dose regimen. The
attribute Tinf is a numeric vector of infusion durations
(hours), one per dose event; it is required for library models of the
"Infusion" type.
The loading dose (400 mg = 2× maintenance) rapidly raises \(C_c\) toward the therapeutic target; subsequent 200 mg doses maintain near-steady-state concentrations.
Sampling times
Six sampling times cover distinct phases of the multi-dose concentration profile:
| Time (h) | Pharmacokinetic phase |
|---|---|
| 1 | End of first infusion — \(C_{\max}\) of dose 1 |
| 12 | Mid-interdose Day 1 — early elimination |
| 24 | Trough before dose 2 — \(C_{\min}\) after dose 1 |
| 44 | Near \(C_{\max}\) of dose 2 (dose at 24 h + ~20 h) |
| 72 | Trough before dose 4 — approaching steady state |
| 120 | 24 h after the last dose — post-steady-state elimination |
Arm and design
A single arm with all 150 subjects on the same regimen (parallel
single-group design). No initialCondition is required: the
library model sets \(C_c(0) = 0\)
internally for IV infusion models.
Define an arm called arm1 of size 150 with
administration administrationRespPK and samplings
samplingTimesRespPK
FIM types: population, individual, and Bayesian
Three Evaluation objects are created, each identical
except for fimType. They differ in which parameter
subvector enters the FIM and how inter-individual variability (IIV) is
handled:
Population FIM
(fimType = "population"). The parameter vector includes
fixed effects (\(\mu_V\), \(\mu_{Cl}\)) and variance components (\(\omega^2_V\), \(\omega^2_{Cl}\)). This is the standard
choice for NLMEM study design.
Individual FIM
(fimType = "individual"). Only fixed effects (\(\mu_V\), \(\mu_{Cl}\)) are estimated; \(\omega^2\) are fixed at their specified
values. This is relevant when IIV is already well-characterized from a
prior study. It tends to give larger (more favorable) FIM values because
it ignores the burden of estimating IIV.
Bayesian FIM (fimType = "Bayesian").
The FIM is augmented by a prior precision matrix derived from \(\omega^2\), reflecting the regularization
imposed by MAP estimation. This is appropriate for therapeutic drug
monitoring or sparse individual PK designs. Shrinkage values are
reported for the Bayesian FIM.
Comparing SE and RSE across the three fimType values for
the same design quantifies the effect of IIV estimation on parameter
precision.
Evaluate the population, individual and Bayesian FIMs
evaluationPop = Evaluation( name = "",
modelFromLibrary = modelFromLibrary,
modelParameters = modelParameters,
modelError = modelError,
outputs = list( "RespPK" ),
designs = list( design1 ),
fimType = "population",
odeSolverParameters = list( atol = 1e-8, rtol = 1e-8 ) )
evaluationFIMPop = run( evaluationPop )
evaluationInd = Evaluation( name = "",
modelFromLibrary = modelFromLibrary,
modelParameters = modelParameters,
modelError = modelError,
outputs = list( "RespPK" ),
designs = list( design1 ),
fimType = "individual",
odeSolverParameters = list( atol = 1e-8, rtol = 1e-8 ) )
evaluationFIMInd = run( evaluationInd )
evaluationBay = Evaluation( name = "",
modelFromLibrary = modelFromLibrary,
modelParameters = modelParameters,
modelError = modelError,
outputs = list( "RespPK" ),
designs = list( design1 ),
fimType = "Bayesian",
odeSolverParameters = list( atol = 1e-8, rtol = 1e-8 ) )
evaluationFIMBay = run( evaluationBay )Display the results of the design evaluations
show( evaluationFIMPop )
show( evaluationFIMInd )
show( evaluationFIMBay )
fisherMatrix = getFisherMatrix( evaluationFIMPop )
getCorrelationMatrix( evaluationFIMPop )
getSE( evaluationFIMPop )
getRSE( evaluationFIMPop )
getShrinkage( evaluationFIMPop )
getDeterminant( evaluationFIMPop )
getDcriterion( evaluationFIMPop )
plotOptions = list( unitTime = c("hour"), unitOutcomes = c("mcg/mL") )
plotsEval2_eval = plotEvaluation( evaluationFIMPop, plotOptions )
plotsEval2_si = plotSensitivityIndices( evaluationFIMPop, plotOptions )
plotOutcomesEvaluationRespPK = plotsEval2_eval$design1$arm1$RespPK
plotSensitivityIndice_RespPK_V = plotsEval2_si$design1$arm1$RespPK$V
plotSensitivityIndice_RespPK_Cl = plotsEval2_si$design1$arm1$RespPK$Cl
plotEval_SE = PFIM::plotSE( evaluationFIMPop )
plotEval_RSE = PFIM::plotRSE( evaluationFIMPop )
***************************************
Population Fisher Matrix
***************************************
μ_V μ_Cl ω²_V ω²_Cl σ_inter_RespPK σ_slope_RespPK
μ_V 0.1393880 -0.2899109 0.00000 0.00000 0.00000 0.0000
μ_Cl -0.2899109 14.3204999 0.00000 0.00000 0.00000 0.0000
ω²_V 0.0000000 0.0000000 404.77110 17.51007 51.06227 259.6216
ω²_Cl 0.0000000 0.0000000 17.51007 427.24316 54.71243 80.1080
σ_inter_RespPK 0.0000000 0.0000000 51.06227 54.71243 1982.21804 1361.3658
σ_slope_RespPK 0.0000000 0.0000000 259.62165 80.10800 1361.36582 1648.5030
***************************************
Fixed effects (μ)
***************************************
μ_V μ_Cl
μ_V 0.1393880 -0.2899109
μ_Cl -0.2899109 14.3204999
***************************************
Variance components (ω², γ², σ)
***************************************
ω²_V ω²_Cl σ_inter_RespPK σ_slope_RespPK
ω²_V 404.77110 17.51007 51.06227 259.6216
ω²_Cl 17.51007 427.24316 54.71243 80.1080
σ_inter_RespPK 51.06227 54.71243 1982.21804 1361.3658
σ_slope_RespPK 259.62165 80.10800 1361.36582 1648.5030
*********************************************
Determinant, condition numbers and D-criterion
***********************************************
Determinant: 380849659428
D-criterion: 85.13842
Condition number (fixed effects): 107.3432
Condition number (variance components): 12.64203
***************************************
Parameters estimation
***************************************
Parameter Value SE RSE(%)
μ_V 50.0000000 2.73670901 5.473418
μ_Cl 5.0000000 0.26999904 5.399981
ω²_V 0.2600000 0.05481857 21.084067
ω²_Cl 0.3400000 0.04861157 14.297522
σ_inter_RespPK 0.5000000 0.03570392 7.140784
σ_slope_RespPK 0.3872983 0.04131312 10.667003
[1] 85.13842
***************************************
Individual Fisher Matrix
***************************************
μ_V μ_Cl σ_inter_RespPK σ_slope_RespPK
μ_V 0.003831386 -0.03578704 0.00000 0.00000
μ_Cl -0.035787035 0.74495930 0.00000 0.00000
σ_inter_RespPK 0.000000000 0.00000000 16.79391 12.15818
σ_slope_RespPK 0.000000000 0.00000000 12.15818 20.61786
***************************************
Fixed effects (μ)
***************************************
μ_V μ_Cl
μ_V 0.003831386 -0.03578704
μ_Cl -0.035787035 0.74495930
***************************************
Variance components (σ)
***************************************
σ_inter_RespPK σ_slope_RespPK
σ_inter_RespPK 16.79391 12.15818
σ_slope_RespPK 12.15818 20.61786
***********************************************
Determinant, condition numbers and D-criterion
***********************************************
Determinant: 0.3122376
D-criterion: 0.7475174
Condition number (fixed effects): 354.3252
Condition number (variance effects): 4.847154
***************************************
Parameters estimation
***************************************
Parameter Value SE RSE(%)
μ_V 50.0000000 21.7585945 43.51719
μ_Cl 5.0000000 1.5604237 31.20847
σ_inter_RespPK 0.5000000 0.3223403 64.46806
σ_slope_RespPK 0.3872983 0.2909168 75.11439
***************************************
Bayesian Fisher Matrix
***************************************
μ_V μ_Cl
μ_V 13.424620 -8.946759
μ_Cl -8.946759 21.565159
***************************************
Fixed effects
***************************************
μ_V μ_Cl
μ_V 13.424620 -8.946759
μ_Cl -8.946759 21.565159
***********************************************
Determinant, condition numbers and D-criterion
***********************************************
Determinant: 209.459560208699
D-criterion: 14.4727177892993
Condition number of the fixed effects: 3.56441811487241
***************************************
Shrinkage
***************************************
μ_V μ_Cl
Shrinkage 39.59854 18.8505
***************************************
Parameters estimation
***************************************
Parameter Value SE RSE(%)
μ_V 50 16.043394 32.08679
μ_Cl 5 1.265817 25.31634
Create and save the report for the design evaluation
plotOptions = list( unitTime = c("hour"), unitOutcomes = c("mcg/mL") )
outputFile = "Example02_EvaluationPopFIM.html"
Report( evaluationFIMPop, paths$reports, outputFile, plotOptions )
outputFile = "Example02_EvaluationIndFIM.html"
Report( evaluationFIMInd, paths$reports, outputFile, plotOptions )
outputFile = "Example02_EvaluationBayFIM.html"
Report( evaluationFIMBay, paths$reports, outputFile, plotOptions )
