Parameter Priors (penalized ML / MAP)

Maturity: experimental — see Feature Maturity for what this means.

A prior pins a parameter near a value you already believe, and lets the data move it only as far as the data can pay for. Declare one next to the parameter:

[parameters]
  theta TVCL(0.2, 0.001, 10.0)  prior(0.15, rse = 25%)
  theta TVV(10.0, 0.1, 500.0)   prior(9.8,  rse = 12%)
  omega ETA_CL ~ 0.09           prior(0.09, rse = 40%)

That is the entire user-facing surface: a central value and a relative uncertainty, both read off the same row of a published parameter table. There is no separate prior problem, no precision or covariance matrix, no inverse-Wishart, and no degrees-of-freedom indexing to keep aligned by hand (compare NONMEM $PRIOR).

When to reach for one

Sparse clinical data — rare disease, pediatrics, TDM-led fits. The likelihood is weakly informative, IIV terms collapse onto boundaries, and structural parameters drift. A prior taken from a published popPK model stabilizes the fit and lets a small local dataset borrow strength from what is already known.

The same mechanism is model updating: take a literature model as the prior, refit on local sparse data, and get a regularized estimate. That is the MAP-Bayesian updating behind much of model-informed precision dosing.

What a prior does to the objective

The penalty is added to the objective the optimizer minimises:

\[ \mathrm{OFV}_{\text{total}} \;=\; \mathrm{OFV}_{\text{data}} \;+\; \sum_k \left(\frac{x_k - m_k}{s_k}\right)^2 \]

where \(x_k\) is the parameter on the scale ferx estimates it, \(m_k\) is the prior centre on that same scale, and \(s_k\) comes from the RSE.

Two consequences worth knowing before you read a result:

  • ofv is the penalized total. The fit report splits it into ofv_data and ofv_prior so you can see how much of the objective the prior is responsible for. Normalisation constants are dropped on both halves (ferx’s −2 log L already omits \(N\log 2\pi\)), so a priored ferx OFV is not absolutely comparable to a priored NONMEM objective — compare a ΔOFV against a null twin instead.
  • AIC and BIC are computed from ofv_data. A penalized objective is not a log-likelihood, and an information criterion built from one would silently reward a tighter prior.

Which scale the prior lives on

The prior’s central value and its spread are always on the scale the parameter was declared on. If you wrote a variance, the prior is a variance; if you wrote (sd), both are standard deviations. Nothing is ever converted by hand.

What the prior becomes internally follows from how ferx estimates the parameter:

Declaration Estimated as Prior family From prior(v, rse = r)
theta X(init, lower ≥ 0, …) \(\ln\theta\) lognormal on \(\theta\) \(m=\ln v\), \(s=\sqrt{\ln(1+r^2)}\)
theta X(init, lower < 0, …) \(\theta\) normal on \(\theta\) \(m=v\), \(s=\lvert v\rvert\,r\)
omega E ~ var \(\ln(\mathrm{SD})\) lognormal on the variance \(m=\tfrac12\ln v\), \(s=\tfrac12\sqrt{\ln(1+r^2)}\)
omega E ~ sd (sd) \(\ln(\mathrm{SD})\) lognormal on the SD \(m=\ln v\), \(s=\sqrt{\ln(1+r^2)}\)
sigma S ~ var / (sd) \(\ln(\mathrm{SD})\) as omega as omega
kappa K ~ var / (sd) \(\ln(\mathrm{SD})\) as omega as omega

Working on the log-variance scale for \(\Omega\) keeps the prior positive and symmetric and avoids any matrix algebra.

The trap: a θ’s prior family follows its bounds

A θ declared with a non-negative lower bound is estimated as \(\ln\theta\), so its prior is lognormal. Widen the bound to allow negative values and ferx estimates \(\theta\) directly, which makes the same prior(...) declaration a normal prior instead. Nothing in the declaration changes; the family does.

