Hello,
We are attempting to model suppression of a biomarker, where a number of
samples (40-60%) are below the quantification limit of the assay and where 2
different assays (with different quantification limits) were used. We are
trying to model these BQL data using the M3 and M4 methods proposed by Ahn et
al (2008).
I would like to know if anyone has any comments or experience implementing the
M3 or M4 methods for biomarker data, where levels are observed at baseline, are
supressed below the LOQ for a given duration, and then return to baseline.
Also please advise if there are other methods to try and incorporate these BQL
data into the model.
I have included the relevant pieces of the control file (for both M3 and M4)
and data from a single subject.
Thanks for your review/suggestions.
Sameer
DATA:
#ID TIME AMT DV CMT EVID TYPE ASSY
1 0 0 65.71 0 0 0 1
1 0 120 0 3 1 0 1
1 168 0 10 0 0 1 1
1 336 0 10 0 0 1 1
1 336 120 0 3 1 0 1
1 504 0 12.21 0 0 0 1
1 672 120 0 3 1 0 1
1 1008 0 10 0 0 1 1
1 1008 120 0 3 1 0 1
1 1344 0 10 0 0 1 1
1 1344 120 0 3 1 0 1
1 1680 0 10 0 0 1 1
1 1680 120 0 3 1 0 1
1 2016 0 10 0 0 0 1
1 2352 0 25.64 0 0 0 1
1 2688 0 59.48 0 0 0 1
MODEL M3:
$DATA data.csv IGNORE=#
$SUB ADVAN8 TRANS1 TOL=6
$MODEL
COMP(central)
COMP(peri)
COMP(depot,DEFDOSE)
COMP(effect)
$DES
DADT(1) = KA*A(3) - K10*A(1) - K12*A(1) + K21*A(2)
DADT(2) = K12*A(1) - K21*A(2)
DADT(3) = -KA*A(3)
CONC = A(1)/V1
DADT(4) = KEO*(CONC-A(4))
$ERROR
CALLFL = 0
LOQ1=10
LOQ2=20
EFF = BL* (1 - IMAX*A(4)**HILL/ (IC50**HILL+A(4)**HILL))
IPRED=EFF
SIGA=THETA(7)
STD=SIGA
IF(TYPE.EQ.0) THEN ; GREATER THAN LOQ
F_FLAG=0
Y=IPRED+SIGA*EPS(1)
IRES =DV-IPRED
IWRES=IRES/STD
ENDIF
IF(TYPE.EQ.1.AND.ASSY.EQ.1) THEN ; BELOW LOQ1
DUM1=(LOQ1-IPRED)/STD
CUM1=PHI(DUM1)
F_FLAG=1
Y=CUM1
IRES = 0
IWRES=0
ENDIF
IF(TYPE.EQ.1.AND.ASSY.EQ.2) THEN ; BELOW LOQ2
DUM2=(LOQ2-IPRED)/STD
CUM2=PHI(DUM2)
F_FLAG=1
Y=CUM2
IRES = 0
IWRES=0
ENDIF
$SIGMA 1 FIX
$ESTIMATION MAXEVAL=9990 NOABORT SIGDIG=3 METHOD=1 INTER LAPLACIAN
POSTHOC PRINT=2 SLOW NUMERICAL
$COVARIANCE PRINT=E SLOW
MODEL M4:
$DATA data.csv IGNORE=#
$SUB ADVAN8 TRANS1 TOL=6
$MODEL
COMP(central)
COMP(peri)
COMP(depot,DEFDOSE)
COMP(effect)
$DES
DADT(1) = KA*A(3) - K10*A(1) - K12*A(1) + K21*A(2)
DADT(2) = K12*A(1) - K21*A(2)
DADT(3) = -KA*A(3)
CONC = A(1)/V1DADT(4) = KEO*(CONC-A(4))
$ERROR
CALLFL = 0
LOQ1=10
LOQ2=20
EFF = BL* (1 - IMX*A(4)**HILL/ (IC50**HILL+A(4)**HILL))
IPRED=EFF
SIGA=THETA(7)
STD=SIGA
IF(TYPE.EQ.0) THEN ; GREATER THAN LOQ
F_FLAG=0
YLO=0
Y=IPRED+SIGA*EPS(1)
IRES =DV-IPRED
IWRES=IRES/STD
ENDIF
IF(TYPE.EQ.1.AND.ASSY.EQ.1) THEN
DUM1=(LOQ1-IPRED)/STD
CUM1=PHI(DUM1)
DUM0=-IPRED/STD
CUMD0=PHI(DUM0)
CCUMD1=(CUM1-CUMD0)/(1-CUMD0)
F_FLAG=1
Y=CCUMD1
IRES = 0
IWRES=0
ENDIF
IF(TYPE.EQ.1.AND.ASSY.EQ.2) THEN
DUM2=(LOQ2-IPRED)/STD
CUM2=PHI(DUM2)
DUM0=-IPRED/STD
CUMD0=PHI(DUM0)
CCUMD2=(CUM2-CUMD0)/(1-CUMD0)
F_FLAG=1
Y=CCUMD2
IRES = 0
IWRES=0
ENDIF
$SIGMA 1 FIX
$ESTIMATION MAXEVAL=9990 NOABORT SIGDIG=3 METHOD=1 INTER LAPLACIAN
POSTHOC PRINT=2 SLOW NUMERICAL
$COVARIANCE PRINT=E SLOW
Sameer Doshi
Pharmacokinetics and Drug Metabolism, Amgen Inc.
(805) 447-6941
Modeling biomarker data below the LOQ
7 messages
5 people
Latest: Nov 20, 2009
Hi Sameer,
Several comments:
--------------------
You did not provide the entire code, but if BL is the observed baseline, it should not be included in the dataset. If you have BL=THETA(*)*EXP(ETA(*)) then the data are fine
--------------------
Additive error is assumed. I would rather use combined error (my guess is that assay STD at DV=65 is larger than STD at DV < LLOQ).
---------------------
M2 method can be implemented using YLO option (BQL observations are included with MDV=1). PRB will give you a model-based probability of
DV > YLO (see YLO EXAMPLE in help).
$ERROR
YLO = LOG1
IF(ASSY.EQ.2) YLO=LOG2
PRB = PR_Y
$EST METH=1 LAPLACIAN SLOW NOABORT
--------------------
I would increase TOL to 9 (if possible). It does not look like a stiff system, so ADVAN6 can be tried
------------------
You have not described the problem: how well these M3 - M4 methods describe your data? If you are not satisfied, could you describe the deficiencies if there are any; these can help to resolve them.
Thanks
Leonid
--------------------------------------
Leonid Gibiansky, Ph.D.
President, QuantPharm LLC
web: www.quantpharm.com
e-mail: LGibiansky at quantpharm.com
tel: (301) 767 5566
Doshi, Sameer wrote:
> Hello,
>
> We are attempting to model suppression of a biomarker, where a number of samples (40-60%) are below the quantification limit of the assay and where 2 different assays (with different quantification limits) were used. We are trying to model these BQL data using the M3 and M4 methods proposed by Ahn et al (2008). I would like to know if anyone has any comments or experience implementing the M3 or M4 methods for biomarker data, where levels are observed at baseline, are supressed below the LOQ for a given duration, and then return to baseline. Also please advise if there are other methods to try and incorporate these BQL data into the model. I have included the relevant pieces of the control file (for both M3 and M4) and data from a single subject. Thanks for your review/suggestions. Sameer DATA:
>
> #ID TIME AMT DV CMT EVID TYPE ASSY
> 1 0 0 65.71 0 0 0 1
> 1 0 120 0 3 1 0 1
> 1 168 0 10 0 0 1 1
> 1 336 0 10 0 0 1 1
> 1 336 120 0 3 1 0 1
> 1 504 0 12.21 0 0 0 1
> 1 672 120 0 3 1 0 1
> 1 1008 0 10 0 0 1 1
> 1 1008 120 0 3 1 0 1
> 1 1344 0 10 0 0 1 1
> 1 1344 120 0 3 1 0 1
> 1 1680 0 10 0 0 1 1
> 1 1680 120 0 3 1 0 1
> 1 2016 0 10 0 0 0 1
> 1 2352 0 25.64 0 0 0 1
> 1 2688 0 59.48 0 0 0 1
>
> MODEL M3:
>
> $DATA data.csv IGNORE=#
> $SUB ADVAN8 TRANS1 TOL=6
> $MODEL
> COMP(central)
> COMP(peri)
> COMP(depot,DEFDOSE)
> COMP(effect)
>
> $DES
>
> DADT(1) = KA*A(3) - K10*A(1) - K12*A(1) + K21*A(2)
> DADT(2) = K12*A(1) - K21*A(2)
> DADT(3) = -KA*A(3)
> CONC = A(1)/V1
> DADT(4) = KEO*(CONC-A(4))
>
> $ERROR
>
> CALLFL = 0
>
> LOQ1=10
>
> LOQ2=20
>
> EFF = BL* (1 - IMAX*A(4)**HILL/ (IC50**HILL+A(4)**HILL))
>
> IPRED=EFF
> SIGA=THETA(7)
> STD=SIGA
> IF(TYPE.EQ.0) THEN ; GREATER THAN LOQ
> F_FLAG=0
> Y=IPRED+SIGA*EPS(1)
> IRES =DV-IPRED
> IWRES=IRES/STD
> ENDIF
> IF(TYPE.EQ.1.AND.ASSY.EQ.1) THEN ; BELOW LOQ1
> DUM1=(LOQ1-IPRED)/STD
> CUM1=PHI(DUM1)
> F_FLAG=1
> Y=CUM1
> IRES = 0
> IWRES=0
> ENDIF
> IF(TYPE.EQ.1.AND.ASSY.EQ.2) THEN ; BELOW LOQ2
> DUM2=(LOQ2-IPRED)/STD
> CUM2=PHI(DUM2)
> F_FLAG=1
> Y=CUM2
> IRES = 0
> IWRES=0
> ENDIF
>
> $SIGMA 1 FIX $ESTIMATION MAXEVAL=9990 NOABORT SIGDIG=3 METHOD=1 INTER LAPLACIAN
>
> POSTHOC PRINT=2 SLOW NUMERICAL
> $COVARIANCE PRINT=E SLOW
>
> MODEL M4:
>
> $DATA data.csv IGNORE=#
> $SUB ADVAN8 TRANS1 TOL=6
> $MODEL
> COMP(central)
> COMP(peri)
> COMP(depot,DEFDOSE)
> COMP(effect)
>
> $DES
>
> DADT(1) = KA*A(3) - K10*A(1) - K12*A(1) + K21*A(2)
> DADT(2) = K12*A(1) - K21*A(2)
> DADT(3) = -KA*A(3)
> CONC = A(1)/V1DADT(4) = KEO*(CONC-A(4))
>
> $ERROR
>
> CALLFL = 0
>
> LOQ1=10
>
> LOQ2=20
>
> EFF = BL* (1 - IMX*A(4)**HILL/ (IC50**HILL+A(4)**HILL))
>
> IPRED=EFF
> SIGA=THETA(7)
> STD=SIGA
> IF(TYPE.EQ.0) THEN ; GREATER THAN LOQ
> F_FLAG=0
> YLO=0
> Y=IPRED+SIGA*EPS(1)
> IRES =DV-IPRED
> IWRES=IRES/STD
> ENDIF
> IF(TYPE.EQ.1.AND.ASSY.EQ.1) THEN
> DUM1=(LOQ1-IPRED)/STD
> CUM1=PHI(DUM1)
> DUM0=-IPRED/STD
> CUMD0=PHI(DUM0)
> CCUMD1=(CUM1-CUMD0)/(1-CUMD0)
> F_FLAG=1
> Y=CCUMD1
> IRES = 0
> IWRES=0
> ENDIF
> IF(TYPE.EQ.1.AND.ASSY.EQ.2) THEN
> DUM2=(LOQ2-IPRED)/STD
> CUM2=PHI(DUM2)
> DUM0=-IPRED/STD
> CUMD0=PHI(DUM0)
> CCUMD2=(CUM2-CUMD0)/(1-CUMD0)
> F_FLAG=1
> Y=CCUMD2
> IRES = 0
> IWRES=0
> ENDIF
>
> $SIGMA 1 FIX $ESTIMATION MAXEVAL=9990 NOABORT SIGDIG=3 METHOD=1 INTER LAPLACIAN
>
> POSTHOC PRINT=2 SLOW NUMERICAL
> $COVARIANCE PRINT=E SLOW
>
> Sameer Doshi
>
> Pharmacokinetics and Drug Metabolism, Amgen Inc.
> (805) 447-6941
Dear Sameer,
We've had this problem with biomarker data and published experiences in
terms of a methodological paper (below). Maybe it can give you some ideas.
Handling data below the limit of quantification in mixed effect models.
Bergstrand M, Karlsson MO.
AAPS J. 2009 Jun;11(2):371-80. Epub 2009 May 19.
Best regards,
Mats
Mats Karlsson, PhD
Professor of Pharmacometrics
Dept of Pharmaceutical Biosciences
Uppsala University
Box 591
751 24 Uppsala Sweden
phone: +46 18 4714105
fax: +46 18 471 4003
Quoted reply history
From: [email protected] [mailto:[email protected]] On
Behalf Of Doshi, Sameer
Sent: Wednesday, November 18, 2009 6:53 PM
To: nmusers
Subject: [NMusers] Modeling biomarker data below the LOQ
Hello,
We are attempting to model suppression of a biomarker, where a number of
samples (40-60%) are below the quantification limit of the assay and where 2
different assays (with different quantification limits) were used. We are
trying to model these BQL data using the M3 and M4 methods proposed by Ahn
et al (2008).
I would like to know if anyone has any comments or experience implementing
the M3 or M4 methods for biomarker data, where levels are observed at
baseline, are supressed below the LOQ for a given duration, and then return
to baseline.
Also please advise if there are other methods to try and incorporate these
BQL data into the model.
I have included the relevant pieces of the control file (for both M3 and M4)
and data from a single subject.
Thanks for your review/suggestions.
Sameer
DATA:
#ID TIME AMT DV CMT EVID TYPE ASSY
1 0 0 65.71 0 0 0 1
1 0 120 0 3 1 0 1
1 168 0 10 0 0 1 1
1 336 0 10 0 0 1 1
1 336 120 0 3 1 0 1
1 504 0 12.21 0 0 0 1
1 672 120 0 3 1 0 1
1 1008 0 10 0 0 1 1
1 1008 120 0 3 1 0 1
1 1344 0 10 0 0 1 1
1 1344 120 0 3 1 0 1
1 1680 0 10 0 0 1 1
1 1680 120 0 3 1 0 1
1 2016 0 10 0 0 0 1
1 2352 0 25.64 0 0 0 1
1 2688 0 59.48 0 0 0 1
MODEL M3:
$DATA data.csv IGNORE=#
$SUB ADVAN8 TRANS1 TOL=6
$MODEL
COMP(central)
COMP(peri)
COMP(depot,DEFDOSE)
COMP(effect)
$DES
DADT(1) = KA*A(3) - K10*A(1) - K12*A(1) + K21*A(2)
DADT(2) = K12*A(1) - K21*A(2)
DADT(3) = -KA*A(3)
CONC = A(1)/V1
DADT(4) = KEO*(CONC-A(4))
$ERROR
CALLFL = 0
LOQ1=10
LOQ2=20
EFF = BL* (1 - IMAX*A(4)**HILL/ (IC50**HILL+A(4)**HILL))
IPRED=EFF
SIGA=THETA(7)
STD=SIGA
IF(TYPE.EQ.0) THEN ; GREATER THAN LOQ
F_FLAG=0
Y=IPRED+SIGA*EPS(1)
IRES =DV-IPRED
IWRES=IRES/STD
ENDIF
IF(TYPE.EQ.1.AND.ASSY.EQ.1) THEN ; BELOW LOQ1
DUM1=(LOQ1-IPRED)/STD
CUM1=PHI(DUM1)
F_FLAG=1
Y=CUM1
IRES = 0
IWRES=0
ENDIF
IF(TYPE.EQ.1.AND.ASSY.EQ.2) THEN ; BELOW LOQ2
DUM2=(LOQ2-IPRED)/STD
CUM2=PHI(DUM2)
F_FLAG=1
Y=CUM2
IRES = 0
IWRES=0
ENDIF
$SIGMA 1 FIX
$ESTIMATION MAXEVAL=9990 NOABORT SIGDIG=3 METHOD=1 INTER LAPLACIAN
POSTHOC PRINT=2 SLOW NUMERICAL
$COVARIANCE PRINT=E SLOW
MODEL M4:
$DATA data.csv IGNORE=#
$SUB ADVAN8 TRANS1 TOL=6
$MODEL
COMP(central)
COMP(peri)
COMP(depot,DEFDOSE)
COMP(effect)
$DES
DADT(1) = KA*A(3) - K10*A(1) - K12*A(1) + K21*A(2)
DADT(2) = K12*A(1) - K21*A(2)
DADT(3) = -KA*A(3)
CONC = A(1)/V1DADT(4) = KEO*(CONC-A(4))
$ERROR
CALLFL = 0
LOQ1=10
LOQ2=20
EFF = BL* (1 - IMX*A(4)**HILL/ (IC50**HILL+A(4)**HILL))
IPRED=EFF
SIGA=THETA(7)
STD=SIGA
IF(TYPE.EQ.0) THEN ; GREATER THAN LOQ
F_FLAG=0
YLO=0
Y=IPRED+SIGA*EPS(1)
IRES =DV-IPRED
IWRES=IRES/STD
ENDIF
IF(TYPE.EQ.1.AND.ASSY.EQ.1) THEN
DUM1=(LOQ1-IPRED)/STD
CUM1=PHI(DUM1)
DUM0=-IPRED/STD
CUMD0=PHI(DUM0)
CCUMD1=(CUM1-CUMD0)/(1-CUMD0)
F_FLAG=1
Y=CCUMD1
IRES = 0
IWRES=0
ENDIF
IF(TYPE.EQ.1.AND.ASSY.EQ.2) THEN
DUM2=(LOQ2-IPRED)/STD
CUM2=PHI(DUM2)
DUM0=-IPRED/STD
CUMD0=PHI(DUM0)
CCUMD2=(CUM2-CUMD0)/(1-CUMD0)
F_FLAG=1
Y=CCUMD2
IRES = 0
IWRES=0
ENDIF
$SIGMA 1 FIX
$ESTIMATION MAXEVAL=9990 NOABORT SIGDIG=3 METHOD=1 INTER LAPLACIAN
POSTHOC PRINT=2 SLOW NUMERICAL
$COVARIANCE PRINT=E SLOW
Sameer Doshi
Pharmacokinetics and Drug Metabolism, Amgen Inc.
(805) 447-6941
Mats, Martin,
In the paper, you mentioned that "an extra additive error
model term was added for samples substituted with LOQ/2". Was it fixed
or estimated? If fixed, how? Have you tried to vary this level?
In many of your examples, LOQ/2 imputations and exclusion of BQL samples
seen to lead to bias in opposite directions; if so, it could be an
optimal value (relative to LOQ) of the fixed extra error term that
provides the least biased parameters.
Laplacian method is often not feasible for receptor (target) models
since they are strongly nonlinear (thus requiring differential
equations) and stiff. Based on the the indirect-response model
simulations considered in your paper, LOQ/2 substitution seems to
provide reasonable results if Laplacian (and thus M2-M3-M4) cannot be used.
Thanks
Leonid
--------------------------------------
Leonid Gibiansky, Ph.D.
President, QuantPharm LLC
web: www.quantpharm.com
e-mail: LGibiansky at quantpharm.com
tel: (301) 767 5566
Mats Karlsson wrote:
> Dear Sameer,
>
>
>
> Weve had this problem with biomarker data and published experiences in
> terms of a methodological paper (below). Maybe it can give you some ideas.
>
>
>
> Handling data below the limit of quantification in mixed effect models.
>
> Bergstrand M, Karlsson MO.
>
> AAPS J. 2009 Jun;11(2):371-80. Epub 2009 May 19.
>
>
>
> Best regards,
>
> Mats
>
>
>
> Mats Karlsson, PhD
>
> Professor of Pharmacometrics
>
> Dept of Pharmaceutical Biosciences
>
> Uppsala University
>
> Box 591
>
> 751 24 Uppsala Sweden
>
> phone: +46 18 4714105
>
> fax: +46 18 471 4003
>
>
>
> *From:* owner-nmusers
> [mailto:owner-nmusers
> *Sent:* Wednesday, November 18, 2009 6:53 PM
> *To:* nmusers
> *Subject:* [NMusers] Modeling biomarker data below the LOQ
>
>
>
> Hello,
>
> We are attempting to model suppression of a biomarker, where a number of
> samples (40-60%) are below the quantification limit of the assay and
> where 2 different assays (with different quantification limits) were
> used. We are trying to model these BQL data using the M3 and M4 methods
> proposed by Ahn et al (2008).
>
>
>
> I would like to know if anyone has any comments or experience
> implementing the M3 or M4 methods for biomarker data, where levels are
> observed at baseline, are supressed below the LOQ for a given duration,
> and then return to baseline.
>
>
>
> Also please advise if there are other methods to try and incorporate
> these BQL data into the model.
>
>
>
> I have included the relevant pieces of the control file (for both M3 and
> M4) and data from a single subject.
>
>
>
> Thanks for your review/suggestions.
>
>
>
> Sameer
>
>
>
> DATA:
>
> #ID TIME AMT DV CMT EVID TYPE ASSY
>
> 1 0 0 65.71 0 0 0 1
>
> 1 0 120 0 3 1 0 1
>
> 1 168 0 10 0 0 1 1
>
> 1 336 0 10 0 0 1 1
>
> 1 336 120 0 3 1 0 1
>
> 1 504 0 12.21 0 0 0 1
>
> 1 672 120 0 3 1 0 1
>
> 1 1008 0 10 0 0 1 1
>
> 1 1008 120 0 3 1 0 1
>
> 1 1344 0 10 0 0 1 1
>
> 1 1344 120 0 3 1 0 1
>
> 1 1680 0 10 0 0 1 1
>
> 1 1680 120 0 3 1 0 1
>
> 1 2016 0 10 0 0 0 1
>
> 1 2352 0 25.64 0 0 0 1
>
> 1 2688 0 59.48 0 0 0 1
>
>
>
> MODEL M3:
>
> $DATA data.csv IGNORE=#
>
> $SUB ADVAN8 TRANS1 TOL=6
>
> $MODEL
>
> COMP(central)
>
> COMP(peri)
>
> COMP(depot,DEFDOSE)
>
> COMP(effect)
>
>
>
> $DES
>
> DADT(1) = KA*A(3) - K10*A(1) - K12*A(1) + K21*A(2)
>
> DADT(2) = K12*A(1) - K21*A(2)
>
> DADT(3) = -KA*A(3)
>
> CONC = A(1)/V1
>
> DADT(4) = KEO*(CONC-A(4))
>
>
>
> $ERROR
>
> CALLFL = 0
>
>
>
> LOQ1
>
> LOQ2
>
>
>
> EFF = BL* (1 - IMAX*A(4)**HILL/ (IC50**HILL+A(4)**HILL))
>
> IPREDF
>
> SIGA=THETA(7)
>
> STD=SIGA
>
> IF(TYPE.EQ.0) THEN ; GREATER THAN LOQ
>
> F_FLAG=0
>
> Y=IPRED+SIGA*EPS(1)
>
> IRES =DV-IPRED
>
> IWRES=IRES/STD
>
> ENDIF
>
> IF(TYPE.EQ.1.AND.ASSY.EQ.1) THEN ; BELOW LOQ1
>
> DUM1=(LOQ1-IPRED)/STD
>
> CUM1=PHI(DUM1)
>
> F_FLAG=1
>
> Y=CUM1
>
> IRES = 0
>
> IWRES=0
>
> ENDIF
>
> IF(TYPE.EQ.1.AND.ASSY.EQ.2) THEN ; BELOW LOQ2
>
> DUM2=(LOQ2-IPRED)/STD
>
> CUM2=PHI(DUM2)
>
> F_FLAG=1
>
> Y=CUM2
>
> IRES = 0
>
> IWRES=0
>
> ENDIF
>
>
>
> $SIGMA 1 FIX
>
>
>
> $ESTIMATION MAXEVAL90 NOABORT SIGDIG=3 METHOD=1 INTER LAPLACIAN
>
> POSTHOC PRINT=2 SLOW NUMERICAL
>
> $COVARIANCE PRINT=E SLOW
>
>
>
> MODEL M4:
>
> $DATA data.csv IGNORE=#
>
> $SUB ADVAN8 TRANS1 TOL=6
>
> $MODEL
>
> COMP(central)
>
> COMP(peri)
>
> COMP(depot,DEFDOSE)
>
> COMP(effect)
>
>
>
> $DES
>
> DADT(1) = KA*A(3) - K10*A(1) - K12*A(1) + K21*A(2)
>
> DADT(2) = K12*A(1) - K21*A(2)
>
> DADT(3) = -KA*A(3)
>
> CONC = A(1)/V1DADT(4) = KEO*(CONC-A(4))
>
>
>
> $ERROR
>
> CALLFL = 0
>
>
>
> LOQ1
>
> LOQ2
>
>
>
> EFF = BL* (1 - IMX*A(4)**HILL/ (IC50**HILL+A(4)**HILL))
>
> IPREDF
>
> SIGA=THETA(7)
>
> STD=SIGA
>
> IF(TYPE.EQ.0) THEN ; GREATER THAN LOQ
>
> F_FLAG=0
>
> YLO=0
>
> Y=IPRED+SIGA*EPS(1)
>
> IRES =DV-IPRED
>
> IWRES=IRES/STD
>
> ENDIF
>
> IF(TYPE.EQ.1.AND.ASSY.EQ.1) THEN
>
> DUM1=(LOQ1-IPRED)/STD
>
> CUM1=PHI(DUM1)
>
> DUM0=-IPRED/STD
>
> CUMD0=PHI(DUM0)
>
> CCUMD1=(CUM1-CUMD0)/(1-CUMD0)
>
> F_FLAG=1
>
> YUMD1
>
> IRES = 0
>
> IWRES=0
>
> ENDIF
>
> IF(TYPE.EQ.1.AND.ASSY.EQ.2) THEN
>
> DUM2=(LOQ2-IPRED)/STD
>
> CUM2=PHI(DUM2)
>
> DUM0=-IPRED/STD
>
> CUMD0=PHI(DUM0)
>
> CCUMD2=(CUM2-CUMD0)/(1-CUMD0)
>
> F_FLAG=1
>
> YUMD2
>
> IRES = 0
>
> IWRES=0
>
> ENDIF
>
>
>
> $SIGMA 1 FIX
>
>
>
> $ESTIMATION MAXEVAL90 NOABORT SIGDIG=3 METHOD=1 INTER LAPLACIAN
>
> POSTHOC PRINT=2 SLOW NUMERICAL
>
> $COVARIANCE PRINT=E SLOW
>
>
>
>
>
>
>
>
>
> Sameer Doshi
>
> Pharmacokinetics and Drug Metabolism, Amgen Inc.
>
> (805) 447-6941
>
>
>
>
>
>
>
>
>
Mats, Martin,
In the paper, you mentioned that "an extra additive error
model term was added for samples substituted with LOQ/2". Was it fixed or estimated? If fixed, how? Have you tried to vary this level?
In many of your examples, LOQ/2 imputations and exclusion of BQL samples seen to lead to bias in opposite directions; if so, it could be an optimal value (relative to LOQ) of the fixed extra error term that provides the least biased parameters.
Laplacian method is often not feasible for receptor (target) models since they are strongly nonlinear (thus requiring differential equations) and stiff. Based on the the indirect-response model simulations considered in your paper, LOQ/2 substitution seems to provide reasonable results if Laplacian (and thus M2-M3-M4) cannot be used.
Thanks
Leonid
--------------------------------------
Leonid Gibiansky, Ph.D.
President, QuantPharm LLC
web: www.quantpharm.com
e-mail: LGibiansky at quantpharm.com
tel: (301) 767 5566
Mats Karlsson wrote:
> Dear Sameer,
>
> We’ve had this problem with biomarker data and published experiences in terms of a methodological paper (below). Maybe it can give you some ideas.
>
> Handling data below the limit of quantification in mixed effect models.
>
> Bergstrand M, Karlsson MO.
>
> AAPS J. 2009 Jun;11(2):371-80. Epub 2009 May 19.
>
> Best regards,
>
> Mats
>
> Mats Karlsson, PhD
>
> Professor of Pharmacometrics
>
> Dept of Pharmaceutical Biosciences
>
> Uppsala University
>
> Box 591
>
> 751 24 Uppsala Sweden
>
> phone: +46 18 4714105
>
> fax: +46 18 471 4003
>
> *From:* [email protected] [ mailto: [email protected] ] *On Behalf Of *Doshi, Sameer
>
> *Sent:* Wednesday, November 18, 2009 6:53 PM
> *To:* nmusers
> *Subject:* [NMusers] Modeling biomarker data below the LOQ
>
> Hello,
>
> We are attempting to model suppression of a biomarker, where a number of samples (40-60%) are below the quantification limit of the assay and where 2 different assays (with different quantification limits) were used. We are trying to model these BQL data using the M3 and M4 methods proposed by Ahn et al (2008). I would like to know if anyone has any comments or experience implementing the M3 or M4 methods for biomarker data, where levels are observed at baseline, are supressed below the LOQ for a given duration, and then return to baseline. Also please advise if there are other methods to try and incorporate these BQL data into the model.
>
> I have included the relevant pieces of the control file (for both M3 and M4) and data from a single subject.
>
> Thanks for your review/suggestions.
>
> Sameer
>
> DATA:
>
> #ID TIME AMT DV CMT EVID TYPE ASSY
>
> 1 0 0 65.71 0 0 0 1
>
> 1 0 120 0 3 1 0 1
>
> 1 168 0 10 0 0 1 1
>
> 1 336 0 10 0 0 1 1
>
> 1 336 120 0 3 1 0 1
>
> 1 504 0 12.21 0 0 0 1
>
> 1 672 120 0 3 1 0 1
>
> 1 1008 0 10 0 0 1 1
>
> 1 1008 120 0 3 1 0 1
>
> 1 1344 0 10 0 0 1 1
>
> 1 1344 120 0 3 1 0 1
>
> 1 1680 0 10 0 0 1 1
>
> 1 1680 120 0 3 1 0 1
>
> 1 2016 0 10 0 0 0 1
>
> 1 2352 0 25.64 0 0 0 1
>
> 1 2688 0 59.48 0 0 0 1
>
> MODEL M3:
>
> $DATA data.csv IGNORE=#
>
> $SUB ADVAN8 TRANS1 TOL=6
>
> $MODEL
>
> COMP(central)
>
> COMP(peri)
>
> COMP(depot,DEFDOSE)
>
> COMP(effect)
>
> $DES
>
> DADT(1) = KA*A(3) - K10*A(1) - K12*A(1) + K21*A(2)
>
> DADT(2) = K12*A(1) - K21*A(2)
>
> DADT(3) = -KA*A(3)
>
> CONC = A(1)/V1
>
> DADT(4) = KEO*(CONC-A(4))
>
> $ERROR
>
> CALLFL = 0
>
> LOQ1=10
>
> LOQ2=20
>
> EFF = BL* (1 - IMAX*A(4)**HILL/ (IC50**HILL+A(4)**HILL))
>
> IPRED=EFF
>
> SIGA=THETA(7)
>
> STD=SIGA
>
> IF(TYPE.EQ.0) THEN ; GREATER THAN LOQ
>
> F_FLAG=0
>
> Y=IPRED+SIGA*EPS(1)
>
> IRES =DV-IPRED
>
> IWRES=IRES/STD
>
> ENDIF
>
> IF(TYPE.EQ.1.AND.ASSY.EQ.1) THEN ; BELOW LOQ1
>
> DUM1=(LOQ1-IPRED)/STD
>
> CUM1=PHI(DUM1)
>
> F_FLAG=1
>
> Y=CUM1
>
> IRES = 0
>
> IWRES=0
>
> ENDIF
>
> IF(TYPE.EQ.1.AND.ASSY.EQ.2) THEN ; BELOW LOQ2
>
> DUM2=(LOQ2-IPRED)/STD
>
> CUM2=PHI(DUM2)
>
> F_FLAG=1
>
> Y=CUM2
>
> IRES = 0
>
> IWRES=0
>
> ENDIF
>
> $SIGMA 1 FIX
>
> $ESTIMATION MAXEVAL=9990 NOABORT SIGDIG=3 METHOD=1 INTER LAPLACIAN
>
> POSTHOC PRINT=2 SLOW NUMERICAL
>
> $COVARIANCE PRINT=E SLOW
>
> MODEL M4:
>
> $DATA data.csv IGNORE=#
>
> $SUB ADVAN8 TRANS1 TOL=6
>
> $MODEL
>
> COMP(central)
>
> COMP(peri)
>
> COMP(depot,DEFDOSE)
>
> COMP(effect)
>
> $DES
>
> DADT(1) = KA*A(3) - K10*A(1) - K12*A(1) + K21*A(2)
>
> DADT(2) = K12*A(1) - K21*A(2)
>
> DADT(3) = -KA*A(3)
>
> CONC = A(1)/V1DADT(4) = KEO*(CONC-A(4))
>
> $ERROR
>
> CALLFL = 0
>
> LOQ1=10
>
> LOQ2=20
>
> EFF = BL* (1 - IMX*A(4)**HILL/ (IC50**HILL+A(4)**HILL))
>
> IPRED=EFF
>
> SIGA=THETA(7)
>
> STD=SIGA
>
> IF(TYPE.EQ.0) THEN ; GREATER THAN LOQ
>
> F_FLAG=0
>
> YLO=0
>
> Y=IPRED+SIGA*EPS(1)
>
> IRES =DV-IPRED
>
> IWRES=IRES/STD
>
> ENDIF
>
> IF(TYPE.EQ.1.AND.ASSY.EQ.1) THEN
>
> DUM1=(LOQ1-IPRED)/STD
>
> CUM1=PHI(DUM1)
>
> DUM0=-IPRED/STD
>
> CUMD0=PHI(DUM0)
>
> CCUMD1=(CUM1-CUMD0)/(1-CUMD0)
>
> F_FLAG=1
>
> Y=CCUMD1
>
> IRES = 0
>
> IWRES=0
>
> ENDIF
>
> IF(TYPE.EQ.1.AND.ASSY.EQ.2) THEN
>
> DUM2=(LOQ2-IPRED)/STD
>
> CUM2=PHI(DUM2)
>
> DUM0=-IPRED/STD
>
> CUMD0=PHI(DUM0)
>
> CCUMD2=(CUM2-CUMD0)/(1-CUMD0)
>
> F_FLAG=1
>
> Y=CCUMD2
>
> IRES = 0
>
> IWRES=0
>
> ENDIF
>
> $SIGMA 1 FIX
>
> $ESTIMATION MAXEVAL=9990 NOABORT SIGDIG=3 METHOD=1 INTER LAPLACIAN
>
> POSTHOC PRINT=2 SLOW NUMERICAL
>
> $COVARIANCE PRINT=E SLOW
>
> Sameer Doshi
>
> Pharmacokinetics and Drug Metabolism, Amgen Inc.
>
> (805) 447-6941
Dear Leonid,
The extra error term for the BQL observations was estimated. This extra
error term was only added for Model A. In the other examples there was
already an additive part in the residual error that sufficiently relaxed the
assumption of the imputation.
I am not a big fan of imputations of BQL samples and think that they should
be avoided in favor of the M3 or M4 method as often as possible. However if
practical issues (probably disappearing as estimation methods and
computational power are improving) hinder you from using the likelihood
based methods I suggest to do a small sensitivity analysis regarding what
imputation to chose. There are no possibility to make general
recommendations on what imputations (LOQ/2 etc.) that will result in
unbiased estimates.
The first thing you should do is to evaluate the simulation properties of
the obtained model parameters. In the article that you refer to we have
described how we think that VPC are best done for datasets with censored
observations such as BQL. If a certain imputation is chosen and the final
model demonstrate good simulation properties both for the contentious
observations (above LOQ) and for the fraction of BQL samples this is a good
indication that the obtained parameter estimates are likely to be fairly
unbiased.
I have also thought of more elaborate approaches with simulation and
re-estimation with different imputations of BQL observations to evaluate the
possible biased introduced by substitution the BQL samples with an assumed
value (similar principle as the Back step method described by Kjellsson MC
et.al. (1)). However this is nothing that I have tested or think is an
attractive solution for many cases.
(1) Kjellsson MC, Jönsson S, Karlsson MO. The back-step method--method for
obtaining unbiased population parameter estimates for ordered categorical
data. AAPS J. 2004 Aug 11;6(3):e19.
By the way I don't agree with your conclusion that LOQ/2 substitution
provides reasonable results in the indirect response model example (C). It
induces a non negligible amount of bias in both the concentration effect
parameter (SLOP) and Kout (5-15%). In my opinion substitution with LOQ/2 are
just as bad as omitting the BQL data (just bad in another way) for that very
example.
Kind regards,
Martin Bergstrand, MSc, PhD student
-----------------------------------------------
Pharmacometrics Research Group,
Department of Pharmaceutical Biosciences,
Uppsala University
-----------------------------------------------
P.O. Box 591
SE-751 24 Uppsala
Sweden
-----------------------------------------------
[email protected]
-----------------------------------------------
Work: +46 18 471 4639
Mobile: +46 709 994 396
Fax: +46 18 471 4003
Quoted reply history
-----Original Message-----
From: Leonid Gibiansky [mailto:[email protected]]
Sent: den 20 november 2009 04:01
To: Mats Karlsson
Cc: 'Doshi, Sameer'; 'nmusers'; Martin Bergstrand
Subject: Re: [NMusers] Modeling biomarker data below the LOQ
Mats, Martin,
In the paper, you mentioned that "an extra additive error
model term was added for samples substituted with LOQ/2". Was it fixed
or estimated? If fixed, how? Have you tried to vary this level?
In many of your examples, LOQ/2 imputations and exclusion of BQL samples
seen to lead to bias in opposite directions; if so, it could be an
optimal value (relative to LOQ) of the fixed extra error term that
provides the least biased parameters.
Laplacian method is often not feasible for receptor (target) models
since they are strongly nonlinear (thus requiring differential
equations) and stiff. Based on the the indirect-response model
simulations considered in your paper, LOQ/2 substitution seems to
provide reasonable results if Laplacian (and thus M2-M3-M4) cannot be used.
Thanks
Leonid
--------------------------------------
Leonid Gibiansky, Ph.D.
President, QuantPharm LLC
web: www.quantpharm.com
e-mail: LGibiansky at quantpharm.com
tel: (301) 767 5566
Mats Karlsson wrote:
> Dear Sameer,
>
>
>
> Weve had this problem with biomarker data and published experiences in
> terms of a methodological paper (below). Maybe it can give you some ideas.
>
>
>
> Handling data below the limit of quantification in mixed effect models.
>
> Bergstrand M, Karlsson MO.
>
> AAPS J. 2009 Jun;11(2):371-80. Epub 2009 May 19.
>
>
>
> Best regards,
>
> Mats
>
>
>
> Mats Karlsson, PhD
>
> Professor of Pharmacometrics
>
> Dept of Pharmaceutical Biosciences
>
> Uppsala University
>
> Box 591
>
> 751 24 Uppsala Sweden
>
> phone: +46 18 4714105
>
> fax: +46 18 471 4003
>
>
>
> *From:* [email protected]
> [mailto:[email protected]] *On Behalf Of *Doshi, Sameer
> *Sent:* Wednesday, November 18, 2009 6:53 PM
> *To:* nmusers
> *Subject:* [NMusers] Modeling biomarker data below the LOQ
>
>
>
> Hello,
>
> We are attempting to model suppression of a biomarker, where a number of
> samples (40-60%) are below the quantification limit of the assay and
> where 2 different assays (with different quantification limits) were
> used. We are trying to model these BQL data using the M3 and M4 methods
> proposed by Ahn et al (2008).
>
>
>
> I would like to know if anyone has any comments or experience
> implementing the M3 or M4 methods for biomarker data, where levels are
> observed at baseline, are supressed below the LOQ for a given duration,
> and then return to baseline.
>
>
>
> Also please advise if there are other methods to try and incorporate
> these BQL data into the model.
>
>
>
> I have included the relevant pieces of the control file (for both M3 and
> M4) and data from a single subject.
>
>
>
> Thanks for your review/suggestions.
>
>
>
> Sameer
>
>
>
> DATA:
>
> #ID TIME AMT DV CMT EVID TYPE ASSY
>
> 1 0 0 65.71 0 0 0 1
>
> 1 0 120 0 3 1 0 1
>
> 1 168 0 10 0 0 1 1
>
> 1 336 0 10 0 0 1 1
>
> 1 336 120 0 3 1 0 1
>
> 1 504 0 12.21 0 0 0 1
>
> 1 672 120 0 3 1 0 1
>
> 1 1008 0 10 0 0 1 1
>
> 1 1008 120 0 3 1 0 1
>
> 1 1344 0 10 0 0 1 1
>
> 1 1344 120 0 3 1 0 1
>
> 1 1680 0 10 0 0 1 1
>
> 1 1680 120 0 3 1 0 1
>
> 1 2016 0 10 0 0 0 1
>
> 1 2352 0 25.64 0 0 0 1
>
> 1 2688 0 59.48 0 0 0 1
>
>
>
> MODEL M3:
>
> $DATA data.csv IGNORE=#
>
> $SUB ADVAN8 TRANS1 TOL=6
>
> $MODEL
>
> COMP(central)
>
> COMP(peri)
>
> COMP(depot,DEFDOSE)
>
> COMP(effect)
>
>
>
> $DES
>
> DADT(1) = KA*A(3) - K10*A(1) - K12*A(1) + K21*A(2)
>
> DADT(2) = K12*A(1) - K21*A(2)
>
> DADT(3) = -KA*A(3)
>
> CONC = A(1)/V1
>
> DADT(4) = KEO*(CONC-A(4))
>
>
>
> $ERROR
>
> CALLFL = 0
>
>
>
> LOQ1=10
>
> LOQ2=20
>
>
>
> EFF = BL* (1 - IMAX*A(4)**HILL/ (IC50**HILL+A(4)**HILL))
>
> IPRED=EFF
>
> SIGA=THETA(7)
>
> STD=SIGA
>
> IF(TYPE.EQ.0) THEN ; GREATER THAN LOQ
>
> F_FLAG=0
>
> Y=IPRED+SIGA*EPS(1)
>
> IRES =DV-IPRED
>
> IWRES=IRES/STD
>
> ENDIF
>
> IF(TYPE.EQ.1.AND.ASSY.EQ.1) THEN ; BELOW LOQ1
>
> DUM1=(LOQ1-IPRED)/STD
>
> CUM1=PHI(DUM1)
>
> F_FLAG=1
>
> Y=CUM1
>
> IRES = 0
>
> IWRES=0
>
> ENDIF
>
> IF(TYPE.EQ.1.AND.ASSY.EQ.2) THEN ; BELOW LOQ2
>
> DUM2=(LOQ2-IPRED)/STD
>
> CUM2=PHI(DUM2)
>
> F_FLAG=1
>
> Y=CUM2
>
> IRES = 0
>
> IWRES=0
>
> ENDIF
>
>
>
> $SIGMA 1 FIX
>
>
>
> $ESTIMATION MAXEVAL=9990 NOABORT SIGDIG=3 METHOD=1 INTER LAPLACIAN
>
> POSTHOC PRINT=2 SLOW NUMERICAL
>
> $COVARIANCE PRINT=E SLOW
>
>
>
> MODEL M4:
>
> $DATA data.csv IGNORE=#
>
> $SUB ADVAN8 TRANS1 TOL=6
>
> $MODEL
>
> COMP(central)
>
> COMP(peri)
>
> COMP(depot,DEFDOSE)
>
> COMP(effect)
>
>
>
> $DES
>
> DADT(1) = KA*A(3) - K10*A(1) - K12*A(1) + K21*A(2)
>
> DADT(2) = K12*A(1) - K21*A(2)
>
> DADT(3) = -KA*A(3)
>
> CONC = A(1)/V1DADT(4) = KEO*(CONC-A(4))
>
>
>
> $ERROR
>
> CALLFL = 0
>
>
>
> LOQ1=10
>
> LOQ2=20
>
>
>
> EFF = BL* (1 - IMX*A(4)**HILL/ (IC50**HILL+A(4)**HILL))
>
> IPRED=EFF
>
> SIGA=THETA(7)
>
> STD=SIGA
>
> IF(TYPE.EQ.0) THEN ; GREATER THAN LOQ
>
> F_FLAG=0
>
> YLO=0
>
> Y=IPRED+SIGA*EPS(1)
>
> IRES =DV-IPRED
>
> IWRES=IRES/STD
>
> ENDIF
>
> IF(TYPE.EQ.1.AND.ASSY.EQ.1) THEN
>
> DUM1=(LOQ1-IPRED)/STD
>
> CUM1=PHI(DUM1)
>
> DUM0=-IPRED/STD
>
> CUMD0=PHI(DUM0)
>
> CCUMD1=(CUM1-CUMD0)/(1-CUMD0)
>
> F_FLAG=1
>
> Y=CCUMD1
>
> IRES = 0
>
> IWRES=0
>
> ENDIF
>
> IF(TYPE.EQ.1.AND.ASSY.EQ.2) THEN
>
> DUM2=(LOQ2-IPRED)/STD
>
> CUM2=PHI(DUM2)
>
> DUM0=-IPRED/STD
>
> CUMD0=PHI(DUM0)
>
> CCUMD2=(CUM2-CUMD0)/(1-CUMD0)
>
> F_FLAG=1
>
> Y=CCUMD2
>
> IRES = 0
>
> IWRES=0
>
> ENDIF
>
>
>
> $SIGMA 1 FIX
>
>
>
> $ESTIMATION MAXEVAL=9990 NOABORT SIGDIG=3 METHOD=1 INTER LAPLACIAN
>
> POSTHOC PRINT=2 SLOW NUMERICAL
>
> $COVARIANCE PRINT=E SLOW
>
>
>
>
>
>
>
>
>
> Sameer Doshi
>
> Pharmacokinetics and Drug Metabolism, Amgen Inc.
>
> (805) 447-6941
>
>
>
>
>
>
>
>
>
Dear Leonid,
Dear Martin,
As Bob, David, Alan, Marc, France and others have implemented the MC-PEM,
EM, and SAEM algorithms in S-ADAPT, ADAPT V, Monolix, and NONMEM 7,
I think we came a lot closer to the "bright future" of not having
to worry as much about computation times. This progress during the
last years is just incredible.
The Laplacian method with interaction to apply the Beal M2, M3, or M4 methods
takes a lot of computation time and tends to cause instability during
estimation.
I worry more about the instability than about time and think that the
Laplacian method is also not as well parallelizable as the MC-PEM algorithm.
I did not realize any run-time difference with or without the Beal M3 option
using MC-PEM in S-ADAPT and guess this is also true for MC-PEM in NONMEM 7,
SAEM in Monolix, and MCMC in WinBugs (Steve/Joy, please correct me, if I am
wrong).
Whenever run-time is an issue and you need to account for LOQs, you might
consider using MC-PEM in parallelized S-ADAPT, SAEM in Monolix, EM in ADAPT V,
or MC-PEM in (soon to be parallelized) NONMEM 7.
Martin, I fully agree with you that imputation methods (such as LOQ/2) are
not the way to go given the availability of the Beal M3/M4 methods. While
imputation may be working in PK, the bias in PD can be substantial and
clinically relevant, at least in infectious diseases. Imputation methods
are probably worse in PD than in PK, since we often do not know our PD
systems as well as PK.
Best wishes
Juergen
Quoted reply history
-----Original Message-----
From: [email protected] [mailto:[email protected]] On
Behalf Of Leonid Gibiansky
Sent: Friday, November 20, 2009 10:21 AM
To: Martin Bergstrand
Cc: 'Mats Karlsson'; 'Doshi, Sameer'; 'nmusers'
Subject: Re: [NMusers] Modeling biomarker data below the LOQ
Dear Martin,
I wish I live in that bright future where "practical issues disappear as
estimation methods and computational power improved"; at this time, many
of my FOCEI models run for hours, even days, and this would make it
difficult to use Laplacian on each and every model.
If you look at the end of this discussion:
http://cognigencorp.com/nonmem/nm/99oct041999.html
you will find Lewis Sheiner description of LOQ/2 imputation
1. Delete all but the first in each continuous series of BQL observations
2. Set the remaining (first) one DV = QL/2
3. Use an additive plus proportional error model with the SD of the
additive part >= QL/2.
Step 3 (additive error with SD > QL/2) is really important. If you
estimated SD for BQL observations without this restriction, this may
invalidate your conclusions concerning this method. In my experience,
this restriction makes a large difference in the results.
The idea is to penalize the model for predictions above LOQ, but not
below LOQ. Ideally, we would like to have a penalty function that is
zero if IPRED < LOQ, and non-zero otherwise. Since we are restricted by
normal residual error, BQL value is imputed as LOQ/2, and SD of the
associated error is fixed at (LOQ/2)^2. Thus, the penalty for
abs(IPRED-LOQ/2)< SD observations is small but increases as IPRED > LOQ.
Note that ignoring observations is equivalent to assigning them a very
large residual error. If you look on the plots in your paper, you will
see that the direction of the bias for the method with ignored
observations (infinite residual error associated with these
observations) is the opposite relative to the bias for LOQ/2 imputation
(presumably, with small associated additive error). It could be that if
you fix an optimal value of the residual error, the bias would get
smaller or disappear completely.
As to the value of the bias, 10% bias is rarely important in real-life
examples: you rarely see parameters estimated with less than 10%
(RSE=0.1) precision (confidence intervals on the parameter estimates are
approximately VALUE*(1 +/- 2*RSE).
With that, I completely agree that if you can use the best method (M3)
you should do it.
Thanks
Leonid
--------------------------------------
Leonid Gibiansky, Ph.D.
President, QuantPharm LLC
web: www.quantpharm.com
e-mail: LGibiansky at quantpharm.com
tel: (301) 767 5566
Martin Bergstrand wrote:
> Dear Leonid,
>
> The extra error term for the BQL observations was estimated. This extra
> error term was only added for Model A. In the other examples there was
> already an additive part in the residual error that sufficiently relaxed the
> assumption of the imputation.
>
> I am not a big fan of imputations of BQL samples and think that they should
> be avoided in favor of the M3 or M4 method as often as possible. However if
> practical issues (probably disappearing as estimation methods and
> computational power are improving) hinder you from using the likelihood
> based methods I suggest to do a small sensitivity analysis regarding what
> imputation to chose. There are no possibility to make general
> recommendations on what imputations (LOQ/2 etc.) that will result in
> unbiased estimates.
>
> The first thing you should do is to evaluate the simulation properties of
> the obtained model parameters. In the article that you refer to we have
> described how we think that VPC are best done for datasets with censored
> observations such as BQL. If a certain imputation is chosen and the final
> model demonstrate good simulation properties both for the contentious
> observations (above LOQ) and for the fraction of BQL samples this is a good
> indication that the obtained parameter estimates are likely to be fairly
> unbiased.
>
> I have also thought of more elaborate approaches with simulation and
> re-estimation with different imputations of BQL observations to evaluate the
> possible biased introduced by substitution the BQL samples with an assumed
> value (similar principle as the Back step method described by Kjellsson MC
> et.al. (1)). However this is nothing that I have tested or think is an
> attractive solution for many cases.
>
> (1) Kjellsson MC, Jönsson S, Karlsson MO. The back-step method--method for
> obtaining unbiased population parameter estimates for ordered categorical
> data. AAPS J. 2004 Aug 11;6(3):e19.
>
> By the way I don't agree with your conclusion that LOQ/2 substitution
> provides reasonable results in the indirect response model example (C). It
> induces a non negligible amount of bias in both the concentration effect
> parameter (SLOP) and Kout (5-15%). In my opinion substitution with LOQ/2 are
> just as bad as omitting the BQL data (just bad in another way) for that very
> example.
>
> Kind regards,
>
> Martin Bergstrand, MSc, PhD student
> -----------------------------------------------
> Pharmacometrics Research Group,
> Department of Pharmaceutical Biosciences,
> Uppsala University
> -----------------------------------------------
> P.O. Box 591
> SE-751 24 Uppsala
> Sweden
> -----------------------------------------------
> [email protected]
> -----------------------------------------------
> Work: +46 18 471 4639
> Mobile: +46 709 994 396
> Fax: +46 18 471 4003
>
>
> -----Original Message-----
> From: Leonid Gibiansky [mailto:[email protected]]
> Sent: den 20 november 2009 04:01
> To: Mats Karlsson
> Cc: 'Doshi, Sameer'; 'nmusers'; Martin Bergstrand
> Subject: Re: [NMusers] Modeling biomarker data below the LOQ
>
> Mats, Martin,
> In the paper, you mentioned that "an extra additive error
> model term was added for samples substituted with LOQ/2". Was it fixed
> or estimated? If fixed, how? Have you tried to vary this level?
>
> In many of your examples, LOQ/2 imputations and exclusion of BQL samples
> seen to lead to bias in opposite directions; if so, it could be an
> optimal value (relative to LOQ) of the fixed extra error term that
> provides the least biased parameters.
>
> Laplacian method is often not feasible for receptor (target) models
> since they are strongly nonlinear (thus requiring differential
> equations) and stiff. Based on the the indirect-response model
> simulations considered in your paper, LOQ/2 substitution seems to
> provide reasonable results if Laplacian (and thus M2-M3-M4) cannot be used.
>
> Thanks
> Leonid
>
> --------------------------------------
> Leonid Gibiansky, Ph.D.
> President, QuantPharm LLC
> web: www.quantpharm.com
> e-mail: LGibiansky at quantpharm.com
> tel: (301) 767 5566
>
>
>
>
> Mats Karlsson wrote:
>> Dear Sameer,
>>
>>
>>
>> We've had this problem with biomarker data and published experiences in
>> terms of a methodological paper (below). Maybe it can give you some ideas.
>>
>>
>>
>> Handling data below the limit of quantification in mixed effect models.
>>
>> Bergstrand M, Karlsson MO.
>>
>> AAPS J. 2009 Jun;11(2):371-80. Epub 2009 May 19.
>>
>>
>>
>> Best regards,
>>
>> Mats
>>
>>
>>
>> Mats Karlsson, PhD
>>
>> Professor of Pharmacometrics
>>
>> Dept of Pharmaceutical Biosciences
>>
>> Uppsala University
>>
>> Box 591
>>
>> 751 24 Uppsala Sweden
>>
>> phone: +46 18 4714105
>>
>> fax: +46 18 471 4003
>>
>>
>>
>> *From:* [email protected]
>> [mailto:[email protected]] *On Behalf Of *Doshi, Sameer
>> *Sent:* Wednesday, November 18, 2009 6:53 PM
>> *To:* nmusers
>> *Subject:* [NMusers] Modeling biomarker data below the LOQ
>>
>>
>>
>> Hello,
>>
>> We are attempting to model suppression of a biomarker, where a number of
>> samples (40-60%) are below the quantification limit of the assay and
>> where 2 different assays (with different quantification limits) were
>> used. We are trying to model these BQL data using the M3 and M4 methods
>> proposed by Ahn et al (2008).
>>
>>
>>
>> I would like to know if anyone has any comments or experience
>> implementing the M3 or M4 methods for biomarker data, where levels are
>> observed at baseline, are supressed below the LOQ for a given duration,
>> and then return to baseline.
>>
>>
>>
>> Also please advise if there are other methods to try and incorporate
>> these BQL data into the model.
>>
>>
>>
>> I have included the relevant pieces of the control file (for both M3 and
>> M4) and data from a single subject.
>>
>>
>>
>> Thanks for your review/suggestions.
>>
>>
>>
>> Sameer
>>
>>
>>
>> DATA:
>>
>> #ID TIME AMT DV CMT EVID TYPE ASSY
>>
>> 1 0 0 65.71 0 0 0 1
>>
>> 1 0 120 0 3 1 0 1
>>
>> 1 168 0 10 0 0 1 1
>>
>> 1 336 0 10 0 0 1 1
>>
>> 1 336 120 0 3 1 0 1
>>
>> 1 504 0 12.21 0 0 0 1
>>
>> 1 672 120 0 3 1 0 1
>>
>> 1 1008 0 10 0 0 1 1
>>
>> 1 1008 120 0 3 1 0 1
>>
>> 1 1344 0 10 0 0 1 1
>>
>> 1 1344 120 0 3 1 0 1
>>
>> 1 1680 0 10 0 0 1 1
>>
>> 1 1680 120 0 3 1 0 1
>>
>> 1 2016 0 10 0 0 0 1
>>
>> 1 2352 0 25.64 0 0 0 1
>>
>> 1 2688 0 59.48 0 0 0 1
>>
>>
>>
>> MODEL M3:
>>
>> $DATA data.csv IGNORE=#
>>
>> $SUB ADVAN8 TRANS1 TOL=6
>>
>> $MODEL
>>
>> COMP(central)
>>
>> COMP(peri)
>>
>> COMP(depot,DEFDOSE)
>>
>> COMP(effect)
>>
>>
>>
>> $DES
>>
>> DADT(1) = KA*A(3) - K10*A(1) - K12*A(1) + K21*A(2)
>>
>> DADT(2) = K12*A(1) - K21*A(2)
>>
>> DADT(3) = -KA*A(3)
>>
>> CONC = A(1)/V1
>>
>> DADT(4) = KEO*(CONC-A(4))
>>
>>
>>
>> $ERROR
>>
>> CALLFL = 0
>>
>>
>>
>> LOQ1=10
>>
>> LOQ2=20
>>
>>
>>
>> EFF = BL* (1 - IMAX*A(4)**HILL/ (IC50**HILL+A(4)**HILL))
>>
>> IPRED=EFF
>>
>> SIGA=THETA(7)
>>
>> STD=SIGA
>>
>> IF(TYPE.EQ.0) THEN ; GREATER THAN LOQ
>>
>> F_FLAG=0
>>
>> Y=IPRED+SIGA*EPS(1)
>>
>> IRES =DV-IPRED
>>
>> IWRES=IRES/STD
>>
>> ENDIF
>>
>> IF(TYPE.EQ.1.AND.ASSY.EQ.1) THEN ; BELOW LOQ1
>>
>> DUM1=(LOQ1-IPRED)/STD
>>
>> CUM1=PHI(DUM1)
>>
>> F_FLAG=1
>>
>> Y=CUM1
>>
>> IRES = 0
>>
>> IWRES=0
>>
>> ENDIF
>>
>> IF(TYPE.EQ.1.AND.ASSY.EQ.2) THEN ; BELOW LOQ2
>>
>> DUM2=(LOQ2-IPRED)/STD
>>
>> CUM2=PHI(DUM2)
>>
>> F_FLAG=1
>>
>> Y=CUM2
>>
>> IRES = 0
>>
>> IWRES=0
>>
>> ENDIF
>>
>>
>>
>> $SIGMA 1 FIX
>>
>>
>>
>> $ESTIMATION MAXEVAL=9990 NOABORT SIGDIG=3 METHOD=1 INTER LAPLACIAN
>>
>> POSTHOC PRINT=2 SLOW NUMERICAL
>>
>> $COVARIANCE PRINT=E SLOW
>>
>>
>>
>> MODEL M4:
>>
>> $DATA data.csv IGNORE=#
>>
>> $SUB ADVAN8 TRANS1 TOL=6
>>
>> $MODEL
>>
>> COMP(central)
>>
>> COMP(peri)
>>
>> COMP(depot,DEFDOSE)
>>
>> COMP(effect)
>>
>>
>>
>> $DES
>>
>> DADT(1) = KA*A(3) - K10*A(1) - K12*A(1) + K21*A(2)
>>
>> DADT(2) = K12*A(1) - K21*A(2)
>>
>> DADT(3) = -KA*A(3)
>>
>> CONC = A(1)/V1DADT(4) = KEO*(CONC-A(4))
>>
>>
>>
>> $ERROR
>>
>> CALLFL = 0
>>
>>
>>
>> LOQ1=10
>>
>> LOQ2=20
>>
>>
>>
>> EFF = BL* (1 - IMX*A(4)**HILL/ (IC50**HILL+A(4)**HILL))
>>
>> IPRED=EFF
>>
>> SIGA=THETA(7)
>>
>> STD=SIGA
>>
>> IF(TYPE.EQ.0) THEN ; GREATER THAN LOQ
>>
>> F_FLAG=0
>>
>> YLO=0
>>
>> Y=IPRED+SIGA*EPS(1)
>>
>> IRES =DV-IPRED
>>
>> IWRES=IRES/STD
>>
>> ENDIF
>>
>> IF(TYPE.EQ.1.AND.ASSY.EQ.1) THEN
>>
>> DUM1=(LOQ1-IPRED)/STD
>>
>> CUM1=PHI(DUM1)
>>
>> DUM0=-IPRED/STD
>>
>> CUMD0=PHI(DUM0)
>>
>> CCUMD1=(CUM1-CUMD0)/(1-CUMD0)
>>
>> F_FLAG=1
>>
>> Y=CCUMD1
>>
>> IRES = 0
>>
>> IWRES=0
>>
>> ENDIF
>>
>> IF(TYPE.EQ.1.AND.ASSY.EQ.2) THEN
>>
>> DUM2=(LOQ2-IPRED)/STD
>>
>> CUM2=PHI(DUM2)
>>
>> DUM0=-IPRED/STD
>>
>> CUMD0=PHI(DUM0)
>>
>> CCUMD2=(CUM2-CUMD0)/(1-CUMD0)
>>
>> F_FLAG=1
>>
>> Y=CCUMD2
>>
>> IRES = 0
>>
>> IWRES=0
>>
>> ENDIF
>>
>>
>>
>> $SIGMA 1 FIX
>>
>>
>>
>> $ESTIMATION MAXEVAL=9990 NOABORT SIGDIG=3 METHOD=1 INTER LAPLACIAN
>>
>> POSTHOC PRINT=2 SLOW NUMERICAL
>>
>> $COVARIANCE PRINT=E SLOW
>>
>>
>>
>>
>>
>>
>>
>>
>>
>> Sameer Doshi
>>
>> Pharmacokinetics and Drug Metabolism, Amgen Inc.
>>
>> (805) 447-6941
>>
>>
>>
>>
>>
>>
>>
>>
>>
>
>