Pumas Docstrings

DataFrames.DataFrameType
DataFrame(events::DosageRegimen, expand::Bool = false)

Create a DataFrame with the information in the dosage regimen. If expand, create a DataFrame with the information in the event list (expanded form).

DataFrames.DataFrameMethod
DataFrame(bts::Union{SimulatedInference, Bootstraps})

Returns a data frame of the parameter estimates.

Distributions.MvNormalMethod
MvNormal(pmi::FittedPumasModelInference)

Returns an MvNormal with a mean equal to the optimal parameter values in the fitted model and a covariance equal to the variance-covariance estimate.

Distributions.MvNormalMethod
MvNormal(model::PumasModel, means::NamedTuple, vcov::NamedTuple)

Returns an MvNormal with mean equal to the vectorized means argument and a covariance equal to the block diagonal matrix of the covariance matrices in the vcov argument. The mean argument should be a named tuple of values in their natural domain. vcov should be a named tuple of numbers (for real parameters) or matrices for vector and matrix parameters. For parameters that are themselves symmetric matrices, e.g:

Ω = [
  Ω₁₁ Ω₁₂ Ω₁₃;
  Ω₂₁ Ω₂₂ Ω₂₃;
  Ω₃₁ Ω₃₂ Ω₃₃;
]

the covariance matrix should be the covariance of the vector: [Ω₁₁, Ω₁₂, Ω₂₂, Ω₁₃, Ω₂₃, Ω₃₃]. Note that because the matrix is symmetric, only the upper triangle is used iterating over the columns first in an outer loop then the over the rows of each column in an inner loop.

Example:

means = (θ1 = 1.0, θ2 = [1.0, 2.0], θ3 = [1.0 0.1; 0.1 1.0])
vcov = (
  θ1 = 0.1,
  θ2 = [0.5 -0.2; -0.2 0.6],
  θ3 = [0.5 -0.2 0.1; -0.2 0.6 -0.3; 0.1 -0.3 0.45],
)
dist = MvNormal(model, means, vcov)
MCMCChains.ChainsMethod
Chains(b::BayesMCMCResults; subject::Union{Nothing, Int} = nothing)

Constructs an MCMCChains.Chains object for the samples of the population or subject-specific parameters. If subject is nothing (the default value) a Chains object of the population parameters' samples is returned. If subject is an integer, the samples of the subject-specific parameters corresponding to the subject with the given index will be returned instead.

Pumas.AnalyticalPKPDProblemType
AnalyticalPKPDProblem(
  f,
  u0,
  tspan,
  events,
  time,
  p,
  bioav)

An analytical PK(PD) problem.

Fields:

  • f: callable that returns the state of the dynamic system
  • u0: initial conditions
  • tspan: time points that define the intervals in which the problem is defined
  • events: events such as dose, reset events, etc
  • time: event times
  • p: a function that returns pre block evaluated at a time point
  • bioav: bioavailability in each compartment
Pumas.AnalyticalPKProblemType

AnalyticalPKProblem(pkprob, prob2, sys2)

A problem that is partially an analytical problem that can be evaluated independently of the rest of the system.

Fields:

  • pkprob: the analytical part of the problem
  • prob2: a problem that represents the rest of the system
  • sys2 (optional): the system of prob2 - used for latexification.
Pumas.BayesMCMCType

BayesMCMC

An instance of the Hamiltonian Monte Carlo (HMC) based No-U-Turn Sampler (NUTS) algorithm. All the options of the algorithm can be passed to the constructor of BayesMCMC as keyword arguments. Example:

alg = BayesMCMC(nsamples = 1000)

The main options that can be set are:

  • target_accept: target acceptance ratio for the NUTS algorithm
  • nsamples: number of Markov Chain Monte Carlo (MCMC) samples to generate
  • nadapts: number of adaptation steps in the NUTS algorithm
  • nchains: number of MCMC chains to sample
  • ess_per_chain: target effective sample size (ESS) per chain, sampling terminates if the target is reached
  • check_every: the number of samples after which the ESS per chain is checked
  • time_limit: a time limit for sampling in seconds, sampling terminates if the time limit is reached
  • ensemblealg: can be set to EnsembleSerial() for serial sampling, EnsembleThreads() for multi-threaded sampling or EnsembleDistributed() for multi-processing (aka distributed parallelism) sampling
  • parallel_chains: can be set to true or false, if set to true the chains will be sampled in parallel using either multi-threading or multi-processing depending on the value of ensemblealg
  • parallel_subjects: can be set to true or false, if set to true the log probability computation will be parallelized over the subjects using either multi-threading or multi-processing depending on the value of ensemblealg
  • rng: the random number generator used
  • diffeq_options: a NamedTuple of all the differential equations solver's options
Pumas.BootstrapMethod
Bootstrap(; rng=Random.default_rng, samples=200, stratify_by=nothing, ensemblealg=EnsembleThreads())

The keyword rng can be used to set the seed for the random resamples of the datasets, samples controls the number of resampled Populations, stratify_by control the number of resampled datasets, and by specifying keyword stratify_by to be a Symbol with the name of a covariate it is possible to stratify by a covariate with a finite number of possible values, ensemblealg controls the parallization of each fit call.

Pumas.CensoredType
Censored(distribution::Distribution, lower::Real, upper::Real)::Censored

Construct a censored distribution based on distribution and the censoring limits lower and upper.

Pumas.Central1Type
Central1()

An analytical model for a one compartment model with dosing into Central. Equivalent to

Central' = -CL/Vc*Central

where clearance, CL, and volume, Vc, are required to be defined in the @pre block.

Pumas.Central1Meta1Type
Central1Meta1()

An analytical model for a model with a central compartment, Central, and a metabolite, Metabolite. Equivalent to Central' = -(CL+CLfm)/VcCentral Metabolite' = CLfm/VcCentral - CLm/Vm*Metabolite where clearances (CL and CLm) and volumes (Vc and Vm) and formation clearance of metabolite CLfm are required to be defined in the @pre block.

Pumas.Central1Periph1Type
Central1Periph1()

An analytical model for a two-compartment model with a central compartment, Central and a peripheral compartment, Peripheral. Equivalent to

Central'    = -(CL+Q)/Vc*Central + Q/Vp*Peripheral
Peripheral' =       Q/Vc*Central - Q/Vp*Peripheral

where clearance, CL, and volumes, Vc and Vp, and distribution clearance, Q, are required to be defined in the @pre block.

Pumas.Central1Periph1Meta1Type
Central1Periph1Meta1()

An analytical model for a two compartment model with a central compartment, Central, with a peripheral compartment, Peripheral, and a metabolite compartment, Metabolite. Equivalent to

Central'     = -(CL+Q+CLfm)/Vc*Central + Q/Vp*CPeripheral
CPeripheral' =            Q/Vc*Central - Q/Vp*CPeripheral
Metabolite'  = -CLm/Vm*Metabolite + CLfm/Vc*Central

where clearances (CL and CLm) and volumes (Vc, Vp and Vm), distribution clearance (Q), and formation clearance of metabolite CLfm are required to be defined in the @pre block.

Pumas.Central1Periph1Meta1Periph1Type
Central1Periph1Meta1Periph1()

An analytical model for a two compartment model with a central compartment, Central, with a peripheral compartment, Peripheral, and a metabolite compartment, Metabolite, with a peripheral compartment, MPeripheral. Equivalent to

Central'     = -(CL+Q+CLfm)/Vc*Central + Q/Vp*CPeripheral
CPeripheral' =            Q/Vc*Central - Q/Vp*CPeripheral
Metabolite'  = -(CLm+Qm)/Vm*Metabolite + Qm/Vmp*MPeripheral + CLfm/Vc*Central
MPeripheral' =        Qm/Vm*Metabolite - Qm/Vmp*MPeripheral

where clearances (CL and CLm) and volumes (Vc, Vp, Vm and Vmp), distribution clearances (Q and Qm) and formation clearance of metabolite CLfm are required to be defined in the @pre block.

Pumas.ConstDomainType
@param x = val
@param x ∈ ConstDomain(val)

Specifies a parameter as a constant.

Pumas.ConstrainedType
Constrained

Constrain a Distribution within a Domain. The most common case is an MvNormal constrained within a VectorDomain. The only supported method for Constrained is logpdf. Notice that the result does not represent a probability distribution since the remaining probability mass is not scaled by the mass excluded by the constraints.

Example

julia> d = Constrained(MvNormal(fill(1.0, 1, 1)), lower=-1, upper=1)
Constrained{ZeroMeanFullNormal{Tuple{Base.OneTo{Int64}}}, VectorDomain{Vector{Int64}, Vector{Int64}, Vector{Float64}}}(ZeroMeanFullNormal(
dim: 1
μ: 1-element Zeros{Float64}
Σ: [1.0;;]
)
, VectorDomain{Vector{Int64}, Vector{Int64}, Vector{Float64}}([-1], [1], [0.0]))


julia> logpdf(d, [ 0])
-0.9189385332046728

julia> logpdf(d, [-2])
-Inf
Pumas.CorrDomainType
CorrDomain(template)

Return a correlation domain.

template provides both the shape and initial values for the resulting correlation parameter.

  • if template is an Int then an identity matrix of size template is returned.
  • if template is a square Matrix then the CorrDomain will match its size and will

have initial values according to the template elements.

Pumas.Depots1Central1Type
Depots1Central1()

An analytical model for a one compartment model with a central compartment, Central, and a depot, Depot. Equivalent to

Depot'   = -Ka*Depot
Central' =  Ka*Depot - CL/Vc*Central

where absoption rate, Ka, clearance, CL, and volume, Vc, are required to be defined in the @pre block.

Pumas.Depots1Central1Periph1Type
Depots1Central1Periph1()

An analytical model for a two-compartment model with a central compartment, Central, a peripheral compartment, Peripheral, and a depot Depot. Equivalent to

Depot'      = -Ka*Depot
Central'    =  Ka*Depot - (CL+Q)/Vc*Central + Q/Vp*Peripheral
Peripheral' =                 Q/Vc*Central - Q/Vp*Peripheral

