Package 'depCensoring'

Description

This function estimates the bootstrap standard errors for the finite-dimensional model parameters and for the non-parametric cumulative hazard function. Parallel computing using foreach has been used to speed up the estimation of standard errors.

Usage

boot.fun(
  init,
  resData,
  X,
  W,
  lhat,
  cumL,
  dist,
  k,
  lb,
  ub,
  Obs.time,
  cop,
  n.boot,
  n.iter,
  ncore,
  eps
)

Arguments

Value

Bootstrap standard errors for parameter estimates and for estimated cumulative hazard function.

Description

This function estimates the bootstrap standard errors for the finite-dimensional model parameters and for the non-parametric cumulative hazard function under the assumption of independent censoring. Parallel computing using foreach has been used to speed up the computation.

Usage

boot.funI(
  init,
  resData,
  X,
  W,
  lhat,
  cumL,
  dist,
  k,
  lb,
  ub,
  Obs.time,
  n.boot,
  n.iter,
  ncore,
  eps
)

Arguments

Value

Bootstrap standard errors for parameter estimates and for estimated cumulative hazard function.

Description

This function estimates the bootstrap standard errors for the finite-dimensional model parameters and for the non-parametric transformation function. Parallel computing using foreach has been used to speed up the estimation of standard errors.

Usage

boot.nonparTrans(init, resData, X, W, n.boot, n.iter, eps)

Arguments

Value

Bootstrap standard errors for parameter estimates and for estimated cumulative hazard function.

Description

This function calculates a bivariate survival probability based on multivariate normal distribution.

Usage

Bvprob(lx, ly, rho)

Arguments

Description

This function transforms the parameters of the Cholesky de- composition to the covariance matrix, represented as a the row-wise con- catenation of the upper-triangular elements.

Usage

chol2par(par.chol1)

Arguments

Value

Covariance matrix corresponding to the provided Cholesky decomposition.

Description

This function transforms the parameters of the Cholesky de- composition to a covariance matrix element. This function is used in chol2par.R.

Usage

chol2par.elem(a, b, par.chol1)

Arguments

Value

Specified element of the covariance matrix resulting from the provided Cholesky decomposition.

Description

This function estimates phi function at fixed time point t

Usage

CompC(theta, t, X, W, ld, cop, dist)

Arguments

Description

Prepare initial values within the control arguments

Usage

control.arguments(maxit = 300, eps = 1e-06, trace = TRUE, ktau.inits = NULL)

Arguments

Description

The distribution function of the Archimedean copula

Usage

copdist.Archimedean(x, copfam, ktau, coppar = NULL)

Arguments

Description

The h-function of the copula

Usage

cophfunc(x, coppar, copfam, condvar = 1)

Arguments

Description

Convert the copula parameter the Kendall's tau

Usage

coppar.to.ktau(coppar, copfam)

Arguments

Description

This function implements the second step likelihood function of the competing risk model defined in Willems et al. (2024+).

Usage

cr.lik(
  n,
  s,
  Y,
  admin,
  cens.inds,
  M,
  Sigma,
  beta.mat,
  sigma.vct,
  rho.mat,
  theta.vct
)

Arguments

Value

Evaluation of the log-likelihood function

References

Willems et al. (2024+). Flexible control function approach under competing risks (in preparation).

Description

This function generates competing risk data that can be used in simulation studies.

Usage

dat.sim.reg.comp.risks(
  n,
  par,
  iseed,
  s,
  conf,
  Zbin,
  Wbin,
  type.cov,
  A.upper = 15
)

Arguments

Value

A generated data set

Description

This function defines the derivative of the transformation function that maps Cholesky parameters to the full covariance matrix.

Usage

dchol2par(par.chol1)

Arguments

Value

Derivative of the function that transforms the cholesky parameters to the full covariance matrix.

Description

This function defines the derivative of the transformation function that maps Cholesky parameters to a covariance matrix element. This function is used in dchol2par.R.

Usage

dchol2par.elem(k, q, a, b, par.chol1)

Arguments

Value

Derivative of the function that transforms the cholesky parameters to the specified element of the covariance matrix, evaluated at the specified arguments.

Description

This function computes distance between two vectors based on L2-norm

This function computes distance between two vectors based on L2-norm

Usage

Distance(b, a)

Distance(b, a)

Arguments

Description

Evaluates the derivative of the Yeo-Johnson transformation at the provided argument.

Usage

DYJtrans(y, theta)

Arguments

Value

The transformed value of y.

Description

This function estimates the control function for the endogenous variable based on the provided covariates. This function is called inside estimate.cmprsk.R.

Usage

estimate.cf(XandW, Z, Zbin, gammaest = NULL)

Arguments

Value

List containing the vector of values for the control function and the regression parameters of the first step.

Description

