STBLE

Estimates survival probabilities and hazard rates for various parametric models.

Required Arguments

XPT — NOBS by NCOL matrix, each row of which contains the covariates for a group for which survival estimates are desired. (Input)

MODEL — Model option parameter. (Input)
MODEL specifies the distribution of the response variable and the relationship of the linear model to a distribution parameter.

MODEL	Distribution
0	Exponential
1	Linear hazard
2	Log‑normal
3	Normal
4	Log‑logistic
5	Logistic
6	Log least extreme value
7	Least extreme value
8	Log extreme value
9	Extreme value
10	Weibull

For further discussion of the models, see the “Description” section.

TIME — Beginning of the time grid for which the survival estimates are desired. (Input)
Survival probabilities and hazard rates are computed for each covariate vector over the grid of time points TIME + i * DELTA for i = 0, 1, …, NPT ‑ 1.

NPT — Number of points on the time grid for which survival probabilities are desired. (Input)

DELTA — Increment between time points on the time grid. (Input)

IFIX — Column number in XPT containing a constant to be added to the linear response. (Input)
The estimated linear response is w + COEF(1) * z(1) + COEF(2) * z(2) + … + COEF(NCOEF) * z(NCOEF), where z is the design vector for the I‑th observation obtained from a row of XPT. w = XPT(I, IFIX) if IFIX is positive, and w = 0 otherwise.

NCOEF — Number of coefficients in the model. (Input)

COEF — Vector of length NCOEF containing the model parameter estimates. (Input)
Usually routine SVGLM is first called to estimate COEF as the first column of matrix COEF in SVGLM. When present in the model, the initial coefficient in COEF is a “nuisance” parameter, and the remaining coefficients are parameters associated with the “linear” model, beginning with the intercept, if present. Nuisance parameters are as follows:

Model	Nuisance Parameter
1	Coefficient of the quadratic term in time, θ
2 – 9	Scale parameter, σ
10	Shape parameter, θ

There is no nuisance parameter for model 0.

SPROB — NPT by 2 * NOBS + 1 matrix. (Output)
SPROB(i, 2) contains the estimated survival probability at time SPROB(i, 1) = TIME + (i ‑ 1) * DELTA for observations with covariates given in row 1 of XPT. SPROB(i, 3) contains the estimate for the hazard rate at this time point. Columns 4 and 5 contain the estimated survival probabilities and hazard rates for observations with covariates given in the second row in XPT, etc., up to columns 2 * NOBS and 2 * NOBS + 1, which contain these statistics for observations with covariates in the last row of XPT.

XBETA — Vector of length NOBS containing the estimated linear response w + COEF(1) * z(1) + … + COEF(NCOEF) * z(NCOEF ) for each row of XPT. (Output)

Optional Arguments

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

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

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

INTCEP — Intercept option. (Input)
Default: INTCEP = 1.

INTCEP	Action
0	No intercept is in the model (unless otherwise provided for by the user).
1	An intercept is automatically included in the model.

NCLVAR — Number of classification variables. (Input)
Dummy or indicator variables are generated for classification variables using the IDUMMY = 2 option of routine GRGLM (see Chapter 2, “Regression”). See also Comment 2.
Default: NCLVAR = 0.

INDCL — Index vector of length NCLVAR containing the column numbers of X that are classification variables. (Input, if NCLVAR is positive, not used otherwise)
If NCLVAR is 0, INDCL is not referenced and can be dimensioned of length 1 in the calling program.

NCLVAL — Vector of length NCLVAR containing the number of values taken on by each classification variable. (Input, if NCLVAR is positive, not referenced otherwise)
NCLVAL(I) is the number of distinct values for the I‑th classification variable. NCLVAL is not referenced and can be dimensioned of length 1 in the calling program if NCLVAR is zero.

CLVAL — Vector of length NCLVAL(1) + NCLVAL(2) + … + NCLVAL(NCLVAR) containing the distinct values of the classification variables. (Input, if NCLVAR is positive; not used otherwise)
The first NCLVAL(1) elements contain the values for the first classification variables, the next NCLVAL(2) elements contain the values for the second classification variable, etc. If NCLVAR is zero, then CLVAL is not referenced and can be dimensioned of length 1 in the calling program.

NEF — Number of effects in the model. (Input)
In addition to effects involving classification variables, simple covariates and the product of simple covariates are also considered effects.
Default: NEF = 0.

NVEF — Vector of length NEF that contains the number of variables associated with each effect. (Input, if NEF is greater than 0; not referenced otherwise)
NVEF is not referenced and can be dimensioned of length 1 in the calling program if NEF is zero.

INDEF — Vector of length NVEF(1) + … + NVEF(NEF) that contains the column numbers in X associated with each effect. (Input, if NEF is greater than 0; not used otherwise)
The first NVEF(1) elements of INDEF contain the column numbers in XPT for the variables in the first effect. The next NVEF(2) elements in INDEF contain the column numbers for the second effect, etc. If NCLVAR is zero, INDEF is not referenced and can be dimensioned of length 1 in the calling program.

