Statistical Analysis of Repeated Measures Data Using SAS (and R)

Download Report

Transcript Statistical Analysis of Repeated Measures Data Using SAS (and R)

Lecture 4
Non-Linear and Generalized
Mixed Effects Models
Ziad Taib
Biostatistics, AZ
MV, CTH
April 2011
1
Date
Part I
Generalized Mixed Effects
Models
2
Date
Outline of part I
1. Generalized Mixed Effects Models
1.
2.
3.
4.
Formulation
Estimation
Inference
Software
2. Non-linear Mixed Effects Models in Pharmacokinetics
1.
2.
3.
4.
Basic Kinetics
Compartmental Models
NONMEM
Software issues
Name, department
3
Date
Various forms of models and relation between them
Classical statistics (Observations are random, parameters are unknown constants)
LM: Assumptions:
1.
independence,
2.
normality,
3.
constant parameters
LMM:
Assumptions 1)
and 3) are modified
GLM: assumption 2)
Exponential family
Repeated measures:
Assumptions 1) and 3)
are modified
GLMM: Assumption 2) Exponential
family and assumptions 1) and 3) are
modified
Longitudinal data
Maximum likelihood
LM - Linear model
Non-linear models
GLM - Generalised linear model
LMM - Linear mixed model
GLMM - Generalised linear mixed model
Name, department
4
Date
Bayesian statistics
Example 1
Toenail Dermatophyte Onychomycosis
 Common toenail infection, difficult to treat, affecting more
than 2% of population. Classical treatments with antifungal
compounds need to be administered until the whole nail
has grown out healthy.
 New compounds have been developed which reduce
treatment to 3 months.
5
Date
Example 1 :
• Randomized, double-blind, parallel group, multicenter study
for the comparison of two such new compounds (A and B)
for oral treatment.
Research question:
Severity relative to treatment of TDO ?
• 2 × 189 patients randomized, 36 centers
• 48 weeks of total follow up (12 months)
• 12 weeks of treatment (3 months)
measurements at months 0, 1, 2, 3, 6, 9, 12.
Name, department
6
Date
Example 2
The Analgesic Trial
 Single-arm trial with 530 patients recruited (491 selected
for analysis).
 Analgesic treatment for pain caused by chronic nonmalignant disease.




Treatment was to be administered for 12 months.
We will focus on Global Satisfaction Assessment (GSA).
GSA scale goes from 1=very good to 5=very bad.
GSA was rated by each subject 4 times during the trial, at
months 3, 6, 9, and 12.
Name, department
7
Date
Questions
 Evolution over time.
 Relation with baseline covariates: age, sex, duration of the pain, type
of pain, disease progression, Pain Control Assessment (PCA), . . .
 Investigation of dropout.
Observed
frequencies
Name, department
8
Date
Generalized linear Models:
Name, department
9
Date
The Bernoulli case
Name, department
10 Date
Name, department
11 Date
Name, department
12 Date
Generalized Linear Models
Name, department
13 Date
Longitudinal Generlized Linear Models
Name, department
14 Date
Generalized Linear Mixed Models
Name, department
15 Date
Name, department
16 Date
Name, department
17 Date
Empirical bayes estimates
Name, department
18 Date
Example 1 (cont’d)
Name, department
19 Date
Name, department
20 Date
Types of inference
Name, department
21 Date
22
Date
Syntax for NLMIXED
http://www.tau.ac.il/cc/pages/docs/sas8/stat/chap46/index.htm













23
PROC NLMIXED options ;
ARRAY array specification ;
BOUNDS boundary constraints ;
BY variables ;
CONTRAST 'label' expression <,expression> ;
ESTIMATE 'label' expression ;
ID expressions ;
MODEL model specification ;
PARMS parameters and starting values ;
PREDICT expression ;
RANDOM random effects specification ;
REPLICATE variable ;
Program statements ; The following sections provide a detailed description of each of
these statements.
Date













24
PROC NLMIXED Statement
ARRAY Statement
BOUNDS Statement
BY Statement
CONTRAST Statement
ESTIMATE Statement
ID Statement
MODEL Statement
PARMS Statement
PREDICT Statement
RANDOM Statement
REPLICATE Statement
Programming Statements
Example
data infection;
input clinic t x n;
datalines;
 This example analyzes the data