This function estimates the parameters in the competing risks model described in Willems et al. (2024+). Note that this model extends the model of Crommen, Beyhum and Van Keilegom (2024) and as such, this function also implements their methodology.

Usage

estimate.cmprsk(
  data,
  admin,
  conf,
  eoi.indicator.names = NULL,
  Zbin = NULL,
  inst = "cf",
  realV = NULL,
  eps = 0.001
)

Arguments

Value

A list of parameter estimates in the second stage of the estimation algorithm (hence omitting the estimates for the control function), as well as an estimate of their variance and confidence intervals.

References

Willems et al. (2024+). Flexible control function approach under competing risks (in preparation).

Crommen, G., Beyhum, J., and Van Keilegom, I. (2024). An instrumental variable approach under dependent censoring. Test, 33(2), 473-495.

Examples

n <- 200

# Set parameters
gamma <- c(1, 2, 1.5, -1)
theta <- c(0.5, 1.5)
eta1 <- c(1, -1, 2, -1.5, 0.5)
eta2 <- c(0.5, 1, 1, 3, 0)

# Generate exogenous covariates
x0 <- rep(1, n)
x1 <- rnorm(n)
x2 <- rbinom(n, 1, 0.5)

# Generate confounder and instrument
w <- rnorm(n)
V <- rnorm(n, 0, 2)
z <- cbind(x0, x1, x2, w) %*% gamma + V
realV <- z - (cbind(x0, x1, x2, w) %*% gamma)

# Generate event times
err <- MASS::mvrnorm(n, mu = c(0, 0), Sigma =
matrix(c(3, 1, 1, 2), nrow = 2, byrow = TRUE))
bn <- cbind(x0, x1, x2, z, realV) %*% cbind(eta1, eta2) + err
Lambda_T1 <- bn[,1]; Lambda_T2 <- bn[,2]
x.ind = (Lambda_T1>0)
y.ind <- (Lambda_T2>0)
T1 <- rep(0,length(Lambda_T1))
T2 <- rep(0,length(Lambda_T2))
T1[x.ind] = ((theta[1]*Lambda_T1[x.ind]+1)^(1/theta[1])-1)
T1[!x.ind] = 1-(1-(2-theta[1])*Lambda_T1[!x.ind])^(1/(2-theta[1]))
T2[y.ind] = ((theta[2]*Lambda_T2[y.ind]+1)^(1/theta[2])-1)
T2[!y.ind] = 1-(1-(2-theta[2])*Lambda_T2[!y.ind])^(1/(2-theta[2]))
# Generate adminstrative censoring time
C <- runif(n, 0, 40)

# Create observed data set
y <- pmin(T1, T2, C)
delta1 <- as.numeric(T1 == y)
delta2 <- as.numeric(T2 == y)
da <- as.numeric(C == y)
data <- data.frame(cbind(y, delta1, delta2, da, x0, x1, x2, z, w))
colnames(data) <- c("y", "delta1", "delta2", "da", "x0", "x1", "x2", "z", "w")

# Estimate the model
admin <- TRUE                # There is administrative censoring in the data.
conf <- TRUE                 # There is confounding in the data (z)
eoi.indicator.names <- NULL  # We will not impose that T1 and T2 are independent
Zbin <- FALSE                # The confounding variable z is not binary
inst <- "cf"                 # Use the control function approach
# Since we don't use the oracle estimator, this argument is ignored anyway
realV <- NULL
estimate.cmprsk(data, admin, conf, eoi.indicator.names, Zbin, inst, realV)

Description

This function allows to estimate the dependency parameter along all other model parameters. First, estimates the cumulative hazard function, and then at the second stage it estimates other model parameters assuming that the cumulative hazard function is known. The details for implementing the dependent censoring methodology can be found in Deresa and Van Keilegom (2024).

Usage

fitDepCens(
  resData,
  X,
  W,
  cop = c("Frank", "Gumbel", "Normal"),
  dist = c("Weibull", "lognormal"),
  start = NULL,
  n.iter = 50,
  bootstrap = TRUE,
  n.boot = 150,
  ncore = 7,
  eps = 1e-04
)

Arguments

Value

This function returns a fit of dependent censoring model; parameter estimates, estimate of the cumulative hazard function, bootstrap standard errors for finite-dimensional parameters, the nonparametric cumulative hazard function, etc.

References

Deresa and Van Keilegom (2024). Copula based Cox proportional hazards models for dependent censoring, Journal of the American Statistical Association, 119:1044-1054.

Examples