where absorption rate, Ka, clearance, CL, and volumes, Vc and Vp, and distribution clearance, Q, are required to be defined in the @pre block.

Pumas.Depots2Central1Type
Depots2Central1()

An analytical model for a one compartment model with a central compartment, Central, and two depots, Depot1 and Depot2. Equivalent to

Depot1'  = -Ka1*Depot1
Depot2'  = -Ka2*Depot2
Central' =  Ka1*Depot1 + Ka2*Depot2 - CL/Vc*Central

where absorption rates, Ka1 and Ka2, clearance, CL, and volume, Vc, are required to be defined in the @pre block.

When using this model during simulation or estimation, it is preferred to have 2 dosing rows for each subject in the dataset, where the first dose goes into cmt =1 (or cmt = Depot1) and the second dose goes into cmt=2 (or cmt=Depot2). Central compartment gets cmt=3 or (cmt = Central). e.g.

ev = DosageRegimen([100,100],cmt=[1,2]) s1 = Subject(id=1, events=ev)

Pumas.DiffEqWrapperType
DiffEqWrapper

Modifies the diffeq system f to add allow for modifying rates

  • f: original diffeq system
  • params: ODE parameters
  • rates_on: should rates be modified?
  • rates: amount to be added to rates
Pumas.DosageRegimenType
DosageRegimen

Lazy representation of a series of Events.

Fields

  • data::DataFrame: The tabular representation of a series of Events.

  • Signature

evts = DosageRegimen(amt::Numeric;
                     time::Numeric = 0,
                     cmt::Union{Numeric,Symbol} = 1,
                     evid::Numeric = 1,
                     ii::Numeric = zero.(time),
                     addl::Numeric = 0,
                     rate::Numeric = zero.(amt)./oneunit.(time),
                     duration::Numeric = zero(amt)./oneunit.(time),
                     ss::Numeric = 0,
                     route::NCA.Route)
  • Examples
julia> DosageRegimen(100, ii = 24, addl = 6)
DosageRegimen
 Row │ time     cmt    amt      evid  ii       addl   rate     duration  ss    route
     │ Float64  Int64  Float64  Int8  Float64  Int64  Float64  Float64   Int8  NCA.Route
─────┼───────────────────────────────────────────────────────────────────────────────────
   1 │     0.0      1    100.0     1     24.0      6      0.0       0.0     0  NullRoute

julia> DosageRegimen(50,  ii = 12, addl = 13)
DosageRegimen
 Row │ time     cmt    amt      evid  ii       addl   rate     duration  ss    route
     │ Float64  Int64  Float64  Int8  Float64  Int64  Float64  Float64   Int8  NCA.Route
─────┼───────────────────────────────────────────────────────────────────────────────────
   1 │     0.0      1     50.0     1     12.0     13      0.0       0.0     0  NullRoute

julia> DosageRegimen(200, ii = 24, addl = 2)
DosageRegimen
 Row │ time     cmt    amt      evid  ii       addl   rate     duration  ss    route
     │ Float64  Int64  Float64  Int8  Float64  Int64  Float64  Float64   Int8  NCA.Route
─────┼───────────────────────────────────────────────────────────────────────────────────
   1 │     0.0      1    200.0     1     24.0      2      0.0       0.0     0  NullRoute

julia> DosageRegimen(200, ii = 24, addl = 2, route = NCA.IVBolus)
DosageRegimen
 Row │ time     cmt    amt      evid  ii       addl   rate     duration  ss    route
     │ Float64  Int64  Float64  Int8  Float64  Int64  Float64  Float64   Int8  NCA.Route
─────┼───────────────────────────────────────────────────────────────────────────────────
   1 │     0.0      1    200.0     1     24.0      2      0.0       0.0     0  IVBolus

From various DosageRegimens

  • Signature

evs = DosageRegimen(regimen1::DosageRegimen, regimen2::DosageRegimen; offset = nothing)

offset specifies if regimen2 should start after an offset following the end of the last event in regimen1.

  • Examples
julia> e1 = DosageRegimen(100, ii = 24, addl = 6)
DosageRegimen
 Row │ time     cmt    amt      evid  ii       addl   rate     duration  ss    route
     │ Float64  Int64  Float64  Int8  Float64  Int64  Float64  Float64   Int8  NCA.Route
─────┼───────────────────────────────────────────────────────────────────────────────────
   1 │     0.0      1    100.0     1     24.0      6      0.0       0.0     0  NullRoute

julia> e2 = DosageRegimen(50, ii = 12, addl = 13)
DosageRegimen
 Row │ time     cmt    amt      evid  ii       addl   rate     duration  ss    route
     │ Float64  Int64  Float64  Int8  Float64  Int64  Float64  Float64   Int8  NCA.Route
─────┼───────────────────────────────────────────────────────────────────────────────────
   1 │     0.0      1     50.0     1     12.0     13      0.0       0.0     0  NullRoute

julia> evs = DosageRegimen(e1, e2)
DosageRegimen
 Row │ time     cmt    amt      evid  ii       addl   rate     duration  ss    route
     │ Float64  Int64  Float64  Int8  Float64  Int64  Float64  Float64   Int8  NCA.Route
─────┼───────────────────────────────────────────────────────────────────────────────────
   1 │     0.0      1    100.0     1     24.0      6      0.0       0.0     0  NullRoute
   2 │     0.0      1     50.0     1     12.0     13      0.0       0.0     0  NullRoute

julia> DosageRegimen(e1, e2, offset = 10)
DosageRegimen
 Row │ time     cmt    amt      evid  ii       addl   rate     duration  ss    route
     │ Float64  Int64  Float64  Int8  Float64  Int64  Float64  Float64   Int8  NCA.Route
─────┼───────────────────────────────────────────────────────────────────────────────────
   1 │     0.0      1    100.0     1     24.0      6      0.0       0.0     0  NullRoute
   2 │   178.0      1     50.0     1     12.0     13      0.0       0.0     0  NullRoute
Pumas.EventType
Event

A single event

Fields:

  • amt: Dose amount (mass units, e.g. g, mg, etc.)
  • time: Time of dose (time unit, e.g. hours, days, etc.)
  • evid: Event identifier, possible values are
    • -1: End of an infusion dose
    • 0: observation (these should have already been removed)
    • 1: dose event
    • 2: other type event
    • 3: reset event (the amounts in each compartment are reset to zero, the on/off status of each compartment is reset to its initial status)
    • 4: reset and dose event
  • cmt: Compartment (an Int or Symbol)
  • rate: The dosing rate:
    • 0: instantaneous bolus dose
    • >0: infusion dose (mass/time units, e.g. mg/hr)
    • -1: infusion rate specifed by model
    • -2: infusion duration specified by model
  • duration: duration of dose, if dose > 0.
  • ss: steady state status
    • 0: dose is not a steady state dose.
    • 1: dose is a steady state dose, and that the compartment amounts are to be reset to the steady-state amounts resulting from the given dose. Compartment amounts resulting from prior dose event records are "zeroed out," and infusions in progress or pending additional doses are cancelled.
    • 2: dose is a steady state dose and that the compartment amounts are to be set to the sum of the steady-state amounts resulting from the given dose plus whatever those amounts would have been at the event time were the steady-state dose not given.
  • ii: interdose interval (units: time). Time between implied or additional doses.
  • base_time: time at which infusion started
  • rate_dir: direction of rate
    • -1 for end-of-infusion
    • +1 for any other doses
  • route: route of administration to be used in NCA analysis
Pumas.FittedPumasEMModelType

FittedPumasEMModel

Contains the results of fitting a PumasEMModel with SAEM.

Fields:

  • model: The PumasEMModel that was fit.
  • data: The population that was fit.
  • μ: Final estimates of @param and @random parameters.
  • ω: Final estimates of ω parameters.
  • σ: Final estimated standard deviations of residuals.
  • opthistory: If approx is SAEM, contains the full history of unconstrained-μ, ω, and, σ estimates across iterations. Otherwise, contains the optimization results.
  • approx: Fitting method used.
  • η: Random effects from the fit. If approx is SAEM, contains a length(pop) × nη × iteration array of η samples from the smoothing phase. Otherwise, contains orthogonal random effects.
  • kwargs: Input keywords stored for use in diagnostics. Includes ensemblealg and diffeq_options.
Pumas.FormulaType

O ~ C | D

For example: CL ~ 1 + wt | Normal -> Formula{:CL,(1,:wt),:Normal}()

Pumas.IDRMethod
IDR(; direct=true, mode=:stimulation_sigmoid)

Get an indirect response model (IDR) for use as the pharmacodynamic part of model macro generation.

  • direct::Bool - toggle whether the drug acts on production (true) or inactivation (false)
  • mode::Symbol - pick the mode of drug action. Valid options are: :stimulationlinear, :stimulationemax, :stimulationsigmoid, :stimulationinfinite, :stimulationlogarithmic, :inhibitionemax, :inhibitionsigmoid, :inhibitionhill,
Pumas.JointMAPType

JointMAP

An instance of the maximum a-posteriori (MAP) algorithm where the mode of the joint posterior of the population and subject-specific parameters is found using the BFGS optimization algorithm. All the options of the algorithm can be passed to the constructor of JointMAP as keyword arguments. Example:

alg = JointMAP(ensemblealg = EnsembleThreads())

The main options that can be set are:

  • ensemblealg: can be set to EnsembleSerial() for serial sampling, EnsembleThreads() for multi-threaded sampling or EnsembleDistributed() for multi-processing (aka distributed parallelism) sampling
  • diffeq_options: a NamedTuple of all the differential equations solver's options
Pumas.NLType
NL(pkmodel; hill=false)

Specify nonlinear clearance of a pkmodel for model generation.

pkmodel::Pumas.Explicitmodel - The pharmacokinetic model, e.g. Depots1Central1(). hill::Bool - toggle the use of a Hill coefficient.

Pumas.ObsAdditiveType
ObsAdditive(obsname::Symbol)

Get an additive error model for automatic model macro generation.

  • obsname::Symbol - choose what the output should be called. If you'll be fitting the model towards data then this symbol should match the identifier of the observation in your data.