from Beitler and Landis (1985),
which represent results from a
multi-center clinical trial
investigating the effectiveness of
two topical cream treatments
(active drug, control) in curing an
infection. For each of eight
clinics, the number of trials and
favorable cures are recorded for
each treatment. The SAS data
set is as follows.
1 1 11 36
1 0 10 37
2 1 16 20
2 0 22 32
3 1 14 19
3 0 7 19
4 1 2 16
4 0 1 17
5 1 6 17
5 0 0 12
6 1 1 11
6 0 0 10
7115
7019
8146
8067
run;
25
Date
 Suppose nij denotes the number of trials for the ith clinic
and the jth treatment (i = 1, ... ,8 j = 0,1), and xij denotes
the corresponding number of favorable cures. Then a
reasonable model for the preceding data is the following
logistic model with random effects:
 The notation tj indicates the jth treatment, and the ui are
assumed to be iid .
26
Date
 The PROC NLMIXED statements to fit this model are as
follows:
proc nlmixed data=infection;
parms beta0=-1 beta1=1 s2u=2;
eta = beta0 + beta1*t + u;
expeta = exp(eta);
p = expeta/(1+expeta);
model x ~ binomial(n,p);
random u ~ normal(0,s2u) subject=clinic;
predict eta out=eta; estimate '1/beta1' 1/beta1; run;
Name, department
27 Date
 The PROC NLMIXED statement invokes the procedure, and the
PARMS statement defines the parameters and their starting values.
The next three statements define pij, and the MODEL statement
defines the conditional distribution of xij to be binomial. The RANDOM
statement defines U to be the random effect with subjects defined by
the CLINIC variable.
 The PREDICT statement constructs predictions for each observation in
the input data set. For this example, predictions of and approximate
standard errors of prediction are output to a SAS data set named ETA.
These predictions include empirical Bayes estimates of the random
effects ui.
 The ESTIMATE statement requests an estimate of the reciprocal of .
28
Date
Parameter Estimates
Paramet
Standar
er
Estimate d Error
DF
t Value Pr > |t|
Alpha
Lower
-2.5123
Upper Gradient
beta0
-1.1974
0.5561
7
-2.15
0.0683
0.05
beta1
0.7385
0.3004
7
2.46
0.0436
0.05 0.02806
1.4488 -2.08E-6
s2u
1.9591
1.1903
7
1.65
0.1438
0.05
-0.8554
4.7736 -2.48E-7
Estimate
Standar
d Error
DF
t Value Pr > |t|
Alpha
Lower
Upper
1.3542
0.5509
7
0.05 0.05146
2.6569
Label
1/beta1
Name, department
29 Date
2.46
0.0436
0.1175
-3.1E-7
Conclusions
 The "Parameter Estimates" table indicates marginal
significance of the two fixed-effects parameters. The
positive value of the estimate of indicates that the
treatment significantly increases the chance of a favorable
cure.
 The "Additional Estimates" table displays results from the
ESTIMATE statement. The estimate of
equals
1/0.7385 = 1.3541 and its standard error equals
0.3004/0.73852 = 0.5509 by the delta method (Billingsley
1986). Note this particular approximation produces a tstatistic identical to that for the estimate of .
30
Date
PROC NLMIXED
Name, department
31 Date
PROC NLMIXED
Name, department
32 Date
Name, department
33 Date
Name, department
34 Date
Name, department
35 Date
Name, department
36 Date
Example 2 (cont’d)
 • We analyze the data using a GLMM, but with different
approximations:
 Integrand approximation: GLIMMIX and MLWIN (PQL1 or PQL2)
 Integral approximation: NLMIXED (adaptive or not) and MIXOR
(non-adaptive)
Results
Name, department
37 Date
Name, department
38 Date
PROC MIXED vs PROC NLMIXED
 The models fit by PROC NLMIXED can be viewed as generalizations of the random
coefficient models fit by the MIXED procedure. This generalization allows the random
coefficients to enter the model nonlinearly, whereas in PROC MIXED they enter linearly.
 With PROC MIXED you can perform both maximum likelihood and restricted maximum
