From Bareiss' algorithm to the stable computation of partial correlations

Abstract

This paper presents a new algorithm for the stable computation of sample partial correlation coefficients. We start by generalizing Bareiss' algorithm for the solution of linear systems of equations with (non-symmetric) Toeplitz coefficient matrix to matrices that are not Toeplitz, and show that it computes their LU and UL factorizations. For symmetric positive-definite matrices B, the normalized version of Bareiss' algorithm is the Hyperbolic Cholesky algorithm, which computes the upper and lower triangular Cholesky factors U and L of B by means of 2×2 hyperbolic rotations. If B is a convariance matrix, certain partial correlation coefficients may be obtained as the geometric means of the multiplier pairs in Bareiss' algorithm or as the hyperbolic tangents of the hyperbolic rotations. The partial correlation coefficients for a given data matrix A may suffer from serious loss of numerical accuracy when computed from the covariance matrix B = A T A. While the data flow graph of Bareiss' algorithm suggests a method, due to Cybenko, for computing partial correlation coefficients directly from the columns of A, the graph of Hyperbolic Cholesky suggests an alternative that requires fewer operations nd exhibits a higher degree of parallelism. In that method one computes the QR decomposition A = QU of the data matrix and applies orthogonal rotations to transform U to L; the sines of the rotation angles constitute a set of partial correlation coefficients. We extend this approach to determine arbitrary partial correlation coefficients at little extra computational cost.