From: "Iyer, Ganesh R [PRI]" <GIyer1@prius.jnj.com>
Subject: RE: Help please...
Date: Mon, 5 Feb 2001 11:53:13 -0500

Hi,

I was wondering if there is any standard format of extracting nonmem output files into S-plus through the command line directly rather than using xpose.

Many Thanks in Advance,
Ganesh


*****


From: R.Port@DKFZ-Heidelberg.de
Subject: Reading NONMEM output into S-Plus
Date: Mon, 5 Feb 2001 19:37:55 +0100

Hi Ganesh:

Re: Reading NONMEM output into S-Plus

I suggest you collecting the desired information in the NONMEM output *table*, then reading it into S-Plus using the read.table command. E.g. if the name of the NONMEM output table is 'NONMEM.output.table' and its first few lines look like this:

TABLE NO. 1 DATA RECS. 1 THROUGH 900
ID PS PAT MALE AGE0 BW
ADM TIME DRCH WDR ETA1 ETA2 ETA3
beta Emax ED50 tau mep DV PRED
RES WRES
1.0000E+00 1.0000E+00 1.0000E+00 1.0000E+00 1.0000E+01 2.5900E+01
1.0000E+00 0.0000E+00 2.0000E+00 7.5000E+03 -5.5092E-01 -2.0011E-01
5.0293E-01 6.6288E-06 0.0000E+00 1.0000E+00 8.4320E+01 6.9725E+00
7.4000E+00 7.4000E+00 0.0000E+00 3.6832E-01
1.0000E+00 1.0000E+00 1.0000E+00 1.0000E+00 1.0000E+01 2.6100E+01
1.0000E+00 2.0000E+00 0.0000E+00 7.5000E+03 -5.5092E-01 -2.0011E-01
5.0293E-01 6.6288E-06 0.0000E+00 1.0000E+00 8.4320E+01 7.0719E+00
8.6000E+00 7.5725E+00 1.0275E+00 1.6684E+00
1.0000E+00 1.0000E+00 1.0000E+00 1.0000E+00 1.0000E+01 2.6000E+01
1.0000E+00 4.0000E+00 0.0000E+00 7.5000E+03 -5.5092E-01 -2.0011E-01
5.0293E-01 6.6288E-06 0.0000E+00 1.0000E+00 8.4320E+01 7.1714E+00
6.5000E+00 7.7450E+00 -1.2450E+00 -9.1380E-01

then, in S-Plus, you would use:

> NONMEM.output.table_read.table('NONMEM.output.table',skip=1,header=T)

Don't forget to use the ONEHEADER option in $TABLE .

Good luck! Ruedi
----------------------------------------------------------------------------
Dr. R.E. Port, German Cancer Research Center, D-0200
P.O. Box 10 19 49, D-69009 Heidelberg, Germany
phone: +49-6221 42-3385
-3347
fax: -3382
e-mail: r.port@dkfz.de


*****


From: "Piotrovskij, Vladimir [JanBe]" <VPIOTROV@janbe.jnj.com>
Subject: RE: Help please...
Date: Tue, 6 Feb 2001 13:45:46 +0100

Ganesh,
I don't think the standard NONMEM output is a suitable sourse file for S-plus. The better way is to save final estimates of THETA, OMEGA and SIGMA is text files that can then be fed to S-plus (scan() function is most suitable for importing in this case). Some time ago I asked Alison how to save these estimates. I give below my question and her answer.

Best regards,
Vladimir
----------------------------------------------------------------------
Vladimir Piotrovsky, Ph.D.
Janssen Research Foundation
Human Pharmacokinetics (ext. 5463)
B-2340 Beerse
Belgium
Email: vpiotrov@janbe.jnj.com
------------------------------------------------------------------
> If one does simulation/fitting and wants to perform, say, 100 simulations,
> he certainly would like to have parameter estimates in one of the tables
> besides a common output file. Is there any way to put THETA, OMEGA and SIGMA
> estimates in a table. Or to get a reduced version of the output file just
> with parameter estimates to facilitate further tabulating/analysis.