likelihood (REML) estimation, whereas PROC NLMIXED only implements maximum
likelihood.
 Finally, PROC MIXED assumes the data to be normally distributed, whereas PROC
NLMIXED enables you to analyze data that are normal, binomial, or Poisson or that have
any likelihood programmable with SAS statements.
 PROC NLMIXED does not implement the same estimation techniques available with the
NLINMIX and GLIMMIX macros. (generalized estimating equations). In contrast, PROC
NLMIXED directly maximizes an approximate integrated likelihood.
39
References
 Beal, S.L. and Sheiner, L.B. (1982), "Estimating Population Kinetics," CRC
Crit. Rev. Biomed. Eng., 8, 195 -222.
 Beal, S.L. and Sheiner, L.B., eds. (1992), NONMEM User's Guide, University
of California, San Francisco, NONMEM Project Group.
 Beitler, P.J. and Landis, J.R. (1985), "A Mixed-effects Model for Categorical
Data," Biometrics, 41, 991 -1000.
 Breslow, N.E. and Clayton, D.G. (1993), "Approximate Inference in
Generalized Linear Mixed Models," Journal of the American Statistical
Association, 88, 9 -25.
 Davidian, M. and Giltinan, D.M. (1995), Nonlinear Models for Repeated
Measurement Data, New York: Chapman & Hall.
 Diggle, P.J., Liang, K.Y., and Zeger, S.L. (1994), Analysis of Longitudinal Data,
Oxford: Clarendon Press.
 Engel, B. and Keen, A. (1992), "A Simple Approach for the Analysis of
Generalized Linear Mixed Models," LWA-92-6, Agricultural Mathematics Group
(GLW-DLO). Wageningen, The Netherlands.
40
Date
 Fahrmeir, L. and Tutz, G. (2002). Multivariate Statistical Modelling Based on
Generalized Linear Models, (2nd edition). Springer Series in Statistics. NewYork: Springer-Verlag.
 Ezzet, F. and Whitehead, J. (1991), "A Random Effects Model for Ordinal
Responses from a Crossover Trial," Statistics in Medicine, 10, 901 -907.
 Galecki, A.T. (1998), "NLMEM: New SAS/IML Macro for Hierarchical Nonlinear
Models," Computer Methods and Programs in Biomedicine, 55, 107 -216.
 Gallant, A.R. (1987), Nonlinear Statistical Models, New York: John Wiley &
Sons, Inc.
 Gilmour, A.R., Anderson, R.D., and Rae, A.L. (1985), "The Analysis of
Binomial Data by Generalized Linear Mixed Model," Biometrika, 72, 593 -599.
 Harville, D.A. and Mee, R.W. (1984), "A Mixed-model Procedure for Analyzing
Ordered Categorical Data," Biometrics, 40, 393 -408.
 Lindstrom, M.J. and Bates, D.M. (1990), "Nonlinear Mixed Effects Models for
Repeated Measures Data," Biometrics, 46, 673 -687.
 Littell, R.C., Milliken, G.A., Stroup, W.W., and Wolfinger, R.D. (1996), SAS
System for Mixed Models, Cary, NC: SAS Institute Inc.
Name, department
41 Date
 Longford, N.T. (1994), "Logistic Regression with Random Coefficients,"
Computational Statistics and Data Analysis, 17, 1 -15.
 McCulloch, C.E. (1994), "Maximum Likelihood Variance Components
Estimation for Binary Data," Journal of the American Statistical Association,
89, 330 -335.
 McGilchrist, C.E. (1994), "Estimation in Generalized Mixed Models," Journal of
the Royal Statistical Society B, 56, 61 -69.
 Pinheiro, J.C. and Bates, D.M. (1995), "Approximations to the Log-likelihood
Function in the Nonlinear Mixed-effects Model," Journal of Computational and
Graphical Statistics, 4, 12 -35.
 Roe, D.J. (1997) "Comparison of Population Pharmacokinetic Modeling
Methods Using Simulated Data: Results from the Population Modeling
Workgroup," Statistics in Medicine, 16, 1241 - 1262.
 Schall, R. (1991). "Estimation in Generalized Linear Models with Random