# Toy data example to illustrate implementation
n = 300
beta = c(0.5)
lambd = 0.35
eta = c(0.9,0.4)
X = cbind(rbinom(n,1,0.5))
W = cbind(rep(1,n),rbinom(n,1,0.5))
# generate dependency structure from Frank
frank.cop <- copula::frankCopula(param = 5,dim = 2)
U = copula::rCopula(n,frank.cop)
T1 = (-log(1-U[,1]))/(lambd*exp(X*beta))           # Survival time
T2 = (-log(1-U[,2]))^(1.1)*exp(W%*%eta)            # Censoring time
A = runif(n,0,15)                                  # administrative censoring time
Z = pmin(T1,T2,A)
d1 = as.numeric(Z==T1)
d2 = as.numeric(Z==T2)
resData = data.frame("Z" = Z,"d1" = d1, "d2" = d2)   # should be data frame
colnames(W) <- c("ones","cov1")
colnames(X) <- "cov.surv"

# Fit dependent censoring model

fit <- fitDepCens(resData = resData, X = X, W = W, bootstrap = FALSE)

# parameter estimates

fit$parameterEstimates

# summary fit results
summary(fit)

# plot cumulative hazard vs time

plot(fit$observedTime, fit$cumhazardFunction, type = "l",xlab = "Time",
ylab = "Estimated cumulative hazard function")

Description

This function allows to estimate all model parameters under the assumption of independent censoring. First, estimates the cumulative hazard function, and then at the second stage it estimates other model parameters assuming that the cumulative hazard is known.

Usage

fitIndepCens(
  resData,
  X,
  W,
  dist = c("Weibull", "lognormal"),
  start = NULL,
  n.iter = 50,
  bootstrap = TRUE,
  n.boot = 150,
  ncore = 7,
  eps = 1e-04
)

Arguments

Value

This function returns a fit of independent censoring model; parameter estimates, estimate of the cumulative hazard function, bootstrap standard errors for finite-dimensional parameters, the nonparametric cumulative hazard function, etc.

Examples

# Toy data example to illustrate implementation
n = 300
beta = c(0.5)
lambd = 0.35
eta = c(0.9,0.4)
X = cbind(rbinom(n,1,0.5))
W = cbind(rep(1,n),rbinom(n,1,0.5))
# generate dependency structure from Frank
frank.cop <- copula::frankCopula(param = 5,dim = 2)
U = copula::rCopula(n,frank.cop)
T1 = (-log(1-U[,1]))/(lambd*exp(X*beta))                  # Survival time'
T2 = (-log(1-U[,2]))^(1.1)*exp(W%*%eta)                   # Censoring time
A = runif(n,0,15)                                         # administrative censoring time
Z = pmin(T1,T2,A)
d1 = as.numeric(Z==T1)
d2 = as.numeric(Z==T2)
resData = data.frame("Z" = Z,"d1" = d1, "d2" = d2)     # should be data frame

colnames(W) <- c("ones","cov1")
colnames(X) <- "cov.surv"

# Fit independent censoring model

fitI <- fitIndepCens(resData = resData, X = X, W = W, bootstrap = FALSE)

# parameter estimates

fitI$parameterEstimates

# summary fit results
summary(fitI)

# plot cumulative hazard vs time

 plot(fitI$observedTime, fitI$cumhazardFunction, type = "l",xlab = "Time",
 ylab = "Estimated cumulative hazard function")

Description

The generator function of the Archimedean copula

Usage

generator.Archimedean(x, coppar, copfam, inverse = FALSE)

Arguments

Description

Computes the inverse Yeo-Johnson transformation of the provided argument.

Usage

IYJtrans(y, theta)

Arguments

Value

The transformed value of y.

Description

Calculate the kernel function

Usage

Kernel(u, name = "Gaussian")

Arguments

Description

Convert the Kendall's tau into the copula parameter

Usage

ktau.to.coppar(ktau, copfam)

Arguments

Description

Loglikehood function under independent censoring

Usage

LikCopInd(theta, resData, X, W, lhat, cumL, dist)

Arguments

Value

Maximized log-likelihood value

Description

Calculate the likelihood function for the fully parametric joint distribution

Usage

Likelihood.Parametric(param, yobs, delta, copfam, margins, cure = FALSE)

Arguments

Description

Calculate the profiled likelihood function with kernel smoothing

Usage

Likelihood.Profile.Kernel(param, yobs, delta, copfam, margins, cure = FALSE)

Arguments

Description

Solve the profiled likelihood function

Usage

Likelihood.Profile.Solve(
  yobs,
  delta,
  copfam,
  margins,
  ktau.init,
  parapar.init,
  cure,
  curerate.init,
  constraints,
  maxit,
  eps
)

Arguments

Description

Calculate the semiparametric version of profiled likelihood function

Usage

Likelihood.Semiparametric(
  param,
  Syobs,
  yobs,
  delta,
  copfam,
  margins,
  cure = FALSE
)

Arguments

Description

This function defines the log-likelihood used to estimate the second step in the competing risks extension of the model described in Willems et al. (2024+).

Usage

LikF.cmprsk(par, data, admin, conf, cf)

