From:Michael.J.Fossler@gsk.com
Subject:[NMusers] Obtaining gradients for the computation of prediction intervals
Date:11:38 AM 9/17/02 -0400
Hi everyone;
Aarons et al. (Br J Clin Pharmacol (1989) 28:305)) describes a method of calculating a prediction
interval for assessing how well a model predicts a validation data set. He first computes a
standardized prediction error for each concentration (SPE) which is
SPE=Observed-Predicted/SD(Observed)
the SD(Observed) is obtained using the following expression (inserted as a graphic)
Can anyone tell us how to obtain the partial derivatives (gradients) which we need
to compute this standard deviation?
Thanks,
Mike Fossler, Clinical PK, Modeling & Simulation
Jagadeesh Aluri, PK Operations
GlaxoSmithKline
Obtaining gradients for the computation of prediction intervals
From:Leonid Gibiansky
Subject:Re: [NMusers] Obtaining gradients for the computation of prediction intervals
Date:Tue, 17 Sep 2002 13:10:16 -0400
Michael,
Unfortunately, I do not have the paper that you cited, and cannot find it on the web.
But I am not sure that you can estimate SD(observed) from anything except
observations. If you have identical dosing and sampling history, SD(observed)
at each time point can be directly computed from the set of observations.
The formula that you mentioned is based on a particular model and looks more
like a SD(predicted). Concentration at any time point can be expressed as a function of
thetas and etas + error term (using previously developed population model ).
Then the formula can be used to compute SD (partial derivatives should be computed at etas=0).
Leonid
From:"Gastonguay, Marc"
Subject:RE: [NMusers] Obtaining gradients for the computation of prediction intervals
Date:Tue, 17 Sep 2002 23:09:14 -0400
Mike,
I recall that the SD you are looking for has been discussed previously in the
context of prediction intervals, although I think the nmusers archive does not have the complete thread (see:
http://www.cognigencorp.com/nonmem/nm/99dec071998.html ). The gradients you are
looking for are stored in the NONMEM commons in the G and HH matrices for etas and epsilons,
respectively. You can access these with verbatim code or presumably by using an INFN subroutine.
A while ago, Tom Ludden and I worked-out some code for extracting the SD of the observation and using it in a 95%
prediction interval (see below). This one works when OMEGA and SIGMA are diagonal (no covariances).
A word of caution... These days, I typically create prediction intervals via simulation, as described in
the archives, so this code has not been used since NONMEM IV. Also recognize that this is my
implementation of the code, and any errors are mine, not Tom's.
Best regards,
Marc Gastonguay
; THIS IS AN NMTRAN CONTROL STREAM FOR GENERATING 95% PREDICTION INTERVALS
; OF THE TYPICAL POPULATION PREDICTION.
; USE FINAL MODEL WITH FINAL PARAMETER ESTIMATES AS INITIAL ESTIMATES AND
; MAXEVAL=0.
; USE DATA FROM A "TYPICAL" INDIVIDUAL - YOU ONLY NEED 1 INDIVIDUAL
; THE PHENOBARBITAL MODEL IS USED AS AN EXAMPLE HERE.
$PROB PHENOBARBITAL POPULATION PK MODEL
$INPUT ID TIME AMT WT APGR DV EVID MDV
$DATA pheno
$SUBROUTINES ADVAN1
$ABBREVIATED COMRES=20
; RESERVES SPACE IN THE GLOBAL COMMON NMPRD4. YOU NEED TO RESERVE AS
; MANY SPACES (OR MORE) AS YOU HAVE ERROR TERMS.
$PK
; BEGIN VERBATIM CODE
" FIRST
" COMMON/ROCM6/THETAF(20), OMEGAF(10,10), SIGMAF(10,10)
" DOUBLE PRECISION THETAF, OMEGAF, SIGMAF
"C YOU MUST DEFINE AS MANY COM() VARIABLES AS YOU HAVE ERROR TERMS.
"C THE INTER- AND INTRA-INDIVIDUAL VARIANCE TERMS ARE THE DIAGONALS OF
"C THE MATRICES OMEGAF AND SIGMAF, RESPECTIVELY.
" COM(1)=OMEGAF(1,1)
" COM(2)=OMEGAF(2,2)
"C .
"C .
"C COM(n)=OMEGAF(n,n)
"C .
" COM(3)=SIGMAF(1,1)
" COM(4)=SIGMAF(2,2)
; BACK TO ABBREVIATED CODE
TVCL=THETA(1)
CL=TVCL*DEXP(ETA(1))
TVV=THETA(2)
V=TVV*DEXP(ETA(2))
K=CL/V
S1=V
$ERROR
Y=F*(1+ERR(1))+ERR(2)
; THE FOLLOWING SECTION CALCULATES 95% PREDICTION INTERVALS OF THE PREDICTION.
" LAST
"C
"C THE DATA FOR THE 95% PREDICTION INTERVAL WILL BE STORED IN FILE=NM.PI
" OPEN(33, FILE = 'NM.PI')
" IF (NEWIND.EQ.0)THEN
" WRITE(33,3900) "ID","TIME","PRED","LOWER","UPPER"
"3900 FORMAT(A7,4(A9))
" END IF
"C
"C THIS NEXT SECTION RETRIEVES VARIANCES WHICH ARE STORED IN COM(1-4).
"C PARTIAL DERIVATIVES ARE ALSO RETRIEVED FROM G(I,J) AND HH(I,J).
"C VETA1 IS THE VARIANCE IN PRED DUE TO ETA1. FOR EXAMPLE:
"C VETA1=(PARTIAL DERIV OF PRED W.R.T. ETA1)**2*OMEGA(1), ETC...
"C
"C ADD VARIABLES BELOW AS NEEDED FOR THE APPROPRIATE # OF ETA & EPS
"C
" IF (ICALL.EQ.2)THEN
" VETA1=G(1,1)**2*COM(1)
" VETA2=G(2,1)**2*COM(2)
"C .
"C .
"C VETAn=G(n,1)**2*COM(n)
"C .
" VEPS1=HH(1,1)**2*COM(3)
" VEPS2=HH(2,1)**2*COM(4)
"C
"C THE TOTAL VARIANCE IS THE SUM OF EACH COMPONENT, SD=SQRT(TOTAL VAR).
"C YOU MUST ADD VETA OR VEPS BELOW AS NEEDED.
"C
" SDPRED=SQRT(VETA1+VETA2+VEPS1+VEPS2)
" WRITE(33,4000) ID,TIME,F,F-(1.96*SDPRED),F+(1.96*SDPRED)
"C
"4000 FORMAT(1X,F6.0,1X,F8.2,1X,3(F8.4,1X))
"
" END IF
$THETA
0.00548
1.4
$OMEGA 0.227 0.146
$SIGMA 0.0001 8
$EST MAXEVAL=0
From:Alice I Nichols [mailto:nichols@bellatlantic.net]
Subject:RE: [NMusers] Obtaining gradients for the computation of prediction intervals
Date:Friday, September 20, 2002 11:24 AM
Mike,
As a followup regarding where to find these derivatives. Here's a way I have used to reliably identify
partials related to this equation.
If you look at the fortran code generated by NMTRAN for your model (fsubs.for) you can find the exact
variable names for both partials w.r.t. etas and partials w.r.t. epsilons. There are comment lines
identifying these equations. Since there are usually several equation lines throughout the code for
each partial, I suggest you follow each equation thru all the code to correctly identify these variables. After
identifing these parameters you can revise your control stream to add these parameters to your output file.
As Marc Gastonguay mentioned however most people will do simulations and calculate SD based on this.
Alice
[Alice I Nichols, PhD
Hawthorne Research and Consulting, INC
132 Hawthorn Rd
King of Prussia, PA 19406
PH: 610-878-9112 FX: 610-878-9113]
From:"Alice I Nichols"
Subject: [NMusers] Obtaining gradients for the computation of prediction intervals
Date:Fri, 20 Sep 2002 11:34:13 -0700
Mike as a followup to my previous note the fortran subroutine to look at is for the error subroutine
written by nmtran for your control stream.
Alice.
[Alice I Nichols, PhD
Hawthorne Research and Consulting, INC
132 Hawthorn Rd
King of Prussia, PA 19406
PH: 610-878-9112 FX: 610-878-9113]
___________________________________