Fit the Additive Trend Filtering Cox Model

Description

This package is called tfCox or trend filtering for Cox model, which is proposed in Jiacheng Wu & Daniela Witten (2019) Flexible and Interpretable Models for Survival Data, Journal of Computational and Graphical Statistics, DOI: 10.1080/10618600.2019.1592758. It provides an approach to fit additive Cox model in which each component function is estimated to be piecewise polynomial with adaptively-chosen knots.

Function tfCox fits the trend filtering Cox model for a range of tuning parameters. Function cv_tfCox returns the optimal tuning parameter selected by K-fold cross validation.

Details

The package includes the following functions: tfCox, cv_tfCox, plot.tfCox, plot.cv_tfCox, predict.tfCox, summary.tfCox, summary.cv_tfCox, sim_dat, plot.sim_dat.

Author(s)

Jiacheng Wu Maintainer: Jiacheng Wu <wujiacheng1992@gmail.com>

References

Jiacheng Wu & Daniela Witten (2019) Flexible and Interpretable Models for Survival Data, Journal of Computational and Graphical Statistics, DOI: 10.1080/10618600.2019.1592758

Fit Trend Filtering Cox model and Choose Tuning Parameter via K-Fold Cross-Validation

Description

Fit additive trend filtering Cox model where each component function is estimated to be piecewise constant or polynomial. Tuning parameter is selected via k-fold cross-validation.

Usage

cv_tfCox(dat, ord=0, alpha=1, discrete=NULL, lambda.seq=NULL,
lambda.min.ratio=0.01, n.lambda=30, n.fold=5, seed=NULL, tol=1e-6,
niter=1000, stepSize=25,  backtracking=0)

Arguments

Details

Note that cv_tfCox does not cross-validate over alpha, and alpha should be provided. However, if the user would like to cross-validate over alpha, then cv_tfCox should be called multiple times for different values of alpha and the same seed. This ensures that the cross-validation folds (fold) remain the same for the different values of alpha. See the example below for details.

Value

An object with S3 class "cv_tfCox".

Author(s)

Jiacheng Wu

References

Jiacheng Wu & Daniela Witten (2019) Flexible and Interpretable Models for Survival Data, Journal of Computational and Graphical Statistics, DOI: 10.1080/10618600.2019.1592758

See Also

summary.cv_tfCox, plot.cv_tfCox, tfCox

Examples

#generate data
set.seed(123)
dat = sim_dat(n=100, zerof=0, scenario=1)

#fit piecewise constant functions
#cross-validation to choose the tuning parameter lambda with fixed alpha=1
cv = cv_tfCox(dat, ord=0, alpha=1, n.fold=2, seed=123)
plot(cv, showSE=TRUE)

Calculate the negative log likelihood from Cox model.

Description

Calculate the negative log likelihood from Cox model from the estimated coefficient matrix theta.

Usage

negloglik(dat, theta)

Arguments

Author(s)

Jiacheng Wu

References

Jiacheng Wu & Daniela Witten (2019) Flexible and Interpretable Models for Survival Data, Journal of Computational and Graphical Statistics, DOI: 10.1080/10618600.2019.1592758

See Also

predict_best_lambda, tfCox_choose_lambda

Examples

#generate training and testing data
dat = sim_dat(n=100, zerof=0, scenario=1)
test_dat = sim_dat(n=100, zerof=0, scenario=1)

#choose the optimal tuning parameter
cv = tfCox_choose_lambda(dat, test_dat, ord=0, alpha=1)
plot(cv$lam_seq, cv$loss)

#optimal tuning parameter
cv$best_lambda

#predict the coefficients of testing covariates from the optimal tuning parameter
#from tfCox_choose_lambda object. 
theta_hat = predict_best_lambda(cv, test_dat$X)

#calculate the loss in the testing data based on the estimated coefficients theta
negloglik(test_dat, theta_hat)

Plots Cross-Validation Curve for Object of Class "cv_tfCox"