Pumas.ObsCombinedType
ObsCombined(obsname::Symbol)

Get a combined (proportional and additive) error model for automatic model macro generation.

  • obsname::Symbol - choose what the output should be called. If you'll be fitting the model towards data then this symbol should match the identifier of the observation in your data.
Pumas.ObsProportionalType
ObsProportional(obsname::Symbol)

Get a proportional error model for automatic model macro generation.

  • obsname::Symbol - choose what the output should be called. If you'll be fitting the model towards data then this symbol should match the identifier of the observation in your data.
Pumas.PDiagDomainType
@param x ∈ PDiagDomain(n::Int; init=ones(n))

Specifies a parameter as a positive diagonal matrix, with initial diagonal elements specified by init.

Pumas.PSDDomainType
@param x ∈ PSDDomain(n::Int; init=Matrix{Float64}(I, n, n))

Specifies a parameter as a symmetric n-by-n positive semi-definite matrix. init sets the initial value of the parameter and should be a positive semi-definite Matrix of Float64s.

Pumas.ParamSetType
ParamSet(params::NamedTuple)

Construct a Pumas.jl parameter set.

params should be a NamedTuple that maps a parameter name to a Domain or a Distribution.

Example

ParamSet((; 
  tvKa = RealDomain(lower=0.),
  tvCL = RealDomain(lower=0.),
  tvω = PDiagDomain(2)
))
Pumas.PumasEMModelType
PumasEMModel{Σ,NRE,C,PRR,DCR,DR,PR}

Parameters: Σ - Covariance structure; (2,1) indicates a matrix with 2x2 and 1x1 diagonal blocks. NRE - The number of formulas without an associated random effect. C - Covariates defined in @covariates block PRR - Parameters defined in @pre block DCR - Parameters defined in @dosecontrol block DR - Parameters defined in @dynamics block POR - Parameters defined in @post block

Currently supported fitting methods are SAEM() and LaplaceI().

Pumas.PumasModelType
PumasModel

A model takes the following arguments (arg [= default]):

  • param: a ParamSet detailing the parameters and their domain
  • random: a mapping from a named tuple of parameters -> ParamSet
  • pre: a mapping from the (params, randeffs, subject) -> ODE parameters.
  • dcp: a mapping from the (params, randeffs, subject) -> dose-control parameters.
  • init: a mapping (col,t0) -> initial conditions
  • prob: a DEProblem describing the dynamics (either exact or analytical)
  • derived: the derived variables and error distributions (param, randeffs, data, ode vals) -> sampling dist
  • observed = (_, _, _, _, _, _, samples) -> samples: simulated values from the error model and post processing: (param, randeffs, data, ode vals, samples) -> vals
  • options::PumasModelOptions = PumasModelOptions(; subject_t0=false): pass additional options, see PumasModelOptions.
  • syms::PumasModelSymbols = PumasModelSymbols(): the symbols of all the variables defined inside the model, see PumasModelSymbols.
  • sys = nothing: the ModelingToolkit.jl System used to generate prob - if applicable.
  • macroexpr::Expr = :(): the expression passed to the @model macro - if applicable.
  • desc::Dict{Symbol, String} = Dict{Symbol, String}(): descriptions of the symbols (variables, parameters, etc.) used in the model.
  • metadata::Dict{Symbol, Any} = Dict{Symbol, Any}(): include a Dict of arbitrary metadata.
Pumas.RealDomainType
@param x ∈ RealDomain(;lower=-∞,upper=∞,init=0)

Specifies a parameter as a real value. lower and upper are the respective bounds, init sets the initial value of the parameter and is used for init_params.

Pumas.SAEMType

SAEM(; iters::Tuple{Int,Int,Int} = (200,100,100), α::Float64 = 0.7, λ::Float64 = 0.975, repeat::Int = 1, verbosity::Int = 0)

  • iters: A tuple indicating the number of iterations to run in the exploratory, middle, and smoothing phases of SAEM.
  • α: The learning rate for the middle phase. On the kth update, α^k weight is placed on the kth iteration, and 1 - α^k on the previous.
  • λ: Maximum standard deviation change rate. A value of 0.95 means that ω_new = clamp(ω_new, ω_old * 0.95, ω_old / 0.95) on each iteration. Slowing the decliine in ωs and σs increases exploration, while limiting the growth prevents excessive flattening of the log likelihood.
  • repeat: Indicates how many times to repeat iters. Defaults to 1.
  • verbosity: How verbose should it be? 0: silent, 1: print on each repeat, 2: print once per iteration.
Pumas.SIRMethod
SIR(; rng=Random.default_rng, samples, resamples)

This is the sampling importance re-sampling (SIR) algorithm. In the SIR algorithm, a number of parameter values are sampled from a truncated multivariate normal distribution using rejection sampling. The mean of the multivariate normal distribution is equal to the maximum likelihood parameter values. The covariance is variance-covariance matrix estimator, and the distribution is truncated to the domains of the parameters. After sampling, the samples are weighed by their respective likelihood values, and a subset of these samples is then selected by weighted re-sampling without replacement. The keyword arguments samples and resamples are the number of samples and re-samples respectively. They do not have default values so they must be input by the user.

Pumas.SubjectType
Subject

The data corresponding to a single subject:

Fields:

  • id: identifier
  • observations: a named tuple of the dependent variables
  • covariates: a named tuple containing the covariates, or nothing.
  • events: a vector of Events.
  • time: a vector of time stamps for the observations

When there are time varying covariates, each covariate is interpolated with a common covariate time support. The interpolated values are then used to build a multi-valued interpolant for the complete time support.

From the multi-valued interpolant, certain discontinuities are flagged in order to use that information for the differential equation solvers and to correctly apply the analytical solution per region as applicable.

Constructor

Subject(;id = "1",
         observations = NamedTuple(),
         events = Event[],
         time = observations isa AbstractDataFrame ? observations.time : nothing,
         event_data = true,
         covariates::Union{Nothing, NamedTuple} = nothing,
         covariates_time = observations isa AbstractDataFrame ? observations.time : nothing,
         covariates_direction = :right)

Subject may be constructed from an <:AbstractDataFrame with the appropriate schema or by providing the arguments directly through separate DataFrames / structures.

Examples:

julia> Subject()
Subject
  ID: 1

julia> Subject(id = 20, events = DosageRegimen(200, ii = 24, addl = 2), covariates = (WT = 14.2, HT = 5.2))
Subject
  ID: 20
  Events: 3
  Covariates: WT, HT

julia> Subject(covariates = (WT = [14.2, 14.7], HT = fill(5.2, 2)), covariates_time = [0, 3])
Subject
  ID: 1
  Covariates: WT, HT
Pumas.SubjectMethod

Subject

Constructor

Subject(simsubject::SimulatedObservations)

Roundtrip the result of simobs, i.e. SimulatedObservations to a Subject/Population

Example:

sims = simobs(model, pop, params)

To convert sims to a Population, broadcast as below

Subject.(sims)
Pumas.SubjectMethod

Subject(subject::Subject; id::String, events)

Construct a Subject from an existing Subject while replacing information according to the keyword arguments. The possible keyword arguments are

  • id, sets the subject identifier to the string id. Defaults to the identifier already in subject.
  • events, sets the subject events to the input DosageRegimen. Defaults to the events already present in subject.
Pumas.TimeToEventType
TimeToEvent{T}

Distribution-like struct to store the hazard and the cumulative hazard in a time-to-event model. The dependent variable in a model that uses TimeToEvent should be a censoring variable that is one when the variable is observed (not censored) zero when the variable is right censored. Currently, no other censoring types are supported.

Example

...
@pre begin
  θeff = θ*DOSE
  λ = λ₀*exp(θeff)
end

@dynamics begin
  Λ' = λ
end

@derived begin
  DV ~ @. TimeToEvent(λ, Λ)
end
...
Pumas.VectorDomainType
@param x ∈ VectorDomain(n::Int; lower=-∞,upper=∞,init=0)

Specifies a parameter as a real vector of length n. lower and upper are the respective bounds, init sets the initial value of the parameter and is used for init_params. The keyword arguments can all take either numbers or vectors. If numbers then the same value will be applied to all n vector elements. If a you specify a vector (of length n) then you can adjust the parameters of the VectorDomain elements individually.

CommonSolve.solveFunction
sol = solve(
  model::AbstractPumasModel,
  population::Union{Population, Subject},
  param::NamedTuple,
  randeffs::NamedTuple=sample_randeffs(rng, model, param, population);
  ensemblealg=EnsembleSerial(),
  diffeq_options=NamedTuple(),
  rng=default_rng()
)

Solve the model applied to the Subject(s) within population using parameters param and random effects randeffs. By default, the times at which to save the solution are taken to be the observation times for each Subject within population. This can be overriden by supplying a vector or range as saveat to diffeq_options - e.g. diffeq_options=(; saveat=[0., 0.5, 1.]).

Arguments

  • model may either be a PumasModel or a PumasEMModel.

  • population may either be a Population of Subjects or a single Subject.

  • param is parameter set in the form of a NamedTuple, e.g. (; tvCL=1., tvKa=2., tvV=1.).

  • randeffs is an optional argument that, if used, should specify the random effects for each subject in population. Such random effects are specified by NamedTuples for PumasModels (e.g. (; tvCL=1., tvKa=2.)) and by Tuples for PumasEMModels (e.g. (1., 2.)). If population is a single Subject (without being enclosed in a vector) then a single such random effect specifier should be passed. If, however, population is a Population of multiple Subjects then randeffs should be a vector containing one such specifier for each Subject. The functions init_randeffs, zero_randeffs, and sample_randeffs are sometimes convenient for generating randeffs:

     randeffs = zero_randeffs(model, param, population)
     solve(model, population, param, randeffs)

If no randeffs is provided, then random ones are generated according to the distribution in the model.

  • ensemblealg is a keyword argument that allows you to choose between different modes of parallelization. Options include EnsembleSerial(), EnsembleThreads() and EnsembleDistributed().

  • diffeq_options is a keyword argument that takes a NamedTuple of options to pass on to the differential equation solver.

  • rng is a keyword argument where you can specify which random number generator to use.