Effects," Biometrika, 78, 719 -727.
 Sheiner L. B. and Beal S. L., "Evaluation of Methods for Estimating Population
Pharmacokinetic Parameters. I. Michaelis-Menten Model: Routine Clinical
Pharmacokinetic Data," Journal of Pharmacokinetics and Biopharmaceutics, 8,
(1980) 553 -571.

42
Date
 Sheiner, L.B. and Beal, S.L. (1985), "Pharmacokinetic Parameter Estimates
from Several Least Squares Procedures: Superiority of Extended Least
Squares," Journal of Pharmacokinetics and Biopharmaceutics, 13, 185 -201.
 Stiratelli, R., Laird, N.M., and Ware, J.H. (1984), "Random Effects Models for
Serial Observations with Binary Response," Biometrics, 40, 961-971.
 Vonesh, E.F., (1992), "Nonlinear Models for the Analysis of Longitudinal Data,"
Statistics in Medicine, 11, 1929 - 1954.
 Vonesh, E.F. and Chinchilli, V.M. (1997), Linear and Nonlinear Models for the
Analysis of Repeated Measurements, New York: Marcel Dekker.
 Wolfinger R.D. (1993), "Laplace's Approximation for Nonlinear Mixed Models,"
Biometrika, 80, 791 -795.
 Wolfinger, R.D. (1997), "Comment: Experiences with the SAS Macro
NLINMIX," Statistics in Medicine, 16, 1258 -1259.
 Wolfinger, R.D. and O'Connell, M. (1993), "Generalized Linear Mixed Models:
a Pseudo-likelihood Approach," Journal of Statistical Computation and
Simulation, 48, 233 -243.
 Yuh, L., Beal, S., Davidian, M., Harrison, F., Hester, A., Kowalski, K., Vonesh,
E., Wolfinger, R. (1994), "Population Pharmacokinetic/Pharmacodynamic
Methodology and Applications: a Bibliography," Biometrics, 50, 566 -575
43
Date
End of Part I
Any Questions
Name, department
44 Date
?
Part II
Introduction to non-linear mixed
models in Pharmakokinetics
Typical data
180
180
Concentration
160
160
140
140
120
120
One curve per patient
100
100
80
80
60
60
40
40
20
20
00
00
55
10
10
15
15
20
20
25
25
30
30
35
35
40
40Time 45
45
Common situation (bio)sciences:
 A continuous response evolves over time (or other condition) within individuals




from a population of interest
Scientific interest focuses on features or mechanisms that underlie individual
time trajectories of the response and how these vary across the population.
A theoretical or empirical model for such individual profiles, typically non-linear
in the parameters that may be interpreted as representing such features or
mechanisms, is available.
Repeated measurements over time are available on each individual in a
sample drawn from the population
Inference on the scientific questions of interest is to be made in the context of
the model and its parameters
Non linear mixed effects models
Nonlinear mixed effects models: or hierarchical non-linear models

A formal statistical framework for this situation

A “hot” methodological research area in the early 1990s

Now widely accepted as a suitable approach to inference, with applications
routinely reported and commercial software available

Many recent extensions, innovations

Have many applications: growth curves, pharmacokinetics, dose-response
etc
PHARMACOKINETICS
 A drugs can administered in many different
ways: orally, by i.v. infusion, by inhalation,
using a plaster etc.
 Pharmacokinetics is the study of the rate
processes that are responsible for the time
course of the level of the drug (or any other
exogenous compound in the body such as
alcohol, toxins etc).
PHARMACOKINETICS
 Pharmacokinetics is about what happens to the drug in the
body. It involves the kinetics of drug absorption, distribution,
and elimination i.e. metabolism and excretion (adme). The
description of drug distribution and elimination is often termed
drug disposition.
 One way to model these processes is to view the body as a
system with a number of compartments through which the
drug is distributed at certain rates. This flow can be described
using constant rates in the cases of absorbtion and
elimination.
Plasma concentration curves (PCC)
 The concentration of a drug in the plasma reflects many of its