Description

This function plots the cross-validation curve for models fitted by a range of tuning parameter lambda using cv_tfCox. The cross-validation error with +/-1 standard error is plotted for each value of lambda. The dotted vertical line indicates the chosen lambda corresponding to the minimum cross-validation error.

Usage

## S3 method for class 'cv_tfCox'
plot(x, showSE=F, ...)

Arguments

Author(s)

Jiacheng Wu

See Also

cv_tfCox

Examples

#generate data
set.seed(123)
dat = sim_dat(n=100, zerof=0, scenario=1)

#fit piecewise constant functions
#cross-validation to choose the tuning parameter lambda with fixed alpha=1
cv = cv_tfCox(dat, ord=0, alpha=1, n.fold=2, seed=123)
plot(cv, showSE=TRUE)

Plot the true covariate effects

Description

This function plots the functional form of covariate effects in four simulation scenarios.

Usage

## S3 method for class 'sim_dat'
plot(x, which.predictor = NULL, n.plot = 4, ...)

Arguments

Author(s)

Jiacheng Wu

See Also

sim_dat

Examples

#generate data
set.seed(123)
dat = sim_dat(n=100, zerof=0, scenario=1)
#plot X versus the true theta
plot.sim_dat(dat)

Plot Fitted Functions from Class "tfCox"

Description

This function plots the fitted functions from a model estimated by tfCox.

Usage

## S3 method for class 'tfCox'
plot(x, which.lambda=1, which.predictor = NULL, n.plot = 4, ...)

Arguments

Author(s)

Jiacheng Wu

See Also

tfCox

Examples

#generate data
set.seed(123)
dat = sim_dat(n=100, zerof=0, scenario=1)

fit = tfCox(dat, ord=0, alpha=1, lambda.seq=0.04)
plot(fit, n.plot=4)

Predict for a New Covariate Matrix and fit from tfCox

Description

This function makes predictions from a specified covariate matrix for a fit of the class "tfCox".

Usage

## S3 method for class 'tfCox'
predict(object, newX, which.lambda=1, ...)

Arguments

Details

Prediction for the new data point is implemented by constant or linear interpolation. 0th order trend filtering will have constant interpolation, and 1th or higher order trend filtering will have linear interpolation.

Value

A n x p matrix containing the fitted theta values.

Author(s)

Jiacheng Wu

See Also

tfCox

Predict from the optimal lambda from tfCox_choose_lambda

Description

Estimate the corresponding theta values from the optimal tuning parameter obtained by tfCox_choose_lambda.

Usage

predict_best_lambda(cv, newX)

Arguments

Value

Estimated theta values.

Author(s)

Jiacheng Wu

References

Jiacheng Wu & Daniela Witten (2019) Flexible and Interpretable Models for Survival Data, Journal of Computational and Graphical Statistics, DOI: 10.1080/10618600.2019.1592758

See Also

tfCox_choose_lambda, negloglik

Examples

#generate training and testing data
dat = sim_dat(n=100, zerof=0, scenario=1)
test_dat = sim_dat(n=100, zerof=0, scenario=1)

#choose the optimal tuning parameter
cv = tfCox_choose_lambda(dat, test_dat, ord=0, alpha=1)
plot(cv$lam_seq, cv$loss)

#optimal tuning parameter
cv$best_lambda

#Estimate the theta values of testing covariates from the optimal tuning parameter
#from tfCox_choose_lambda object. 
theta_hat = predict_best_lambda(cv, test_dat$X)

Simulate Data from a Variety of Functional Scenarios

Description

This function generates survival data according to the simulation scenarios considered in Section 4 of Wu, J., and Witten, D. (2019) Flexible and interpretable models for survival data. Cox model has the form

\lambda(t|x) = \lambda_0(t) exp(\sum_{j=1}^p f_j(x))

. Failure time is generated by Weibull distribution with baseline hazard

\lambda_0(t) = scale * shape * t ^ {shape-1}