Distributions.logpdfMethod
logpdf(d::Pumas.Constrained, x)

Evaluate the logarithm of the probability density of the constrained distribution, d, at x.

  • If d is a constrained univariate distribution then x should be a scalar.
  • If d is a constrained multivariate distribution then x should be a vector.

Evaluations of x outside of d's constraints returns -Inf.

Note that d itself is not a true distribution since its probability mass is not rescaled to 1.

GlobalSensitivity.gsaFunction

GlobalSensitivity.gsa(model, population, params, method, vars, prangelow, prangehigh; kwargs...)

Function to perform global sensitivty analysis

The arguments are:

  • model: a PumasModel, either defined by the @model DSL or the function-based interface.
  • population: a Population.
  • params: a named tuple of parameters.
  • method: one of the GSAMethods from GlobalSensitivity.jl, Sobol(), Morris(), eFAST(), RegressionGSA().
  • vars: a list of the derived variables to run GSA on.
  • p_range_low & p_range_high: the lower and upper bounds for the parameters you want to run the GSA on.

For method specific arguments that are passed with the method constructor you can refer to the GlobalSensitivity.jl documentation.

LinearAlgebra.condMethod
cond(pmi::FittedPumasModelInference)

Return the condition number of the variance-covariance matrix stored in pmi. Throws an error if pmi is the result of a call to infer with Pumas.Bootstrap or if the

variance-covariance calculation failed.

Pumas._analytical_solveMethod
_analytical_solve(
  m::ExplicitModel,
  t::Real,
  t0::Real,
  amounts::AbstractVector,
  doses::AbstractVector,
  pre,
  rates::AbstractVector)

Solve the model at a given point in time t. The reference time t0 is associated with the initial system state amounts, bolus doses and the start of an infusion with a rates going into compartments given by rates. For each evaluation, pre will be evaluated at t.

Pumas._convergencedataMethod
_convergencedata(obj; metakey="x")

Returns the "timeseries" of optimization as a matrix, with series as columns.

Warn

This must return parameter data in the same order that _paramnames returns names.

Pumas._fixed_to_constant_paramsetMethod
_fixed_to_constant_paramset(paramset::ParamSet, param::NamedTuple, fixed::NamedTuple, omegas::Tuple)

Replace individual parameter Domains in paramset with ConstantTranform if the parameter has an entry in fixed. Return a new parameter ParamSet with the values in fixed in place of the values in input param.

Pumas._lpdfMethod
_lpdf(d,x)

The log pdf: this differs from Distributions.logdpf definintion in a couple of ways:

  • if d is a non-distribution it assumes the Dirac distribution.
  • if x is NaN or Missing, it returns 0.
  • if d is a NamedTuple of distributions, and x is a NamedTuple of observations, it computes the sum of the observed variables.
Pumas._objectivefunctionvaluesMethod
_objectivefunctionvalues(obj)

Returns the objective function values during optimization. Must return a Vector{Number}.

Pumas._paramnamesMethod
_paramnames(obj)

Returns the names of the parameters which convergence is being checked for.

Warn

This must return parameter names in the same order that _convergencedata returns data.

Pumas.add_tMethod

add_t(ex, vars)

Append (t) to all the vars symbols in the expression ex.

Pumas.apply_formattingMethod

applyformatting(inputex, vars; show_t=true, italicize=true)

Apply formatting to the vars of an input expression before latexification.

Pumas.build_model_strMethod
build_model_str(pk, pd, obs, dcp)

Generate a string that fully specifies an @model macro call.

  • pk - a pharmacokinetic model, for example Depots1Central1()
  • pd::Union{Nothing, AbstractPDModule} - a pharmacodynamic model, for example IDR(; direct=true, mode=:stimulation_emax)
  • obs::AbstractObservationModule - an observational error model, for example ObsCombined(:MyOutputName)
  • dcp::DoseControl - a dose control specification, for example DoseControl(; bioav=true)
Pumas.conditional_nllMethod
conditional_nll(m::AbstractPumasModel, subject::Subject, param, randeffs; diffeq_options)

Compute the conditional negative log-likelihood of model m for subject with parameters param and random effects randeffs. diffeq_options is a NamedTuple of options passed to the ODE solver. Requires that the derived produces distributions.

Pumas.cor2covMethod

cor2cov(ρ::AbstractVector, ω::AbstractVector)

cor2cov(ω::AbstractMatrix)

cor2cov(ω::AbstractVector)


Converts correlations and standard deviations to a covariance matrix. There are 3 ways to use the cor2cov function.


The first method has 2 inputs:

  • ρ: a vector of correlations which is a vector of the strict lower triangular elements of the correlation matrix of the random variables.
  • ω: a vector of standard deviations of the random variables

Example:

ω11, ρ12, ω22 = 0.2, 0.1, 0.3

ρ = [ρ12]
ω = [ω11, ω22]
Pumas.cor2cov(ρ, ω)

gives the covariance matrix:

[ ω11^2        ρ12*ω11*ω22;
  ρ12*ω11*ω22  ω22^2      ]

which is equal to:

2×2 Matrix{Float64}:
 0.04   0.006
 0.006  0.09

The second method has a single input matrix ω which is a symmetric matrix whose diagonal is the standard deviations of the random variables and off-diagonals are their correlations. Example:

ω = [
  ω11  ρ12;
  ρ12  ω22
]
Pumas.cor2cov(ω)

gives the covariance matrix:

[ ω11^2        ρ12*ω11*ω22;
  ρ12*ω11*ω22  ω22^2      ]

which is equal to:

2×2 Matrix{Float64}:
 0.04   0.006
 0.006  0.09

The third method has a single input vector ω which is a vector of the lower triangular elements of the symmetric matrix whose diagonal is the standard deviations of the random variables and off-diagonals are their correlations.

Example:

ω = [ω11, ρ12, ω22]
Pumas.cor2cov(ω)

gives the covariance matrix:

[ ω11^2        ρ12*ω11*ω22;
  ρ12*ω11*ω22  ω22^2      ]

which is equal to:

2×2 Matrix{Float64}:
 0.04   0.006
 0.006  0.09
Pumas.covariates_interpolantMethod
covariates_interpolant(
  covariates_nt::NamedTuple,
  covariates_times::NamedTuple,
  id::String;
  [covariates_direction::Symbol]) -> time_grid, interpolant

Creates an interpolation of the elements of covariates_nt and covariates_times using the a piecewise constant interpolation scheme. The two NamedTuples should have the same keys. The values of covariates_nt and covariates_times should contain the interpolation values and time point respectively. The argument id is a string indicating the subject ID and the argument is only used for providing useful error messages.

Returns a callable interpolation struct as well as the common time grid for the covariate observations. This is safe for values which are not time-varying and it also allows for mixing subjects with multiple measurements and subjects with a single measurement. The default is left-sided constant interpolation.

Pumas.covariates_interpolantMethod
covariates_interpolant(
  covariates::Vector{Symbol},
  data::AbstractDataFrame,
  time::Symbol,
  id::String;
  [covariates_direction::Symbol]) -> time_grid, interpolant

Creates an interpolation of the columns of data listed in the covariates vector at the time points in the time column of data using the a piecewise constant interpolation scheme. The argument id is a string indicating the subject ID and the argument is only used for providing useful error messages.

Returns a callable interpolation struct as well as the common time grid for the covariate observations. This is safe for values which are not time-varying and it also allows for mixing subjects with multiple measurements and subjects with a single measurement. The default is left-sided constant interpolation.

Pumas.dosecontrolMethod
dosecontrol(fpm::AbstractFittedPumasModel)::Vector{ConstantInterpolationStructArray}

Return the dosecontrol parameters from fpm. The dosecontrol parameters are the variables bioav, lags, duration, rate if they're defined in the @dosecontrol block in the @model macro. The function returns a vector of covariate interpolation objects. Each of them can be evaluated at an arbitrary time point, but the interpolant is piecewise constant between the event times. Note, that in the model solution, subject simulation, and parameter fitting the dose control parameters can be time-varying, but they only affect the dosing according to their value at the beginning of each dose event. This means that the dose control parameters are only relevevant at these time points. Each of the covariate interpolation objects can be converted to a DataFrame by calling the DataFrame constructor on each element, or a DataFrame with all subjects can be obtained by calling the DataFrame constructor on the vector output by this function. This DataFrame will only containt the times that correspond to the start of each dose event for the reason discussed above.

Pumas.eiwresMethod
eiwres(model::PumasModel, subject::Subject, param::NamedTuple, nsim::Integer)

Calculate the Expected Simulation based Individual Weighted Residuals (EIWRES).

Pumas.em_to_m_paramsMethod
em_to_m_params(m::PumasEMModel, params::NamedTuple)

Given that params is a NamedTuple of coefficients used as inits for fitting m or returned by coef(fpm::FittedPumasEMModel), then em_to_m_params will return a NamedTuple of coefficients usable as inits for fitting or performing diagnostics on a PumasModel created by calling convert(PumasModel, m::PumasEMModel).

Pumas.empirical_bayesMethod
empirical_bayes(fpm::Pumas.FittedPumasModel)
empirical_bayes(fpm::Pumas.FittedPumasEMModel)
empirical_bayes(insp::Pumas.FittedPumasModelInspection)

Return sampled random effects or empirical bayes estimates from a fit or model inspection. If the model was estimated with the Pumas.FO likelihood approximation methods the empirical bayes estimates will be obtained using the Pumas.LaplaceI approximation. If either Pumas.FOCE or Pumas.LaplaceI was used the final empirical bayes estimates will be returned. If Pumas.SAEM was used to fit the empirical bayes estimates are obtained using the Pumas.LaplaceI approximation.

Pumas.empirical_bayes_distMethod
empirical_bayes_dist(fpm::AbstractFittedPumasModel)

Estimate the distribution of random effects (typically a Normal approximation of the empirical Bayes posterior) for the subjects in a FittedPumasModel. The result is returned as a vector of MvNormals.

