---
title: "Censored Quantile Regression & Monotonicity-Respecting Restoring"
author: "Yijian Huang (yhuang5@emory.edu)"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Censored Quantile Regression & Monotonicity-Respecting Restoring}
  %\VignetteEngine{knitr::rmarkdown}
  \usepackage[utf8]{inputenc}
---

```{r setup, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
```

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.

## Installation
`cequre` is available on CRAN:
```{r install, eval=FALSE, message=FALSE, warning=FALSE}
install.packages("cequre")
```

## Censored quantile regression of Huang (2010)
This procedure is illustrated with the Mayo primary biliary cholangitis dataset
as given in package `survival`.
```{r cequre, eval=TRUE, message=FALSE, warning=FALSE}
## 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:
```{r cequreplot, eval=TRUE, message=FALSE, warning=FALSE, fig.height = 5, fig.width = 7, dpi = 300, out.width = "80%"}
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)
```

## Restoration of monotonicity respecting using adaptive interpolation of Huang (2017)
The monotonicity-respecting restoration is illustrated with the preceding
censored quantile regression using the Mayo primary biliary cholangitis dataset:
```{r monodr, eval=TRUE, message=FALSE, warning=FALSE, fig.height = 5, fig.width = 7, dpi = 300, out.width = "80%"}
# 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)
```

## References
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.