. In the paper, however, failure time is generated by a simplied weibull distribution: exponential(1) baseline hazard corresponding to shape=1 and scale=1. Censoring time is generated independently by exponential distribution with intensity censoring.rate. Thus the observed time is the minimum of failure time and censoring time. Each scenario has four covariates that have some non-linear association with the outcome. There is the option to also generate a user-specified number of covariates that have no association with the outcome.

Usage

sim_dat(n, zerof=0, scenario=1, scale=1, shape=1, censoring.rate=0.01, n.discrete=0)

Arguments

Value

Author(s)

Jiacheng Wu

References

Jiacheng Wu & Daniela Witten (2019) Flexible and Interpretable Models for Survival Data, Journal of Computational and Graphical Statistics, DOI: 10.1080/10618600.2019.1592758

See Also

plot.sim_dat

Examples

#generate data
set.seed(123)
dat = sim_dat(n=100, zerof=0, scenario=1)
#plot X versus the true theta
plot.sim_dat(dat)

Summarize cv_tfCox object

Description

This function summarizes cv_tfCox object and identifies the tuning parameter chosen by cross-validation.

Usage

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

Arguments

Author(s)

Jiacheng Wu

See Also

cv_tfCox, plot.cv_tfCox

Examples

#generate data
set.seed(1234)
dat = sim_dat(n=100, zerof=0, scenario=1)

#cross-validation to choose the tuning parameter lambda with fixed alpha=1
cv = cv_tfCox(dat, ord=0, alpha=1, n.fold=2)
#summarize the cross-validation
summary(cv)
#plot the cross-validation curve
plot(cv)

Summarize tfCox object

Description

This function summarizes tfCox object

Usage

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

Arguments

Details

Summarize the fit by the number of knots and percent sparsity achieved. Percent sparsity is the percentage of features estimated to have no relationship with the outcome.

Author(s)

Jiacheng Wu

See Also

tfCox, plot.tfCox

Examples

#generate data
set.seed(1234)
dat = sim_dat(n=100, zerof=0, scenario=1)

#fit piecewise constant for alpha=1 and a range of lambda
fit = tfCox(dat, ord=0, alpha=1)

#summarize the fit by the number of knots and percent sparsity achieved.
#Percent sparsity is the percentage of features estimated to have
#no relationship with outcome
summary(fit)

Fit the additive trend filtering Cox model with a range of tuning parameters

Description

Fit additive trend filtering Cox model where each component function is estimated to be piecewise constant or polynomial.

Usage

tfCox(dat, ord=0, alpha=1, lambda.seq=NULL, discrete=NULL, n.lambda=30,
lambda.min.ratio = 0.01, tol=1e-6, niter=1000, stepSize=25, backtracking=0)

Arguments

Details

The optimization problem has the form

l(\theta)+\alpha\lambda\sum_{j=1}^p |D_jP_j\theta_j|_1+(1-\alpha)\lambda\sum_{j=1}^p|\theta_j|_2

where l is the loss function defined as the negative log partial likelihood divided by n, and \alpha provides a trade-off between trend filtering penalty and group lasso penalty. Covariate matrix X is not standardized before solving the optimization problem.

Value

An object with S3 class "tfCox".

Author(s)

Jiacheng Wu

References

Jiacheng Wu & Daniela Witten (2019) Flexible and Interpretable Models for Survival Data, Journal of Computational and Graphical Statistics, DOI: 10.1080/10618600.2019.1592758

See Also

summary.tfCox, predict.tfCox, plot.tfCox, cv_tfCox

Examples

###################################################################
#constant trend filtering (fused lasso) with adaptively chosen knots
#generate data from simulation scenario 1 with piecewise constant functions
set.seed(1234)
dat = sim_dat(n=100, zerof=0, scenario=1)

#fit piecewise constant for alpha=1 and a range of lambda
fit = tfCox(dat, ord=0, alpha=1)
summary(fit)
#plot the fit of lambda index 15 and the first predictor
plot(fit, which.lambda=15, which.predictor=1)