properties. A PCC gives a hint as to how the ADME processes
interact. If we draw a PCC in a logarithmic scale after an i.v. dose, we
expect to get a straight line since we assume the concentration of the
drug in plasma to decrease exponentially. This is first order- or linear
kinetics. The elimination rate is then proportional to the concentration
in plasma. This model is approximately true for most drugs.
Plasma concentration curve
Concentrati
on
Tim
e
Pharmacokinetic models
Various types of
models
One-compartment model with rapid intravenous
administration: The pharmacokinetics parameters




Half life
Distribution volume
AUC
Tmax and Cmax
i.v.
k
D, V
D
•D: Dose
•VD: Volume
•k: Elimination rate
•Cl: Clearance
One compartment model
 General model
 Tablet
Dose ka
C (t )  F
(e  k e t  e  k a t )
V k a  ke
dC
 v in  v out
dt
 IV
C(t) , V
Vin
Ve
ka
ke
dC
  kC 0
dt
D
 Cl 
Ct  exp 
t
V
 V 
Typical example in kinetics
A typical kinetics experiment is performed on a number, m, of
groups of h patients.
Individuals in different groups receive the same formulation of an
active principle, and different groups receive different formulations.
The formulations are given by IV route at time t=0.
The dose, D, is the same for all formulations.
For all formulations, the plasma concentration is measured at certain
sampling times.
Random or fixed ?
The formulation
Fixed
Dose
Fixed
The sampling times
Fixed
The concentrations
Random
Analytical error
Departure to kinetic model
The patients
Random
Population kinetics
Fixed
Classical kinetics
An example
180
One PCC per patients
Concentration
160
140
120
100
80
60
40
20
0
0
5
10
15
20
25
30
35
40
45
Time
Step 1 : Write a (PK/PD) model
A statistical model
Mean model :
functional relationship
Variance model :
Assumptions on the residuals
Step 1 : Write a deterministic (mean)
model to describe the individual kinetics
140
120
100
80
60
40
20
0
0
10
20
30
40
50
60
70
One compartment model with constant
intravenous infusion rate
D
C (t )  C0 exp  kt ; C0  ; Cl  kV
V
D
 Cl 
 C (t )  exp   t 
V
 V 
t
D
 Cl 
C  exp  t 
V
 V 
Step 1 : Write a deterministic (mean)
model to describe the individual kinetics
140
D
 Cl 
C (t )  exp   t 
V
 V 
120
100
80
60
40
20
0
0
10
20
30
40
50
60
70
Step 1 : Write a deterministic (mean)
model to describe the individual kinetics
140
120
100
residual
80
60
40
20
0
0
10
20
30
40
50
60
70
Step 1 : Write a model (variance) to
describe the magnitude of departure to
the kinetics
25
20
15
Residual
10
5
0
0
10
20
30
40
50
60
70
-5
-10
-15
-20
-25
Time
Step 1 : Write a model (variance) to
describe the magnitude of departure to
the kinetics
25
20
15
Residual
10
5
0
0
10
20
30
40
50
60
70
-5
-10
-15
-20
-25
Time
Step 1 : Describe the shape of departure
to the kinetics
Residual
0
10
20
30
40
50
60
70
Time
Step 1 :Write an "individual" model
Yi , j jth concentration measured on the ith patient
ti , j
jth sample time of the ith patient
residual
 Cli 
 Cli 
D
D
Yi , j  exp  
ti , j    exp  
ti , j  i , j
Vi
Vi
 Vi

 Vi

Gaussian residual with unit variance
Step 2 : Describe variation between
individual parameters
0
Population of patients
0.1
0.2
0.3
0.4
Distribution of clearances
Clearance
Step 2 : Our view through a sample of
patients
Sample of patients
Sample of clearances
Step 2 : Two main approaches:parametric
and semi-parametric
Sample of clearances
Semi-parametric approach
Step 2 : Two main approaches
Sample of clearances
Semi-parametric approach
(e.g. kernel estimate)
Step 2 : Semi-parametric approach
• Does require a large sample size to provide
results
• Difficult to implement
• Is implemented on “commercial” PK software
Bias?
Step 2 : Two main approaches
0
Sample of clearances
0.1
0.2
0.3
0.4
Parametric approach
Step 2 : Parametric approach
• Easier to understand
• Does not require a large sample size to provide
(good or poor) results
• Easy to implement
• Is implemented on the most popular pop PK
software
(NONMEM, S+, SAS,…)
Step 2 : Parametric approach
 Cli 
 Cli 
