A Hierarchical Bayesian Cox Model with Variational Inference for High-Dimensional Survival Analysis

Large-scale molecular datasets have become increasingly important in bioinformatics, providing a foundation for prognostic modeling and the identification of biologically relevant genes. In survival analysis, however, censored observations and the high-dimensional nature of genomic data pose significant challenges for reliable prediction and variable selection. To address these challenges, we develop a Penalized Hierarchical Bayesian Cox (PHBC) modeling framework that extends hierarchical shrinkage and global–local prior formulations using Adaptive Gaussian Priors (AGP) and a modified Global–Local Gaussian–Inverse-Gamma Hyperprior (Gaussian–IGH) structure. The framework consists of two hierarchical components: (i) Adaptive Gaussian priors assigned to regression coefficients, with coefficient-specific variances updated from the data through an empirical-Bayes rule, and (ii) Inverse-Gamma hyperpriors on variance components to provide flexible and adaptive shrinkage. Inference is based on the Cox partial likelihood, yielding a partial-likelihood (pseudo-Bayesian) posterior. We derive the corresponding posterior kernels and construct scalable variational inference approximations, resulting in two computationally efficient procedures: AGP-VI and Global–Local Gaussian–IGH-VI. The methods are evaluated through extensive simulation studies and applications to six low- and high-dimensional genomic survival datasets. Performance is assessed using the concordance index (C-index), Integrated Brier Score, time-dependent AUC, and true/false positive rates for variable selection. A dedicated simulation study further examines sensitivity to multicollinearity under weak, moderate, and high predictor correlations, with the proposed methods remaining competitive across the different correlation levels. Across the real-data applications, no single method uniformly dominated all predictive criteria. The proposed variational methods showed strong but dataset-dependent performance, with Global–Local Gaussian–IGH-VI achieving the best overall predictive performance on one dataset and the best C-index on another, while AGP-VI achieved the best discrimination on two others. Classical penalized and Bayesian Cox approaches remained competitive, particularly in calibration. These results demonstrate the complementary strengths and trade-offs of the proposed methods in high-dimensional genomic survival analysis.