BivPPL: Estimation of a Bivariate Correlated Frailty Model for...

In lilywang1988/BivPPL: Bivariate frailty model estimation

Description

Fit the bivariate frailty model of survival data using penalized partial log-likelihood. The estimation procedure is based on a Laplace approximation.

Usage

BivPPL(data, max.iter = 1000, eps = 1e-05, NP.max.iter = 500, NP.eps = 1e-06, huge = FALSE, independence = FALSE, D.initial = 0.5, raw = TRUE, check.log.lik = FALSE)

Arguments

data	Data input is a list object. If the input dataset is generated by gen.data, set raw=TRUE), or the input can be a transformed dataset by tocoxme to fit in the R package coxme for comparison (raw=FALSE). The function also permits any external data lists with variables: z1 and z2 for covariate matrices, x and y for follow-up times, and delta1 and delta2 for event status, for each event type respectively (where we still set raw=FALSE). Note that, correspondingly ordered and identical numbers of events are required for both event types. Thus one will need to sort the observations, trim or add trivial rows (follow-up time to be 0) to make them balanced.
max.iter	The maximum total number of iterations (inner and outer loops)
eps	Convergene criterion either for the parameter (check.log.lik=FALSE) or for the apprixmate marginal likelihood (check.log.lik=TRUE)
NP.max.iter	The maximum iterations of the Newton-Raphson in the inner loop
NP.eps	The covergence criteria of the Newton-Raphson in the inner loop
huge	logical; if TRUE it will sparsen the Hessian matrix to reduce computation time and memory usage; recommended when there are many cluster (e.g., n>500).
independence	logical; if TRUE it will restrict the off-digonal or the covariance entries of the estimated D matrix (variance-covariance matrix of the estimated bivariate frailty vector) to be 0. Mainly used to conduct likelihood ratio tests.
D.initial	The initial value for the diagonal entries of the D matrix, we let the initial D to be diagonal and the two variances are identical.
raw	logical; if TRUE, it will be generated from gen.data and thus every variables are recorded in a hierarchical list structure; if FALSE, it will detect whether it has
check.log.lik	logical; if FALSE, the convergence of the parameter estimation will be checked; if TRUE, the convergence of the approximate marginal likelihood will be checked.

Value

beta_hat	Regression parameter estimates
beta_ASE	Asymptotic standard error estimates for the regression parameter estimators
D_hat	Variance-covariance estimate for the bivariate frailties
r_hat	Estimated frailties
LogMargProb	Approximate marginal log-likelihood
num.iter	Number of iterations
dist	Distance of the last two iterations for convergence check

See Also

coxme

Examples

set.seed(100)
N<-100 # number of clusters
beta1 <- c(0.5,-0.3,0.5) # regression parameters for event type 1
beta2 <- c(0.8,-0.2,0.3) # regression parameters for event type 2
beta  <- c(beta1,beta2)
theta <- c(0.25,0.25,-0.125) # variance-covariance matrix for the bivariate frailty (denoted as D), it is a vector (D[1,1],D[2,2], D[1,2])
lambda01 <- 1
lambda02 <- 1
cen <- 10 # maximum censoring time
centype <- TRUE # fixed censoring time at cen

data <- gen.data(N,beta1,beta2,theta,lambda01,lambda02,c=cen,Ctype=centype)
ptm<-proc.time()
res <- BivPPL(data)
proc.time() - ptm
res$beta_hat
res$beta_ASE
res$D_hat

# fit the model assuming the independence between the two frailties
ptm<-proc.time()
res2 <- BivPPL(data,independence=T)
proc.time() - ptm
res2$beta_hat
res2$beta_ASE
res2$D_hat

# A likelihood ratio test
LRT<-2*(res$LogMargProb-res2$LogMargProb)
print(round(pchisq(abs(LRT),df=1,lower.tail = F),3))


# sparsen the Hessian matrix
ptm<-proc.time()
res3 <- BivPPL(data,huge=TRUE)
proc.time() - ptm
res3$beta_hat
res3$beta_ASE
res3$D_hat

lilywang1988/BivPPL documentation built on Aug. 9, 2019, 6:14 p.m.