From: drfreedman@drfreedmaninc.com
Subject: [NMusers] AR(1) residuals
Date: 12/17/2003 6:42 PM

Colleagues
Viewing the presentation 
NONMEM Alphabet Soup: L2, AR(1), and INFN
by Marc Gastonguay posted at 
;http://www.globomaxservice.com/presentations/presentations.h
tml , I tried the following code (all drug-specific) 
information removed. The code including the "FIRST statement 
appears to be inserted into FSUBS.FOR in an executable 
section not permitting specification statements.

Can you suggest any changes? Can this approach to serial 
correlation work with the ADBVAN subroutines?

Thank-you in advance

Dr. Immanuel Freedman

;				run001.mod

$DATA input.csv IGNORE=C
;
$SUBROUTINE ADVAN2 TRANS2
;
$PK
;PK MODEL		
 F1=value;			Bioavailability 
 CL=THETA(1)*EXP(ETA(1));       Clearance (mL/kg/d)
 V=THETA(2)*EXP(ETA(2));	Central volume (mL/kg)
 KA=THETA(3)*EXP(ETA(3));       Absorption rate (/d)
 S2=V;				Scale factor
;							
$ERROR
"	FIRST
"	COMMON /NMPRD5/ C(50,29)
"    	DIMENSION T(50)
"	IF (NEWIND .NE. 2) I=0
"	IF (MDV .EQ. 0) THEN
"		I=I+1
"		T(I)=TIME
"		DO 5 J=1,I
" 5		C(J,I) = EXP(-THETA(4)*(TIME-T(J)))	
"	ENDIF
CP=A(2)/V;
Y=CP*EXP(ERR(1))+ERR(2);				
			
; 
$THETA 
;values
$OMEGA
;values
$SIGMA
;values

$ESTIMATION NOABORT MAXEVAL=9999 PRINT=10 METHOD=0  POSTHOC 
MSFO=run001.msf
$COVARIANCE
;
$TABLE NOPRINT ONEHEADER FILE=sdtab001
C ID TIME AMT DV CMT TYPE ROUT
F1 CL V  KA;
;
_______________________________________________________

From: "Gastonguay, Marc" - marc.gastonguay@snet.net
Subject: RE: [NMusers] AR(1) residuals
Date: 12/17/2003 8:01 PM

Dear Dr. Freedman,

It is not clear what problem you are running into, but this approach can be
used with the ADVAN subroutines. Below, I've listed an example control
stream using ADVAN1 and part of the corresponding data set (the phenobarb
example included in the NONMEM distribution). The FIRST statement ensures
that the verbatim code is inserted immediately after variable declarations
in the error subroutine of fsubs.for.

You'll find a discussion and example of this AR(1) residual error model
described in:
Karlsson MO, Beal SL, Sheiner LB. Three new residual error models for
population PK/PD analyses. J Pharmacokinet Biopharm 1995; 23(6):651-672.

To be honest, I rarely find this residual error model to be useful. I'd
prefer to address the serial correlation in residuals by improving the
structural or covariate model whenever possible.

Hope this helps,
Marc

(By the way, I should point out that the presentation you referred to was
the result of a collaboration between myself, Tom Ludden and Stuart Beal.)

---------Control stream----------------------
$PROB 201, PHENOBARBITAL POPULATION PK MODEL AR1 RESIDUAL
$INPUT X ID TIME AMT WT APGR DV EVID=DROP MDV=DROP
$DATA 002.CSV IGNORE=C
$SUBROUTINES ADVAN1
$PK
   TVCL=THETA(1)
   CL=TVCL*EXP(ETA(1))
   TVV=THETA(2)
   V=TVV*EXP(ETA(2))
   K=CL/V
   S1=V

$ERROR
"      FIRST
"      COMMON /NMPRD5/ C(250,29)
"      DIMENSION T(50)
"      IF (NEWIND.NE.2) I=0
"      IF (MDV.EQ.0) THEN
"         I=I+1
"         T(I)=TIME
"         DO 5 J=1,I
"  5      C(J,1)=EXP(-THETA(3)*(TIME-T(J)))
"      ENDIF
   E1=EPS(1)
   E2=EPS(2)
   IPRD=F
   IRES=DV-IPRD

   Y=F + F*E1 + E2