#cross-validation to choose the tuning parameter lambda with fixed alpha=1
cv = cv_tfCox(dat, ord=0, alpha=1, n.fold=2)
summary(cv)
cv$best.lambda
#plot the cross-validation curve
plot(cv)

#fit the model with the best tuning parameter chosen by cross-validation
one.fit = tfCox(dat, ord=0, alpha=1, lambda.seq=cv$best.lambda)
#predict theta from the fitted tfCox object
theta_hat = predict(one.fit, newX=dat$X, which.lambda=1)

#plot the fitted theta_hat (line) with the true theta (dot)
for (i in 1:4) {
  ordi = order(dat$X[,i])
  plot(dat$X[ordi,i], dat$true_theta[ordi,i],
    xlab=paste("predictor",i), ylab="theta" )
  lines(dat$X[ordi,i], theta_hat[ordi,i], type="s")
}


#################################################################
#linear trend filtering with adaptively chosen knots
#generate data from simulation scenario 3 with piecewise linear functions
set.seed(1234)
dat = sim_dat(n=100, zerof=0, scenario=3)

#fit piecewise constant for alpha=1 and a range of lambda
fit = tfCox(dat, ord=1, alpha=1)
summary(fit)
#plot the fit of lambda index 15 and the first predictor
plot(fit, which.lambda=15, which.predictor=1)

#cross-validation to choose the tuning parameter lambda with fixed alpha=1
cv = cv_tfCox(dat, ord=1, alpha=1, n.fold=2)
summary(cv)
#plot the cross-validation curve
plot(cv)

#fit the model with the best tuning parameter chosen by cross-validation
one.fit = tfCox(dat, ord=1, alpha=1, lambda.seq=cv$best.lambda)
#predict theta from the fitted tfCox object
theta_hat = predict(one.fit, newX=dat$X, which.lambda=1)

#plot the fitted theta_hat (line) with the true theta (dot)
for (i in 1:4) {
  ordi = order(dat$X[,i])
  plot(dat$X[ordi,i], dat$true_theta[ordi,i],
       xlab=paste("predictor",i), ylab="theta" )
  lines(dat$X[ordi,i], theta_hat[ordi,i], type="l")
}

Choose the tuning parameter lambda using training and testing dataset

Description

Fit additive trend filtering Cox model where each component function is estimated to be piecewise constant or polynomial. Tuning parameter is selected via training and testing dataset described in Wu and Witten (2019). Training data is used to build the model, and testing data is used for selecting tuning parameter based on log likelihood. It is a convenience function to replicate the simulation results in Wu and Witten (2019).

Usage

tfCox_choose_lambda(dat, test_dat, ord = 0, alpha = 1, discrete = NULL, 
lam_seq = NULL, nlambda = 30, c = NULL, tol = 1e-06, niter=1000, 
stepSize=25, backtracking=0)

Arguments

Value

Author(s)

Jiacheng Wu

References

Jiacheng Wu & Daniela Witten (2019) Flexible and Interpretable Models for Survival Data, Journal of Computational and Graphical Statistics, DOI: 10.1080/10618600.2019.1592758

See Also

predict_best_lambda, negloglik

Examples

#generate training and testing data
dat = sim_dat(n=100, zerof=0, scenario=1)
test_dat = sim_dat(n=100, zerof=0, scenario=1)

#choose the optimal tuning parameter
cv = tfCox_choose_lambda(dat, test_dat, ord=0, alpha=1)
plot(cv$lam_seq, cv$loss)

#optimal tuning parameter
cv$best_lambda

#predict the coefficients of testing covariates from the optimal tuning parameter
#from tfCox_choose_lambda object. 
theta_hat = predict_best_lambda(cv, test_dat$X)

#calculate the loss in the testing data based on the estimated coefficients theta
negloglik(test_dat, theta_hat)