ParNPHCox: ParNPHCox

In AB5103/doubleCoxr: Parametric Gamma-Frailty Models with Non-Proportional Hazard Functions

Description

Provides parameter estimates, their standard errors, loglikelihood and other attributes for the parametric gamma-frailty model with non-proportional hazard functions

Usage

ParNPHCox(formula.scale, formula.shape, cluster, dist, data, expr, strata)

Arguments

formula.scale	A formula object, with the response on the left of a ~ operator, and the terms on the right. The response must be a survival object as returned by the Surv() function. The status indicator 'event' in the Surv object must be 0=censored, 1=non-censored. This object describes the effect of several factors in the Cox proportional hazards model.
formula.shape	A formula object, with the response on the left of a ~ operator, and the terms on the right. The response must be a survival object as returned by the Surv() function. The status indicator 'event' in the Surv object must be 0=censored, 1=non-censored. This object describes the effect of several factors on the shape parameter.
cluster	The name of a cluster variable in data (is equal to NULL for the fixed-effect model)
dist	Baseline hazard function ('Weibull' or 'Gompertz').
data	A data.frame in which to interpret the variables named in the formula. The data set includes the following fields: time-to-failure and censoring in the case without left truncation or time-of-start, time-of-failure, and censoring in the case with left truncation at the time of begin (censoring must be either 0 for no event or 1 for event); covariates (continuous or categorical) used in a study (can be empty set).
expr	The vector of expressions for contrasts. NULL, otherwise.
strata	Full list of the factors used in the study and their levels in the form: strata=list(factor1=level of the factor1,...,last factor=level of the last factor), NULL, otherwise.

Details

Two kinds of the baseline hazards are used in this function:

1. Weibull with baseline cumulative hazard function

   H_{0}(t;a,b)=(t/a)^b;

2. Gompertz with baseline cumulative hazard function

   H_{0}(t;a,b)=(a/b)(exp(bt)-1)).

The cumulative hazard function for the vector of covariates u is defined by

H_{cum}(t;a,b,β _{scale},β _{shape},u)=exp(β _{scale}u)H_{0}(t;a,exp(β _{shape}u)b).

The marginal survival function is defined by

S(t;a,b,β _{scale},β _{shape},σ ^2,u)=Eexp(-ZH_{cum}(t;a,b,β _{scale},β _{shape},u))=(1+σ ^2H_{cum}(t;a,b,β _{scale},β _{shape},u))^{-1/σ ^2}

for gamma-distributed frailty Z with mean 0 and variance σ ^2. The 'formula.scale' and 'formula.shape' are formula objects used in the R-package 'survival' and have the form

Surv(time, Cens) ~ factor_1+...+factor_k).

Interactions between factors are allowed. If no factors are used in the Cox regression the value of 1 stands on the right side of the formula. The records with NA value for at least one factor in both formulas are excluded from the analysis.

Remark 1. The concordance is evaluated using the function 'survConcordance' from the R-package 'survival'.

Remark 2. The mean contrasts, CIs, and p-values are calculated on the basis of the empirical distribution for contrasts constructed using 10^6 of the generations of the vector of parameters.

Remark 3. In some cases (flat likelihood function, multimodality, etc.) the hessian cannot be correctly estimated and the function issues the following error message:

Error in ParNPHCox(formula.scale, formula.shape, cluster, dist, data = ...) : hessian cannot be correctly calculated.

Change the model and try again. It can occur if, for example, the data set is not large enough the estimates of the shape and scale parameters are poorly calculated. In this case it is recommended to simplify the model excluding some factors from the formulas for shape or scale and try to calculate the estimates again.

Remark 4. The number of factors can increase after conversion of any categorical factor in binary ones. For example, if factor "type" has three levels "A", "B", and "C" we will get after conversion two binary factors - "typeB" and "typeC" (factor "typeA" is the baseline one and does not appear in the list of binary factors). The command "strata=list(...)" must include at most one non-zero level for each non-baseline level of the converted categorical factor. For example, it is correct to write "strata=list(typeB=1,typeC=0,...)" or briefly "strata=list(typeB=1,...)" and it is not correct to write "strata=list(typeB=1,typeC=1,...)".

Remark 5. If the covariate is a numerical one but we want to consider it as a categorical variable we can convert it to categorical one using the function "as.factor(variable)".

Value

List containing the following components:

• par: parameter estimates

• se: standard errors for parameter estimates

• LogLik: the value of the loglikelihood

• Tabs: the table of parameter estimates, their standard errors, and p-values in the Latex format

• Names: the names of estimated parameters

• Conc: concordance and its standard error

• pval: the vector of p-values for parameter estimates. For null hypothesis the values of parameters are equal to zero

• p.contrast: the data frame for means, CIs, and p-values for contrasts

• pstrata: data frame for times, means, and 95 cumulative hazards and instant hazards. The column names are: "Time", "Mean.survival", "Low.survival", "Upper.survival", "Mean.cumulative.hazard", "Low.cumulative.hazard", "Upper.cumulative.hazard", "Mean.hazard", "Low.hazard", "Upper.hazard". The duplicates are excluded

Examples

########## Example with 'Weibull baseline hazard function ############
require(survival)
require(numDeriv)
data("lung",package="survival")      #download data set "lung" from the r-package "survival"
lung$sex=lung$sex-1                  #converting numerical variable 'sex' with values 1 and 2
#to a binary one with values 0 and 1
lung$status=lung$status-1            #status variable must be 0 or 1
lung$ph.ecog=as.factor(lung$ph.ecog) #convering categorical variable in binary ones
formula.scale=as.formula('Surv(time, status)  ~ age + sex*ph.ecog')
formula.shape=as.formula('Surv(time, status) ~ 1') # no dependency
cluster='ph.karno'
dist='Weibull'
expr=expression(a -500* b, log(msurv(19)/msurv(14))-1, age.scale,sex.scale)
NamFact(lung,formula.scale,formula.shape)
#[1] "age" "ph.ecog1" "ph.ecog2" "ph.ecog3" "sex" "sex:ph.ecog1" "sex:ph.ecog2" "sex:ph.ecog3"
strata=list(age=0,sex=1,ph.ecog1=1)
dcWeibull=ParNPHCox(formula.scale,formula.shape,cluster,dist,data=lung,expr,strata)

########## Example with 'Gompertz baseline hazard function ############
formula.scale=as.formula('Surv(time, status) ~ sex')
formula.shape=as.formula('Surv(time, status) ~ sex')
cluster='pat.karno'
dist='Gompertz'
expr=expression(a - b, msurv(19) - msurv(14), sex.scale + sex.shape)
NamFact(lung,formula.scale,formula.shape)
#[1] "sex"
strata=list(sex=1)
dcGompertz=ParNPHCox(formula.scale,formula.shape,cluster,dist,data=lung,expr,strata)

AB5103/doubleCoxr documentation built on Feb. 20, 2022, 2:20 p.m.