Pumas.epredictMethod
epredict(fpm::AbstractFittedPumasModel; nsim::Integer, rng::AbstractRNG=default_rng())

Calculate the expected simulation based population predictions.

The keyword nsim controls the number of times a subject is simulated and the keyword rng can be used to provide a custom random number generator.

Pumas.eventnumMethod
eventnum(t, events) -> # of doses for each time

Creates an array that matches t in length which denotes the number of events that have occurred prior to the current point in time. If t is a scalar, outputs the number of events before that time point.

Pumas.expectationMethod
expectation(g, ::MonteCarloExpectation, model::PumasModel, subject::Subject, dist::MvNormal; imaxiters = 10000, diffeq_options::NamedTuple = NamedTuple(), rng::AbstractRNG = Pumas.default_rng())

Computes the expected predictions for subject with respect to the population parameters and the subject-specific random effects using Monte Carlo integration. dist is a multivariate normal distribution used to sample the population parameters e.g. dist = MvNormal(pmi) where pmi is the output of infer. Sampling is done using rejection sampling to ensure the parameter values are in the parameters' domains. imaxiters is the number of Monte Carlo samples used. diffeq_options can be used to customize the options of the differential equation solver used. rng is the random number generator used during the sampling, which defaults to Pumas.default_rng().

Pumas.expectationMethod
expectation(g, ::MonteCarloExpectation, model::PumasModel, subject::Subject, param_dists::NamedTuple; imaxiters = 10000, diffeq_options::NamedTuple = NamedTuple(), rng::AbstractRNG = Pumas.default_rng())

Computes the expected predictions for subject with respect to the population parameters and the subject-specific random effects using Monte Carlo integration. param_dists should be a named tuple of distributions for the population parameters. imaxiters is the number of Monte Carlo samples to use. diffeq_options can be used to customize the options of the differential equation solver used. rng is the random number generator used in the sampling, which defaults to Pumas.default_rng().

Pumas.expectationMethod
expectation(g, quant::QuadratureExpectation, model::PumasModel, subject::Subject, param_dists::NamedTuple; ireltol=1e-3, iabstol=1e-3, imaxiters=10000, diffeq_options::NamedTuple = NamedTuple())

Computes the expected predictions for subject with respect to the population parameters using quadrature methods. param_dists should be a named tuple of distributions for the population parameters. Currently, random effects are not supported. ireltol and iabstol are the tolerances for the integration. The integration algorithm will terminate if the tolerance or imaxiters is reached, whichever is first. diffeq_options can be used to customize the options of the differential equation solver.

Pumas.findinfluentialFunction
findinfluential(fpm::AbstractFittedPumasModel, k::Integer=5)

Return a vector of the k most influencial observations based on the value of (minus) the log-likelihood function.

Pumas.get_model_combinationsMethod
get_model_combinations()

Get all valid arguments to build_model_str or build_model_expr.

Returns a tuple (pks, pds, obs, dcps) where pks is a vector of all possible PK specifications, and so on, for the model generation.

Example usage:

for args in Iterators.product(get_model_combinations())
  obs = make_obs(args[3], args[1], args[2])
  build_model_str(args[1], args[2], obs, args[4])
end

to create all model strings (~15k different ones, 1M lines of code).

Useful for testing and possibly for automatic building, etc.

Pumas.getblockMethod
getblock(blockname, model::AbstractPumasModel) -> array of expressions for a PumasModel macro block.

The blockname should match the name of the blocks used for model creation with the @model or @emmodel macro. Only the blocks that were explicitly used in model creation can be used. Valid such blocknames are :covariates, :derived, :dynamics, :dosecontrol, :init, :observations, :options, :param, :post, :pre, :random, :reactions and :vars. Only models that are defined using the @model or @emmodel macro are supported.

Pumas.gewekediagMethod
gewekediag(b::BayesMCMCResults; subject::Union{Nothing,Int} = nothing, first::Real = 0.1, last::Real = 0.5)

Computes the Geweke diagnostic [Geweke1991] for each chain outputting a p-value per parameter. If subject is nothing (default), the chains diagnosed are those of the population parameters. If subject is set to an integer index, the chains diagnosed are those of the subject-specific parameters corresponding to the subject with the input index.

The Geweke diagnostic compares the sample means of two disjoint sub-chains X₁ and X₂ of the entire chain using a normal difference of means hypothesis test where the null and alternative hypotheses are defined as:

  • H₀: μ₁ = μ₂
  • H₁: μ₁ ≂̸ μ₂

where μ₁ and μ₂ are the population means. The first sub-chain X₁ is taken as the first (first * 100)% of the samples in the chain, where first is a keyword argument defaulting to 0.1. The second sub-chain X₂ is taken as the last (last * 100)% of the samples in the chain, where last is a keyword argument defaulting to 0.5.

The test statistic used is: z₀ = x̅₁ - x̅₂ / √(s₁² + s₂²) where x̅₁ and x̅₂ are the sample means of X₁ and X₂ respectively, and s₁ and s₂ are the Markov Chain standard error (MCSE) estimates of X₁ and X₂ respectively. Auto-correlation is assumed within the samples of each individual sub-chain, but the samples in X₁ are assumed to be independent of the samples in X₂. The p-value output is an estimate of P(|z| > |z₀|), where z is a standard normally distributed random variable.

Low p-values indicate one of the following:

  • The first and last parts of the chain are sampled from distributions with different means, i.e. non-convergence,
  • The need to discard some initial samples as burn-in, or
  • The need to run the sampling for longer due to lack of samples or high auto-correlation.

High p-values indicate the inability to conclude that the means of the first and last parts of the chain are different with statistical significance. However, this alone does not guarantee convergence to a fixed posterior distribution because:

  • Either the standard deviations or higher moments of X₁ and X₂ may be different, or
  • The independence assumption between X₁ and X₂ may not be satisfied when high auto-correlation exists.
Pumas.heideldiagMethod
heideldiag(b::BayesMCMCResults; subject::Union{Nothing,Int} = nothing, alpha::Real = 0.05, eps::Real = 0.1, start::Integer = 1)

Compute the Heidelberger and Welch diagnostic [Heidelberger1983] for each chain. If subject is nothing (default), the chains diagnosed are those of the population parameters. If subject is set to an integer index, the chains diagnosed are those of the subject-specific parameters corresponding to the subject with the input index.

The Heidelberger diagnostic attempts to:

  • Identify a cutoff point for the initial transient phase for each parameter, after which the samples can be assumed to come from a steady-state posterior distribution. The initial transient phase can be removed as burn-in. The cutoff point for each parameter is given in the burnin column of the output dataframe.
  • Estimate the relative confidence interval (confidence interval divided by the mean) for the mean of the steady-state posterior distribution of each parameter, assuming such steady-state distribution exists in the samples. A large confidence interval implies either the lack of convergence to a stationary distribution or lack of samples. Half the relative confidence interval is given in the halfwidth column of the output dataframe. The test column will be true (1) if the halfwidth is less than the input target eps (default is 0.1) and false (0) otherwise. Note that parameters with mean value close to 0 can have erroneous relative confidence intervals because of the division by the mean. The test value can therefore be expected to be false (0) for those parameters without concluding lack of convergence.
  • Quantify the extent to which the distribution of the samples is stationary using statistical testing. The returned p-value, shown in the pvalue column of the output dataframe, can be considered a measure of mean stationarity. A p-value lower than the input threshold alpha (default is 0.05) implies lack of stationarity of the mean, i.e. the posterior samples did not converge to a steady-state distribution with a fixed mean.

The Heidelberger diagnostic only tests for the mean of the distribution. Therefore, it can only be used to detect lack of convergence and not to prove convergence. In other words, even if all the numbers seem normal, one cannot conclude that the chain converged to a stationary distribution or that it converged to the true posterior.

Pumas.icoefMethod
icoef(fpm::AbstractFittedPumasModel; obstimes)::Vector{ConstantInterpolationStructArray}

Return the individual coefficients from fpm. The individual coefficients are the variables defined in the @pre block in the @model macro. The function returns a vector of covariate interpolation objects defined on obstimes. If no obstimes are given, the time points used to define the covariates will be used. Each of them can be evaluated at any time point, but the interpolant is piecewise constant between the obstimes. Each of the covariate interpolation objects can be converted to a DataFrame by calling the DataFrame constructor on each element, or a DataFrame with all subjects can be obtained by calling the DataFrame constructor on the vector output by this function.

Pumas.inferMethod
infer(fpm::FittedPumasModel; level=0.95, rethrow_error::Bool=false, sandwich_estimator::Bool=true) -> FittedPumasModelInference

Compute the vcov matrix and return a struct used for inference based on the fitted model fpm. The confidence intervals are calculated as the (1-level)/2 and (1+level)/2 quantiles of the estimated parameters. sandwich_estimator is a boolean that switches on or off the sandwich estimator. If rethrow_error is false (the default value), no error will be thrown if the variance-covariance matrix estimator fails. If it is true, an error will be thrown if the estimator fails.

Pumas.inferMethod

infer(fpm::FittedPumasModel, bts::Pumas.Bootstrap; level=0.95)

Perform bootstrapping by resampling the Subjects from the Population stored in fpm. The confidence intervals are calculated as the (1-level)/2 and (1+level)/2 quantiles of the estimated parameters. The number of samples used in the bootstrapping is bts.samples. bts.ensemblealg specifies the ensemblealg used here. If ensemblealg is EnsembleSerial(), a single thread will be used. If it is EnsembleThreads() (the default value), multiple threads will be used. See the documentation of Bootstrap for more details on constructing an instance of Bootstrap.

Pumas.inferMethod
infer(fpm::FittedPumasModel, sir::SIR; level=0.95, ensemblealg = EnsembleThreads()) -> FittedPumasModelInference

