BigDataStatMeth (version 1.0.3)

bdInvCholesky_hdf5: Matrix Inversion using Cholesky Decomposition for HDF5-Stored Matrices

Description

Computes the inverse of a symmetric positive-definite matrix stored in an HDF5 file using the Cholesky decomposition method. This approach is more efficient and numerically stable than general matrix inversion methods for symmetric positive-definite matrices.

Usage

bdInvCholesky_hdf5(
  filename,
  group,
  dataset,
  outdataset,
  outgroup = NULL,
  fullMatrix = NULL,
  overwrite = NULL,
  threads = 2L,
  elementsBlock = 1000000L
)

Value

List with components:

fn

Character string with the HDF5 filename

ds

Character string with the full dataset path to the inverse Cholesky decomposition A^(-1) result (group/dataset)

Arguments

filename

Character string. Path to the HDF5 file containing the input matrix.

group

Character string. Path to the group containing the input dataset.

dataset

Character string. Name of the input dataset to invert.

outdataset

Character string. Name for the output dataset.

outgroup

Character string. Optional output group path. If not provided, results are stored in the input group.

fullMatrix

Logical. If TRUE, stores the complete inverse matrix. If FALSE (default), stores only the lower triangular part to save space.

overwrite

Logical. If TRUE, allows overwriting existing results.

threads

Integer. Number of threads for parallel computation (default = 2).

elementsBlock

Integer. Maximum number of elements to process in each block (default = 1,000,000). For matrices larger than 5000x5000, automatically adjusted to number of rows or columns * 2.

Details

This function implements an efficient matrix inversion algorithm that leverages the special properties of symmetric positive-definite matrices. Key features:

• Uses Cholesky decomposition for improved numerical stability

• Block-based computation for large matrices

• Optional storage formats (full or triangular)

• Parallel processing support

• Memory-efficient block algorithm

The algorithm proceeds in two main steps:

1. Compute the Cholesky decomposition A = LL'

2. Solve the system LL'X = I for X = A^(-1)

Advantages of this method:

• More efficient than general matrix inversion

• Better numerical stability

• Preserves matrix symmetry

• Exploits positive-definiteness for efficiency

References

• Golub, G. H., & Van Loan, C. F. (2013). Matrix Computations, 4th Edition. Johns Hopkins University Press.

• Higham, N. J. (2002). Accuracy and Stability of Numerical Algorithms, 2nd Edition. SIAM.

See Also

• bdCholesky_hdf5 for the underlying Cholesky decomposition

• bdSolve_hdf5 for solving linear systems

Examples

if (FALSE) {
library(rhdf5)

# Create a symmetric positive-definite matrix
set.seed(1234)
X <- matrix(rnorm(100), 10, 10)
A <- crossprod(X)  # A = X'X is symmetric positive-definite

# Save to HDF5
h5createFile("matrix.h5")
h5write(A, "matrix.h5", "data/matrix")

# Compute inverse using Cholesky decomposition
bdInvCholesky_hdf5("matrix.h5", "data", "matrix",
                   outdataset = "inverse",
                   outgroup = "results",
                   fullMatrix = TRUE,
                   threads = 4)

# Verify the inverse
Ainv <- h5read("matrix.h5", "results/inverse")
max(abs(A %*% Ainv - diag(nrow(A))))  # Should be very small
}