IPRINT — Printing option. (Input)
Default: IPRINT = 0.

IPRINT	Action
0	No printing is performed.
1	Printing is performed.

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

FORTRAN 90 Interface

Generic: CALL STBLE (XPT, MODEL, TIME, NPT, DELTA, IFIX, NCOEF, COEF, SPROB, XBETA [, …])

Specific: The specific interface names are S_STBLE and D_STBLE.

FORTRAN 77 Interface

Single: CALL STBLE (NOBS, NCOL, XPT, LDXPT, MODEL, TIME, NPT, DELTA, IFIX, INTCEP, NCLVAR, INDCL, NCLVAL, CLVAL, NEF, NVEF, INDEF, NCOEF, COEF, IPRINT, SPROB, LDSPRO, XBETA)

Double: The double precision name is DSTBLE.

Description

Routine STBLE computes estimates of survival probabilities and hazard rates for the parametric survival/reliability models fit by routine SVGLM for one or more vectors of covariate values. Because estimates for the parameters of the model must be given, routine SVGLM is usually invoked to obtain these estimates prior to invoking STBLE.

Let η = xTβ be the linear parameterization, where x is a design vector obtained in STBLE via routine GRGLM (see Chapter 2, “Regression”) from a row of XPT, and β is a vector of parameters associated with the linear model. Let T denote the random response variable and S(t) denote the probability that T > t. All models considered also allow a fixed parameter w (input in column IFIX of XPT). Use of this parameter is discussed in the document for routine SVGLM. There may also be nuisance parameters θ > 0, or σ > 0. Let Φ denote the cumulative normal distribution. The survival models available in STBLE are

Model	Name	S(t)
0	Exponential
1	Linear hazard
2	Log‑normal
3	Normal
4	Log‑logistic
5	Logistic
6	Log least extreme value
7	Least extreme value
8	Log extreme value
9	Extreme value
10	Weibull

Let λ(t) denote the hazard rate at time t. Then λ(t) and S(t) are related as

Models 0, 1, 2, 4, 6, 8, and 10 require that T > 0 (in which case, we assume λ(s) = 0 for s < 0), while the remaining models allow arbitrary values for T, ‑∞ < T < ∞. The computations proceed in routine STBLE as follows:

1. The input arguments are checked for consistency and validity.

2. For each row of XPT, the explanatory variables are generated from the classification and variables and the covariates using routine GRGLM with the IDUMMY = 2 option. (When IDUMMY is two, GRGLM assigns an indicator variable the value 1.0 when the observation is in the class, assigns the value 0.0 otherwise, and omits the last indicator variable from the design vector. See the manual documentation for GRGLM.) Given the explanatary variables x, η is computed as η = xTβ, where β is input in COEF.

3. For each time point requested in the time grid, the survival probabilities and hazard rates are computed.

Comments

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

CALL S2BLE (NOBS, NCOL, XPT, LDXPT, MODEL, TIME, NPT, DELTA, IFIX, INTCEP, NCLVAR, INDCL, NCLVAL, CLVAL, NEF, NVEF, INDEF, NCOEF, COEF, IPRINT, SPROB, LDSPRO, XBETA, CHWK, Z, RWK)

The additional arguments are as follows:

CHWK — CHARACTER * 10 work vector of length NCOL.

Z — Work vector of length NCOEF.

RWK — Work vector of length MAX(7, NCOL) if IPRINT = 1, or of length 1 if IPRINT = 0.

2. Dummy variables are generated for the classification variables as follows: The list of all distinct values of each classification variable is as stored in CLVAL. Dummy variables are generated for each but the last of these distinct values. Each dummy variable is zero unless the classification variable equals the list value corresponding to the dummy variable, in which case the dummy variable is one. See argument IDUMMY for IDUMMY = 2 in routine GRGLM (see Chapter 2, “Regression”).

3. Informational errors

Type	Code	Description
3	1	Some survival probabilities are less than or equal to zero. The corresponding hazard values cannot be computed.
4	2	The specified number of coefficients, NCOEF, is incorrect.
4	3	The model specified is not defined for negative time.

Example

The example is a continuation of the first example given for routine SVGLM. Prior to calling STBLE, SVGLM is invoked to compute the parameter estimates. The example is taken from Lawless (1982, page 287) and involves the mortality of patients suffering from lung cancer.

USE SVGLM_INT

USE SCOPY_INT

USE STBLE_INT

IMPLICIT NONE

INTEGER ICEN, IFIX, ILT, INFIN, IPRINT, IRT, LDCASE, &

LDCOEF, LDCOV, LDSPRO, LDX, LDXPT, MAXCL, &