Arguments

Value

Log-likelihood evaluation of the second step.

References

Willems et al. (2024+). Flexible control function approach under competing risks (in preparation).

Description

This function parametrizes the covariance matrix using its Cholesky decomposition, so that optimization of the likelihood can be done based on this parametrization, and positive-definiteness of the covariance matrix is guaranteed at each step of the optimization algorithm.

Usage

likF.cmprsk.Cholesky(par.chol, data, admin, conf, cf, eps = 0.001)

Arguments

Value

Log-likelihood evaluation of the second step.

Description

This function defines the maximum likelihood used to estimate the control function in the case of a continuous endogenous variable.

Usage

LikGamma1(gamma, Z, M)

Arguments

Value

Returns the log-likelihood function corresponding to the data, evaluated at the point gamma.

Description

This function defines the maximum likelihood used to estimate the control function in the case of a binary endogenous variable.

Usage

LikGamma2(gamma, Z, M)

Arguments

Value

Returns the log-likelihood function corresponding to the data, evaluated at the point gamma.

Description

This function defines the log-likelihood used in estimating the second step in the competing risks extension of the model described in Willems et al. (2024+). The results of this function will serve as starting values for subsequent optimizations (LikI.comprsk.R and LikF.cmprsk.R)

Usage

LikI.bis(par, data, admin, conf, cf)

Arguments

Value

Starting values for subsequent optimization function used in the second step of the estimation procedure.

References

Willems et al. (2024+). Flexible control function approach under competing risks (in preparation).

Description

This function defines the log-likelihood used to estimate the second step in the competing risks extension assuming independence of some of the competing risks in the model described in Willems et al. (2024+).

Usage

LikI.cmprsk(par, data, eoi.indicator.names, admin, conf, cf)

Arguments

Value

Log-likelihood evaluation for the second step in the esimation procedure.

References

Willems et al. (2024+). Flexible control function approach under competing risks (in preparation).

Description

This function does the same as LikI.cmprsk (in fact, it even calls said function), but it parametrizes the covariance matrix using its Cholesky decomposition in order to guarantee positive definiteness. This function is never used, might not work and could be deleted.

Usage

LikI.cmprsk.Cholesky(
  par.chol,
  data,
  eoi.indicator.names,
  admin,
  conf,
  cf,
  eps = 0.001
)

Arguments

Value

Log-likelihood evaluation for the second step in the estimation procedure.

Description

This function defines the 'full' likelihood of the model. Specifically, it includes the estimation of the control function in the computation of the likelihood. This function is used in the estimation of the variance of the estimates (variance.cmprsk.R).

Usage

likIFG.cmprsk.Cholesky(
  parhatG,
  data,
  eoi.indicator.names,
  admin,
  conf,
  Zbin,
  inst
)

Arguments

Value

Full model log-likelihood evaluation.

Description

Computes the logarithm of a number.

Usage

log_transform(y)

Arguments

Value

This function returns the logarithm of the provided argument y if it is greater than zero. If y is smaller than zero, it will return 0.

Description

This likelihood function is maximized to estimate the model parameters under the Clayton copula.

Usage

loglike.clayton.unconstrained(para, Y, Delta, Dist.T, Dist.C)

Arguments

Value

Maximized log-likelihood value.

Description

This likelihood function is maximized to estimate the model parameters under the Frank copula.

Usage

loglike.frank.unconstrained(para, Y, Delta, Dist.T, Dist.C)

Arguments

Value

Maximized log-likelihood value.

Description

This likelihood function is maximized to estimate the model parameters under the Gaussian copula.

Usage

loglike.gaussian.unconstrained(para, Y, Delta, Dist.T, Dist.C)

Arguments

Value

Maximized log-likelihood value.

Description

This likelihood function is maximized to estimate the model parameters under the Gumbel copula.

Usage

loglike.gumbel.unconstrained(para, Y, Delta, Dist.T, Dist.C)

Arguments

Value

Maximized log-likelihood value.

Description

This likelihood function is maximized to estimate the model parameters under the independence copula.

Usage

loglike.indep.unconstrained(para, Y, Delta, Dist.T, Dist.C)

Arguments

Value

Maximized log-likelihood value.

Description

Change hazard and cumulative hazard to long format

Usage

Longfun(Z, T1, lhat, Lhat)

Arguments

Description

Change a nonparametric transformation function to long format

Usage

LongNPT(Z, T1, H)

Arguments

Description

This function allows to estimate the dependency parameter along all other model parameters. First, estimates a non-parametric transformation function, and then at the second stage it estimates other model parameters assuming that the non-parametric function is known. The details for implementing the dependent censoring methodology can be found in Deresa and Van Keilegom (2021).

Usage

