Thanks to Matthew Albrecht for responding to Joyce’s email about this! Below is a passage from the supplement of his paper in Conservation Biology about this very issue. If I understand correctly, they used a cox regression model with mixed effects in the coxme package in R, but used a phylogenetic variance-covariance matrix generated through MCMCglmm.
Incorporating phylogenetic information into random effects models can change whether an effect size is significantly different from zero (Chamberlain et al. 2012). Thus, we use a correlated random effects model as a more explicit way to account for unobserved heterogeneity in time to recruitment due to variation in phylogenetic relatedness. Correlated random effects models define the covariance structure of random terms in a mixed model using a variance- covariance matrix (Therne 2015). In this case, we defined the covariance structure of the random effect of species by the degree of phylogenetic relatedness among species. Hence, a phylogenetic variance-covariance matrix was used to characterize the variance family for the random effect of species. To implement this model, we generated a phylogeny of the 27 species represented in the dataset on the online Phylomatic web-service, using an already dated plant phylogeny (R20070305.bl.new) available from Phylomatic’s local store of mega-trees (Webb & Donoghue 2005). The final topologycontained 27 tips and 21 internal nodes (Appendix S3). We derived a phylogenetic covariance matrix from this phylogeny using the inverseA function of the MCMCglmm package (Hadfield 2016), and forced it to a symmetric matrix with the forcenSymmetric function of the Matrix package (Bates & Maechler 2016) in R. This phylogenetic covariance matrix served as the variance-covariance structure in a correlated random effects model with species as a random intercept.