Here is an answer that I sent to another user with a very similar question.


Hope this helps. Alison
===============

There are 2 ways to go about this.

1) The final parameter estimates, standard errors, and minimum value of the objective function can easily be extracted from a NONMEM output report using the little program nmsee that I have mentioned in email to this group before. Usage on Unix for standard errors:

nmsee prefix out.sim | grep of

Nmsee can be obtained at the NONMEM repository.
Its in the directory NMSEE.DIR, under NONMEM.DIR.
The full address is /public/nonmem.dir/nmsee.dir
This is a reminder on how to access the repository:
ftp 204.161.113.34
ftp pkpd.icon.palo-alto.med.va.gov
World Wide Web at URL:
http://pkpd.icon.palo-alto.med.va.gov/
<http://pkpd.icon.palo-alto.med.va.gov/>

2) The final parameter estimates, standard errors, and other quantities are available in NONMEM read-only commons ROCM6, ROCM7, ROCM8, and ROCM9. These are discussed in Guide VIII, or try "nmhelp rocm6" etc.

If you use $PRED, you can obtain them in a block of code that tests for ICALL.EQ.3.

Or if you use PREDPP, you can use code such as the following in an INFN routine. E.g., if the code is in a file named myinfn, then specify
$SUBROUTINES ... INFN=myinfn

Here's the code. It is suitable for diagonal omega and sigma, but needs some modification for full omega or sigma. (Constants 40 and 30 in the declarations for THETAF etc. should be the same as the values given to LTH and LVR in file NSIZES, in case you changed those values.) It also performs a pass thru the data file to obtain individual etas, which you may delete if you do not want them.

SUBROUTINE INFN(ICALL,THETA,DATREC,INDXS,NEWIND)
DIMENSION THETA(*),DATREC(*),INDXS(*)
DOUBLE PRECISION THETA
COMMON /ROCM6/ THETAF(40), OMEGAF(30,30), SIGMAF(30,30)
COMMON /ROCM7/ SETH(40),SEOM(30,30),SESIG(30,30)
COMMON /ROCM8/ OBJECT
COMMON /ROCM9/ IERE,IERC
DOUBLE PRECISION THETAF, OMEGAF, SIGMAF
DOUBLE PRECISION OBJECT
REAL SETH,SEOM,SESIG
DOUBLE PRECISION ETA(10)
INTEGER J,I
INTEGER IERE,IERC
INTEGER MODE
INTEGER NTH,NETA,NEPS
C THE NEXT 3 LINES SHOULD BE EDITED FOR THE NUMBERS OF ETAs, THETAs
C AND EPSILONs IN YOUR PROBLEM
NTH=8
NETA=6
NEPS=2
IF (ICALL.EQ.0) THEN
OPEN(49,FILE='ETA')
OPEN(50,FILE='OBJECT')
OPEN(51,FILE='THETA')
OPEN(52,FILE='SETHETA')
OPEN(53,FILE='OMEGA')
OPEN(54,FILE='SEOMEGA')
OPEN(55,FILE='SIGMA')
OPEN(56,FILE='SESIGMA')
OPEN(57,FILE='IERE_C')
ENDIF
C THE NEXT 11 LINES WRITE INDIVIDUAL ETAS IF 'POSTHOC'
C OPTION CHOSEN IN THE ESTIMATION STEP
IF (ICALL.EQ.3) THEN
MODE=0
CALL PASS(MODE)
MODE=1
20 CALL PASS(MODE)
IF (MODE.EQ.0) GO TO 30
IF (NEWIND.NE.2) THEN
CALL GETETA(ETA)
WRITE (49,97) (ETA(I),I=1,NETA)
ENDIF
GO TO 20
30 CONTINUE
WRITE (50,99) OBJECT
WRITE (51,99) (THETAF(J),J=1,NTH)
WRITE (52,99) (SETH(J),J=1,NTH)
7000 WRITE (53,99) (OMEGAF(J,J),J=1,NETA)
WRITE (54,99) (SEOM(J,J),J=1,NETA)
7999 WRITE (55,99) (SIGMAF(J,J),J=1,NEPS)
WRITE (56,99) (SESIG(J,J), J=1,NEPS)
WRITE (57,98) IERE,IERC
ENDIF
99 FORMAT (20F15.8)
98 FORMAT (2I8)
97 FORMAT (F8.0, 10E15.8)
RETURN
END