NonParTrans(
  resData,
  X,
  W,
  start = NULL,
  n.iter = 15,
  bootstrap = FALSE,
  n.boot = 50,
  eps = 0.001
)

Arguments

Value

This function returns a fit of a semiparametric transformation model; parameter estimates, estimate of the non-parametric transformation function, bootstrap standard errors for finite-dimensional parameters, the nonparametric cumulative hazard function, etc.

References

Deresa, N. and Van Keilegom, I. (2021). On semiparametric modelling, estimation and inference for survival data subject to dependent censoring, Biometrika, 108, 965–979.

Examples

# Toy data example to illustrate implementation
n = 300
beta = c(0.5, 1); eta = c(1,1.5); rho = 0.70
sigma = matrix(c(1,rho,rho,1),ncol=2)
err = MASS::mvrnorm(n, mu = c(0,0) , Sigma=sigma)
err1 = err[,1]; err2 = err[,2]
x1 = rbinom(n,1,0.5); x2 = runif(n,-1,1)
X = matrix(c(x1,x2),ncol=2,nrow=n); W = X   # data matrix
T1 = X%*%beta+err1
C =  W%*%eta+err2
T1 = exp(T1); C = exp(C)
A = runif(n,0,8); Y = pmin(T1,C,A)
d1 = as.numeric(Y==T1)
d2 = as.numeric(Y==C)
resData = data.frame("Z" = Y,"d1" = d1, "d2" = d2)   # should be data frame
colnames(X) = c("X1", "X2")
colnames(W) = c("W1","W2")

#  Bootstrap is false by default
output = NonParTrans(resData = resData, X = X, W = W, n.iter = 2)
output$parameterEstimates

Description

Estimates the model parameters by maximizing the log-likelihood.

Usage

optimlikelihood(Y, Delta, Copula, Dist.T, Dist.C, start)

Arguments

Value

A list containing the minimized negative log-likelihood using the independence copula model, the estimated parameter values for the model with the independence copula, the minimized negative log-likelihood using the specified copula model and the estimated parameter values for the model with the specified copula.

Description

Obtain the value of the density function

Usage

parafam.d(x, parameter, distribution, truncation = NULL)

Arguments

Description

Obtain the value of the distribution function

Usage

parafam.p(x, parameter, distribution, truncation = NULL)

Arguments

Description

Obtain the adjustment value of truncation

Usage

parafam.trunc(truncation, parameter, distribution)

Arguments

Description

Note that it is not assumed that the association parameter of the copula function is known, unlike most other papers in the literature. The details for implementing the methodology can be found in Czado and Van Keilegom (2023).

Usage

ParamCop(Y, Delta, Copula, Dist.T, Dist.C, start = c(1, 1, 1, 1))

Arguments

Value

A table containing the minimized negative log-likelihood using the independence copula model, the estimated parameter values for the model with the independence copula, the minimized negative log-likelihood using the specified copula model and the estimated parameter values for the model with the specified copula.

References

Czado and Van Keilegom (2023). Dependent censoring based on parametric copulas. Biometrika, 110(3), 721-738.

Examples

tau = 0.75
Copula = "frank"
Dist.T = "weibull"
Dist.C = "lnorm"
par.T = c(2,1)
par.C = c(1,2)
n=1000

simdata<-TCsim(tau,Copula,Dist.T,Dist.C,par.T,par.C,n)

Y = simdata[[1]]
Delta = simdata[[2]]

ParamCop(Y,Delta,Copula,Dist.T,Dist.C)

Description

Generate constraints of parameters

Usage

Parameters.Constraints(copfam, margins, cure)

Arguments

Description

Computes a given power of a number.

Usage

power_transform(y, pw)

Arguments

Value

This function returns the result of raising y to the power pw when y > 0. Otherwise, it will return 1.

Description

The PseudoL function is maximized in order to estimate the finite dimensional model parameters, including the dependency parameter. This function assumes that the cumulative hazard function is known.

Usage

PseudoL(theta, resData, X, W, lhat, cumL, cop, dist)

Arguments

Value

maximized log-likelihood value

Description

This function computes the score vectors and the Jacobean matrix for finite model parameters.

Usage

ScoreEqn(theta, resData, X, W, H)

Arguments

Description

Function to indicate position of t in observed survival time

Usage

SearchIndicate(t, T1)

Arguments

Description

This function estimates the nonparametric transformation function H when the survival time and censoring time are dependent given covariates. The estimating equation of H was derived based on the martingale ideas. More details about the derivation of a nonparmaetric estimator of H and its estimation algorithm can be found in Deresa and Van Keilegom (2021).

Usage

SolveH(theta, resData, X, W)

Arguments

Value

Returns the estimated transformation function H for a fixed value of parameters theta.

References

Deresa, N. and Van Keilegom, I. (2021). On semiparametric modelling, estimation and inference for survival data subject to dependent censoring, Biometrika, 108, 965–979.

