Bayesian semiparametric modeling of survival data based on mixtures of $B$-spline distributions

Abstract

The nonparametric part of a semiparametric regression model usually involves prior specification for an infinite-dimensional parameter $F$. This paper introduces a class of finite mixture models based on $B$-spline distributions as an approximation to priors on the set of cumulative distribution functions. This class includes the mixture of beta distributions of Diaconis and Ylvisaker (1985) and the mixtures of triangular distributions of Perron and Mengersen (2001) as special cases. We describe how this approach can be used to model the baseline hazards in a Bayesian stratified proportional hazards model. A numerical illustration is given using survival data from a multicenter clinical AIDS trial, thus generalizing the approach by Carlin and Hodges (1999). Using conditional predictive ordinates and the deviance information criterion, we compare the fit of hierarchical proportional hazards regression models based on mixtures of $B$-spline distributions of various degrees.

Introduction

During the past two decades, discrete and continuous mixture models have been increasingly employed as flexible nonparametric methods to model unknown distributional shapes. The extent and the potential of their applications have widened considerably due to advances in simulation-based model fitting. In the Bayesian framework, the focus has been on estimating the unknown number of components in finite mixture models (Richardson and Green, 1997, Robert and Mengersen, 1999, Celeux et al., 2000, Stephens, 2000, Jasra et al., 2005, Nobile and Fearnside, 2007, among others), developing inference techniques for the Dirichlet process (DP) mixture setting that specifies an infinite dimensional mixture by assuming a random mixture distribution which is not restricted to a specified parametric family (e.g. MacEachern and Müller, 1998, Neal, 2000, Müller and Quintana, 2004) and extending the DP to dependent distributions (Gelfand et al., 2005, Teh et al., 2006, Omiros and Roberts, 2008).

Bayesian nonparametric methods to analyzing survival data allow for a flexible modeling not restricted to parametric assumptions, handling of censored and truncated data, and inference that does not rely on asymptotic approximations for large samples. Bayesian semiparametric approaches that incorporate covariate information follow either the framework of proportional hazards or accelerated failure time regression models and require a prior random process that generates the nonparametric part, the baseline hazard function or (log-) failure time distribution, respectively. As reviewed by Sinha and Dey (1997), Ibrahim et al. (2001) and Kottas (2006), Dirichlet, gamma, beta and correlated process priors, DP mixture and Pólya tree priors have been utilized as well as random finite mixture process priors for that purpose.

The approach taken in this paper falls within the class of finite mixture process priors that have been developed by Gelfand and Mallick (1995) to model monotone functions. Their work is based on previous research on the approximation of a general prior by a more tractable mixture of conjugate priors by Diaconis and Ylvisaker (1985) who showed that discrete mixtures of beta densities provide a continuous dense class of models for densities on [0, 1]. Using a suitable transformation of the cumulative hazard to a cumulative distribution function (cdf) on [0, 1], they modeled this cdf as a finite but random mixture of beta cdfs. This strategy was used by Carlin and Hodges (1999) for proportional hazards regression of highly stratified data. The general approach of using finite mixtures to model monotone functions was modified by Perron and Mengersen (2001) who suggested mixtures of triangular distributions and showed the superiority of a mixture of triangular over beta distributions in approximating monotone functions.

So far, the two nonparametric techniques using mixtures of beta and triangular distributions have been presented as two distinct approaches. In this article, we show that these two methods can be integrated into the conceptual framework of mixtures of $B$-spline distributions. We point out that the mixture of beta and triangular distributions are two special cases of the mixture of $B$-spline distributions. More precisely, when the degree of the $B$-spline is one, the mixture of $B$-spline distributions reduces to the mixture of triangular distributions. When there are no interior knots, the $B$-spline basis becomes the Bernstein polynomial basis of the mixture of beta distributions. Thus, we present a generalized finite mixture model approach based on $B$-spline distributions, which contains beta mixtures, triangular mixtures and mixtures of higher degree $B$-spline distributions. But we also note that other choices of basis functions, such as the truncated power series basis (cf. Nott and Li, 2010, Ruppert et al., 2003), or Bayesian $P$-splines (Lang and Brezger, 2004) are feasible.

We illustrate the flexibility of this new class of mixture of $B$-spline distributions in the context of semiparametric survival analysis of highly stratified data using the dataset of a multi-center controlled clinical AIDS trial as in Carlin and Hodges (1999) and compare results to those reported in that paper. In particular, the strange and anomalous behavior of the baseline hazard when using the mixture of beta distributions as in Carlin and Hodges (1999) vanishes when using mixture of triangular and higher order $B$-spline distributions. The conditional predictive ordinate (CPO) statistic (Geisser and Eddy, 1979)  and the deviance information criterion (DIC, Spiegelhalter et al., 2002) are adopted for model comparisons.

The remainder of the article is organized as follows. Section 2 starts with an illustrative example of a multi-center clinical trial where a flexible nonparametric modeling of unit-specific hazard functions is called for. Section 3 defines the class of mixtures of $B$-spline distributions which includes mixtures of beta and triangular distributions as special cases. Section 4 presents a class of semiparametric proportional hazards models where the nonparametric part is based on mixtures of $B$-spline distributions. This section includes a discussion of prior specification and model comparison issues. In Section 5, we use these semiparametric models for analyzing the dataset introduced in Section 2. Finally, Section 6 concludes with a summary and discussion of the results.