D
D
Yi , j  exp  
ti , j    exp  
ti , j  i , j
Vi
Vi
 Vi

 Vi

A simple model :
ln Cli  Cl   i

V
 ln Vi  V   i
Cl
ln V
ln Cl
Step 2 : Population parameters
ln V
V
 Cl,V 
V
 Cl
Cl
Mean parameters
Cl V
ln Cl
  Cl2
 Cl V   Variance parameters :

  
2
 measure inter-individual




V
 Cl V

variability
Step 2 : Parametric approach
 Cli 
 Cli 
D
D
Yi , j  exp  
ti , j    exp  
ti , j  i , j
Vi
Vi
 Vi

 Vi

A model including covariates
ln Cli   Cl  θ1 X 1i  θ2 X 2i   i

V
 ln Vi  V   i
Cl
Step 2 : A model including covariates
ln Cli   Cl  1 X 1i   2 X 2i   iCl
ln Cl
 Cl
i
 Cl  1 X 1   2 X 2
X2i
Age
X1i
BMI
Step 3 :Estimate the parameters of the
current model
Several methods with different properties
1. Naive pooled data
2. Two-stages
3. Likelihood approximations
1.
Laplacian expansion based metho
2.
Gaussian quadratures
4. Simulations methods
1. Naive pooled data : a single patient
Naïve Pooled Data combines all the data as if they came
from a single reference individual and fit into a model using
classical fitting procedures. It is simple, but can not
investigate fixed effect sources of variability, distinguish
between variability within and between individuals.
 Cl 
 Cl 
D
D
Y j  exp  
t j    exp  
t j  j
 V

 V

V
V




Concentration
180
The naïve approach does not
allow to estimate interindividual variation.
160
140
120
100
80
60
40
20
0
0
5
10
15
20
25
30
35
40
Time
45
Concentration
2. Two stages method: stage 1
Within individual variability
 Cli 
 Cli 
D
D
Yi , j  exp  
ti , j    exp  
ti , j  i , j
Vi
Vi
180
 Vi

 Vi

160
Cˆ l1 , Vˆ1
140
Cˆ l2 , Vˆ2
120
Cˆ l ,Vˆ



100
3
3



.
.
.
80
60
Cˆl ,Vˆ 
40
n
20
n
0
00
55
10
15
20
25
30
35
40
45
Time
Two stages method : stage 2
Between individual variability
ln Cˆ li   Cl   Cl
i
 ˆ
V
ln
V





i
V
i
• Does not require a specific software
• Does not use information about the distribution
• Leads to an overestimation of  which tends
to zero when the number of observations per
animal increases.
• Cannot be used with sparse data
3. The Maximum Likelihood Estimator
Let
   Cl , V ,  ,  ,  , 
   ,
i
l  y ,  
Cl
V
i
i

2
Cl
2
V
2

N
 ln  exp  h  , y , d
i
i 1
ˆ  Arg inf 
i
i
i
N
 ln  exp  h  , y , d
i 1
i
i
i
i
The Maximum Likelihood Estimator
ˆ
•Is the best estimator that can be obtained
among the consistent estimators
•It is efficient (it has the smallest variance)
•Unfortunately, l(y,) cannot be computed exactly
•Several approximations of l(y,) are used.
3.1 Laplacian expansion based
methods
First Order (FO) (Beal, Sheiner 1982) NONMEM
Linearisation about 0
 Cli 
 Cli 
D
D
Yi , j  exp  
ti , j    exp  
ti , j  i , j
Vi
Vi
 Vi

 Vi

 exp  Cl  
D

exp  
ti , j   Z1  iCl  Z 2  iV  Z 3  iViCl
exp V 
 exp V  
 exp  Cl  
D

exp  
ti , j  i , j
exp V 
 exp V  