Perform sampling importance re-sampling for the Subjects from the Population stored in fpm. The confidence intervals are calculated as the (1-level)/2 and (1+level)/2quantiles of the estimated parameters.ensemblealgcan beEnsembleThreads()(the default value) to use multi-threading orEnsembleSerial()` to use a single thread.

Pumas.init_paramsFunction
init_params(model::PumasModel)

Create a parameter set populated with the initial parameters of PumasModel.

Pumas.init_randeffsFunction
init_randeffs(model::AbstractPumasModel, param::NamedTuple, [, pop::Population])

Create an object with random effects locked to their mean values.

The optional argument pop takes a Population and changes the output to a vector of such random effects with one element for each subject within the population.

Pumas.inspectMethod
inspect(fpm::AbstractFittedPumasModel; wres_approx::LikelihoodApproximation, nsim::Int, rng)

Output a summary of the model predictions, residuals, Empirical Bayes estimates, and NPDEs (when requested).

Called on a fit output and allows the keyword argument wres_approx for approximation method to be used in residual calculation. The default value is the approximation method used for the marginal likelihood calculation in the fit that produced fpm. The keyword nsim controls the number of times each subject is simulated for npde computation. A FittedPumasModelInspection object with pred, wres, ebes, and npdes is output.

Pumas.lrtestMethod
lrtest(fpm_0::AbstractFittedPumasModel, fpm_A::AbstractFittedPumasModel)::LikelihoodRatioTest

Compute the likelihood ratio test statistic of the null hypothesis defined by fpm_0 against the the alternative hypothesis defined by fpm_A. The pvalue function can be used for extracting the p-value based on the asymptotic Χ²(k) distribution of the test statistic.

Pumas.marginal_nllFunction
marginal_nll(model, subject, param, approx, ...)
marginal_nll(model, population, param, approx, ...)

Compute the marginal negative loglikelihood of a subject or dataset, using the integral approximation approx.

See also loglikelihood.

Pumas.npdeMethod
npde(fpm::AbstractFittedPumasModel; nsim::Integer)

Calculate the normalised prediction distribution errors (NPDE).

The keyword nsim controls the number of times each subject is simulated.

Pumas.penalized_conditional_nllMethod
penalized_conditional_nll(m::AbstractPumasModel, subject::Subject, param::NamedTuple, randeffs::NamedTuple; diffeq_options)

Compute the penalized conditional negative log-likelihood. This is the same as conditional_nll, except that it incorporates the penalty from the prior distribution of the random effects.

Here param can be either a NamedTuple or a vector (representing a transformation of the random effects to Cartesian space).

Pumas.probstableMethod

probstable(fpm::FittedPumasModel)

Return a DataFrame with outcome probabilities of all discrete dependent variables.

Pumas.pvalueMethod
pvalue(t::LikelihoodRatioTest)::Real

Compute the p-value of the likelihood ratio test t based on the asymptotic Χ²(k) distribution of the test statistic.

Pumas.rafterydiagMethod
rafterydiag(b::BayesMCMCResults; subject::Union{Nothing,Int} = nothing, q = 0.025, r = 0.005, s = 0.95, eps = 0.001)

Compute the Raftery and Lewis diagnostic [Raftery1992]. This diagnostic is used to determine the number of iterations required to estimate a specified quantile q within a desired degree of accuracy. The diagnostic is designed to determine the number of autocorrelated samples required to estimate a specified quantile $\theta_q$, such that $\Pr(\theta \le \theta_q) = q$, within a desired degree of accuracy. In particular, if $\hat{\theta}_q$ is the estimand and $\Pr(\theta \le \hat{\theta}_q) = \hat{P}_q$ the estimated cumulative probability, then accuracy is specified in terms of r and s, where $\Pr(q - r < \hat{P}_q < q + r) = s$. Thinning may be employed in the calculation of the diagnostic to satisfy its underlying assumptions. However, users may not want to apply the same (or any) thinning when estimating posterior summary statistics because doing so results in a loss of information. Accordingly, sample sizes estimated by the diagnostic tend to be conservative (too large).

Furthermore, the argument r specifies the margin of error for estimated cumulative probabilities and s the probability for the margin of error. eps specifies the tolerance within which the probabilities of transitioning from initial to retained iterations are within the equilibrium probabilities for the chain. This argument determines the number of samples to discard as a burn-in sequence and is typically left at its default value.

Pumas.read_pumasMethod
read_pumas(filepath::String, args...; kwargs...)
read_pumas(
  df::AbstractDataFrame;
  observations=Symbol[:dv],
  covariates=Symbol[],
  id::Symbol=:id,
  time::Symbol=:time,
  evid::Union{Nothing,Symbol}=nothing,
  amt::Symbol=:amt,
  addl::Symbol=:addl,
  ii::Symbol=:ii,
  cmt::Symbol=:cmt,
  rate::Symbol=:rate,
  ss::Symbol=:ss,
  route::Symbol=:route,
  mdv::Symbol=nothing,
  event_data::Bool=true,
  covariates_direction::Symbol=:right,
  check::Bool=event_data,
  adjust_evid34::Bool=true)

Import NMTRAN-formatted data.

  • df : DataFrame contaning the data to be converted to a Vector{<:Subject}
  • observations : dependent variables specified by a vector of column names
  • covariates : covariates specified by a vector of column names
  • id : specifies the ID column of the dataframe
  • time : specifies the time column of the dataframe
  • evid : specifies the event ID column of the dataframe. See ?Pumas.Event for more details.
  • amt : specifies the dose amount column of the dataframe. See ?Pumas.Event for more details.
  • addl : specifies the column of the dataframe that indicated the number of repeated dose events. If not specified then the value is zero.
  • ii : specifies the dose interval column of the dataframe. See ?Pumas.Event for more details.
  • cmt : specifies the compartment column of the dataframe. See ?Pumas.Event for more details.
  • rate : specifies the infusion rate column of the dataframe. See ?Pumas.Event for more details.
  • ss : specifies the steady state column of the dataframe. See ?Pumas.Event for more details.
  • route : specifies the route of administration column of the dataframe. See ?Pumas.Event for more details.
  • mdv : specifies the the column of the dataframe indicating if observations are missing.
  • event_data : toggles assertions applicable to event data
  • covariates_direction : specifies direction of covariate interpolation. Either :left or :right (default)
  • check : toggles NMTRAN compliance check of the input data
  • adjust_evid34 : toggles adjustment of time vector for reset events (evid=3 and evid=4). If true (the default) then the time of the previous event is added to the time on record to ensure that the time vector is monotonically increasing.
Pumas.remove_italicsMethod

remove_italics(ex, vars)

Ensure that the variables specified in vars are not italicized during latexification of ex.

Pumas.sample_randeffsFunction
sample_randeffs([rng::AbstractRNG=Random.default_rng(),] model::AbstractPumasModel, param::NamedTuple [, pop::Population])

Generate a random set of random effects for model model, using parameters param. Optionally, a random number generator object rng can be passed as the first argument.

The optional argument pop takes a Population and changes the output to a vector of such random effects with one element for each subject within the population.

Pumas.simobsFunction
simobs(pmi::FittedPumasModelInference, population::Population = pmi.fpm.data, randeffs::Union{Nothing, AbstractVector{<:NamedTuple}} = nothing; samples::Int, rng = default_rng(), kwargs...,)

Simulates observations from the fitted model using a truncated multivariate normal distribution for the parameter values when a covariance estimate is available and using non-parametric sampling for simulated inference and bootstrapping. When a covariance estimate is available, the optimal parameter values are used as the mean and the variance-covariance estimate is used as the covariance matrix. Rejection sampling is used to avoid parameter values outside the parameter domains. Each sample uses a different parameter value. samples is the number of samples to sample. population is the population of subjects which defaults to the population associated the fitted model. randeffs can be set to a vector of named tuples, one for each sample. If randeffs is not specified (the default behaviour), it will be sampled from its distribution.

Pumas.simobsFunction
simobs(
  model::AbstractPumasModel,
  population::Union{Subject, Population}
  param,
  randeffs=sample_randeffs(model, param, population);
  obstimes=nothing,
  ensemblealg=EnsembleSerial(),
  diffeq_options=NamedTuple(),
  rng=Random.default_rng(),
)

Simulate random observations from model for population with parameters param at obstimes (by default, use the times of the existing observations for each subject). If no randeffs is provided, then random ones are generated according to the distribution in the model.

Arguments

  • model may either be a PumasModel or a PumasEMModel.

  • population may either be a Population of Subjects or a single Subject.

  • param should be either a single parameter set, in the form of a NamedTuple, or a vector of such parameter sets. If a vector then each of the parameter sets in that vector will be applied in turn. Example: (; tvCL=1., tvKa=2., tvV=1.)

  • randeffs is an optional argument that, if used, should specify the random effects for each subject in population. Such random effects are specified by NamedTuples for PumasModels (e.g. (; tvCL=1., tvKa=2.)) and by Tuples for PumasEMModels (e.g. (1., 2.)). If population is a single Subject (without being enclosed in a vector) then a single such random effect specifier should be passed. If, however, population is a Population of multiple Subjects then randeffs should be a vector containing one such specifier for each Subject. The functions init_randeffs, zero_randeffs, and sample_randeffs are sometimes convenient for generating randeffs:
     randeffs = zero_randeffs(model, param, population)
     solve(model, population, param, randeffs)

If no randeffs is provided, then random ones are generated according to the distribution in the model.

  • obstimes is a keyword argument where you can pass a vector of times at which to simulate observations. The default, nothing, ensures the use of the existing observation times for each Subject.

  • ensemblealg is a keyword argument that allows you to choose between different modes of parallelization. Options include EnsembleSerial(), EnsembleThreads() and EnsembleDistributed().

  • diffeq_options is a keyword argument that takes a NamedTuple of options to pass on to the differential equation solver.

  • rng is a keyword argument where you can specify which random number generator to use.

Pumas.simobsFunction
simobs(fpm::FittedPumasModel, [population::Population,] vcov::AbstractMatrix, randeffs::Union{Nothing, AbstractVector{<:NamedTuple}} = nothing; samples::Int, rng = default_rng(), kwargs...)

Simulates observations from the fitted model using a truncated multi-variate normal distribution for the parameter values. The optimal parameter values are used for the mean and the user supplied variance-covariance (vcov) is used as the covariance matrix. Rejection sampling is used to avoid parameter values outside the parameter domains. Each sample uses a different parameter value. samples is the number of samples to sample. population is the population of subjects which defaults to the population associated the fitted model so it's optional to pass. randeffs can be set to a vector of named tuples, one for each sample. If randeffs is not specified (the default behaviour), it will be sampled from its distribution.

Pumas.simobsMethod
simobs(b::BayesMCMCResults; subject = nothing, samples::Int = 10, simulate_error = false, rng::AbstractRNG = default_rng())

Simulates model predictions using the posterior samples for the population and subject-specific parameters. The predictions are either for the subject with index subject or the whole population if no index is input. samples is the number of simulated predictions returned per subject. If simulate_error is false, the mean of the predictive distribution's error model will be returned, otherwise a sample from the predictive distribution's error model will be returned for each subject. rng is the pseudo-random number generator used in the sampling.

Pumas.simobsMethod
simobs(model::PumasModel, population::Population; samples::Int = 10, simulate_error = false, rng::AbstractRNG = default_rng())

Simulates model predictions for population using the prior distributions for the population and subject-specific parameters. samples is the number of simulated predictions returned per subject. If simulate_error is false, the mean of the predictive distribution's error model will be returned, otherwise a sample from the predictive distribution's error model will be returned for each subject. rng is the pseudo-random number generator used in the sampling.

Pumas.simobsMethod
simobs(model::PumasModel, subject::Subject; samples::Int = 10, simulate_error = false, rng::AbstractRNG = default_rng())

Simulates model predictions for subject using the prior distributions for the population and subject-specific parameters. samples is the number of simulated predictions returned. If simulate_error is false, the mean of the predictive distribution's error model will be returned, otherwise a sample from the predictive distribution's error model will be returned for each subject. rng is the pseudo-random number generator used in the sampling.

Pumas.simobstteFunction
simobstte(
  model::PumasModel,
  subject::Union{Subject,Population},
  param::NamedTuple,
  randeffs::Union{Vector{<:NamedTuple}, NamedTuple, Nothing}=nothing;
  minT=0.0,
  maxT=nothing,
  nT=10,
  repeated=false,
  rng = default_rng())

Simulate observations from a time-to-event model, i.e. one with a TimeToEvent dependent variable and return either a new Subject or Population with random event time stamps in the time vector. The observations vector contains the censoring information (see TimeToEvent for more details on the data coding).

The simobstte function first computes nT values of the survival probability from t=minT to maxT and then interpolates with a cubic spline to get a smooth survival funtion. Given a survival funtion, it is possible to simulate from the distribution by using inverse cdf sampling. Instead of sampling a uniform variate to be compared with the survival function, we sample an exponential variate and compare with cumulative hazard. The two approaches are equivalent. The Roots package is then used for computing the root. If repeated=true then new event times are drawn until maxT has been reached.

Pumas.strip_tMethod

strip_t(ex, vars)

Strip the (t) from the vars in the expression ex.

Works both when (t) is appended using expressions, :(x(t)), or symbols, Symbol("x(t)").

Pumas.tadMethod
tad(t, events) -> time after most recent dose

Converts absolute time t (scalar or array) to relative time after the most recent dose. If t is earlier than all dose times, then the (negative) difference between t and the first dose time is returned instead. If no dose events exist, t is returned unmodified.

Pumas.truncateMethod
truncate(b::BayesMCMCResults; burnin = 0, ratio = 1.0)

Removes burnin samples from the beginning of each chain and then keeps only ratio * 100% of the remaining samples in each chain (i.e. thinning).

Pumas.variablesFunction

variables(dynobj)

Get a list of variable names for a model or dynamics object.

Pumas.vechMethod

vech(A; diag = true)

Converts the symmetric matrix A to a vector of lower triangular elements. If diag is set to true, the vector will include the diagonal elements of A. If diag is false, the vector will only include the strict lower triangle and the diagonal elements of A will be ignored. For example when diag is true:

A = [
  1.0 0.1;
  0.1 2.0
]
Pumas.vech(A, diag = true)

gives

3-element Vector{Float64}:
 1.0
 0.1
 2.0

And when diag is false:

Pumas.vech(A, diag = false)

gives

1-element Vector{Float64}:
 0.1
Pumas.vechinvMethod

vechinv(x; diag = true)

Converts a vector of lower triangular elements to a symmetric matrix. If diag is set to true, the vector is assumed to include diagonal elements. If diag is false, the vector is assumed to only include the strict lower triangle and the diagonal elements of the output matrix will be 1. For example when diag is true:

Pumas.vechinv([1.0, -0.1, 2.0], diag = true)

gives

2×2 Matrix{Float64}:
  1.0  -0.1
 -0.1   2.0

And when diag is false:

Pumas.vechinv([-0.1], diag = false)

gives

2×2 Matrix{Float64}:
  1.0  -0.1
 -0.1   1.0
Pumas.vpcMethod
vpc(fpm::AbstractFittedPumasModel;
       samples::Integer = 499
       qreg_method = IP(),
       observations::Array{Symbol} = [keys(fpm.data[1].observations)[1]],
       stratify_by::Array{Symbol} = Symbol[],
       quantiles::NTuple{3,Float64}=(0.1, 0.5, 0.9),
       level::Real=0.95,
       ensemblealg=EnsembleSerial(),
       bandwidth=nothing,
       maxnumstrats=[6 for i in 1:length(stratify_by)],
       covariates::Array{Symbol} = [:time],
       smooth::Bool = true,
       rng::AbstractRNG=default_rng(),
       obstimes::AbstractVector = [],
       prediction_correction = false)

Computes the quantiles for VPC for a FittedPumasModel or FittedPumasEMModel with simulated confidence intervals around the empirical quantiles based on samples simulated populations.

The following keyword arguments are supported:

  • samples: The number of simulated populations to generate, defaults to 499
  • quantiles::NTuple{3,Float64}: A three-tuple of the quantiles for which the quantiles will be computed. The default is (0.1, 0.5, 0.9) which computes the 10th, 50th and 90th percentile.
  • level::Real: Probability level to use for the simulated confidence intervals. The default is 0.95.
  • observations::Array{Symbol}: The name of the dependent variable to use for the VPCs. The default is the first dependent variable in the Population.
  • stratify_by: The covariates to be used for stratification. Takes an array of the Symbols of the stratification covariates.
  • ensemblealg: This is passed to the simobs call for the samples simulations. For more description check the docs for simobs.
  • bandwidth: The kernel bandwidth in the quantile regression. If you are seeing NaNs or an error, increasing the bandwidth should help in most cases. With higher values of the bandwidth you will get more smoothened plots of the quantiles so it's a good idea to check with your data to pick the right bandwidth. For a DiscreteVPC, the default value is nothing to indicate automatic bandwidth tuning. For a CensoredVPC, the bandwidth should be given as a tuple where the first value applies to the ContinuousVPC and the second to the DiscreteVPC for the censoring variable. For example, (4.0, nothing) applies a bandwidth of 4.0 to the ContinuousVPC and automatically tunes the DiscreteVPC. Automatic tuning is only applicable to DiscreteVPCs at this moment.
  • maxnumstrats: The maximum number of strata for each covariate. Takes an array with the number of strata for the corresponding covariate, passed in stratify_by. It defaults to 6 for each of the covariates.
  • covariates: The independent variable for VPC, defaults to [:time].
  • smooth: In case of discrete VPC used to determine whether to interpolate the dependent variable if independent variable is continuous, defaults to true.
  • rng: A random number generator, uses the default_rng from Random as default.
  • obstimes: The times for simulation in case of continuous VPC, same as obstimes in simobs. Defaults to union of all subject's unique times in the data where observation is not missing.
  • qreg_method: Defaults to IP(). For most users the method used in quantile regression is not going to be of concern, but if you see large run times switching qreg_method to IP(true) should help in improving the performance with a tradeoff in the accuracy of the fitting.
  • prediction_correction: Default to false. Set to true to enable prediction correction that multiplies all observations with the ratio between the mean prediction and each individuals population prediction.
Pumas.wresidualsFunction
wresiduals(fpm::AbstractFittedPumasModel, approx::LikelihoodApproximation; nsim=nothing)

Calculate the individual and population weighted residual.

Takes a fit result, an approximation method for the marginal likelihood calculation which defaults to the method used in the fit and the number of simulations with the keyword argument nsim. If nsim is specified only the Expected Simulation based Individual Weighted Residuals (EIWRES) is included in the output as individual residual and population residual is not computed. Using the FO approximation method corresponds to the WRES and while FOCE(I) corresponds to CWRES. The output is a SubjectResidual object that stores the population (wres) and individual (iwres) residuals along with the subject and approximation method (approx).

Pumas.zero_randeffsFunction
zero_randeffs(model::AbstractPumasModel, param::NamedTuple [, pop::Population])

Create an object to signify that the random effects should be zero.

The optional argument pop takes a Population and changes the output to a vector of such random effects with one element for each subject within the population.

Pumas.ηshrinkageMethod
ηshrinkage(fpm::AbstractFittedPumasModel)

Calculate the η-shrinkage.

Takes the result of a fit as the only input argument. A named tuple of the random effects and corresponding η-shrinkage values is output.

Pumas.ϵshrinkageMethod
ϵshrinkage(fpm::AbstractFittedPumasModel)

Calculate the ϵ-shrinkage.

Takes the result of a fit as the only input argument. A named tuple of derived variables and corresponding ϵ-shrinkage values is output.

StatsAPI.aicMethod
aic(fpm::AbstractFittedPumasModel)

Calculate the Akaike information criterion (AIC) of the fitted Pumas model fpm.

StatsAPI.bicMethod
bic(fpm::AbstractFittedPumasModel)

Calculate the Bayesian information criterion (BIC) of the fitted Pumas model fpm.

StatsAPI.coeftableMethod
coeftable(pmi::FittedPumasModelInference; level = pmi.level) -> DataFrame

Construct a DataFrame of parameter names, estimates, standard error, and confidence interval with confidence level level. The default value of level is the one that was passed in infer.

StatsAPI.coeftableMethod
coeftable(fpm::FittedPumasModel) -> DataFrame

Construct a DataFrame of parameter names and estimates from fpm.

StatsAPI.coeftableMethod
coeftable(cfpm::Vector{<:FittedPumasModel}) -> DataFrame

Construct a DataFrame of parameter names and estimates and their standard deviation from vector of fitted single-subject models vfpm.

StatsAPI.fitMethod

Distributions.fit( model::PumasModel, data::Population, param::NamedTuple, alg::BayesMCMC )

Samples from the posterior distribution of the model's population and subject-specific parameters given the observations in data using the Hamiltonian Monte Carlo (HMC) based No-U-Turn Sampler (NUTS) algorithm. The MCMC sampler is initialized starting from param for the population parameters. The subject-specific parameters are also apprioriately initialized. The Bayesian inference algorithm's options are passed in the alg struct. Run ?BayesMCMC for more details on the options that can be set.

StatsAPI.fitMethod
fit(
  model::PumasModel,
  population::Population,
  param::NamedTuple,
  approx::Union{LikelihoodApproximation, MAP};
  optimize_fn = DefaultOptimizeFN(),
  constantcoef::NamedTuple = NamedTuple(),
  omegas::Tuple = tuple(),
  ensemblealg::DiffEqBase.EnsembleAlgorithm = EnsembleThreads(),
  checkidentification=true,
  diffeq_options))

Fit the Pumas model model to the dataset population with starting values param using the estimation method approx. Currently supported values for the approx argument are FO, FOCE, LaplaceI, NaivePooled, and BayesMCMC. See the online documentation for more details about the different methods.

The argument optimize_fn is used for optimizing the objective function for all approx methods except BayesMCMC. The default optimization function uses the quasi-Newton routine BFGS method from the Optim package. Optimization specific arguments can be passed to DefaultOptimizeFN, e.g. the optimization trace can be disabled and the algorithm can be changed to L-BFGS by passing optimize_fn=DefaultOptimizeFN(Optim.LBFGS(); show_trace=false) to fit. The positional argument is a zero or first or method from Optim and the keywords are used to the available options in Optim. See Optim for more defails.

It is possible to fix one or more parameters of the fit by passing a NamedTuple as the constantcoef argument with keys and values corresponding to the names and values of the fixed parameters, e.g. constantcoef=(σ=0.1,).

When models include an @random block and fitting with NaivePooled is requested, it is required that the user supplies the names of the parameters of the random effects as the omegas argument such that these can be ignored in the optimization, e.g. omegas=(Ω,).

Parallelization of the optimization is supported for most estimation methods via the ensemble interface of DifferentialEquations.jl. Currently, the only supported options are:

  • EnsembleThreads(): the default. Accelerate fits by using multiple threads.
  • EnsembleSerial(): fit using a single thread.
  • EnsembleDistributed(): fit by using multiple worker processes.

The fit function will check if any gradients and throw an exception if any of the elements are exactly zero unless checkidentification is set to false.

Further keyword arguments can be passed via the diffeq_options argument. This allows for passing arguments to the differential equations solver such as alg, abstol, and reltol. The default values for these are AutoVern7(Rodas5(autodiff=true)), 1e-12, and 1e-8 respectively. See the DifferentialEquations.jl documentation for more details.

StatsAPI.fitMethod
fit(
  model::PumasEMModel,
  population::Population,
  param::NamedTuple,
  approx::Union{SAEM,LaplaceI};
  ensemblealg = EnsembleThreads(),
  rng = default_rng()
)

Fit the PumasEMModel to the dataset population using the initial values specified in param. The supported methods for the approx argument are currently SAEM and LaplaceI. See the online documentation for more details about these two methods.

When fitting with SAEM the variance terms of the random effects (Ω) and the dispersion parameters for the error model (σ) are initialized to the identity matrix or 1 as appropriate. They may also be specified. With SAEM, it is reccomended to choose an init larger than the true values to facilitate exploration of the parameter space and avoid getting trapped in local optima early. Currently LaplaceI fits with a PumasEMModel require Ω to be diagonal.

Options for ensemblealg are EnsembleSerial() and EnsembleThreads() (the default). The fit will be parallel if ensemblealg == EnsembleThreads(). SAEM imposes the additional requirement for threaded fits that either the rng used supports Future.randjump or that rng be a vector of rngs, one per thread. The default_rng() supports randjump().

StatsAPI.loglikelihoodMethod
loglikelihood(fpm::AbstractFittedPumasModel)

Compute the loglikelihood of a fitted Pumas model.

StatsAPI.predictMethod
predict(
  fpm::AbstractFittedPumasModel,
  [population::Union{Subject,Population};
  [obstimes::AbstractVector]
)::Union{SubjectPrediction,Vector{SubjectPrediction}}

Compute population and individual predictions for the fitted model fpm. By default, the predictions are computed for the estimation data but the predictions can also be computed for user supplied data by passing either a single subject or a vector of subjects (Population) as the population argument.

If the optional obstimes argument is passed then the time points in obstimes are used for the predictions. Otherwise, the time points of the observations for each subject in the population are used for the predictions.

Any optional keyword arguments used when fitting fpm are reused when computing the predictions.

StatsAPI.predictMethod
predict(
  model::AbstractPumasModel,
  population::Union{Subject,Population},
  param::NamedTuple;
  [obstimes::AbstractVector,
  diffeq_options::NamedTuple]
)::Union{SubjectPrediction,Vector{SubjectPrediction}}

Compute population and individual predictions for either the single subject or vector of subjects (Population) population based on model and the population parameters param.

If the optional obstimes argument is passed then the time points in obstimes are used for the predictions. Otherwise, the time points of the observations for each subject in the population are used for the predictions.

The function allows for extra keyword arguments to be passed on to the differential equations solver through the diffeq_options keyword. See the online documentation for more details.

StatsAPI.predictMethod
predict(
  model::AbstractPumasModel,
  subject::Subject,
  param::NamedTuple,
  [randeffs,];
  [obstimes::AbstractVector,
  diffeq_options::NamedTuple]
)::Union{SubjectPrediction,Vector{SubjectPrediction}}

Compute population and individual predictions for the single subject based on model and the population parameters param. A NamedTuple of random effects, randeffs, can be omitted or provided by the user. If they are omitted, they will be estimated from the data in the subject.

If the optional obstimes argument is passed then the time points in obstimes are used for the predictions. Otherwise, the time points of the observations of the subject are used for the predictions.

The function allows for extra keyword arguments to be passed on to the differential equations solver through the diffeq_options keyword. See the online documentation for more details.

StatsAPI.stderrorMethod
stderror(f::AbstractFittedPumasModel) -> NamedTuple

Compute the standard errors of the population parameters and return the result as a NamedTuple matching the NamedTuple of population parameters.

StatsAPI.vcovMethod
vcov(f::AbstractFittedPumasModel) -> Matrix

Compute the covariance matrix of the population parameters

Pumas.@emmodelMacro
@emmodel

Define a PumasEMModel. It may have the following blocks:

  • @param and @random: Defines fixed and random effects, e.g.
@random begin
    CL     ~ 1 + wt |   LogNormal
    θbioav ~ 1      | LogitNormal
end

Distributions specify the parameter's support. Only Normal (-∞,∞) LogNormal (0,∞), and LogitNormal (0,1) are currently supported.

These define 1 unconstrained parameter per term in the formula, for example CL = exp(μ + wt * β + η) where μ, β ∈ (-∞,∞) and wt is a covariate constant with respect to time. η is a random effect. The @param block is equivalent, with the exception that there is no random effect η. These variables are accessible in @pre, @dosecontrol, @dynamics, @post, and @error blocks.

  • @covariates: Covariates available in the @pre, @dosecontrol, @dynamics, and @post blocks.
  • @pre: Block evaluated before the @dynamics.
  • @dosecontrol: Block specifying dose control parameters. Function of @param and @random variables only.
  • @init: return initial conditions for dynamical variables.
  • @dynamics: Dynamics, equivalent to the functionality in the @model macro.
  • @post: Block evalauted after the @dynamics and before the @error.
  • @error: Block describing the error model. Dispersion parameters are implicit. Y ~ ProportionalNormal(μ) indicates that Y has a Normal distribution with mean μ. Y must be observed subject data, while μ can can be defined in the @param, @random, @pre, @dynamics, or @post blocks.
Pumas.@modelMacro
@model

Defines a Pumas model, with the following possible blocks:

  • @param: defines the model's fixed effects and corresponding domains, e.g. tvcl ∈ RealDomain(lower = 0).
  • @random: defines the model's random effects and corresponding distribution, e.g. η ~ MvNormal(Ω).
  • @covariates: Names of covariates used in the model.
  • @pre: pre-processes variables for the dynamic system and statistical specification.
  • @dosecontrol: defines dose control parameters as a function of fixed and random effects.
  • @vars: define variables usable in other blocks. It is recomended to define them in @pre instead if not a function of dynamic variables.
  • @init: defines initial conditions of the dynamic system.
  • @dynamics: defines the dynamic system.
  • @derived: defines the statistical model of the dependent variables.
  • @observed: lists model information to be stored in the solution.
  • @options specifies model options, including checklinear, inplace, and subject_t0.
Pumas.@ncaMacro
 @nca

Defines a macro that can be used in the observed block to compute NCA metrics of the computed model solutions.

Example

...
@derived begin
  cp := Central/Vc
  DV ~ @. Normal(cp, σ*cp)
end

@observed begin
  nca = @nca cp
  auc = NCA.auc(nca, interval = (0,4))
  cmax = NCA.cmax(nca)
end
...
  • Geweke1991Geweke, J. F. (1991). Evaluating the accuracy of sampling-based approaches to the calculation of posterior moments (No. 148). Federal Reserve Bank of Minneapolis.
  • Heidelberger1983Heidelberger, P., & Welch, P. D. (1983). Simulation run length control in the presence of an initial transient. Operations Research, 31(6), 1109-1144.
  • Raftery1992A L Raftery and S Lewis. Bayesian Statistics, chapter How Many Iterations in the Gibbs Sampler? Volume 4. Oxford University Press, New York, 1992.