This is visible rather than hidden: the fit report prints the realised family and the implied 95% prior interval for every priored parameter, so check that column rather than inferring the family from the model file.

Declaring a prior

prior(<value>, rse = <r>%)      # relative uncertainty, as a percent
prior(<value>, rse = <r>)       # …or as a fraction below 1
prior(<value>, sd  = <s>)       # absolute SD, on the declared scale

value = may be written explicitly (prior(value = 0.15, rse = 25%)). Exactly one of rse / sd is required.

RSE units are validated, not guessed. rse = 25% and rse = 0.25 both mean 25%. A bare number greater than 1 is a hard error: rse = 25 is overwhelmingly a percent written without its sign, and reading it as 2500% would produce a prior so flat the fit is indistinguishable from an unpriored one — a mistake you would never see, because the fit still converges.

Use sd = when a relative uncertainty is meaningless: a θ that can be negative (a covariate exponent, a drug-effect slope) whose prior centre sits at or near zero. On a log-estimated parameter sd is read as the relative spread it implies (sd = 0.03 on a value of 0.15 is rse = 20%), so the two spellings never disagree.

Model updating: a prior from a previous fit

Transcribing a whole parameter table by hand is the part of model updating that goes wrong. [priors] from_fit reads a previous ferx run and builds the priors from it:

[priors]
  from_fit = "parent-model-fit.yaml"

Every parameter the source run reports with a usable standard error, and whose name and family match a [parameters] declaration in this model, gets prior(estimate, rse = SE/estimate). That is exactly the prior you would have typed off the source run’s own output — no new convention, and the report reads the same for an imported prior as for a typed one.

A relative path is resolved against the model file’s directory, like [data] path, so a model and the fit it updates from can sit together in a run directory.

The source fit is read when the model is parsed, not once per estimation stage. Two things follow that are worth knowing:

  • a from_fit that cannot be honoured (missing file, wrong format, nothing to import) fails at parse time, so ferx check refuses the model and the error names the path rather than surfacing halfway through a fit;
  • the prior the optimizer minimised, the prior the covariance step used and the prior the report prints are the same one, even if the source file changes while the fit is running.

(A model that needs data-derived bindings — a symbolic [covariate_model] centring statistic, or a theta NAME[COL] level block — is re-parsed once the dataset is known, and reads the source fit again on that pass. Both reads happen before any estimation, so the guarantee above still holds.)

Which file, and from where

File Written by Precision
{model}-fit.yaml every run 6 decimals
{model}-fit.json --output-format json exact
run.fitrx --output run.fitrx exact

The YAML is the one every run writes, so it is the default answer. It is formatted to six decimal places, which matters only when a parameter is smaller than about 1e-4 — an Ω variance near a boundary, an additive σ in small units. Use .json or .fitrx when it is.

ImportantRun the parent from the directory the update lives in

The CLI writes {model}-fit.yaml into the working directory, while from_fit resolves against the model file’s directory. Those are the same place only when you run from where the models are:

cd my-run
ferx parent.ferx --data ../data/study.csv   # writes ./parent-fit.yaml
ferx update.ferx --data ../data/local.csv   # reads  ./parent-fit.yaml

Running the parent from one directory and the update from another gives failed to read …: No such file or directory, naming the path it looked for.

What is imported, and what is skipped

Imported: θ, Ω diagonals, Σ, and κ (Ω_IOV) diagonals.

Skipped and reported in the fit’s warnings, so a prior you expected is visible by name rather than merely absent:

  • a parameter with no standard error in the source — a FIXed one, or a run whose covariance step did not execute. Without an SE there is no spread;
  • a parameter that is FIXed in this model, which a prior could never move;
  • a block_omega / block_sigma element or a [mixture] override, which v1 does not support.

Skipped silently, because both are what you asked for:

  • a parameter that already carries an inline prior(...) here. The inline one wins, so you can tighten or loosen one imported prior without giving up the rest — and the result is listed in the fit’s prior table either way;
  • a source parameter this model simply does not have. A published model is usually larger than the one being updated, so one warning per absent parameter would bury the ones that matter.

