Maximum approximate Bernstein/Beta likelihood estimation in Cox's proportional hazards regression model based on interal censored event time data.
mable.ph(
formula,
data,
M,
g = NULL,
p = NULL,
pi0 = NULL,
tau = Inf,
x0 = NULL,
controls = mable.ctrl(),
progress = TRUE
)A list with components
m the selected/preselected optimal degree m
p the estimate of \(p = (p_0, \dots, p_m, p_{m+1})\), the coefficients of Bernstein polynomial of degree m
coefficients the estimated regression coefficients of the PH model
SE the standard errors of the estimated regression coefficients
z the z-scores of the estimated regression coefficients
mloglik the maximum log-likelihood at an optimal degree m
tau.n maximum observed time \(\tau_n\)
tau right endpoint of support \([0, \tau)\)
x0 the working baseline covariates
egx0 the value of \(e^{\gamma'x_0}\)
convergence an integer code, 1 indicates either the EM-like
iteration for finding maximum likelihood reached the maximum iteration for at least one m
or the search of an optimal degree using change-point method reached the maximum candidate degree,
2 indicates both occured, and 0 indicates a successful completion.
delta the final delta if m0 = m1 or the final pval of the change-point
for searching the optimal degree m;
and, if m0<m1,
M the vector (m0, m1), where m1 is the last candidate degree when the search stoped
lk log-likelihoods evaluated at \(m \in \{m_0,\ldots, m_1\}\)
lr likelihood ratios for change-points evaluated at \(m \in \{m_0+1, \ldots, m_1\}\)
pval the p-values of the change-point tests for choosing optimal model degree
chpts the change-points chosen with the given candidate model degrees
regression formula. Response must be cbind. See 'Details'.
a data frame containing variables in formula.
a positive integer or a vector (m0, m1). If M = m or m0 = m1,
then m0 is a preselected degree. If m0<m1 it specifies the set of
consective candidate model degrees m0:m1 for searching an optimal degree,
where m1-m0>3.
initial guess of \(d\)-vector of regression coefficients. See 'Details'.
an initial coefficients of Bernstein polynomial model of degree m0,
default is the uniform initial.
Initial guess of \(\pi(x_0) = F(\tau_n|x_0)\). Without right censored data, pi0 = 1. See 'Details'.
right endpoint of support \([0, \tau)\) must be greater than or equal to the maximum observed time
a data frame specifying working baseline covariates on the right-hand-side of formula. See 'Details'.
Object of class mable.ctrl() specifying iteration limit
and other control options. Default is mable.ctrl.
if TRUE a text progressbar is displayed
Zhong Guan <zguan@iu.edu>
Consider Cox's PH model with covariate for interval-censored failure time data:
\(S(t|x) = S(t|x_0)^{\exp(\gamma^T(x-x_0))}\), where \(x_0\) satisfies \(\gamma^T(x-x_0)\ge 0\),
where \(x\) and \(x_0\) may
contain dummy variables and interaction terms. The working baseline x0 in arguments
contains only the values of terms excluding dummy variables and interaction terms
in the right-hand-side of formula. Thus g is the initial guess of
the coefficients \(\gamma\) of \(x-x_0\) and could be longer than x0.
Let \(f(t|x)\) and \(F(t|x) = 1-S(t|x)\) be the density and cumulative distribution
functions of the event time given \(X = x\), respectively.
Then \(f(t|x_0)\) on \([0, \tau_n]\) can be approximated by
\(f_m(t|x_0, p) = \tau_n^{-1}\sum_{i=0}^m p_i\beta_{mi}(t/\tau_n)\),
where \(p_i \ge 0\), \(i = 0, \ldots, m\), \(\sum_{i=0}^mp_i = 1-p_{m+1}\),
\(\beta_{mi}(u)\) is the beta denity with shapes \(i+1\) and \(m-i+1\), and
\(\tau_n\) is the largest observed time, either uncensored time, or right endpoint of interval/left censored,
or left endpoint of right censored time. So we can approximate \(S(t|x_0)\) on \([0, \tau_n]\) by
\(S_m(t|x_0; p) = \sum_{i=0}^{m+1} p_i \bar B_{mi}(t/\tau_n)\), where
\(\bar B_{mi}(u)\), \(i = 0, \ldots, m\), is the beta survival function with shapes
\(i+1\) and \(m-i+1\), \(\bar B_{m,m+1}(t) = 1\), \(p_{m+1} = 1-\pi(x_0)\), and
\(\pi(x_0) = F(\tau_n|x_0)\). For data without right-censored time, \(p_{m+1} = 1-\pi(x_0) = 0\).
Response variable should be of the form cbind(l, u), where (l, u) is the interval
containing the event time. Data is uncensored if l = u, right censored
if u = Inf or u = NA, and left censored data if l = 0.
The associated covariate contains \(d\) columns. The baseline x0 should chosen so that
\(\gamma'(x-x_0)\) is nonnegative for all the observed \(x\) and
all \(\gamma\) in a neighborhood of its true value.
A missing initial value of g is imputed by ic_sp() of package icenReg.
The search for optimal degree m stops if either m1 is reached or the test
for change-point results in a p-value pval smaller than sig.level.
This process takes longer than maple.ph to select an optimal degree.
Guan, Z. Maximum Approximate Bernstein Likelihood Estimation in Proportional Hazard Model for Interval-Censored Data, Statistics in Medicine. 2020; 1–21. https://doi.org/10.1002/sim.8801.
maple.ph
# \donttest{
# Ovarian Cancer Survival Data
require(survival)
futime2<-ovarian$futime
futime2[ovarian$fustat==0]<-Inf
ovarian2<-data.frame(age=ovarian$age, futime1=ovarian$futime,
futime2=futime2)
ova<-mable.ph(cbind(futime1, futime2) ~ age, data = ovarian2,
M=c(2,35), g=.16, x0=data.frame(age=35))
op<-par(mfrow=c(2,2))
plot(ova, which = "likelihood")
plot(ova, which = "change-point")
plot(ova, y=data.frame(age=60), which="survival", add=FALSE, type="l",
xlab="Days", main="Age = 60")
plot(ova, y=data.frame(age=65), which="survival", add=FALSE, type="l",
xlab="Days", main="Age = 65")
par(op)
# }
Run the code above in your browser using DataLab