Description

This function obtains an estimating equation of H at the first observed survival time t1.

Usage

SolveHt1(Ht1, Z, nu, t, X, W, theta)

Arguments

Description

This function estimates the cumulative hazard function of survival time (T) under dependent censoring (C). The estimation makes use of the estimating equations derived based on martingale ideas.

Usage

SolveL(
  theta,
  resData,
  X,
  W,
  cop = c("Frank", "Gumbel", "Normal"),
  dist = c("Weibull", "lognormal")
)

Arguments

Value

This function returns an estimated hazard function, cumulative hazard function and distinct observed survival times;

Examples

n = 200
beta = c(0.5)
lambd = 0.35
eta = c(0.9,0.4)
X = cbind(rbinom(n,1,0.5))
W = cbind(rep(1,n),rbinom(n,1,0.5))
frank.cop <- copula::frankCopula(param = 5,dim = 2)
U = copula::rCopula(n,frank.cop)
T1 = (-log(1-U[,1]))/(lambd*exp(X*beta))         # Survival time'
T2 = (-log(1-U[,2]))^(1.1)*exp(W%*%eta)          # Censoring time
A = runif(n,0,15)                                # administrative censoring time
Z = pmin(T1,T2,A)
d1 = as.numeric(Z==T1)
d2 = as.numeric(Z==T2)
resData = data.frame("Z" = Z,"d1" = d1, "d2" = d2)
theta = c(0.3,1,0.3,1,2)

# Estimate cumulative hazard function
cumFit <- SolveL(theta, resData,X,W)
cumhaz = cumFit$cumhaz
time = cumFit$times

# plot hazard vs time

plot(time, cumhaz, type = "l",xlab = "Time",
ylab = "Estimated cumulative hazard function")

Description

This function estimates the cumulative hazard function of survival time (T) under the assumption of independent censoring. The estimating equation is derived based on martingale ideas.

Usage

SolveLI(theta, resData, X)

Arguments

Value

This function returns an estimated hazard function, cumulative hazard function and distinct observed survival times;

Examples

n = 200
beta = c(0.5)
lambd = 0.35
eta = c(0.9,0.4)
X = cbind(rbinom(n,1,0.5))
W = cbind(rep(1,n),rbinom(n,1,0.5))
frank.cop <- copula::frankCopula(param = 5,dim = 2)
U = copula::rCopula(n,frank.cop)
T1 = (-log(1-U[,1]))/(lambd*exp(X*beta))           # Survival time'
T2 = (-log(1-U[,2]))^(1.1)*exp(W%*%eta)            # Censoring time
A = runif(n,0,15)                                  # administrative censoring time
Z = pmin(T1,T2,A)
d1 = as.numeric(Z==T1)
d2 = as.numeric(Z==T2)
resData = data.frame("Z" = Z,"d1" = d1, "d2" = d2)
theta = c(0.3,1,0.3,1)

# Estimate cumulative hazard function

cumFit_ind <- SolveLI(theta, resData,X)

cumhaz = cumFit_ind$cumhaz
time = cumFit_ind$times

# plot hazard vs time

plot(time, cumhaz, type = "l",xlab = "Time",
ylab = "Estimated cumulative hazard function")

Description

This function estimates the model parameters

Usage

SolveScore(theta, resData, X, W, H, eps = 0.001)

Arguments

Description

Summary of depCensoringFit object

Usage

## S3 method for class 'depFit'
summary(object, ...)

Arguments

Value

Summary of dependent censoring model fit in the form of table

Description

Summary of indepCensoringFit object

Usage

## S3 method for class 'indepFit'
summary(object, ...)

Arguments

Value

Summary of independent censoring model fit in the form of table

Description

Provide semiparametric approaches that can be used to model right-censored survival data under dependent censoring (without covariates). The copula-based approach is adopted and there is no need to explicitly specify the association parameter. One of the margins can be modeled nonparametrically. As a byproduct, both marginal distributions of survival and censoring times can be considered as fully parametric. The existence of a cured fraction concerning survival time can also be taken into consideration.

Usage

SurvDC(
  yobs,
  delta,
  tm = NULL,
  copfam = "frank",
  margins = list(survfam = NULL, censfam = "lnorm"),
  cure = FALSE,
  Var = list(do = TRUE, nboot = 200, level = 0.05),
  control = list(maxit = 300, eps = 1e-06, trace = TRUE, ktau.inits = NULL)
)

Arguments

Details

This unified function provide approaches that can be used to model right-censored survival data under dependent censoring (without covariates). Various specifications of marginal distributions can be considered by choosing different combinations of the provided arguments. Generally speaking, the following two scenarios are what we mainly focused on:

nonparametric survival margin and parametric censoring margin (without cure)

survfam = NULL, censfam is not NULL and cure = FALSE.

