Pumas
Docstrings
DataFrames.DataFrame
— TypeDataFrame(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.DataFrame
— MethodDataFrame(bts::Union{SimulatedInference, Bootstraps})
Returns a data frame of the parameter estimates.
Distributions.MvNormal
— MethodMvNormal(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.MvNormal
— MethodMvNormal(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)
Pumas.AnalyticalPKPDProblem
— TypeAnalyticalPKPDProblem(
f,
u0,
tspan,
events,
time,
p,
bioav)
An analytical PK(PD) problem.
Fields:
f
: callable that returns the state of the dynamic systemu0
: initial conditionstspan
: time points that define the intervals in which the problem is definedevents
: events such as dose, reset events, etctime
: event timesp
: a function that returns pre block evaluated at a time pointbioav
: bioavailability in each compartment
Pumas.AnalyticalPKProblem
— TypeAnalyticalPKProblem(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 problemprob2
: a problem that represents the rest of the systemsys2
(optional): the system ofprob2
- used for latexification.
Pumas.Censored
— TypeCensored(distribution::Distribution, lower::Real, upper::Real)::Censored
Construct a censored distribution based on distribution
and the censoring limits lower
and upper
.
Pumas.Central1
— TypeCentral1()
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.Central1Periph1
— TypeCentral1Periph1()
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.Central1Periph1Meta1
— TypeCentral1Periph1Meta1()
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.Central1Periph1Meta1Periph1
— TypeCentral1Periph1Meta1Periph1()
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.Constrained
— TypeConstrained
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.CorrDomain
— TypeCorrDomain(template)
Return a correlation domain.
template
provides both the shape and initial values for the resulting correlation parameter.
- if
template
is anInt
then an identity matrix of sizetemplate
is returned. - if
template
is a squareMatrix
then theCorrDomain
will match its size and will
have initial values according to the template
elements.
Pumas.Depots1Central1
— TypeDepots1Central1()
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.Depots1Central1Periph1
— TypeDepots1Central1Periph1()
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.Depots2Central1
— TypeDepots2Central1()
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.DosageRegimen
— TypeDosageRegimen
Lazy representation of a series of Events.
Fields
data::DataFrame
: The tabular representation of a series ofEvent
s.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 DosageRegimen
s
- 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.PDiagDomain
— Type@param x ∈ PDiagDomain(n::Int; init=ones(n))
Specifies a parameter as a positive diagonal matrix, with initial diagonal elements specified by init
.
Pumas.PSDDomain
— Type@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 Float64
s.
Pumas.ParamSet
— TypeParamSet(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.Population
— TypeA Population
is an AbstractVector
of Subject
s.
Pumas.PumasModel
— TypePumasModel
A model takes the following arguments (arg [= default]
):
param
: aParamSet
detailing the parameters and their domainrandom
: 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 conditionsprob
: a DEProblem describing the dynamics (either exact or analytical)derived
: the derived variables and error distributions (param, randeffs, data, ode vals) -> sampling distobserved = (_, _, _, _, _, _, samples) -> samples
: simulated values from the error model and post processing: (param, randeffs, data, ode vals, samples) -> valsoptions::PumasModelOptions = PumasModelOptions(; subject_t0=false)
: pass additional options, seePumasModelOptions
.syms::PumasModelSymbols = PumasModelSymbols()
: the symbols of all the variables defined inside the model, seePumasModelSymbols
.sys = nothing
: the ModelingToolkit.jlSystem
used to generateprob
- 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 aDict
of arbitrary metadata.
Pumas.RealDomain
— Type@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.Subject
— TypeSubject
The data corresponding to a single subject:
Fields:
id
: identifierobservations
: a named tuple of the dependent variablescovariates
: a named tuple containing the covariates, ornothing
.events
: a vector ofEvent
s.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 = nothing,
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.Subject
— MethodSubject
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.TimeToEvent
— TypeTimeToEvent{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 zero if the variable isn't censored and one if 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.VectorDomain
— Type@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.solve
— Functionsol = 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 aPumasModel
or aPumasEMModel
.population
may either be aPopulation
ofSubject
s or a singleSubject
.param
is parameter set in the form of aNamedTuple
, e.g.(; tvCL=1., tvKa=2., tvV=1.)
.randeffs
is an optional argument that, if used, should specify the random effects for each subject inpopulation
. Such random effects are specified byNamedTuple
s forPumasModels
(e.g.(; tvCL=1., tvKa=2.)
) and byTuples
forPumasEMModel
s (e.g.(1., 2.)
). Ifpopulation
is a singleSubject
(without being enclosed in a vector) then a single such random effect specifier should be passed. If, however,population
is aPopulation
of multipleSubject
s thenrandeffs
should be a vector containing one such specifier for eachSubject
. The functionsinit_randeffs
,zero_randeffs
, andsample_randeffs
are sometimes convenient for generatingrandeffs
:
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 includeEnsembleSerial()
,EnsembleThreads()
andEnsembleDistributed()
.diffeq_options
is a keyword argument that takes aNamedTuple
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.logpdf
— Methodlogpdf(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 thenx
should be a scalar. - If
d
is a constrained multivariate distribution thenx
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.gsa
— FunctionGlobalSensitivity.gsa(model, population, params, method, vars, prangelow, prangehigh; kwargs...)
Function to perform global sensitivty analysis
The arguments are:
model
: aPumasModel
, either defined by the@model
DSL or the function-based interface.population
: aPopulation
.params
: a named tuple of parameters.method
: one of theGSAMethod
s 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.cond
— Methodcond(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.conditional_nll
— Methodconditional_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.empirical_bayes
— Methodempirical_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.eventnum
— Methodeventnum(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.expectation
— Methodexpectation(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.expectation
— Methodexpectation(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.expectation
— Methodexpectation(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.findinfluential
— Functionfindinfluential(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.icoef
— Methodicoef(fpm::AbstractFittedPumasModel)::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 return a vector of covariate interpolation objects. Each of which can be evaluated at any time point. Each of the covariate interpolation objects can be converted to a DataFrame
by calling the DataFrame
constructor. Hence, a complete DataFrame
of the individual coefficients can be obtained by calling reduce(vcat, DataFrame.(icoef(fpm)))
.
Pumas.infer
— Methodinfer(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.infer
— Methodinfer(fpm::FittedPumasModel, bts::Pumas.Bootstrap; level=0.95)
Perform bootstrapping by resampling the Subject
s 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.infer
— Methodinfer(fpm::FittedPumasModel, sir::SIR; level=0.95, ensemblealg = EnsembleThreads()) -> FittedPumasModelInference
Perform sampling importance re-sampling for the Subject
s 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 be
EnsembleThreads()(the default value) to use multi-threading or
EnsembleSerial()` to use a single thread.
Pumas.init_params
— Functioninit_params(model::PumasModel)
Create a parameter set populated with the initial parameters of PumasModel
.
Pumas.init_randeffs
— Functioninit_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.inspect
— Methodinspect(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.lrtest
— Methodlrtest(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.probstable
— Methodprobstable(fpm::FittedPumasModel)
Return a DataFrame with outcome probabilities of all discrete dependent variables.
Pumas.pvalue
— Methodpvalue(t::LikelihoodRatioTest)::Real
Compute the p-value of the likelihood ratio test t
based on the asymptotic Χ²(k)
distribution of the test statistic.
Pumas.read_pumas
— Methodread_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 aVector{<:Subject}
observations
: dependent variables specified by a vector of column namescovariates
: covariates specified by a vector of column namesid
: specifies the ID column of the dataframetime
: specifies the time column of the dataframeevid
: 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 datacovariates_direction
: specifies direction of covariate interpolation. Either:left
or:right
(default)check
: toggles NMTRAN compliance check of the input dataadjust_evid34
: toggles adjustment of time vector for reset events (evid=3
andevid=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.sample_randeffs
— Functionsample_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.simobs
— Functionsimobs(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.simobs
— Functionsimobs(
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 aPumasModel
or aPumasEMModel
.population
may either be aPopulation
ofSubject
s or a singleSubject
.param
should be either a single parameter set, in the form of aNamedTuple
, 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 inpopulation
. Such random effects are specified byNamedTuple
s forPumasModels
(e.g.(; tvCL=1., tvKa=2.)
) and byTuples
forPumasEMModel
s (e.g.(1., 2.)
). Ifpopulation
is a singleSubject
(without being enclosed in a vector) then a single such random effect specifier should be passed. If, however,population
is aPopulation
of multipleSubject
s thenrandeffs
should be a vector containing one such specifier for eachSubject
. The functionsinit_randeffs
,zero_randeffs
, andsample_randeffs
are sometimes convenient for generatingrandeffs
:
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 eachSubject
.ensemblealg
is a keyword argument that allows you to choose between different modes of parallelization. Options includeEnsembleSerial()
,EnsembleThreads()
andEnsembleDistributed()
.diffeq_options
is a keyword argument that takes aNamedTuple
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.simobs
— Functionsimobs(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 multi-variate normal distribution for the parameter values. The optimal parameter values are used 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.simobstte
— Functionsimobstte(
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 function first computes nT
values of survival probability from t=minT
to maxT
and then interpolate 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 use with the survival probability, we sample an exponential and compare to the cumulative hazard which is 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.tad
— Methodtad(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.vpc
— Methodvpc(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=2,
maxnumstrats=[6 for i in 1:length(stratify_by)],
covariates::Array{Symbol} = [:time],
smooth::Bool = true,
rng::AbstractRNG=default_rng(),
obstimes::AbstractVector = [])
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 to499
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 is0.95
.observations::Array{Symbol}
: The name of the dependent variable to use for the VPCs. The default is the first dependent variable in thePopulation
.stratify_by
: The covariates to be used for stratification. Takes an array of theSymbol
s of the stratification covariates.ensemblealg
: This is passed to thesimobs
call for thesamples
simulations. For more description check the docs forsimobs
.bandwidth
: The kernel bandwidth in the quantile regression. If you are seeingNaN
s or an error, increasing the bandwidth should help in most cases. With higher values of thebandwidth
you will get more smoothened plots of the quantiles so it's a good idea to check with your data the rightbandwidth
.maxnumstrats
: The maximum number of strata for each covariate. Takes an array with the number of strata for the corresponding covariate, passed instratify_by
. It defaults to6
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 totrue
.rng
: A random number generator, uses thedefault_rng
fromRandom
as default.obstimes
: The times for simulation in case of continuous VPC, same asobstimes
insimobs
. Defaults to union of all subject's unique times in the data where observation is not missing.qreg_method
: Defaults toIP()
. For most users the method used in quantile regression is not going to be of concern, but if you see large run times switchingqreg_method
toIP(true)
should help in improving the performance with a tradeoff in the accuracy of the fitting.
Pumas.wresiduals
— Functionwresiduals(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_randeffs
— Functionzero_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.ηshrinkage
— Methodη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.ϵshrinkage
— Methodϵ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.
StatsBase.aic
— Methodaic(fpm::AbstractFittedPumasModel)
Calculate the Akaike information criterion (AIC) of the fitted Pumas model fpm
.
StatsBase.bic
— Methodbic(fpm::AbstractFittedPumasModel)
Calculate the Bayesian information criterion (BIC) of the fitted Pumas model fpm
.
StatsBase.coeftable
— Methodcoeftable(pmi::FittedPumasModelInference) -> DataFrame
Construct a DataFrame of parameter names, estimates, standard error, and confidence interval from a pmi
.
StatsBase.coeftable
— Methodcoeftable(fpm::FittedPumasModel) -> DataFrame
Construct a DataFrame of parameter names and estimates from fpm
.
StatsBase.coeftable
— Methodcoeftable(cfpm::Vector{<:FittedPumasModel}) -> DataFrame
Construct a DataFrame of parameter names and estimates and their standard deviation from vector of fitted single-subject models vfpm
.
StatsBase.fit
— Methodfit(
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()
.
StatsBase.fit
— Methodfit(
model::PumasModel,
population::Population,
param::NamedTuple,
approx::Union{LikelihoodApproximation, MAP};
optimize_fn = DefaultOptimizeFN(),
constantcoef::NamedTuple = NamedTuple(),
omegas::Tuple = tuple(),
ensemblealg::DiffEqBase.EnsembleAlgorithm = EnsembleSerial(),
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. The default is EnsembleSerial()
. Currently, the only supported parallelization for model fitting is EnsembleThreads()
.
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.
StatsBase.loglikelihood
— Methodloglikelihood(fpm::AbstractFittedPumasModel)
Compute the loglikelihood of a fitted Pumas model.
StatsBase.predict
— Methodpredict(
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.
StatsBase.predict
— Methodpredict(
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.
StatsBase.predict
— Methodpredict(
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.
StatsBase.stderror
— Methodstderror(f::AbstractFittedPumasModel) -> NamedTuple
Compute the standard errors of the population parameters and return the result as a NamedTuple
matching the NamedTuple
of population parameters.
StatsBase.vcov
— Methodvcov(f::AbstractFittedPumasModel) -> Matrix
Compute the covariance matrix of the population parameters
Pumas.@emmodel
— Macro@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 thatY
has aNormal
distribution with meanμ
.Y
must be observed subject data, whileμ
can can be defined in the@param
,@random
,@pre
,@dynamics
, or@post
blocks.
Pumas.@model
— Macro@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, includingchecklinear
,inplace
, andsubject_t0
.
Pumas.@nca
— Macro @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
...