MODEL, NCLVAR, NCOL, NEF, NOBS, NPT

REAL DELTA, TIME, XPWR

PARAMETER (DELTA=20.0, ICEN=2, IFIX=0, ILT=0, INFIN=0, IPRINT=1, &

IRT=1, LDCASE=40, LDCOEF=9, LDCOV=9, LDX=40, LDXPT=2, &

MAXCL=6, MODEL=0, NCLVAR=2, NCOL=7, NEF=5, NOBS=40, &

NPT=10, TIME=10.0, XPWR=0.0, LDSPRO=NPT)

INTEGER IADDS(NOBS), INDCL(NCLVAR), INDEF(5), NCLVAL(NCLVAR), &

NCOEF, NRMISS, NVEF(NEF)

REAL ALGL, CASE(LDCASE,5), CLVAL(MAXCL), COEF(LDCOEF,4), &

COV(LDCOV,LDCOV), GR(LDCOV), SPROB(LDSPRO,2*NOBS+1), &

X(LDX,NCOL), XBETA(NOBS), XMEAN(LDCOV), XPT(LDXPT,NCOL)

DATA X/411, 126, 118, 92, 8, 25, 11, 54, 153, 16, 56, 21, 287, &

10, 8, 12, 177, 12, 200, 250, 100, 999, 231, 991, 1, 201, &

44, 15, 103, 2, 20, 51, 18, 90, 84, 164, 19, 43, 340, 231, &

5*0, 1, 16*0, 1, 5*0, 1, 11*0, 7, 6, 7, 4, 4, 7, 7, 8, 6, &

3, 8, 4, 6, 4, 2, 5, 5, 4, 8, 7, 6, 9, 5, 7, 2, 8, 6, 5, 7, &

4, 3, 3, 4, 6, 8, 7, 3, 6, 8, 7, 64, 63, 65, 69, 63, 48, &

48, 63, 63, 53, 43, 55, 66, 67, 61, 63, 66, 68, 41, 53, 37, &

54, 52, 50, 65, 52, 70, 40, 36, 44, 54, 59, 69, 50, 62, 68, &

39, 49, 64, 67, 5, 9, 11, 10, 58, 9, 11, 4, 14, 4, 12, 2, &

25, 23, 19, 4, 16, 12, 12, 8, 13, 12, 8, 7, 21, 28, 13, 13, &

22, 36, 9, 87, 5, 22, 4, 15, 4, 11, 10, 18, 7*1, 7*2, 2*3, &

5*4, 7*1, 4*2, 3*3, 5*4, 21*0, 19*1/

DATA NVEF/1, 1, 1, 1, 1/, INDEF/3, 4, 5, 6, 7/, INDCL/6, 7/

CALL SVGLM (X, MODEL, ILT, IRT, MAXCL, NCOEF, COEF, ALGL, &

COV, XMEAN, CASE, GR, IADDS, ICEN=ICEN, &

NCLVAR=NCLVAR, INDCL=INDCL, NEF=NEF, NVEF=NVEF, &

INDEF=INDEF, NCLVAL=NCLVAL, CLVAL=CLVAL)

CALL SCOPY (NCOL, X(1:,1), LDX, XPT(1:,1), LDXPT)

CALL SCOPY (NCOL, X(2:,1), LDX, XPT(2:,1), LDXPT)

CALL STBLE (XPT, MODEL, TIME, NPT, DELTA, IFIX, &

NCOEF, COEF, SPROB, XBETA, NCLVAR=NCLVAR, INDCL=INDCL,&

NCLVAL=NCLVAL, CLVAL=CLVAL, NEF=NEF, NVEF=NVEF,&

INDEF=INDEF, IPRINT=IPRINT)

END

Output

group 1

xpt

1 2 3 4 5 6

411 0 7 64 5 1

design vector

1 2 3 4 5 6

1 7 64 5 1 0

7 8

0 1

xbeta = -5.57097

group 2

xpt

1 2 3 4 5 6

126 0 6 63 9 1

design vector

1 2 3 4 5 6

1 6 63 9 1 0

7 8

0 1

xbeta = -5.03551

survival and hazard estimates

(sprob)

time s1 h1 s2 h2

10.00 0.9626 0.003807 0.9370 0.006503

30.00 0.8921 0.003807 0.8228 0.006503

50.00 0.8267 0.003807 0.7224 0.006503

70.00 0.7661 0.003807 0.6343 0.006503

90.00 0.7099 0.003807 0.5570 0.006503

110.00 0.6579 0.003807 0.4890 0.006503

130.00 0.6096 0.003807 0.4294 0.006503

150.00 0.5649 0.003807 0.3770 0.006503

170.00 0.5235 0.003807 0.3310 0.006503

190.00 0.4852 0.003807 0.2907 0.006503

Note that in simple exponential models the hazard rate is constant over time.