An import that lands nothing at all is a hard error, for the same reason an unresolvable typed prior is: the alternative is an unpenalized fit that looks exactly like a penalized one.

Matching is on name and family, not name alone. θ names and η names are separate namespaces, so theta CL and omega CL can both exist in one model; a source θ called CL lands on the θ and never on the Ω.

Scales are converted for you

A fit report always states Ω and κ as variances and Σ as an SD, whatever the source model declared them as — so a source written (sd) imports identically to one written on the variance scale. What decides the prior is this model’s declaration, and the relative SE is carried across with it:

Source reports Declared here as Prior centre Prior RSE
Ω variance v omega X ~ v v SE/v
Ω variance v omega X ~ s (sd) √v SE/(2v)
Σ SD s sigma X ~ v s² 2·SE/s
Σ SD s sigma X ~ s (sd) s SE/s

(The factors of two are the delta method: the relative SE of xᵖ is |p| times the relative SE of x.)

A θ estimated on the identity scale — one whose declared lower bound is negative — imports its absolute SE instead, since a relative one is meaningless where the estimate may be zero.

Reading the output

--- Objective Function ---
OFV:  1584.2213
    data:  1583.9501
    prior: 0.2712
AIC:  1595.9501
BIC:  1612.1530

--- Parameter Priors (penalized ML / MAP) ---
  PARAMETER               PRIOR     ESTIMATE    SHIFT    PENALTY  PRIOR 95%
  TVCL                  0.15000      0.13195    -0.52      0.271  [0.0926, 0.2430] (lognormal)
  SHIFT is (estimate − prior) in prior SDs. Standard errors above are the curvature of
  the penalized objective — MAP standard errors, not posterior SDs.

SHIFT is the column that earns its place. A raw difference only says the estimate moved; the standardized one says whether the data and the prior actually disagree — which is the question a MAP fit is asked. A shift beyond about ±2 means the data is fighting the prior, and is worth investigating before the result is used.

The same values are written to the fit YAML under objective_function (as ofv_data / ofv_prior) and parameter_priors.

Standard errors are MAP standard errors

The reported SEs are the curvature of the penalized objective, not posterior standard deviations. This is deliberate and it is the only defensible choice: a sparse fit carries a prior precisely because the data alone does not identify the direction, and that is exactly where the unpenalized Hessian goes flat or non-positive-definite. Reporting the data-only curvature would give a standard error for a quantity the data cannot estimate.

It also means a tighter prior gives a smaller SE, mechanically. Do not read a narrow MAP interval as evidence the data was informative — read SHIFT and ofv_prior for that.

Where it applies

Priors are applied by the FOCE family: foce, focei, laplace, gn, gn_hybrid. They enter the objective, the gradient, and the covariance step’s Hessian.

If the last estimating stage of a chain does not apply priors, the fit is refused rather than run unpenalized. An unapplied prior is invisible in the output — the fit converges, the estimates look reasonable, and nothing says the prior was dropped — so a silent fallback would be worse than no feature. Chaining is fine as long as that stage qualifies: method = [saem, focei] works, because the FOCEI stage is the one whose estimates are reported.

“Last estimating” is deliberate: a trailing evaluation-only stage (imp with imp_eval_only, laplace with agq_eval_only) only scores the preceding stage’s parameters, so it neither applies a prior nor undoes one. method = [focei, imp] with imp_eval_only = true is accepted — FOCEI applied the prior — while [saem, laplace] with agq_eval_only = true is refused, because Laplace there is only evaluating SAEM’s unpriored estimates.

Covariance estimators

covariance_method = r (the default) is the only estimator accepted under a prior; it is the one whose Hessian carries the prior’s curvature, so its R⁻¹ is the penalized (MAP) covariance.

