Chapter 13: Survival Analysis, Life Testing, and Reliability

KAPMR

Computes Kaplan-Meier estimates of survival probabilities in stratified samples.

Required Arguments

XNOBS by NCOL matrix containing the data.   (Input)

IRT — Column number in X containing the response variable.   (Input)
For the i-th right-censored observation, X(i, IRT) contains the right-censoring time. Otherwise, X(i, IRT) contains the failure time. (See ICEN.)

SPROBNOBS by 2 matrix.   (Output)
SPROB(i, 1) contains the estimated survival probability at time X(i, IRT) in the i-th observation’s stratum, while SPROB(i, 2) contains Greenwood’s estimate of the standard deviation of this estimated probability. If the i-th observation contains censor codes out of range or if a variable is missing, then the corresponding elements of SPROB are set to missing (NaN, not a number). Similarly, if an element in SPROB is not defined, then it is set to missing.

Optional Arguments

NOBS — Number of observations.   (Input)
Default: NOBS = size (X,1).

NCOL — Number of columns in X.   (Input)
Default: NCOL = size (X,2).

LDX — Leading dimension of X exactly as specified in the dimension statement in the calling program.   (Input)
Default: LDX = size (X,1).

IFRQ — Column number in X containing the frequency of response for this observation.   (Input)
If IFRQ = 0, a response frequency of 1 for each observation is assumed.
Default: IFRQ = 0.

ICEN — Column number in X containing the censoring code for this observation.   (Input)
Default: ICEN = 0.
If ICEN = 0, a censoring code of 0 is assumed. Valid censoring codes are:

            Code                       Meaning

            0                              Exact failure at X(i, IRT).

            1                              Right censored. The response is greater than X(i, IRT).

            If X(i, ICEN) is not 0 or 1, then the i-th observation is omitted from the analysis.

IGRP — Column number in X containing the stratum number for this observation.   (Input)
If IGRP = 0, the data is assumed to be from one stratum. Otherwise, column IGRP of X contains a unique value for each stratum in the data. Kaplan-Meier estimates are computed within each stratum.
Default: IGRP = 0.

ISRT — Sorting option.   (Input)
If ISRT = 1, column IRT of X is assumed to be sorted in ascending order within each stratum. Otherwise, a detached sort will be performed by KAPMR . If sorting is performed by KAPMR, all censored individuals are assumed to follow tied failures.
Default: ISRT = 0.

LDSPRO — Leading dimension of SPROB exactly as specified in the dimension statement in the calling program.   (Input)
Default: LDSPRO = size (SPROB,1).

NRMISS — Number of rows of data in X that contain any missing values.   (Output)

FORTRAN 90 Interface

Generic:                              CALL KAPMR (X, IRT, SPROB [,…])

Specific:                             The specific interface names are S_KAPMR and D_KAPMR.

FORTRAN 77 Interface

Single:                                CALL KAPMR (NOBS, NCOL, X, LDX, IRT, IFRQ, ICEN, IGRP, ISRT, SPROB, LDSPRO, NRMISS)

Double:                              The double precision name is DKAPMR.

Description

Routine KAPMR computes Kaplan-Meier (or product-limit) estimates of survival probabilities for a sample of failure times that possibly contain right censoring. A survival probability S(t) is defined as 1 F(t), where F(t) is the cumulative distribution function of the failure times (t). Greenwood’s estimate of the standard errors of the survival probability estimates are also computed. (See Kalbfleisch and Prentice, 1980, pages 13 and 14.)

Let (ti, δi), for i = 1,…, n denote the failure/censoring times and the censoring codes for the n observations in a single sample. Here, ti = X(i, IRT) is a failure time if δi is 0, where
δi = X(i, ICEN). Also, ti is a censoring time if δi is 1. Rows in X containing values other than 0 or 1 for δi are ignored. Let the number of observations in the sample that have not failed by time s(ι) be denoted by n(ι), where s(ι) is an ordered (from smallest to largest) listing of the distinct failure times (censoring times are omitted). Then the Kaplan-Meier estimate of the survival probabilities is a step function, which in the interval from s(ι) to s(i+1) (including the lower endpoint) is given by

where d(j) denotes the number of failures occurring at time s(j). Note that one row of X may correspond to more than one failed (or censored) observation when the frequency option is in effect (IFRQ is not zero). The Kaplan-Meier estimate of the survival probability prior to time s(1) is 1.0, while the Kaplan-Meier estimate of the survival probability after the last failure time is not defined.

Greenwood’s estimate of the variance of

in the interval from s(i) to s(i+1) is given as

Routine KAPMR computes the single sample estimates of the survival probabilities for all samples of data included in X during a single call. This is accomplished through the IGRP column of X, which if present, must contain a distinct code for each sample of observations. If IGRP = 0, there is no grouping column, and all observations are assumed to be from the same sample.

When failures and right-censored observations are tied and the data are to be sorted by KAPMR (ISRT is not 1), KAPMR assumes that the time of censoring for the tied-censored observations is immediately after the tied failure (within the same sample). When the ISRT = 1 option is in effect, the data are assumed to be sorted from smallest to largest according to column IRT of X within each stratum. Furthermore, a small increment of time is assumed (theoretically) to elapse between the failed and censored observations that are tied (in the same sample). Thus, when the ISRT = 1 option is in effect, the user must sort all of the data in X from smallest to largest according to column IRT (and column IGRP, if present). By appropriate sorting of the observations, the user can handle censored and failed observations that are tied in any manner desired.

Comments

1.         Workspace may be explicitly provided, if desired, by use of K2PMR/DK2PMR. The reference is:

