ESCMID Postgraduate Technical Workshop
Advanced Antimicrobial Pharmacokinetic and Pharmacodynamic Modelling & Simulation
Population PK/PD modelling:
Software Choices
Monday, October 6th, 2014
Liverpool, UK
Jürgen B Bulitta, PhD
Senior Research Fellow
Monash University, Melbourne, Australia
Adjunct Assistant Professor
SUNY Buffalo, NY, USA
[email protected]
Author’s Copyright © 2014. All rights reserved.
Outline
1. “The right task”
2. Software tools for
•
•
•
•
•
•
•
•
Non-compartmental analysis
Exposure-response & exposure-toxicity relationships
Simulation of PK and PK/PD profiles
Estimation of PK and PK/PD model parameters
Optimal dosing of a patient population (Monte Carlo simulation)
Optimal dosing of individual patients
Optimal clinical trial design
Clinical trial simulation
3. Concluding remarks
2
Author’s Copyright © 2014. All rights reserved.
Most appropriate tool for the required task
Define and discuss the objective(s)
Evaluate the time available
Gain a thorough understanding of the data
Discuss / select relevant approaches
 May need to revisit the approach / time
Select the most suitable
software tool and data
analytical approach.
Secondary, as long as
one of suitable tools
is chosen and used
by a skilled person.
Plausibility checks!
A model without a task is not very useful.
Often, several tools can yield suitable results.
There is no single tool which does it all!
Communication of approach & results is critical!
Author’s Copyright © 2014. All rights reserved.
3
Relevant tasks for anti-infective PK/PD
1. Determine Cmax, AUC, T>MIC, clearance (CL), volume (Vd), and other
parameters via non-compartmental analysis (NCA).
2. Fit an exposure-response (or exposure-toxicity) relationship
correlating fAUC/MIC with bacterial counts at 24 h, for example.
3. Simulate antibiotic concentrations for dosage regimens in mice or humans.
4. Develop a population PK or PK/PD model to:
a) Propose optimized empiric dosage regimens for a patient population
via Monte Carlo simulations.
b) Identify covariates that affect dosing (e.g. renal function, body size).
c) Evaluate proposed mechanisms of action, resistance, and synergy.
 To design mechanistically optimized (combination) dosage regimens.
