frailtypack (version 2.0-2)

frailtyPenal: Fit Shared Gamma Frailty model using penalized likelihood estimation

Description

Fit a shared gamma frailty model using a Penalized Likelihood on the hazard function. Left truncated and censored data and strata (max=2) are allowed. It allows to obtain a non-parametric smooth hazard of survival function. This approach is different from the partial penalized likelihood approach of Therneau et al.

Usage

frailtyPenal(formula, data, Frailty = TRUE, recurrentAG=FALSE, 
             cross.validation=FALSE, n.knots, kappa1, kappa2, maxit=350)

Arguments

formula

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 like in survival package.

data

a data.frame in which to interpret the variables named in the 'formula'.

Frailty

Logical value. Is model with frailties fitted? If so variance of frailty parameter is estimated. If not, Cox proportional hazards model is estimated using Penalized Likelihood on the hazard function

recurrentAG

Logical value. Is Andersen-Gill model fitted? If so indicates that recurrent event times with the counting process approach of Andersen and Gill is used. This formulation can be used for dealing with time-dependent cov

cross.validation

Logical value. Is cross validation procedure used for estimating smoothing parameter? If so a search of the smoothing parameter using cross validation is done, with kappa1 as the seed. The cross validation is not implem

n.knots

integer giving the number of knots to use. Value required. It corresponds to the (n.knots+2) splines functions for the approximation of the hazard or the survival functions. Number of knots must be between 4 and 20.

kappa1

positive smoothing parameter. The coefficient kappa of the integral of the squared second derivative of hazard function in the fit (penalized log likelihood). We advise the user to identify several possible tuning parameters, note t

kappa2

positive smoothing parameter for the second stratum, when data are stratified. See kappa1.

maxit

maximum number of iterations for the Marquardt algorithm. Default is 350

Value

• an object of class '"frailtyPenal"'. Methods defined for 'frailtyPenal' objects are provided for print and plot. The following components are included in a 'frailtyPenal'object.
• n
  the number of observations used in the fit.
• groups
  the maximum number of groups used in the fit
• n.events
  the number of events observed in the fit
• logVerComPenal
  the complete marginal penalized log-likelihood
• theta
  variance of frailty parameter
• coef
  the coefficients of the linear predictor, which multiply the columns of the model matrix.
• varH
  the variance matrix of theta and of the coefficients.
• varHIH
  the robust estimation of the variance matrix of theta and of the coefficients.
• x1
  vector of times where both survival and hazard function are estimated. By default seq(0,max(time),length=99), where time is the vector of survival times.
• lam
  matrix of hazard estimates at x1 times and confidence bands.
• surv
  matrix of baseline survival estimates at x1 times and confidence bands.
• x2
  see x1 value for the second stratum
• lam2
  the same value as lam for the second stratum
• surv2
  the same value as surv for the second stratum

Details

The estimated parameter are obtained using the robust Marquardt algorithm (Marquardt, 1963) which is a combination between a Newton-Raphson algorithm and a steepest descent algorithm. When frailty parameter is small, numerical problems may arise. To solve this problem, an alternative formula of the penalized log-likelihood is used (see Rondeau, 2003 for further details). Cubic M-splines of order 4 are used for the hazard function, and I-splines (integrated M-splines) are used for the cumulative hazard function. PARAMETERS As frailtypack is written in Fortran 77 some parameters had to be hard coded in. The default values of these parameters are maximum number of observations: 60000 maximum number of groups: 5000 maximum number of subjects: 30000 If these parameters are not large enough (an error message will let you know this), you need to reset them in frailtypack.f and recompile. In particular, the statements defining these parameters are PARAMETER (ndatemax = 60000) PARAMETER (ngmax = 5000) PARAMETER (nsujetmax = 30000)

References

D. Marquardt (1963). An algorithm for least-squares estimation of nonlinear parameters. SIAM Journal of Applied Mathematics, 431-441. V. Rondeau, D Commenges, and P. Joly (2003). Maximum penalized likelihood estimation in a gamma-frailty model. Lifetime Data Analysis 9, 139-153. McGilchrist CA and Aisbett CW (1991). Regression with frailty in survival analysis. Biometrics 47, 461-466.

See Also

print.frailtyPenal, summary.frailtyPenal

Examples

data(kidney)
#Shared frailty model
frailtyPenal(Surv(time,status)~sex+age+cluster(id),
             n.knots=12,kappa1=1000,data=kidney)


#model without frailties (e.g., Cox proportional hazards 
#                          estimated via penalized likelihood)
frailtyPenal(Surv(time,status)~sex+age+cluster(id),
             n.knots=12,kappa1=1000,data=kidney,Frailty=FALSE)


# truncated data

# first, we create a hypothetical truncated data
kidney$tt0<-rep(0,nrow(kidney))
kidney$tt0[1:3]<-c(2,9,13)

# then, we fit the model
frailtyPenal(Surv(tt0,time,status)~sex+age+cluster(id),
             n.knots=12,kappa1=1000,data=kidney)


#stratified data. Let's use another dataset
data(readmission)
frailtyPenal(Surv(time,event)~as.factor(dukes)+cluster(id)+strata(sex),
             n.knots=10,kappa1=10000,kappa2=10000,data=readmission)

#Andersen-Gill counting-process approach with time-dependent covariate
frailtyPenal(Surv(t.start,t.stop,event)~as.factor(sex)+as.factor(dukes)+
             as.factor(charlson)+cluster(id),data=readmission, Frail=TRUE,
             n.knots=6,kappa1=100000,recurrentAG=TRUE)

# with the use of the cross validation approach, to find the smoothing parameter 
frailtyPenal(Surv(t.start,t.stop,event)~as.factor(sex)+as.factor(dukes)+
             as.factor(charlson)+cluster(id),data=readmission, Frail=TRUE,
             n.knots=6,kappa1=5000,recurrentAG=TRUE,cross.validation=TRUE)