This repo provides implementation for the proposed estimation and inference procedure in the following paper:

Survival Analysis via Ordinary Differential Equation.

Weijing Tang, Kevin He, Gongjun Xu, and Ji Zhu.

The implementation of the proposed method is build on top of the following toolboxes in Matlab.

• Optimization Toolbox for solving constrained and unconstrained optimization problems
• Curve Fitting Toolbox for constructing the sieve space using B-splines

Fitting the Cox model: with unspecified

• The data generation: run cox/generator.m and set the setting as 1. The generated data will be saved in "data/survode/".
• The proposed ODE-Cox approach where is approximated by the cubic B-splines: run cox/main.m to obtain the point estimate and run cox/inference.m to obtain the standard error estimate. See cox/summary.m for example.
• The comparison with the existing maximum partial likelihood estimation: the implementation baseline/cox_mple.R is based on the "coxph" function in the R package "survival".
• The comparison with the the parametric method in Royston and Parmar (2002), where the log-transformed baseline cumulative hazard is modeled as a natural cubic spline function of the log-transformed time: the implementation baseline/cox_flexsurv.R is based on the "flexsurvspline" function in the R package "flexsurv".

Fitting the time-varying Cox model: with unspecified and

• The data generation: run baseline/cox_time_varying_generate_data.R with a given sample size. The generated data will be saved in "data/cox_time_varying_model/".
• The proposed ODE-Cox-tv approach where and are approximated by the cubic B-splines: see cox_time_varying/main.m for example.
• The comparison with the existing maximum partial likelihood estimation: the implementation baseline/cox_time_varying_effects.R is based on the "coxph" function the “tt” argument set as the same cubic B-spline transformation of time in the R package "survival".

Fitting the semi-parametric accelerated failure time (AFT) model: with unspecified

• The data generation: run aft/generator.m and set the setting as 2. The generated data will be saved in "data/survode/".
• The proposed ODE-AFT approach where is approximated by the cubic B-splines: see aft/main.m for example.
• I have implemented both forward and adjoint methods to compute the gradients: in aft/main.m, set forward=true for using forward method and forward=false for using adjoint method along with parallel computing. The default is using the forward method if the memory permits.
• The comparison with the existing rank-based method: the implementation baseline/aft_rank.R is based on the "aftsrr" function in R package "aftgee".

Fitting the semi-parametric linear transformation model: with specified and unspecified . In particular, the solution of the initial value problem with corresponds to the G-transformation function in Zeng and Lin(2007), which covers the class of Box-Cox transformations and the class of logarithmic transformations.

• The data generation: run Gtransform/generator.m and set the setting as 3. The generated data will be saved in "data/survode/".
• The proposed ODE-LT approach where is approximated by the cubic B-splines: see Gtransform/main.m for example.
• The comparison with the existing nonparametric maximum likelihood estimation in Zeng and Lin(2007): the implementation baseline/NPMLE_Trasnform/main.m is modified based on the original implementation using the EM algorithm in Zeng and Lin(2007).

Fitting the general (or nonparametric) linear transformation model: with both and unspecified

• The data generation: run survode_cgd/generator.m and set the setting between 1 and 4. In particular, the Cox, the AFT, and the logarithmic transformation model are mis-specified under setting 4. The generated data will be saved in "data/survode/".
• The proposed ODE-Flex approach: see survode_cgd/main.m and survode_cgd/mle.m for example.
• The comparison with the existing smoothed partial rank method in Song et al. (2006): the implementation baseline/SPR_LT.py is able to reproduce the simulation results in Song et al. (2006). Note that SPR introduces an additional parameter c in the objective function to improve the estimation accuracy. We evaluated SPR with various values of the parameter c and the sample size N under our data settings 1)-4) in the manuscript.

I have implemented two placements of knots for cubic splines. The argument knots_setting in "main.m" can be set as:

• knots_setting = "quantile": the interior knots arelocated at the quantiles of the distinct observation time points.
• knots_setting = "equal": the interior knots equally separate the time interval from 0 to the maximum of observed times.

The default placement is knots_setting = "quantile".

Weijing Tang (University of Michigan) - weijtang@umich.edu