Uhaz: U-shaped Hazard Function Estimation
In npsurv: Nonparametric Survival Analysis

Description Usage Arguments Details Value Author(s) References See Also Examples

Uhaz computes the nonparametric maximum likelihood esimate (NPMLE) of a U-shaped hazard function from exact or interval-censored data, or a mix of the two types of data.

1	Uhaz(data, w = 1, deg = 1, maxit = 100, tol = 1e-06, verb = 0)

`data`	vector or matrix, or an object of class `icendata`.
`w`	weights or multiplicities of the observations.
`deg`	nonnegative real number for spline degree (i.e., p in the formula below).
`maxit`	maximum number of iterations.
`tol`	tolerance level for stopping the algorithm. It is used as the threshold on the increase of the log-likelihood after each iteration.
`verb`	verbosity level for printing intermediate results in each iteration.

If data is a vector, it contains only exact observations, with weights given in w.

If data is a matrix with two columns, it contains interval-censored observations, with the two columns storing their left and right end-points, respectively. If the left and right end-points are equal, then the observation is exact. Weights are provided by w.

If data is a matrix with three columns, it contains interval-censored observations, with the first two columns storing their left and right end-points, respectively. The weight of each observation is the third-column value multiplied by the corresponding weight value in w.

The algorithm used for the computing the NPMLE of a hazard function under the U-shape restriction is is proposed by Wang and Fani (2015). Such a hazard function is given by

A U-shaped hazard function is given by

h(t) = alpha + sum_{j=1}^k nu_j (tau_j - t)_+^p + sum_{j=1}^m mu_j (t - eta_j)_+^p,

where alpha, nu_j, mu_j ≥ 0, tau_1 < ... < tau_k <= eta_1 < ... < eta_m, and p >= 0 is the the spline degree which determines the smoothness of the U-shaped hazard. As p increases, the family of hazard functions becomes increasingly smoother, but at the time, smaller. When p = 0, the hazard function is U-shaped, as studied by Bray et al. (1967). When p = 1, the hazard function is convex, as studied by Jankowski and Wellner (2009a,b).

Note that deg (i.e., p in the above mathematical display) can take on any nonnegative real value.

An object of class Uhaz, which is a list with components:

`convergence`	= `TRUE`, converged successfully; = `FALSE`, maximum number of iterations reached.
`grad`	gradient values at the knots.
`numiter`	number of iterations used.
`ll`	log-likelihood value of the NPMLE `h`.
`h`	NPMLE of the U-shaped hazard function, an object of class `uh`.

Yong Wang <yongwang@auckland.ac.nz>

Bray, T. A., Crawford, G. B., and Proschan, F. (1967). Maximum Likelihood Estimation of a U-shaped Failure Rate Function. Defense Technical Information Center.

Jankowski, H. K. and Wellner, J. A. (2009a). Computation of nonparametric convex hazard estimators via profile methods. Journal of Nonparametric Statistics, 21, 505-518.

Jankowski, H. K. and Wellner, J. A. (2009b). Nonparametric estimation of a convex bathtub-shaped hazard function. Bernoulli, 15, 1010-1035.

Wang, Y. and Fani, S. (2018). Nonparametric maximum likelihood computation of a U-shaped hazard function. Statistics and Computing, 28, 187-200.

icendata, nzmort.

## Interval-censored observations
data(ap)
(r = Uhaz(ap, deg=0))
plot(r, ylim=c(0,.3), col=1)
for(i in 1:6) plot(Uhaz(ap, deg=i/2), add=TRUE, col=i+1)
legend(15, 0.01, paste0("deg = ", 0:6/2), lwd=2, col=1:7, xjust=1, yjust=0)

## Exact observations
data(nzmort)
x = with(nzmort, nzmort[ethnic=="maori",])[,1:2]   # Maori mortality
(h0 = Uhaz(x[,1]+0.5, x[,2], deg=0)$h)     # U-shaped hazard
(h1 = Uhaz(x[,1]+0.5, x[,2], deg=1)$h)     # convex hazard
(h2 <- Uhaz(x[,1]+0.5, x[,2], deg=2)$h)    # smooth U-shaped hazard

plot(h0, pch=2)                            # plot hazard functions
plot(h1, add=TRUE, col="green3", pch=1)
plot(h2, add=TRUE, col="red3", pch=19)

age = 0:max(x[,1])                         # plot densities
count = integer(length(age))
count[x[,"age"]+1] = x[,"deaths"]
barplot(count/sum(count), space=0, col="lightgrey")
axis(1, pos=NA, at=0:10*10)
plot(h0, fn="d", add=TRUE, pch=2)
plot(h1, fn="d", add=TRUE, col="green3", pch=1)
plot(h2, fn="d", add=TRUE, col="red3", pch=19)

plot(h0, fn="s", pch=2)                    # plot survival functions
plot(h1, fn="s", add=TRUE, col="green3", pch=1)
plot(h2, fn="s", add=TRUE, col="red3", pch=19)

## Exact and right-censored observations
data(gastric)
plot(h0<-Uhaz(gastric, deg=0)$h)           # plot hazard functions
plot(h1<-Uhaz(gastric, deg=1)$h, add=TRUE, col="green3")
plot(h2<-Uhaz(gastric, deg=2)$h, add=TRUE, col="red3")

plot(npsurv(gastric), fn="s", col="grey")  # plot survival functions
plot(h0, fn="s", add=TRUE)
plot(h1, fn="s", add=TRUE, col="green3")
plot(h2, fn="s", add=TRUE, col="red3")