With (say) $OMEGA BLOCK(3) .... the change is as follows:

Locate:
NETA=6
7000 WRITE (53,99) (OMEGAF(J,J),J=1,NETA)
WRITE (54,99) (SEOM(J,J),J=1,NETA)

Change to:
NETA=3
7000 WRITE (53,99) ((OMEGAF(J,I),J=1,I),I=1,NETA)
WRITE (54,99) ((SEOM(J,I),J=1,I),I=1,NETA)


*****


From: R.Port@DKFZ-Heidelberg.de
Subject: Reading NONMEM output into S-Plus
Date: Tue, 6 Feb 2001 17:42:54 +0100

Hi Ganesh:

a few more suggestions...

If this is your $TABLE record:

$TABLE STUD ID CL V2 Q V3 KA TIME IPRED NOPRINT ONEHEADER FILE=BLUESKY2.TAB

Run NONMEM as usually to create BLUESKY2.TAB . Then, in S-Plus, issue:

> a_read.table('BLUESKY2.TAB',skip=1,header=T)

# If you want to plot residuals versus time:

> a1_a[a[,'DV']!=0,] # Extract lines where DV!=0
> plot(a1[,'TIME'],a1[,'RES'])

If you want to use your final THETA estimates in S-Plus (e.g. to create smooth prediction curves) you can also make them a data item of the output table:

$ERROR ...

IF(ICALL.EQ.3) THEN
TH1 = THETA(1)
TH2 = THETA(2)
TH3 = THETA(3)
TH4 = THETA(4)
ENDIF

$TABLE STUD ID TH1 TH2 TH3 TH4 CL V2 Q V3 KA TIME IPRED
NOPRINT ONEHEADER FILE=BLUESKY2.TAB

Then, in S-Plus:

> a_read.table('BLUESKY2.TAB',skip=1,header=T)
> theta1_a[1,'TH1'] # theta1 is mean population k_a

If you want to see final OMEGA^2 estimates in S-Plus, or objective function values, things become more complicated, and you will need something like what Vladimir Piotrovsky suggests, i.e. extracting lines from the NONMEM report file and scan them into S-Plus' using the scan () function, or making use of NONMEM commons (ROCM6...) as described by Alison.

Dr. Nicholas Holford has written very useful awk procedures to be used under UNIX that extract the more important information from the NONMEM report file (I guess they can be found in the NONMEM repository or on his own Web site).

By the way, if you are doing simulations, you might like to invoke NONMEM from within an S-Plus command file so that everything is done with one command on the S-Plus command line. E.g. if you are working under UNIX and the script that invokes NONMEM is named "start.NONMEM":

# S-Plus command file, e.g. 'graph1': ----------------------------------------

!start.NONMEM control.file report.file

a_read.table('BLUESKY2.TAB',skip=1,header=T)
a1_a[a[,'DV']!=0,]
plot(a1[,'TIME'],a1[,'RES'])

# ----------------------------------------------------------------------------

On the S-Plus command line: > source('graph1')

Good luck! Ruedi


*****


From: Nick Holford <n.holford@auckland.ac.nz>
Subject: Re: Reading NONMEM output into S-Plus
Date: Wed, 07 Feb 2001 09:50:03 +1300

R.Port@DKFZ-Heidelberg.de wrote:

> Dr. Nicholas Holford has written very useful awk procedures to be used under
> UNIX that extract the more important information from the NONMEM report file
> (I guess they can be found in the NONMEM repository or on his own Web site).

