Package cequre
performs censored quantile regression of Huang (2010), and restores
monotonicity respecting via adaptive interpolation for dynamic
regression of Huang (2017). The monotonicity-respecting restoration
applies to general dynamic regression models including (uncensored or
censored) quantile regression model, additive hazards model, and dynamic
survival models of Peng and Huang (2007), among others.
This procedure is illustrated with the Mayo primary biliary
cholangitis dataset as given in package survival
.
## Mayo PBC data
library(survival)
pbc_analy <- as.matrix(na.omit(pbc[,c("time","status","age","edema","bili","albumin","protime")]))
# log transformation for time, bili, albumin, and protime
pbc_analy[,c(1,5:7)] <- log(pbc_analy[,c(1,5:7)])
colnames(pbc_analy)[c(1,5:7)] <- paste("log",colnames(pbc_analy)[c(1,5:7)])
# convert status to censoring indicator
pbc_analy[,2] <- pbc_analy[,2]>1
## Censored quantile regression
library(cequre)
taus <- .1*(1:9)
fit <- cequre(pbc_analy[,1],pbc_analy[,2],pbc_analy[,-c(1,2)],
taus=taus,res=200)
Plot the estimated regression coefficient processes and associated pointwise 95% CI’s:
var.names <- c("baseline",colnames(pbc_analy)[-(1:2)])
par(mfrow=c(3,2),mai=c(.5,.5,.1,.1))
for(k in 1:6) {
plot(fit$curve[7,],fit$curve[k,],type="s",lwd=2,xlab="",ylab="",
ylim=c(-1,1)*max(abs(fit$curve[k,])))
lines(c(0,1),c(0,0),lty=2)
text(1,max(abs(fit$curve[k,])),var.names[k],adj=c(1,1))
# pointwise 95% CI
for(i in 1:length(taus)) lines(c(1,1)*taus[i],fit$bt[k,i]+
c(-1,1)*qnorm(.975)*sqrt(fit$va[k,k,i]))
}
mtext("Estimated regression coefficient",side=2,line=-1,outer=T,cex=.8)
mtext("Probability",side=1,line=-1,outer=T,cex=.8)
The monotonicity-respecting restoration is illustrated with the preceding censored quantile regression using the Mayo primary biliary cholangitis dataset:
# zch needs to include constant 1 as being the covariate for the intercept
mfit <- monodr(fit$curve,
zch = cbind(rep(1,dim(pbc_analy)[1]),pbc_analy[,-(1:2)]),
initau=fit$tau.bnd/2)
# plot the original regression coefficient (in black) and the one after the
# monotonicity-respecting restoration (in orange)
par(mfrow=c(3,2),mai=c(.5,.5,.1,.1))
for(k in 1:6) {
plot(fit$curve[7,],fit$curve[k,],type="s",xlab="",ylab="",
ylim=c(-1,1)*max(abs(fit$curve[k,])))
lines(c(0,1),c(0,0),lty=2)
text(1,max(abs(fit$curve[k,])),var.names[k],adj=c(1,1))
lines(mfit$airc[7,],mfit$airc[k,],lwd=2,col="dark orange")
}
mtext("Estimated regression coefficient",side=2,line=-1,outer=T,cex=.8)
mtext("Probability",side=1,line=-1,outer=T,cex=.8)
Huang, Y. (2010) Quantile calculus and censored regression, The Annals of Statistics 38, 1607–1637.
Huang, Y. (2017) Restoration of monotonicity respecting in dynamic regression. Journal of the American Statistical Association 112, 613–622.