Why a Parameter Transfer Needs More Than a List of Estimates
A NONMEM run can provide estimates for typical parameters, interindividual variability, and residual error. Those estimates do not all enter a simulation in the same way. A typical clearance estimate is a fixed-effect value; an OMEGA entry describes the distribution of an individual random effect; and a SIGMA entry concerns residual variability. Treating all three as ordinary parameter values can produce a simulation that runs but represents a different model.
The aim here is to carry a fitted NONMEM parameter set into an mrgsolve model while retaining its meanings. This assumes the structural model and parameterization have already been mapped correctly, including units, covariate relationships, and event definitions. The earlier material on mapping parameter names and values is useful for that preliminary alignment; the focus here is on extracting estimates and using them for typical-subject and population simulation.
Read the Estimate in Its Model Context
Suppose a one-compartment IV bolus model uses lognormally distributed interindividual variability in clearance and volume:
Here \(TVCL\) and \(TVV\) are fixed effects, while \(\Omega\) is a variance-covariance matrix on the log-random-effect scale. With this parameterization, an individual with both random effects equal to zero has \(CL_i=TVCL\) and \(V_i=TVV\). For the lognormal parameters, those fixed effects are medians across the population, not generally arithmetic means.
A NONMEM report may show estimates as THETA values and OMEGA or SIGMA entries. Read the model code and the output labels together: THETA index alone does not tell you the parameter's units or whether a transformation is applied. Likewise, an OMEGA entry is not necessarily a standard deviation. In the model above, a diagonal OMEGA entry is a variance of an \(\eta\), and an off-diagonal entry is a covariance. If the output or a secondary export instead presents standard deviations or correlations, convert them back to the intended covariance matrix before entering that structure in the simulation model.
Worked Example: Extracting and Checking an Estimate Set
Consider a hypothetical NONMEM run for the model above. The final estimates are \(TVCL=4.00\ \mathrm{L/h}\), \(TVV=20.0\ \mathrm{L}\), and \[ \Omega = \begin{pmatrix} 0.0900 & 0.0120\\ 0.0120 & 0.0400 \end{pmatrix}. \] The rows and columns are ordered as \(\eta_{CL}\), then \(\eta_V\). Identify what belongs in the fixed parameters and check the covariance entries.
The fixed parameters are \(TVCL=4.00\ \mathrm{L/h}\) and \(TVV=20.0\ \mathrm{L}\). The diagonal entries are \(\operatorname{Var}(\eta_{CL})=0.0900\) and \(\operatorname{Var}(\eta_V)=0.0400\); the off-diagonal entries are \(\operatorname{Cov}(\eta_{CL},\eta_V)=0.0120\). The corresponding standard deviations are \(\sqrt{0.0900}=0.300\) and \(\sqrt{0.0400}=0.200\), but these standard deviations are not the values to enter when the model expects OMEGA variances.
The implied correlation is \[ \rho=\frac{0.0120}{\sqrt{0.0900}\sqrt{0.0400}} =\frac{0.0120}{0.300(0.200)}=0.200. \] The determinant is \(0.0900(0.0400)-(0.0120)^2=0.003456>0\), and both diagonal entries are positive. Thus this symmetric two-dimensional matrix is positive definite and is a valid covariance matrix for the stated random effects.
Map the Estimates into an mrgsolve Model
For this model, the fixed-effect values become named mrgsolve parameters, and the OMEGA estimates define the random-effect distribution. A compact model specification can be written as follows; the model code is illustrative and assumes the state is the amount in a central compartment.
The order is essential: the first random effect in the model must correspond to the first row and column of \(\Omega\). If the OMEGA order is reversed without changing the parameter equations, the marginal variances may appear plausible while the covariance is attached to the wrong pair of parameters. Keep a written index map beside the model code, especially when the NONMEM run has multiple OMEGA blocks.
Only include a SIGMA estimate when the simulation is intended to include the observation-error process. For example, a proportional residual model has a different interpretation from an additive residual model, and its implementation must match the specified observation equation. Do not place a residual variance in OMEGA or use an OMEGA variance as if it were a residual variance. A deterministic concentration trajectory can be simulated without residual error when the target is the structural prediction rather than simulated observations.
Typical-Subject Simulation
A typical-subject prediction sets the individual random effects to their reference values, usually \(\eta_{CL}=\eta_V=0\). For lognormal parameter models, this makes the individual parameters equal to their typical values. It does not represent the average concentration over a heterogeneous population: nonlinear transformations mean the average of individual predictions need not equal the prediction at typical parameters.
For an IV bolus dose \(D\) into the central compartment, the amount after dose is \(D\exp[-(CL_i/V_i)t]\), and the concentration is that amount divided by \(V_i\). This follows by solving the first-order differential equation \(dCENT/dt=-(CL_i/V_i)CENT\) with initial amount \(CENT(0)=D\). The solution is \(CENT(t)=D\exp[-(CL_i/V_i)t]\), so \(C_i(t)=D/V_i\exp[-(CL_i/V_i)t]\).
Worked Example: Predicting a Typical Subject
Use the hypothetical estimates above, set both \(\eta\) values to zero, and give a \(100\ \mathrm{mg}\) IV bolus at time zero. Calculate the concentration at 8 h.
The typical parameters are \(CL=4.00\ \mathrm{L/h}\) and \(V=20.0\ \mathrm{L}\), so the elimination rate constant is \(CL/V=4.00/20.0=0.200\ \mathrm{h^{-1}}\). Substituting into the concentration equation gives \[ C(8)=\frac{100\ \mathrm{mg}}{20.0\ \mathrm{L}} \exp[-(0.200\ \mathrm{h^{-1}})(8\ \mathrm{h})] =5.00\ \mathrm{mg/L}\times\exp(-1.600) \approx 1.009\ \mathrm{mg/L}. \]
The result is the structural concentration for the typical subject under the stated bolus model. It does not include residual error, and it is not a population-average concentration.
Population Simulation and a Numerical Precision Check
To simulate a population, retain the fixed effects and draw each subject's random-effect vector from the transferred distribution, \(\eta_i\sim N(0,\Omega)\). In mrgsolve, this means the OMEGA block must contain the intended covariance matrix and the simulation must use it as the source of individual random effects. A typical-subject run instead uses zero random effects; do not assume that supplying OMEGA alone makes every subject typical. Check the simulation settings and whether random effects are drawn or supplied as data.
The effect of a draw is applied through the parameter equations. For example, if a sampled subject has \(\eta_{CL}=-0.200\) and \(\eta_V=0.100\), then \(CL_i=4.00\exp(-0.200)=3.275\ \mathrm{L/h}\), and \(V_i=20.0\exp(0.100)=22.10\ \mathrm{L}\). Its concentration after a \(100\ \mathrm{mg}\) bolus at 8 h is \[ C_i(8)=\frac{100}{22.10}\exp\left[-\left(\frac{3.275}{22.10}\right)8\right] \approx 1.383\ \mathrm{mg/L}. \] This is one hypothetical subject conditional on the stated random-effect values, not an estimate of a population summary.
Calculation precision matters when reporting simulated concentrations. For instance, if the unrounded parameter formulas give \(D/V_i\approx2.7634\ \mathrm{mg/L}\) and \(\exp[-(CL_i/V_i)8]\approx0.25925\), then their product is \(2.7634(0.25925)=0.716411\ldots\), which rounds to \(C_i(8)\approx0.7164\ \mathrm{mg/L}\). If the factors are instead rounded to \(2.763\) and \(0.2592\) before multiplication, their product is \(0.7161696\), or \(0.7162\ \mathrm{mg/L}\). Use full precision for the calculation and round only the final reported result.
Derivation and Checks for the Random-Effect Transfer
A covariance matrix must be symmetric and positive semidefinite: for every vector \(a\), \(a^\mathsf{T}\Omega a\geq0\). This condition ensures that \(\Omega\) can describe a multivariate distribution and can be used to generate random effects. For a positive-definite matrix, a Cholesky factorization gives \(\Omega=LL^\mathsf{T}\). If \(z\) is a vector of independent standard normal draws, then \(\eta=Lz\) has mean zero and covariance \[ \operatorname{Cov}(\eta)=\operatorname{Cov}(Lz) =L\operatorname{Cov}(z)L^\mathsf{T} =LI L^\mathsf{T} =\Omega. \] This is why preserving the complete matrix, including its off-diagonal entries and ordering, matters for population simulation.
For the lognormal model, \(\eta_{CL}\sim N(0,\omega^2_{CL})\) implies \(\log(CL_i)=\log(TVCL)+\eta_{CL}\). Its median is \(TVCL\), because the median of a zero-mean normal random variable is zero. Its arithmetic mean is \(TVCL\exp(\omega^2_{CL}/2)\). With \(TVCL=4.00\ \mathrm{L/h}\) and \(\omega^2_{CL}=0.0900\), the arithmetic mean is \(4.00\exp(0.0450)\approx4.184\ \mathrm{L/h}\). Thus the typical value and the arithmetic mean are distinct, even though both can be called population summaries in different contexts.
Diagnostics, Alternatives, and Common Failure Modes
Before relying on a transferred simulation, inspect the inputs and outputs. Confirm that the final estimates—not initial values or an intermediate run—were copied; that the THETA mapping reflects the actual model equations; that units agree; and that OMEGA order, diagonals, and covariances are preserved. Check that the matrix is symmetric and positive semidefinite. Also verify the random-effect source: a simulation that inadvertently uses zero \(\eta\) values will reproduce a typical subject even if OMEGA is present.
A useful numerical diagnostic is to compare a small hand calculation with the simulation for a known subject. Set the random effects to zero, administer a single bolus, and compare the concentration at a selected time with the analytical expression. Then test a supplied nonzero random-effect vector and verify that the resulting clearance and volume match their exponential equations. These checks reveal many scale, index, and event-setup errors before a large population run.
The main alternatives are to transfer the estimated OMEGA covariance matrix and let mrgsolve generate random effects, or to generate and supply \(\eta\) values externally. The first is concise and keeps the distribution in the model; the second can be useful for controlled draws, reproducible comparisons, or nonstandard sampling schemes. In either case, the random effects must follow the intended distribution and use the same ordering and transformation. Neither alternative repairs a mismatch in model parameterization.
- Confusing variance with standard deviation. If the model expects an OMEGA variance of \(0.0900\), entering \(0.300\) substitutes its standard deviation and changes the variability.
- Dropping covariance. Replacing the example's \(0.0120\) covariance with zero changes the joint distribution, even though each marginal variance stays the same.
- Calling the zero-random-effect prediction the population mean. It is the typical-subject prediction under this parameterization; nonlinearity generally makes it different from the average over simulated subjects.
- Rounding intermediate factors. Premature rounding can alter the last reported digit, as the concentration calculation above shows.
These transfers are only as reliable as the correspondence between the fitted and simulation models. If the structural equations, transformations, covariates, or residual model differ, copied estimates alone cannot guarantee that the simulation represents the fitted model. The task here is faithful parameter transfer, not a substitute for checking the model definitions or the assumptions behind the estimation.
Check Your Understanding
Answer each question, then compare with the solutions below.
- For lognormal clearance with \(TVCL=5.00\ \mathrm{L/h}\) and \(\operatorname{Var}(\eta_{CL})=0.160\), what is the clearance at \(\eta_{CL}=0\), and what is the population arithmetic mean under the model?
- A report gives an OMEGA standard deviation of \(0.250\) for one random effect. What variance belongs in a model specification that expects OMEGA variance?
- Why does reversing the order of two random effects without changing their parameter equations cause a transfer error when OMEGA has a nonzero covariance?
- For the example's matrix, calculate the correlation between clearance and volume random effects.
Solutions. (1) At \(\eta=0\), clearance is \(5.00\ \mathrm{L/h}\). The arithmetic mean is \(5.00\exp(0.160/2)=5.00\exp(0.0800)\approx5.416\ \mathrm{L/h}\). (2) The variance is \(0.250^2=0.0625\). (3) The covariance is associated with matrix row and column order; reversing indices assigns the joint dependence to the wrong random-effect pair unless the equations and matrix are permuted together. (4) The correlation is \(0.0120/(\sqrt{0.0900}\sqrt{0.0400})=0.0120/(0.300\times0.200)=0.200\), dimensionless.