Exact methods for exponentiation of matrices of dimension $N$ can be computationally expensive in terms of execution time ($N^{3}$) and memory requirements ($N^{2}$) not to mention numerical precision issues. A type of matrix often exponentiated in the sciences is the rate matrix. Here we explore five methods to exponentiate rate matrices some of which apply even more broadly to other matrix types. Three of the methods leverage a mathematical analogy between computing matrix elements of a matrix exponential and computing transition probabilities of a dynamical processes (technically a Markov jump process, MJP, typically simulated using Gillespie). In doing so, we identify a novel MJP-based method relying on restricting the number of ``trajectory" jumps based on the magnitude of the matrix elements with favorable computational scaling. We then discuss this method's downstream implications on mixing properties of Monte Carlo posterior samplers. We also benchmark two other methods of matrix exponentiation valid for any matrix (beyond rate matrices and, more generally, positive definite matrices) related to solving differential equations: Runge-Kutta integrators and Krylov subspace methods. Under conditions where both the largest matrix element and the number of non-vanishing elements scale linearly with $N$ -- reasonable conditions for rate matrices often exponentiated -- computational time scaling with the most competitive methods (Krylov and one of the MJP-based methods) reduces to $N^2$ with total memory requirements of $N$.
翻译:对于维度为$N$的矩阵,精确的指数计算方法在执行时间(复杂度为$N^{3}$)和内存需求(复杂度为$N^{2}$)上计算成本高昂,更不用说数值精度问题。科学领域中常被指数化的一种矩阵是速率矩阵。本文探讨了五种指数化速率矩阵的方法,其中部分方法甚至可更广泛地应用于其他矩阵类型。三种方法利用了计算矩阵指数的矩阵元素与计算动态过程(技术上称为马尔可夫跳跃过程MJP,通常使用Gillespie算法模拟)的转移概率之间的数学类比。在此过程中,我们提出了一种基于MJP的新方法,该方法通过矩阵元素的量级限制“轨迹”跳跃次数,具有优越的计算扩展性。随后,我们讨论了该方法对蒙特卡洛后验采样器混合性质的潜在影响。我们还对另外两种适用于任意矩阵(不限于速率矩阵,更一般地,正定矩阵)的矩阵指数求值方法进行了基准测试,它们与微分方程求解有关:龙格-库塔积分器和Krylov子空间方法。当最大矩阵元素和非零元素数量均与$N$呈线性关系时(这通常是需指数化速率矩阵的合理条件),最具竞争力的方法(Krylov方法及一种基于MJP的方法)的计算时间可降至$N^2$,总内存需求为$N$。