$THETA
   (0, 1)
   (0, 5)
   (0, 0.01)

$OMEGA
 0.16
 0.16
$SIGMA
 6
 0.04
;$MSFI 001.MSF
$EST MAXEVAL=9999 PRINT=20 NOABORT POSTHOC MSF=201.MSF
;$COV
$TABLE ID TIME IPRD CL WT APGR E1 E2 NOPRINT ONEHEADER FILE=201.TAB


---------Data (first 2 individuals of pheno data
set)---------------------------

C	ID	TIME		AMT	WT	APGR	DV	EVID	MDV
	1	0		25	1.4	7	.	1	1
	1	2		.	1.4	7	17.3	0	0
	1	12.5		3.5	1.4	7	.	1	1
	1	24.5		3.5	1.4	7	.	1	1
	1	37		3.5	1.4	7	.	1	1
	1	48		3.5	1.4	7	.	1	1
	1	60.5		3.5	1.4	7	.	1	1
	1	72.5		3.5	1.4	7	.	1	1
	1	85.3		3.5	1.4	7	.	1	1
	1	96.5		3.5	1.4	7	.	1	1
	1	108.5	3.5	1.4	7	.	1	1
	1	112.5	.	1.4	7	31	0	0
	2	0		15	1.5	9	.	1	1
	2	2		.	1.5	9	9.7	0	0
	2	4		3.8	1.5	9	.	1	1
	2	16		3.8	1.5	9	.	1	1
	2	27.8		3.8	1.5	9	.	1	1
	2	40		3.8	1.5	9	.	1	1
	2	52		3.8	1.5	9	.	1	1
	2	63.5		.	1.5	9	24.6	0	0
	2	64		3.8	1.5	9	.	1	1
	2	76		3.8	1.5	9	.	1	1
	2	88		3.8	1.5	9	.	1	1
	2	100		3.8	1.5	9	.	1	1
	2	112		3.8	1.5	9	.	1	1
	2	124		3.8	1.5	9	.	1	1
	2	135.5	.	1.5	9	33	0	0
_______________________________________________________

From: "Gastonguay, Marc" - marc.gastonguay@snet.net
Subject: RE: [NMusers] AR(1) residuals
Date: 12/19/2003 11:14 AM

Dear Dr Freedman,
A follow-up to my previous message:

I was reminded by Tom Ludden that a recent code-fix was announced on nmusers
(Subject:[NMusers] Epsilons correlated across L2 records, Date: 4Dec2003;
not yet in the archives), which directly applies to the model you've
inquired about when using conditional estimation methods. Implementation of
the example I supplied under conditional estimation methods assumes that you
have made the recommended code change as described in the Dec 4
announcement.

Marc

Webmaster Update:
The post referenced above can be found here:
http://www.cognigencorp.com/nonmem/nm/99dec042003.html
Thank you.

_______________________________________________________

From: drfreedman@drfreedmaninc.com
Subject: [NMusers] NONMEM AR(1) residual bug?
Date: 12/23/2003 8:59 PM

Colleagues

As a follow-up to my previous message on this topic, the 
residual serial correlation control stream below yields the 
following FSUBS.FOR fragment which does not compile with the 
Compaq Fortran 6.6a compiler.

Is this a bug in NONMEM?

Thank-you in advance
Immanuel Freedman, Ph.D., MIEEE
(619) 884-1347 


Control Stream:
$ERROR
"      FIRST
"      COMMON /NMPRD5/ C(250,29)
"      DIMENSION T(50)
"      IF (NEWIND.NE.2) I=0
"      IF (MDV.EQ.0) THEN
"         I=I+1
"         T(I)=TIME
"         DO 5 J=1,I
"  5      C(J,1)=EXP(-THETA(11)*(TIME-T(J)))
"      ENDIF

FSUBS.FOR fragment:

     COMMON /NMPRD5/ C(250,29)
      DIMENSION T(50)
      IF (NEWIND.NE.2) I=0
      IF (EVTREC(N,14).EQ.0) THEN
         I=I+1
         T(I)=EVTREC(N,03)
         DO 5 J=1,I
  5      EVTREC(N,01)(J,1)=EXP(-THETA(11)*(EVTREC(N,03)-T
(J)))
      ENDIF

_______________________________________________________