Rybicki Press algorithm

The Rybicki–Press algorithm is a fast direct algorithm for inverting a matrix, whose entries are given by ${\displaystyle A(i,j)=\exp(-a\vert t_{i}-t_{j}\vert )}$, where ${\displaystyle a\in \mathbb {R} }$.[1] It is a computational optimization of a general set of statistical methods developed to determine whether two noisy, irregularly sampled data sets are, in fact, dimensionally shifted representations of the same underlying function.[2][3] The most common use of the algorithm is in the detection of periodicity in astronomical observations.[3]

Extended Sparse Matrix arising from a ${\displaystyle 10\times 10}$ semi-separable matrix whose semi-separable rank is ${\displaystyle 4}$.

Recently, this method has been extended (Generalized Rybicki Press algorithm) for inverting matrices whose entries of the form ${\displaystyle A(i,j)=\sum _{k=1}^{p}a_{k}\exp(-\beta _{k}\vert t_{i}-t_{j}\vert )}$.[4] The key observation in the Generalized Rybicki Press (GPP) algorithm is that the matrix ${\displaystyle A}$ is a semi-separable matrix with rank ${\displaystyle p}$. More precisely, if the matrix ${\displaystyle A\in \mathbb {R} ^{n\times n}}$ has a semi-separable rank is ${\displaystyle p}$, the cost for solving the linear system ${\displaystyle Ax=b}$ and obtaining the determinant of the matrix scales as ${\displaystyle {\mathcal {O}}\left(p^{2}n\right)}$, thereby making it extremely attractive for large matrices. This implementation of the GPP algorithm can be found here.[5] The key idea is that the dense matrix ${\displaystyle A}$ can be converted into a sparser matrix of a larger size (see figure on the right), whose sparsity structure can be leveraged to reduce the computational complexity.

The fact that matrix ${\displaystyle A}$ is a semi-separable matrix also forms the basis for celerite[6] library, which is a library for fast and scalable Gaussian Process (GP) Regression in one dimension[7] with implementations in C++, Python, and Julia. The celerite method[7] also provides an algorithm for generating samples from a high-dimensional distribution. The method has found attractive applications in a wide range of fields, especially in astronomical data analysis.[8][9]