Joint lifetime modeling with matrix distributions

Abstract Acyclic phase-type (PH) distributions have been a popular tool in survival analysis, thanks to their natural interpretation in terms of aging toward its inevitable absorption. In this article, we consider an extension to the bivariate setting for the modeling of joint lifetimes. In contrast to previous models in the literature that were based on a separate estimation of the marginal behavior and the dependence structure through a copula, we propose a new time-inhomogeneous version of a multivariate PH (mIPH) class that leads to a model for joint lifetimes without that separation. We study properties of mIPH class members and provide an adapted estimation procedure that allows for right-censoring and covariate information. We show that initial distribution vectors in our construction can be tailored to reflect the dependence of the random variables and use multinomial regression to determine the influence of covariates on starting probabilities. Moreover, we highlight the flexibility and parsimony, in terms of needed phases, introduced by the time inhomogeneity. Numerical illustrations are given for the data set of joint lifetimes of Frees et al., where 10 phases turn out to be sufficient for a reasonable fitting performance. As a by-product, the proposed approach enables a natural causal interpretation of the association in the aging mechanism of joint lifetimes that goes beyond a statistical fit.