nonparametric survival margin and parametric censoring margin (with cure)

survfam = NULL, censfam is not NULL and cure = TRUE.

As byproducts, several other scenarios (the distribution of the underlying survival time is not nonparametric but fully parametric) can also be considered by this R function:

parametric survival and censoring margins (without cure)

both survfam and censfam are not NULL and cure = FALSE.

parametric survival and censoring margins (with cure)

both survfam and censfam are not NULL and cure = TRUE.

parametric survival margin and nonparametric censoring margin (without cure)

survfam is not NULL, censfam = NULL and cure = FALSE.

Furthermore, one might expect that a scenario with "parametric survival margin and nonparametric censoring margin (with cure)" can also be included. Indeed, it can be done based on: survfam is not NULL, censfam = NULL and cure = TRUE. However, from a theoretical perspective of view, whether this type of modeling is reasonable or not still needs further investigations.

We emphasize that the first scenario (in byproducts) has also be considered in another function of this package. Specifically, the scenario of "parametric survival margin and nonparametric censoring margin (without cure)" can be fitted based on ParamCop(). However, the default joint modeling of survival and censoring times are based on their joint survival function in line with the semiparametric case (instead of modeling joint distribution function directly as in Czado and Van Keilegom (2023) <doi:10.1093/biomet/asac067>), but the idea of estimation methodology are exactly the same.

@references Czado and Van Keilegom (2023). Dependent censoring based on parametric copulas. Biometrika, 110(3), 721-738. @references Delhelle and Van Keilegom (2024). Copula based dependent censoring in cure models. TEST (to appear). @references Ding and Van Keilegom (2024). Semiparametric estimation of the survival function under dependent censoring (in preparation).

Value

A list of fitted results is returned. Within this outputted list, the following elements can be found:

probs

survival probabilities of the survial margin at tm.

ktau

Kendall's tau.

parapar

