A QUASI–MONTE CARLO EXPECTATION–MAXIMIZATION FRAMEWORK FOR MULTIVARIATE JOINT MODELLING WITH PENALIZED COX REGRESSION

Abstract

Abstract The joint modelling of longitudinal and time-to-event data provides an effective framework for analysing the association between biomarker trajectories and survival outcomes. However, when multiple longitudinal processes are modelled simultaneously, the resulting high-dimensional random effects structure leads to substantial computational challenges due to complex likelihood integrals. This paper proposes an efficient estimation framework for multivariate joint models using Quasi Monte Carlo (QMC) integration within a Monte Carlo Expectation Maximization (MCEM) algorithm. The approach employs deterministic low-discrepancy sequences, specifically Sobol and Halton sequences to approximate the intractable integrals in the E-step, thereby improving convergence speed and estimation accuracy compared to classical Monte Carlo (MC) methods. To enhance model flexibility and stability in the survival component, a ridge-penalized Cox regression is incorporated into the joint modelling framework. The penalty term regularizes the parameter estimates, mitigates overfitting, and improves generalization in the presence of correlated covariates or high-dimensional predictors. Extensive simulation studies demonstrate that QMC methods yield comparable or superior estimation accuracy relative to standard MC integration, with substantially reduced computation time. Among the QMC variants, the Halton sequence exhibits the most stable convergence and lowest estimation variance, particularly in small-sample and high-dimensional scenarios. The proposed approach is further illustrated using PBC dataset, where the multivariate joint model simultaneously links two or three longitudinal biomarkers with survival outcomes. The empirical results confirm that QMC integration significantly accelerates computation while maintaining robust and precise parameter estimates. This study highlights the potential of QMC-based numerical integration as a practical and efficient alternative for fitting high-dimensional joint models in biomedical applications..