Just to clarify:
A set of awk scripts are included in the Wings for NONMEM utilities. WFN is primarily intended for use on Windoze (i.e. DOS) platforms. The awk scripts are essentially platform independent and can be used under Unix but the "glue" for WFN comes from DOS batch command files rather than unix shell scripts. WFN is free and can be downloaded from:

http://www.geocities.com/wfn2k

WFN provides DOS users a means to run NONMEM with output extracted from the standard NONMEM listing file in tab delimited files suitable for reading into Excel, Splus, etc. No special code is required e.g. INFN subroutines. Support is provided for non-parametric bootstrap, multiple problems within one control stream, and simulation with multiple subproblems.

--
Nick Holford, Divn Pharmacology & Clinical Pharmacology
University of Auckland, 85 Park Rd, Private Bag 92019, Auckland, New Zealand
email:n.holford@auckland.ac.nz tel:+64(9)373-7599x6730 fax:373-7556
http://www.phm.auckland.ac.nz/Staff/NHolford/nholford.htm


*****


From: "Niclas Jonsson" <Niclas.Jonsson@farmbio.uu.se>
Subject: Re: Reading NONMEM output into S-Plus
Date: Wed, 7 Feb 2001 08:37:58 +0000 (GMT)

If you want to get hold of the parameter estimates in a NONMEM output file you could use teh read.lst function that comes with Xpose, e.g:

> tmp <- read.lst("run13.lst")
> tmp
$term:
[1] "MINIMIZATION SUCCESSFUL"
[2] "NO. OF FUNCTION EVALUATIONS USED: 208"
[3] "NO. OF SIG. DIGITS IN FINAL EST.: 3.9"

$ofv:
[1] 2116.521

$thetas:
[1] 20.30000 84.90000 1.77000 0.00566

$omega:
$omega[[1]]:
[1] 0.136

$omega[[2]]:
[1] 0.000 0.182

$omega[[3]]:
[1] 0.00 0.00 1.39


$sigma:
$sigma[[1]]:
[1] 0.0238


$sethetas:
[1] 1.07000 4.51000 0.39000 0.00317

$seomegas:
$seomegas[[1]]:
[1] 0.0297

$seomegas[[2]]:
[1] 0.0000 0.0299

$seomegas[[3]]:
[1] 0.000 0.000 0.663


$sesigmas:
$sesigmas[[1]]:
[1] 0.00469

In other words, you will get back a list with various components like ofv, thetas and so on. If you want to access, for example, the thetas from the tmp object created above, you would type:

> tmp$thetas

and that would give you a vector of the theta estimates, i.e.

[1] 20.30000 84.90000 1.77000 0.00566

Best,
Niclas
--
Department of Pharmacy
Uppsala University
Box 580
SE-751 23 Uppsala
Sweden
Phone: +46 18 471 43 85
Fax: +46 18 471 40 03
Mobile: +46 70 485 61 98
E-mail: niclas.jonsson@biof.uu.se


*****


From: "Piotrovskij, Vladimir [JanBe]" <VPIOTROV@janbe.jnj.com>
Subject: RE: Reading NONMEM output into S-Plus
Date: Wed, 7 Feb 2001 12:44:56 +0100

I have just checked the read.lst() function of Xpose3. I have not installed the whole library, but read.lst() only (plus a couple of functions invoked by it). I find this function very nice and I do not see any reasons for not using it. It exctracts necessary information directly from the standard NONMEM output file. Consequently, no additional files with THETAs, OMEGAs, etc. are needed. INFN subroutine is not needed anymore.

Thanks to Niclas and Mats!

Best regards,
Vladimir
----------------------------------------------------------------------
Vladimir Piotrovsky, Ph.D.
Janssen Research Foundation
Human Pharmacokinetics (ext. 5463)
B-2340 Beerse
Belgium
Email: vpiotrov@janbe.jnj.com