Both estimators built on S are refused. S is a sum of per-subject score cross-products, and a prior contributes one score for the whole population rather than one per subject, so nothing assembled from S can represent it:

  • covariance_method = s is S⁻¹, the unpenalized information outright;
  • covariance_method = rsr is the R⁻¹ S R⁻¹ sandwich, which wraps a prior-free S in two prior-shrunk R⁻¹ factors. The result is neither the penalized estimator nor the unpenalized one — it under-states the standard error on every priored coordinate — so it is rejected rather than reported with a caveat.

For the same reason, a large-parameter fit that would normally be auto-routed from r to s stays on r when a prior is declared, and says so.

Under SIR (sir = true), the importance weights score the penalized objective, so the resampled intervals reflect the prior rather than drifting back toward the unpenalized MLE.

Not supported in v1

Each of these is an error, never a silently ignored declaration:

  • a prior on a block_omega / block_sigma / block_kappa element — those are estimated as Cholesky entries rather than as variances, and a correlated block needs the inverse-Wishart machinery this feature exists to avoid. Declare the diagonal variances separately to prior them;
  • a prior on a FIXed parameter (it could never move it, so it is a typo);
  • a prior on a θ level block (theta PLACEBO[STUDY, TIME](…)) — one declaration standing for many parameters would make one prior mean “the same prior on every level”, which is a useful ridge but a different feature;
  • a prior on a per-class [mixture] Ω/Σ override.

Also out of scope by design: joint priors across parameters (no precision matrices anywhere), and full Bayesian inference — this is MAP / penalized ML, and method = bayes carries its own priors.

Comparison with NONMEM $PRIOR

NONMEM’s $PRIOR NWPRI works, and is a known source of friction: a separate prior problem, a variance block for the thetas, an inverse-Wishart on \(\Omega\) via a prior matrix plus degrees of freedom, and all the matrix indexing kept aligned by hand.

NONMEM $PRIOR NWPRI ferx prior(...)
Where it is written a second $PROBLEM inline, next to the parameter
θ prior $THETAP + $THETAPV (a variance-covariance matrix) prior(value, rse = …)
Ω prior inverse-Wishart: $OMEGAP + $OMEGAPD degrees of freedom prior(value, rse = …) on the declared variance
Correlated priors supported (full matrices) not in v1
Scale conversions yours internal

The θ side is directly equivalent — a normal prior, differing only in whether it sits on \(\theta\) or \(\ln\theta\), which in ferx follows the declared bounds. The \(\Omega\) side is genuinely a different prior: lognormal on the variance rather than inverse-Wishart. For a diagonal \(\Omega\) the two are close in practice and neither is more “correct”; they differ in tail behaviour and in how a very small variance is penalized.

The θ prior is anchored against NWPRI

Not just equivalent on paper. With both engines writing the model on the log scale and evaluating at the same fixed point, ferx’s penalty reproduces NONMEM’s:

NONMEM 7.5.1 ferx
null (no prior) 177.82632978040530 177.826330
prior centred on θ 177.826 (== null) —
prior offset from θ 177.94393305151641 177.943933
prior contribution 0.11760327111 0.117603

The centred run — prior mean placed exactly on the parameter — reproducing the null is what establishes that NWPRI drops the prior’s normalisation constant just as ferx does, so these are comparable absolutely and not only as a difference. Control streams and the measured .lst files live in nonmem_anchor/; the comparison runs as a test (the_prior_penalty_matches_nonmem_nwpri).

The Ω prior has no NWPRI counterpart to difference against, since NONMEM’s is an inverse-Wishart. It is pinned instead by closed-form unit tests, including a differential pair that catches the variance-vs-SD scale.

References

  • Gisleskog PO, Karlsson MO, Beal SL. Use of prior information to stabilize a population data analysis. J Pharmacokinet Pharmacodyn 29(5):473–505, 2002. The canonical frequentist-prior reference, and exactly the sparse-data clinical use case above.
  • NONMEM 7 guides, $PRIOR (NWPRI, TNPRI).