CALL K2PMR (NOBS, NCOL, X, LDX, IRT, IFRQ, ICEN, IGRP, ISRT, SPROB, LDSPRO, NRMISS, IGP, IPERM, INDDR, IWK, WK, IPER)

The additional arguments are as follows:

IGP — Work vector of length NOBS.

IPERM — Work vector of length NOBS + NCOL.

INDDR — Work vector of length NOBS.

IWK — Work vector of length max(NOBS, NCOL).

WK — Work vector of length 2 * max(NOBS, NCOL).

IPER — Work vector of length NOBS.

2.         Missing values may occur in any of the columns of X. Any row of X that contains missing values in the IRT, ICEN, or IFRQ columns (when the ICEN and IFRQ columns are present) is omitted from the analysis. Missing values in the IGRP column, if present, are classified into an additional “missing” group.

Example

The following example is taken from Kalbfleisch and Prentice (1980, page 1). The first column in X contains the death/censoring times for rats suffering from vaginal cancer. The second column contains information as to which of two forms of treatment were provided, while the third column contains the censoring code. Finally, the fourth column contains the frequency of each observation. The product-limit estimates of the survival probabilities are computed for both groups with one call to KAPMR. In this example, the output in SPROB has been equivalenced with columns 5 and 6 of X so that the input and output matrices could be printed together. Routine KAPMR could have been called with the ISRT = 1 option in effect if the censored observations had been sorted with respect to the failure time variable.

 

      USE KAPMR_INT

      USE WRRRL_INT

      USE UMACH_INT

 

      IMPLICIT   NONE

      INTEGER    ICEN, IFRQ, IGRP, IRT, ISRT, LDSPRO, LDX, NCOL, NOBS

      PARAMETER  (ICEN=3, IFRQ=4, IGRP=2, IRT=1, ISRT=0, LDSPRO=33, &

                 LDX=33, NCOL=6, NOBS=33)

!

      INTEGER    NOUT, NRMISS

      REAL       SPROB(LDSPRO,2), X(LDX,NCOL)

      CHARACTER  XLABEL(7)*6, YLABEL(1)*6

!

      EQUIVALENCE (X(1,5), SPROB)

!

      DATA XLABEL/'OBS', 'TIME', 'GROUP', 'CENSOR', 'FREQ', 'S-HAT', &

          'SE'/

      DATA YLABEL/'NUMBER'/

      DATA X/143, 164, 188, 190, 192, 206, 209, 213, 216, 220, 227, &

          230, 234, 246, 265, 304, 216, 244, 142, 156, 163, 198, 205, &

          232, 233, 239, 240, 261, 280, 296, 323, 204, 344, 18*5, &

          15*7, 16*0, 2*1, 13*0, 4*1, 2, 20*1, 2, 4, 3*1, 2*2, 3*1, &

          66*0/

!

      CALL KAPMR (X, IRT, SPROB, IFRQ=IFRQ, ICEN=ICEN, IGRP=IGRP, &

                  NRMISS=NRMISS)

!

      CALL WRRRL ('X/SPROB', X, YLABEL, XLABEL, FMT='(W10.6)')

      CALL UMACH (2, NOUT)

      WRITE (NOUT,'(//'' NRMISS = '', I5)') NRMISS

      END

Output

 

                                 X/SPROB
 OBS     TIME       GROUP      CENSOR        FREQ       S-HAT          SE
 1    143.000       5.000       0.000       1.000       0.947       0.051
 2    164.000       5.000       0.000       1.000       0.895       0.070
 3    188.000       5.000       0.000       2.000       0.789       0.094
 4    190.000       5.000       0.000       1.000       0.737       0.101
 5    192.000       5.000       0.000       1.000       0.684       0.107
 6    206.000       5.000       0.000       1.000       0.632       0.111
 7    209.000       5.000       0.000       1.000       0.579       0.113
 8    213.000       5.000       0.000       1.000       0.526       0.115
 9    216.000       5.000       0.000       1.000       0.474       0.115
10    220.000       5.000       0.000       1.000       0.414       0.115
11    227.000       5.000       0.000       1.000       0.355       0.112
12    230.000       5.000       0.000       1.000       0.296       0.108
13    234.000       5.000       0.000       1.000       0.237       0.101
14    246.000       5.000       0.000       1.000       0.158       0.093
15    265.000       5.000       0.000       1.000       0.079       0.073
16    304.000       5.000       0.000       1.000       0.000         NaN
17    216.000       5.000       1.000       1.000       0.474       0.115
18    244.000       5.000       1.000       1.000       0.237       0.101
19    142.000       7.000       0.000       1.000       0.952       0.046
20    156.000       7.000       0.000       1.000       0.905       0.064
21    163.000       7.000       0.000       1.000       0.857       0.076
22    198.000       7.000       0.000       1.000       0.810       0.086
23    205.000       7.000       0.000       1.000       0.759       0.094
24    232.000       7.000       0.000       2.000       0.658       0.105
25    233.000       7.000       0.000       4.000       0.455       0.111
26    239.000       7.000       0.000       1.000       0.405       0.110
27    240.000       7.000       0.000       1.000       0.354       0.107
28    261.000       7.000       0.000       1.000       0.304       0.103
29    280.000       7.000       0.000       2.000       0.202       0.090
30    296.000       7.000       0.000       2.000       0.101       0.068
31    323.000       7.000       0.000       1.000       0.051       0.049
32    204.000       7.000       1.000       1.000       0.810       0.086
33    344.000       7.000       1.000       1.000         NaN         NaN

NRMISS =     0



http://www.vni.com/
PHONE: 713.784.3131
FAX:713.781.9260