estimation of all parameters (except Kendall's tau) contained in the parametric part.

GoF

goodness-of-test results.

curerate

cure rate. If cure = FALSE, it is NULL.

Examples

#----------------------------------------------------------#
# Basic preparations before running subsequent examples ####
#----------------------------------------------------------#

# library necessary packages

#------------------------------------------------------------------------#
# simulated data from Frank copula log-Normal margins (without cure)
#------------------------------------------------------------------------#

# generate the simulated data

# - the sample size of the generated data
n <- 1000

# information on the used copula
copfam.true <- "frank"
ktau.true <- 0.5
coppar.true <- 5.74

# parameters of the underlying log-normal marginal distributions
survpar.true <- c(2.20,1.00)
censpar.true <- c(2.20,0.25)

# - true underlying survival and censoring times
set.seed(1)
u.TC <- copula::rCopula(
  n        = n,
  copula   = copula::archmCopula(
    family = copfam.true,
    param  = coppar.true,
    dim    = 2
  )
)
yobs.T <- qlnorm(1-u.TC[,1],survpar.true[1],survpar.true[2])
yobs.C <- qlnorm(1-u.TC[,2],censpar.true[1],censpar.true[2])

# observations
yobs  <- pmin(yobs.T,yobs.C)
delta <- as.numeric(yobs.T<=yobs.C)
cat("censoring rate is", mean(1-delta))

# model the data under different scenarios

# scenario 1: nonparametric survival margin and parametric censoring margin
set.seed(1)
sol.scenario1 <- SurvDC(
  yobs    = yobs,
  delta   = delta,
  tm      = quantile(yobs, c(0.25,0.50,0.75)),
  copfam  = copfam.true,
  margins = list(survfam = NULL, censfam = "lnorm"),
  Var     = list(do = FALSE, nboot = 50)
)
sol.scenario1$probs
sol.scenario1$ktau
sol.scenario1$parapar

# scenario 2: parametric survival and censoring margins
set.seed(1)
sol.scenario2 <- SurvDC(
  yobs    = yobs,
  delta   = delta,
  tm      = quantile(yobs, c(0.25,0.50,0.75)),
  copfam  = copfam.true,
  margins = list(survfam = "lnorm", censfam = "lnorm"),
  Var     = list(do = FALSE, nboot = 50)
)
sol.scenario2$probs
sol.scenario2$ktau
sol.scenario2$parapar

# scenario 3: parametric survival margin and nonparametric censoring margin
set.seed(1)
sol.scenario3 <- SurvDC(
  yobs    = yobs,
  delta   = delta,
  tm      = quantile(yobs, c(0.25,0.50,0.75)),
  copfam  = copfam.true,
  margins = list(survfam = "lnorm", censfam = NULL),
  Var     = list(do = FALSE, nboot = 50)
)
sol.scenario3$probs
sol.scenario3$ktau
sol.scenario3$parapar

#------------------------------------------------------------------------
# simulated data from Frank copula log-Normal margins (with cure)
#------------------------------------------------------------------------

# generate the simulated data

#  true underlying cure rate
curerate.true <- 0.2

# true underlying survival and censoring times
set.seed(1)
u.TC <- copula::rCopula(
  n        = n,
  copula   = copula::archmCopula(
    family = copfam.true,
    param  = coppar.true,
    dim    = 2
  )
)
yobs.T <- sapply(u.TC[,1],function(uT){
  if(uT<=curerate.true){ val <- Inf }else{
    val <- EnvStats::qlnormTrunc((1-uT)/(1-curerate.true),survpar.true[1],survpar.true[2],0,15)
  }
  return(val)
})
yobs.C <- qlnorm(1-u.TC[,2],censpar.true[1],censpar.true[2])
cat("cure rate is",mean(yobs.T==Inf))

#  observations
yobs  <- pmin(yobs.T,yobs.C)
delta <- as.numeric(yobs.T<=yobs.C)
cat("censoring rate is",mean(1-delta))

# model the data under different scenarios (with cure)

# scenario 4: parametric survival and censoring margins
set.seed(1)
sol.scenario4 <- SurvDC(
  yobs    = yobs,
  delta   = delta,
  tm      = quantile(yobs, c(0.25,0.50,0.75)),
  copfam  = copfam.true,
  margins = list(survfam = "lnorm", censfam = "lnorm"),
  Var     = list(do = FALSE, nboot = 50),
  cure    = TRUE
)
sol.scenario4$probs
sol.scenario4$ktau
sol.scenario4$parapar
sol.scenario4$curerate

# scenario 5: nonparametric survival margin and parametric censoring margin
set.seed(1)
sol.scenario5 <- SurvDC(
  yobs    = yobs,
  delta   = delta,
  tm      = quantile(yobs, c(0.25,0.50,0.75)),
  copfam  = copfam.true,
  margins = list(survfam = NULL, censfam = "lnorm"),
  Var     = list(do = FALSE, nboot = 50),
  cure    = TRUE
)
sol.scenario5$probs
sol.scenario5$ktau
sol.scenario5$parapar
sol.scenario5$curerate

Description

Calculate the goodness-of-fit test statistic

Usage

SurvDC.GoF(
  yobs,
  delta,
  copfam,
  margins,
  ktau,
  parapar,
  cure = FALSE,
  curerate = NULL
)

Arguments

Description

Estimated survival function based on copula-graphic estimator (Archimedean copula only)

Usage

SurvFunc.CG(tm = NULL, yobs, delta, copfam, ktau, coppar = NULL)

Arguments

Description

Estimated survival function based on Kaplan-Meier estimator

Usage

SurvFunc.KM(tm = NULL, yobs, delta, type = "right")

Arguments

Description

Maximum likelihood estimator for a given parametric distribution

Usage

SurvMLE(
  yobs,
  delta,
  distribution,
  truncation = NULL,
  cure = FALSE,
  maxit = 300
)

Arguments

Description

Likelihood for a given parametric distribution

Usage

SurvMLE.Likelihood(
  param,
  yobs,
  delta,
  distribution,
  truncation = NULL,
  cure = FALSE
)

Arguments

Description

Generates the follow-up time and censoring indicator according to the specified model.

Usage

TCsim(
  tau = 0,
  Copula = "frank",
  Dist.T = "lnorm",
  Dist.C = "lnorm",
  par.T = c(0, 1),
  par.C = c(0, 1),
  n = 10000
)

Arguments

Value

A list containing the generated follow-up times and censoring indicators.

Examples

tau = 0.5
Copula = "gaussian"
Dist.T = "lnorm"
Dist.C = "lnorm"
par.T = c(1,1)
par.C = c(2,2)
n=1000

simdata <- TCsim(tau,Copula,Dist.T,Dist.C,par.T,par.C,n)
Y = simdata[[1]]
Delta = simdata[[2]]
hist(Y)
mean(Delta)

Description

Checks the required preconditions of the data and possibly restructures the data.

Usage

uniformize.data(
  data,
  admin = FALSE,
  conf = FALSE,
  comp.risks = FALSE,
  Zbin = NULL,
  Wbin = NULL
)

Arguments

Value

Returns the uniformized data set.

Description

This function computes the variance of the estimates computed by the 'estimate.cmprsk.R' function.

Usage

variance.cmprsk(
  parhatc,
  gammaest,
  data,
  admin,
  conf,
  inst,
  cf,
  eoi.indicator.names,
  Zbin,
  use.chol,
  n.trans,
  totparl
)

Arguments

Value

Variance estimates of the provided vector of estimated parameters.

Description

Computes the Yeo-Johnson transformation of the provided argument.

Usage

YJtrans(y, theta)

Arguments

Value

The transformed value of y.