d) Achieve therapeutic target goals most precisely in an individual patient.
Account for MIC / pathogen, renal function, other diseases, etc.
e) Propose a clinical trial design that provides safe and effective regimens
and optimally informs a PK or PK/PD model for patients.
f) Evaluate the robustness of a clinical trial design.
4
Author’s Copyright © 2014. All rights reserved.
NCA (Cmax,
AUC, T>MIC
Exposure
response
Safety model
/ relationship
Simulation of drug
concentrations
Optimal
dosing
Individual
patient
Patient
population
Tasks and Tools
(Population)
PK model
(Population)
PD model
Predict effect and
safety for in vitro and
animal models, and
ultimately patients
Optimal trial design
 PK/PD model
Robust Clinical
Trial design
5
Author’s Copyright © 2014. All rights reserved.
NCA (Cmax,
AUC, T>MIC
Any (!)
estimation
tool
Exposure
response
Skill level:
Basic
Intermediate
Advanced
NONMEM, ADAPT, S-ADAPT,
Pmetrics/USC*PACK, Monolix,
Phoenix/NLME, WinBUGS, SAS,
Matlab, R-packages, others
WinNonlin,
Kinetica,
others
Safety model
/ relationship
(Population)
PK model
(Population)
PD model
Berkeley Madonna
Any tool for sim. / est.
Simulation of drug
concentrations
Predict effect and
safety for in vitro and
animal models, and
ultimately patients
Optimal
dosing
BestDose
ID-ODS
DoseMe
TCI Works
MWPharm
First-dose Individual
patient
CADDy
WinAUIC
Berkeley Madonna, MicLab
Any population tool
Patient
population
Optimal trial design
 PK/PD model
Robust Clinical
Trial design
PFIM,
WinPOPT,
POPT,
PopDes,
POPED,
etc.
Clinical Trial
Simulator
6
Author’s Copyright © 2014. All rights reserved.
Concept of population modeling
– fit one subject in the perspective
of all the other subjects –
7
Least squares method
Observations
Fitted line
Offsets
y
o
o
o
o
 Minimize sum of
squares of the offsets
o
o
o
Option to weight residuals
(e.g. by 1/x)
 weighted least squares
o
x
This method is NOT used for population PK/PD
estimation.
8
Author’s Copyright © 2014. All rights reserved.
Naïve pooling
4
2
1
0.5
0
2
Time
(h)6
4
16
8
4
2
1
0.5
8
Time (h)
CLAv = 9.5 L/h
VAv = 39 L
0
2
Time
4 (h)
6
Time (h)
CLPooling = 9.7 L/h
VPooling = 40 L
Naïve averaging and naïve pooling
- All observations considered to come from one subject
- No individual CL and V estimates
- Ignore between subject variability (BSV)
- cannot use actual sampling / dosing times
Naïve averaging
- Average or median is
computed
- actively reduces available
observations
Naïve pooling
- One ID per dose level
- Preferable to averaging
- Only yields fairly unbiased
estimates, if BSV is small
(maybe: %CV <15%)
Directly estimated parameters
Calculated from parameter estimates.
8
16
16
4
CL = 12 L/h8
V = 36 L 4
2
2
Concentration (mg/L)
8
Standard-two-stage
method
(mg/L) (mg/L)
Concentration
(mg/L) Concentration
Concentration
16
Concentration (mg/L)
Concentration (mg/L)
Naïve averaging
8
1
0.5
16 0
2
8
1
0.5
CL4 = 86 L/h8
V = 41 L
Time (h)
4
CL = 10 L/h
V = 40 L
2
1
0
2
4
6
8
Time (h)
“Fit each
subject
separately”
0.5
0
2
4
6(h)
Time
8
Time (h)
Descriptive statistics, Average  SD
CL = 10.0  2.00 L/h
V = 39.0  2.65 L
Standard-two-state method
- usually yields reasonable
estimates for average PK
parameters, if data are rich
- usually yields biased (too
These are not
large) estimates for the
population methods! between subject variability
Author’s Copyright © 2014. All rights reserved.
9
Borrowing of information
Sparse data
IV bolus dose, 1-compartment model
Drug concentration
100
10


Li  - 2  log  lyi | θ, σ   hθ | μ, Ω dθ
 -

1
lyi | θ, σ  How well does the curve
0.1
100 0
Drug concentration
Considering all relevant values of
CL and V can be done by sampling e.g.
1000 combinations of CL and V for
patient i and then calculating the
individual subject objective function (Li).
4
8
12
16
20
24
Time (h)
generated from one set
of θ fit the observations?
hθ| μ, Ω How likely is it to obtain
10
the current set of
parameters (θ) given
the parameter variability
model.
1
0.1
0
4
8
12
Time (h)
16
20
24
Author’s Copyright © 2014. All rights reserved.
10
Applications of population PK/PD models
Model
structure
Observations
(“data”)
Optimizing
individual
doses
Estimation
Accept or
Reject or
Revise
Parameter
estimates
Qualify
(validate)
Re-estimate
New
experimental
data
Expert
knowledge
Author’s Copyright © 2014. All rights reserved.
11
Tasks of Population PK/PD Software
Simulation
Optimizing
individual
doses
Estimation
Data processing & Plotting
Optimal
design
12
Author’s Copyright © 2014. All rights reserved.
Tasks of Population PK/PD Software
Optimizing indiv.
Doses: BestDose,
Simulation
Berkeley Madonna,
Model Maker,
acslXtreme,
Stella, Gepasi,
and all estimation
programs
ID-ODS, DoseMe
TCI Works, MWPharm
First-dose, CADDy
WinAUIC
Data processing &
Any programming language,
Perl, Python, Visual Basic,
Matlab, S-Plus, R, SAS, etc.
Optimal
design
Review of population optimal design
software: Mentré F et al. PAGE 2007.
Estimation
NONMEM,
ADAPT V,
Plotting
S-ADAPT,
Pmetrics/NPAG,
R, S-Plus,
WinNonlin,
WinNonlin,
Matlab, EXCEL, Phoenix/NLME,
SigmaPlot, etc. PDx-MC-PEM,
SAAM II,
many more
ADAPT, S-ADAPT,
PFIM, WinPOPT, POPT,
PKStaMP, PopDes,
POPED,
Clinical Trial Simulator
13
Author’s Copyright © 2014. All rights reserved.
Estimation approaches
No good
reason to use
these algorithms (Sheiner
(Naïve averaging
or Naïve pooling)
& Beal, JPB
1981; 9:635-51).
(Standard-twostage approach)
FO method (e.g.
in NONMEM)
is now outdated
by at least
2 generations
of algorithms
or >15 years.
Estimation
Two Stage Hierarchical
Population Methods
(Adapt V, S-Adapt, NONMEM,
Pmetrics/NPAG, Monolix, PEM,
Phoenix/NLME, SAS, R)
Not recommended to get initial
estimates via naïve methods or an
outdated population algorithm:
 risk of local minima.
Three Stage Hierarchical
Population Methods
(BUGS, S-Adapt, Monolix,
NONMEM, Pmetrics, others)
Author’s Copyright © 2014. All rights reserved.
14
Simulate the antibiotic
concentrations for dosage
regimens in mice or humans
15
Simulations without between patient variability
(also called deterministic simulations)
Concentration (mg/L)
80
70
Central
compartment
60
50
Peripheral
compartment
40
30
20
10
10
•
400 mg
800 mg
0.1
0.01
0
12
24
36
12
24
48
48
MichaelisMenten
elimination
Concentration (mg/L)
100
Almost any PK/PD and other software
package (incl. Excel®) will do this task.
Time to implement the model (i.e. coding
and debugging) guaranteed to be longer
than computation time (usually ~msec).
Robust differential equation solver critical.
Ability to plot multiple / all variables highly
important to better understand a model.
Berkeley Madonna very powerful for this
task (and easily usable by beginners).
36
Time (h)
Time (h)
•
•
200 mg
0.001
0
•
100 mg
1
0
•
Linear 2-compartment model at
different doses
100
Concentration (mg/L)
Rise to steady-state after
multiple dosing
10
100 mg
200 mg
400 mg
1
800 mg
0.1
Km =
0.1 mg/L
0.01
0.001
0
2
Author’s Copyright © 2014. All rights reserved.
4
6
Time (h)
8
10
12
16
Propose optimized empiric dosage
regimens for a patient population
via Monte Carlo simulations.
Identify covariates that affect dosing
(e.g. renal function, body size, etc.)
17
The average patient?
Sörgel F. Chemotherapie Journal (2003)
Between Subject Variability
PK Model for 1 patient
PK model for a
population of patients
Variability in PK parameters
One set of PK parameters
for Clearance (CL), volume
of distribution (V), and
absorption rate constant (ka).
Bulitta JB, PhD Thesis, 2006.
Author’s Copyright © 2014. All rights reserved.
Absorption rate
constant ka
not shown.
Derivation of PKPD susceptibility breakpoints
Calculate fT>MIC for each
patient and for all relevant MICs
Population PK model
V2
V3
10
Simulate 10,000
virtual patients
8
TABS
6
4
MIC = 2 mg/L
2
0
0
2
4
6
8 10 12 14 16 18 20 22 24
Time (h)
Bulitta JB, PhD Thesis, 2006.
Fraction
below the
40%
target:
30%
Target:
9.6h
62%
38%
20%
10%
0%
0 2 4 6 8 10 12 14 16 18
fT>MIC
(h) MIC
for MIC
2 mg/L
fT>MIC for
= 2ofmg/L
Target: fT>MIC TARGET=0.4
40% of 24h = 9.6h
Probability of
Probability of
target
targetattainment
attainment
V1
50%
Fraction of samples
Fraction
Fractionofofpatients
subjects
CL
100%
90%
80%
70%
60%
50%
40%
TLAG
30%
20%
10%
0%
PKPD
breakpoint:
1 mg/L
0.25
Author’s Copyright © 2014. All rights reserved.
0.5
1
2
MIC (mg/L)
38%
4
8
20
Types of parameter variability models
Mixture of parametric
distributions
Parametric
distribution
Non-parametric
every patient can
be its own subpopulation
one sub-population
representing several
comprising the entire
sub-populations
population
0.25
0.25
0.25
0.2
0.15
0.1
0.2
Probability
Probability
Probability
0.2
0.15
0.1
0.15
0.1
0.05
0.05
0.05
0
0
0
0
2
4
6
8
10 12 14
Clearance (L/h)
0
2
4
6
8
10 12 14
Clearance (L/h)
NONMEM, ADAPT, S-ADAPT,
Number of sub-populations
Monolix, Phoenix/NLME,
needs to be user-specified.
WinBUGS, SAS, Matlab,
Software: ADAPT V, S-ADAPT,
R-packages, Pmetrics /
NONMEM, Monolix, others.
USC*PACK, others.
0
2
4
6
8
10 12 14
Clearance (L/h)
Available in Pmetrics /
USC*PACK (NPEM,
NPAG) and in NPML.
21
Author’s Copyright © 2014. All rights reserved.
Software for population PK modelling
and Monte Carlo simulations
Monte Carlo
Simulation
Berkeley Madonna,
acslXtreme, Crystal Ball,
others
Specialized tool in
antimicrobial PK/PD:
MicLab (Medimatics)
Simulation and Estimation
NONMEM, ADAPT V, S-ADAPT, BUGS,
Pmetrics/NPAG, Phoenix/NLME, R, others
Important software characteristics:
• Can the software handle correlation
between PK parameters (e.g. CL and V)?
• Availability of post-processing tools to
summarize the results (e.g. Perl or R scripts).
• Is the simulation tool robust (e.g. differential
equation solver; capabilities to assure / enforce
positive semi-definite covariance matrix).
• Ability to simulate with a prior distribution
(Full Bayesian approach, 3-stage method).
22
Author’s Copyright © 2014. All rights reserved.
Mechanism-based modeling /
systems pharmacology
Evaluating proposed mechanisms of
action, resistance and synergy
 Design mechanistically optimized
monotherapy and combination
dosage regimens.
23
Mechanism-based modeling of antibiotic action
and resistance
Resistance
often limits
access to
target site.
Time course
& mechanisms
of activity and
resistance.
Efflux pumps,
Beta-lactamase
activity
Error-prone
replication
Bulitta JB et al. Curr Pharm Biotechnol, 2011: 2044-2061
Rationally Optimizing Combination Chemotherapy
DRUG A
Receptor 1
Killing via
pathway A
DRUG B
Receptor 2
Killing via
pathway B
Receptor 3
Enhances
killing
Bacterial killing
Eradicate
Inhibit
(upregulation)
Phenotypic
resistance
mechanism(s)
Genotypically Genotypically Mutation Genotypically Phenotypically
‘resistant’ cells resistant nonsusceptible
intermediate
replicating
Pre-existing vs.
bacteria
bacteria
persisters
de novo formation
Spontaneous or error-prone mutation
Mechanism-based modeling integrates time course & probabilities
Authors’ copyright © 2011, all rights reserved.
25
Systems pharmacology
Mechanism-based PK/PD modeling
Estimation and simulation
Software:
S-ADAPT, ADAPT V,
NONMEM, Pmetrics,
Monolix
Important software characteristics:
• Robustness and efficiency of estimation algorithm  models with many
(often >20) parameters to be simultaneously estimated.
• Pre- and post-processing tools (e.g. SADAPT-TRAN, NM-TRAN,
AMGET [for ADAPT V]) extremely important.
• Automated code enhancing & debugging enables beginner and
intermediate users to perform such modeling and often accelerates
model coding >10-fold.
• Robust computational tool (differential equation solver, other features)
• Customization of result plots highly important.
• Parallelized estimation.
26
Author’s Copyright © 2014. All rights reserved.
Exact calculation of the integral of one simplified
function approximating the log-likelihood
True probability distribution
characterizing the uncertainty
of clearance for the ith patient
This scenario with one
support point applies to the
LAPLACIAN method (FOCE
relies on more approximations
than LAPLACIAN; FO method
far worse)
Density
Function with a
precisely known integral
BUT: It is unknown
how well the simplified
function approximates
the log-likelihood!
mode
Model Parameter (e.g. clearance)
Author’s Copyright © 2014. All rights reserved.
27
Approximating the true log-likelihood
as precisely as needed
Density
True probability distribution
for clearance of patient i
Use a sampling function to
randomly draw points on the
x-axis. The probability of a
point drawn at a position is
determined by the sampling
function (also called
proposal density).
Model Parameter (e.g. clearance)
Author’s Copyright © 2014. All rights reserved.
28
Approximating the true log-likelihood
as precisely as needed
Density
True probability distribution
for clearance of patient i
~zero
density
at this tail
For each randomly drawn point of the
sampling function, calculate the exact
value of the true density.
 Interpolate between support points.
Important sampling
algorithm approximates
the true log-likelihood as
closely as needed.
Model Parameter (e.g. clearance)
Author’s Copyright © 2014. All rights reserved.
29
Model predictions
at 106 CFU/mL
inoculum
Bulitta JB, et al.
Current Pharmaceutical
Biotechnology 2011;
12: 2044-2061.
Log10 (CFU/mL)
Prospective ‘validation’ based on external in vivo data
fT>MIC
40%
60%
93%
A: 30 min
75%
B: 5 h
100%
Time (h)
Targets in
patients:
CFU: Colony-forming unit.
Time > MIC
Log10 CFU per
lung at 24 h
Ambrose PG
et al. CID
2007, 44:
79-86.
Neutropenic mouse lung infection model (at 24 h)
Bacteriostasis
target ~35% fT>MIC
Craig WA. Clin Infect Dis 1998, 26:1-12.
Andes D & Craig WA. IJAA 2002; 19:261-8.
 The model
quantitatively
predicted the
PKPD target
values for
cephalosporins.
Near-maximal cell
killing target ~65% fT>MIC
30
Achieve therapeutic target goals most
precisely in an individual patient
Account for MIC / pathogen,
renal function, other diseases, etc.
31
Different Bayesian updating methods to individualize
PK parameters in an unstable critically ill patient
Available in
Pmetrics,
Best Dose
Michael Neely,
Roger Jelliffe
et al.
Bulitta JB et al. Curr Pharm Biotechnol, 2011: 2044-2061
Optimized dosing of individual patients
Simulated concentrations
Target
2
4
5
1.0
6
0.1
0
4
8
12
16
20
24
Time (h)
Time
(h)
(mg/L)
Concentration
Concentration
100.0
Observation +/- SE
12
0.0
13
5
6
7
12
0.1
4
8
12
16
Time (h)
Time (h)
20
0 2 4 6 8 10 12 14
24
0.5
2
4
1.0
0.2
0.1
3
13
Typical
patient
Loading dose + Continuous inf
Dose optimization based
on “typical patient”
0.3
7
Clearance
(L/h)
10.0
0
Probability
3
10.0
Clearance distribution
0.4
MAP
Bayesian
Clearance
(L/h)
individualization
0.4
Probability
(mg/L)
Concentration
Concentration
100.0
0.5
0.3
Concentration (mg/L)
Clearance
(L/h)
Dose optimization based on full
non-parametric distribution
(Multiple-Model Dosage Design)
Hit target most precisely!
0.2
0.1
0.0
0 2 4 6 8 10 12 14
Clearance (L/h)
Author’s Copyright © 2014. All rights reserved.
Time (h)
33
Real data from Dr Roger Jelliffe.
Optimal Dose Selection programs for antibiotics
Roberts JA et al., Lancet Infect Dis 2014; 14:498-509
34
Software choices for antimicrobial PK/PD
1. Carefully defining the objective should always be the first step.
2. A variety of powerful software tools are available and accessible
also to beginner and intermediate users. No single tool does it all.
3. Very significant improvements in software usability, efficiency and
robustness of algorithms were achieved over the last 10-15 years.
4. Model estimation time is usually no longer a real limitation, even for
complex models with >30 parameters. (Parallelized estimation!)
In the future, a semi-automated code generator will be very helpful.
5. Performing a Monte Carlo simulation to optimize empiric dosage
regimens is very helpful. However, this is NOT the same as
selecting an optimal dosage regimen for an individual patient.
6. Softwares for optimal dosing of individual patients are available and
are being enhanced for different devices (incl. smart phones).
7. Communication / explanation of results by a skilled modeler is critical.
35
Backup slides
36
Exposure response and
exposure toxicity
relationships
37
Non-compartmental analysis
Many software packages:
Commercial software:
Phoenix / WinNonlin
Kinetica
EquivTest
Matlab/SimBiology
PKSolutions (Excel based)
Topfit
Calculate / obtain:
Cmax, AUC, T>MIC, CL, Vd, etc.
Free software:
Bear
PK for R
S-ADAPT
Pharm PK
archives &
Wikipedia
WinNonlin industry
standard (documentation
excellent; FDA CFR Part 11)
38
Author’s Copyright © 2014. All rights reserved.
Exposure response – continuous outcome data
Log10 CFU per lung at 24 h
Cefotaxime vs K. pneumoniae in neutropenic lung infection model (after 24 hours of therapy)
Cmax / MIC
Modeling approach:
Amount of data
Type of output data
Signal
Time-course data
Between subject variability
Recommended algorithm
AUC / MIC
Time > MIC
Usually significant
Continuous
Often strong
No (or not used)
Yes (but not used in analysis)
Maximum likelihood or
(Weighted least squares)
Craig WA. Clin Infect Dis 1998, 26:1-12.
Pictures from: Drusano GL. Nat Rev Microbiol 2004; 2:289-300.
Software tools:
Many nonlinear regression
tools.
Suggestions:
ADAPT (maximum likelihood,
free, user-friendly).
WinNonlin (as commercial
package).
Many other tools equally
capable.
39
Probability of cure or toxicity
Exposure response – Dichotomous (Yes/No data)
MIC: 0.25 mg/L
Afibrile
MIC: 1 mg/L
on day 7
MIC: 4 mg/L
Nephotoxicity for Q24h dosing
Nephotoxicity for Q12h dosing
AUC
Modeling approach:
Amount of data
Type of output data
Signal
Time-course of risk
Between subject variability
Less (especially for tox.)
Dichotomous
(e.g. Live/Dead, Yes/No)
Often weaker
No (or usually not used)
Yes (but not used in analysis)
Recommended algorithm
Maximum likelihood or
Bayesian algorithms
Drusano & Louie. Antimicrob Agents Chemother 2011; 55:2528-31.
Software tools:
Statistical packages for logistic
regression and CART analysis.
Parametric hazard models to
describe time-dependent risks
(eg. of death or adverse events).
Suggestions:
Systat for logistic regression.
(Other tools equally capable.)
Population modelling tools for
parametric hazard models:
NONMEM, S-ADAPT, Monolix,
etc.
40
Optimize
empiric dosage
regimens via
Monte Carlo
simulation
Individual patient’s PK
and MIC are unknown.
Can incorporate
influential covariates
(e.g. renal impairment).
Bulitta JB et al.
Curr Pharm Biotechnol,
2011: 2044-2061
Achieve Target Goals
Sequential MultipleModel (MM) Bayesian
updating
Pmetrics
www.lapk.org
Individual
patient data
Population
PK model
Individual
patient data
Interacting MultipleModel (MM) approach.
Unstable
patients!
Here, PK parameters
can change over time
(unstable patients!)
Slide kindly provided
by Dr. Roger Jelliffe.
42
Authors’ Disclosure
I have never received any funding or any other financial incentive from
a software company and am not associated with any commercial
software package.
I am the creator and developer of the free, open-source
SADAPT-TRAN tool and am a passionate user of Berkeley Madonna.
I have been an active user of a Pharsight Academic Excellence
License for several years since 1999.
I have received collaborative research grants from Pfizer, Trius,
Cempra, CSL, Cubist, Novartis, and Boryung.
None of this collaborative work is related to this presentation.
43
Author’s Copyright © 2014. All rights reserved.