Laplacian expansion based methods
First Order Conditional Estimation (FOCE) (Beal, Sheiner) NONMEM
Non Linear Mixed Effects models (NLME) (Pinheiro, Bates)S+, SAS
(Wolfinger)
Linearisation about the current prediction of the individual parameter
 Cli 
 Cli 
D
D
Yi , j  exp  
ti , j    exp  
ti , j  i , j
Vi
Vi
 Vi

 Vi

 Cˆ li 
D
 exp  
ti , j   Z1  ,ˆi  iCl  ˆiCl  Z 2  ,ˆi  iV  ˆiV
ˆ
Vˆi
V
i


 Cˆ li 
D
Cl
Cl
V
V
 Z 3  ,ˆi  i  ˆi i  ˆi   exp  
ti , j  i , j
ˆ
Vˆi
 Vi








Gaussian quadratures
Approximation of the integrals by discrete sums
N
l  y,    ln  exp  hi i , yi , di
i 1
N
P
i 1
k 1
 
  ln  exp  hi ik , yi ,

4. Simulations methods
Simulated Pseudo Maximum Likelihood (SPML)
Minimize
1
2
i  2 yi  i  , D  Vi1  ,D,   ln Vi  , D, 
Cl
K


exp



1
D
Cl


i , j    
exp
t
i
,
j
 exp V   V

K k 1 exp V   Cl



i ,K



Vi   simulated variance
i ,K
i ,K


Properties
Criterion
When
Advantages
Drawbacks
Naive pooled data
Never
Easy to use
Does not provide
consistent estimate
Two stages
Rich data/
initial estimates
Does not require
a specific software
Overestimation of
variance components
FO
Initial estimate
quick computation
Gives quickly a result
Does not provide
consistent estimate
FOCE/NLME
Rich data/ small
Give quickly a result.
intra individual available on specific
variance
softwares
Biased estimates when
sparse data and/or
large intra
Gaussian
quadrature
Always
consistent and
efficient estimates
provided P is large
The computation is long
when P is large
SMPL
Always
consistent estimates
The computation is long
when K is large
Model check: Graphical analysis
Predicted concentrations
ln Cli   Cl   Cli

V
ln
V





i
V
i
ln Cli   Cl  1 BWi   2 agei   Cli

V
ln
V





i
V
i
180
160
160
140
140
Variance reduction
120
120
100
100
80
80
60
60
40
40
20
20
0
0
0
20
40
60
80
100
120
140
0
20
40
Observed concentrations
60
80
100
120
140
Graphical analysis
ˆi, j
3
3
2
2
1
1
0
0
10
20
30
40
50
0
0
5
10
15
20
25
30
35
40
45
-1
-1
-2
-2
-3
-4
The PK model seems good
Time
-3
The PK model is inappropriate
Graphical analysis
ˆV
i
ˆ Cl
i
under gaussian assumption
ˆ Cl
ˆV
i
i
Normality should be questioned
Normality acceptable
add other covariates
or try semi-parametric model
The Theophylline example
 An alkaloid derived from tea or produced synthetically; it is a smooth
muscle relaxant used chiefly for its bronchodilator effect in the
treatment of chronic obstructive pulmonary emphysema, bronchial
asthma, chronic bronchitis and bronchospastic distress. It also has
myocardial stimulant, coronary vasodilator, diuretic and respiratory
center stimulant effects.
http://www.tau.ac.il/cc/pages/docs/sas8/stat/chap46/sect38.htm
References
 Davidian, M. and Giltinan, D.M. (1995). Nonlinear Models for Repeated
Measurement Data. Chapman & Hall/CRC Press.
 Davidian, M. and Giltinan, D.M. (2003). Nonlinear models for repeated
measurement data: An overview and update. Journal of Agricultural,
Biological, and Environmental Statistics 8, 387–419.
 Davidian, M. (2009). Non-linear mixed-effects models. In Longitudinal Data
Analysis, G. Fitzmaurice, M. Davidian, G. Verbeke, and G. Molenberghs
(eds). Chapman & Hall/CRC Press, ch. 5, 107–141.
 (An outstanding overview ) “Pharmacokinetics and pharmaco- dynamics ,”
by D.M. Giltinan, in Encyclopedia of Biostatistics, 2nd edition.
?
Any Questions