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
7
0
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
7
0
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.
