Title: | Riemannian Metrics for Symmetric Positive Definite Matrices |
Version: | 0.1.0 |
Description: | Implements various Riemannian metrics for symmetric positive definite matrices, including AIRM (Affine Invariant Riemannian Metric, see Pennec, Fillard, and Ayache (2006) <doi:10.1007/s11263-005-3222-z>), Log-Euclidean (see Arsigny, Fillard, Pennec, and Ayache (2006) <doi:10.1002/mrm.20965>), Euclidean, Log-Cholesky (see Lin (2019) <doi:10.1137/18M1221084>), and Bures-Wasserstein metrics (see Bhatia, Jain, and Lim (2019) <doi:10.1016/j.exmath.2018.01.002>). Provides functions for computing logarithmic and exponential maps, vectorization, and statistical operations on the manifold of positive definite matrices. |
License: | MIT + file LICENSE |
Encoding: | UTF-8 |
LazyData: | true |
RoxygenNote: | 7.3.2 |
Depends: | R (≥ 4.3.0), Matrix |
Imports: | methods, expm, R6, purrr, MASS, furrr |
Suggests: | testthat (≥ 3.0.0), knitr, rmarkdown |
VignetteBuilder: | knitr |
URL: | https://nicoesve.github.io/riemtan/ |
BugReports: | https://github.com/nicoesve/riemtan/issues |
Config/testthat/edition: | 3 |
Maintainer: | Nicolas Escobar <nescoba@iu.edu> |
NeedsCompilation: | no |
Packaged: | 2025-04-19 21:30:24 UTC; nicolasescobar |
Author: | Nicolas Escobar |
Repository: | CRAN |
Date/Publication: | 2025-04-23 10:10:02 UTC |
riemtan: Riemannian Metrics for Symmetric Positive Definite Matrices
Description
Implements various Riemannian metrics for symmetric positive definite matrices, including AIRM (Affine Invariant Riemannian Metric, see Pennec, Fillard, and Ayache (2006) doi:10.1007/s11263-005-3222-z), Log-Euclidean (see Arsigny, Fillard, Pennec, and Ayache (2006) doi:10.1002/mrm.20965), Euclidean, Log-Cholesky (see Lin (2019) doi:10.1137/18M1221084), and Bures-Wasserstein metrics (see Bhatia, Jain, and Lim (2019) doi:10.1016/j.exmath.2018.01.002). Provides functions for computing logarithmic and exponential maps, vectorization, and statistical operations on the manifold of positive definite matrices.
Author(s)
Maintainer: Nicolas Escobar nescoba@iu.edu (ORCID)
Other contributors:
Jaroslaw Harezlak [thesis advisor]
See Also
Useful links:
CSample Class
Description
This class represents a sample of connectomes, with various properties and methods to handle their tangent and vectorized images.
Active bindings
connectomes
Connectomes data
tangent_images
Tangent images data
vector_images
Vector images data
sample_size
Sample size
matrix_size
Matrix size
mfd_dim
Manifold dimension
is_centered
Centering status
frechet_mean
Frechet mean
riem_metric
Riemannian Metric used
variation
Variation of the sample
sample_cov
Sample covariance
ref_point
Reference point for tangent or vectorized images
Methods
Public methods
Method new()
Initialize a CSample object
Usage
CSample$new( conns = NULL, tan_imgs = NULL, vec_imgs = NULL, centered = NULL, ref_pt = NULL, metric_obj )
Arguments
conns
A list of connectomes (default is NULL).
tan_imgs
A list of tangent images (default is NULL).
vec_imgs
A matrix whose rows are vectorized images (default is NULL).
centered
Boolean indicating whether tangent or vectorized images are centered (default is NULL).
ref_pt
A connectome (default is identity)
metric_obj
Object of class
rmetric
representing the Riemannian metric used.
Returns
A new CSample
object.
Method compute_tangents()
This function computes the tangent images from the connectomes.
Usage
CSample$compute_tangents(ref_pt = default_ref_pt(private$p))
Arguments
ref_pt
A reference point, which must be a
dppMatrix
object (default isdefault_ref_pt
).
Details
Error if ref_pt
is not a dppMatrix
object or if conns
is not specified.
Returns
None
Method compute_conns()
This function computes the connectomes from the tangent images.
Usage
CSample$compute_conns()
Details
Error if tangent images are not specified.
Returns
None
Method compute_vecs()
This function computes the vectorized tangent images from the tangent images.
Usage
CSample$compute_vecs()
Details
Error if tangent images are not specified.
Returns
None
Method compute_unvecs()
This function computes the tangent images from the vector images.
Usage
CSample$compute_unvecs()
Details
Error if vec_imgs
is not specified.
Returns
None
Method compute_fmean()
This function computes the Frechet mean of the sample.
Usage
CSample$compute_fmean(tol = 0.05, max_iter = 20, lr = 0.2)
Arguments
tol
Tolerance for the convergence of the mean (default is 0.05).
max_iter
Maximum number of iterations for the computation (default is 20).
lr
Learning rate for the optimization algorithm (default is 0.2).
Returns
None
Method change_ref_pt()
This function changes the reference point for the tangent images.
Usage
CSample$change_ref_pt(new_ref_pt)
Arguments
new_ref_pt
A new reference point, which must be a
dppMatrix
object.
Details
Error if tangent images have not been computed or if new_ref_pt
is not a dppMatrix
object.
Returns
None
Method center()
Center the sample
Usage
CSample$center()
Details
This function centers the sample by computing the Frechet mean if it is not already computed, and then changing the reference point to the computed Frechet mean. Error if tangent images are not specified. Error if the sample is already centered.
Returns
None. This function is called for its side effects.
Method compute_variation()
Compute Variation
Usage
CSample$compute_variation()
Details
This function computes the variation of the sample. It first checks if the vector images are null, and if so, it computes the vectors, computing first the tangent images if necessary. If the sample is not centered, it centers the sample and recomputes the vectors. Finally, it calculates the variation as the mean of the sum of squares of the vector images. Error if vec_imgs
is not specified.
Returns
None. This function is called for its side effects.
Method compute_sample_cov()
Compute Sample Covariance
Usage
CSample$compute_sample_cov()
Details
This function computes the sample covariance matrix for the vector images. It first checks if the vector images are null, and if so, it computes the vectors, computing first the tangent images if necessary.
Returns
None. This function is called for its side effects.
Method clone()
The objects of this class are cloneable with this method.
Usage
CSample$clone(deep = FALSE)
Arguments
deep
Whether to make a deep clone.
TangentImageHandler Class
Description
TangentImageHandler Class
TangentImageHandler Class
Details
This class handles tangent images on a manifold. It provides methods to set a reference point, compute tangents, and perform various operations using a provided metric. Initialize the TangentImageHandler
Error if the tangent images have not been specified
Error if the reference point is not an object of class dppMatrix
Error if the matrix is not of type dspMatrix Tangent images getter
Active bindings
ref_point
A matrix of type dppMatrix
tangent_images
A list of dspMatrix objects
Methods
Public methods
Method new()
Usage
TangentImageHandler$new(metric_obj, reference_point = NULL)
Arguments
metric_obj
An rmetric object for operations.
reference_point
An optional reference point on the manifold.
Returns
A new instance of TangentImageHandler. Set a new reference point.
Method set_reference_point()
If tangent images have been created, it recomputes them by mapping to the manifold and then to the new tangent space.
Usage
TangentImageHandler$set_reference_point(new_ref_pt)
Arguments
new_ref_pt
A new reference point of class dppMatrix.
Returns
None. Computes the tangent images from the points in the manifold
Method compute_tangents()
Usage
TangentImageHandler$compute_tangents(manifold_points)
Arguments
manifold_points
A list of connectomes
Returns
None Computes vectorizations from tangent images
Method compute_vecs()
Usage
TangentImageHandler$compute_vecs()
Returns
A matrix, each row of which is a vectorization Computes connectomes from tangent images
Method compute_conns()
Usage
TangentImageHandler$compute_conns()
Returns
A list of connectomes Setter for the tangent images
Method set_tangent_images()
Usage
TangentImageHandler$set_tangent_images(reference_point, tangent_images)
Arguments
reference_point
A connectome
tangent_images
A list of tangent images
Returns
None Appends a matrix to the list of tangent images
Method add_tangent_image()
Usage
TangentImageHandler$add_tangent_image(image)
Arguments
image
Matrix to be added
Method get_tangent_images()
Usage
TangentImageHandler$get_tangent_images()
Returns
list of tangent matrices Wrapper for set_reference_point
Method relocate_tangents()
Usage
TangentImageHandler$relocate_tangents(new_ref_pt)
Arguments
new_ref_pt
The new reference point
Returns
None
Method clone()
The objects of this class are cloneable with this method.
Usage
TangentImageHandler$clone(deep = FALSE)
Arguments
deep
Whether to make a deep clone.
Compute the AIRM Exponential
Description
This function computes the Riemannian exponential map for the Affine-Invariant Riemannian Metric (AIRM).
Usage
airm_exp(sigma, v)
Arguments
sigma |
A symmetric positive-definite matrix of class |
v |
A tangent vector of class |
Value
A symmetric positive-definite matrix of class dppMatrix
.
Examples
if (requireNamespace("Matrix", quietly = TRUE)) {
library(Matrix)
sigma <- diag(2) |>
Matrix::nearPD() |>
_$mat |>
Matrix::pack()
v <- diag(c(1, 0.5)) |>
Matrix::symmpart() |>
Matrix::pack()
airm_exp(sigma, v)
}
Compute the AIRM Logarithm
Description
This function computes the Riemannian logarithmic map for the Affine-Invariant Riemannian Metric (AIRM).
Usage
airm_log(sigma, lambda)
Arguments
sigma |
A symmetric positive-definite matrix of class |
lambda |
A symmetric positive-definite matrix of class |
Value
A symmetric matrix of class dspMatrix
, representing the tangent space image of lambda
at sigma
.
Examples
if (requireNamespace("Matrix", quietly = TRUE)) {
library(Matrix)
sigma <- diag(2) |>
Matrix::nearPD() |>
_$mat |>
Matrix::pack()
lambda <- diag(c(2, 3)) |>
Matrix::nearPD() |>
_$mat |>
Matrix::pack()
airm_log(sigma, lambda)
}
Compute the Inverse Vectorization (AIRM)
Description
Converts a vector back into a tangent matrix relative to a reference point using AIRM.
Usage
airm_unvec(sigma, w)
Arguments
sigma |
A symmetric positive-definite matrix of class |
w |
A numeric vector, representing the vectorized tangent image. |
Value
A symmetric matrix of class dspMatrix
, representing the tangent vector.
Examples
if (requireNamespace("Matrix", quietly = TRUE)) {
library(Matrix)
sigma <- diag(2) |>
Matrix::nearPD() |>
_$mat |>
Matrix::pack()
w <- c(1, sqrt(2), 2)
airm_unvec(sigma, w)
}
Compute the AIRM Vectorization of Tangent Space
Description
Vectorizes a tangent matrix into a vector in Euclidean space using AIRM.
Usage
airm_vec(sigma, v)
Arguments
sigma |
A symmetric positive-definite matrix of class |
v |
A symmetric matrix of class |
Value
A numeric vector, representing the vectorized tangent image.
Examples
if (requireNamespace("Matrix", quietly = TRUE)) {
library(Matrix)
sigma <- diag(2) |>
Matrix::nearPD() |>
_$mat |>
Matrix::pack()
v <- diag(c(1, 0.5)) |>
Matrix::symmpart() |>
Matrix::pack()
airm_vec(sigma, v)
}
Compute the Bures-Wasserstein Exponential
Description
This function computes the Riemannian exponential map using the Bures-Wasserstein metric for symmetric positive-definite matrices. The map operates by solving a Lyapunov equation and then constructing the exponential.
Usage
bures_wasserstein_exp(sigma, v)
Arguments
sigma |
A symmetric positive-definite matrix of class |
v |
A symmetric matrix of class |
Value
A symmetric positive-definite matrix of class dppMatrix
, representing the point on the manifold.
Compute the Bures-Wasserstein Logarithm
Description
This function computes the Riemannian logarithmic map using the Bures-Wasserstein metric for symmetric positive-definite matrices.
Usage
bures_wasserstein_log(sigma, lambda)
Arguments
sigma |
A symmetric positive-definite matrix of class |
lambda |
A symmetric positive-definite matrix of class |
Value
A symmetric matrix of class dspMatrix
, representing the tangent space image of lambda
at sigma
.
Compute the Bures-Wasserstein Inverse Vectorization
Description
Compute the Bures-Wasserstein Inverse Vectorization
Usage
bures_wasserstein_unvec(sigma, w)
Arguments
sigma |
A symmetric positive-definite matrix of class |
w |
A numeric vector representing the vectorized tangent image |
Value
A symmetric matrix of class dspMatrix
Compute the Bures-Wasserstein Vectorization
Description
Compute the Bures-Wasserstein Vectorization
Usage
bures_wasserstein_vec(sigma, v)
Arguments
sigma |
A symmetric positive-definite matrix of class |
v |
A symmetric matrix of class |
Value
A numeric vector representing the vectorized tangent image
Compute the Frechet Mean
Description
This function computes the Frechet mean of a sample using an iterative algorithm.
Usage
compute_frechet_mean(sample, tol = 0.05, max_iter = 20, lr = 0.2)
Arguments
sample |
An object of class |
tol |
A numeric value specifying the tolerance for convergence. Default is 0.05. |
max_iter |
An integer specifying the maximum number of iterations. Default is 20. |
lr |
A numeric value specifying the learning rate. Default is 0.2. |
Details
The function iteratively updates the reference point of the sample until the change in the reference point is less than the specified tolerance or the maximum number of iterations is reached. If the tangent images are not already computed, they will be computed before starting the iterations.
Value
The computed Frechet mean.
Examples
if (requireNamespace("Matrix", quietly = TRUE)) {
library(Matrix)
# Load the AIRM metric object
data(airm)
# Create a CSample object with example data
conns <- list(
diag(2) |> Matrix::nearPD() |> _$mat |> Matrix::pack(),
diag(c(2, 3)) |> Matrix::nearPD() |> _$mat |> Matrix::pack()
)
sample <- CSample$new(conns = conns, metric_obj = airm)
# Compute the Frechet mean
compute_frechet_mean(sample, tol = 0.01, max_iter = 50, lr = 0.1)
}
Default reference point
Description
Default reference point
Usage
default_ref_pt(p)
Arguments
p |
the dimension |
Value
A diagonal matrix of the desired dimension
Differential of Matrix Exponential Map
Description
Computes the differential of the matrix exponential map located at a point a, evaluated at x
Usage
dexp(a, x)
Arguments
a |
A symmetric matrix of class dspMatrix |
x |
A symmetric matrix representing tangent vector of class dspMatrix |
Value
A positive definite symmetric matrix representing the differential located at a and evaluated at x, of class dppMatrix
Differential of Matrix Logarithm Map
Description
Computes the differential of the matrix logarithm map at a point Sigma, evaluated at H
Usage
dlog(sigma, h)
Arguments
sigma |
A symmetric positive definite matrix of class dspMatrix |
h |
A symmetric matrix representing tangent vector of class dsyMatrix |
Value
A symmetric matrix representing the differential evaluated at H of class dsyMatrix
Compute the Euclidean Exponential
Description
Compute the Euclidean Exponential
Usage
euclidean_exp(sigma, v)
Arguments
sigma |
A reference point. |
v |
A tangent vector to be mapped back to the manifold at |
Value
The point on the manifold corresponding to the tangent vector at sigma
.
Compute the Euclidean Logarithm
Description
Compute the Euclidean Logarithm
Usage
euclidean_log(sigma, lambda)
Arguments
sigma |
A reference point. |
lambda |
A point on the manifold. |
Value
The tangent space image of lambda
at sigma
.
Compute the Inverse Vectorization (Euclidean)
Description
Converts a vector back into a tangent matrix relative to a reference point using Euclidean metric.
Usage
euclidean_unvec(sigma, w)
Arguments
sigma |
A symmetric matrix. |
w |
A numeric vector, representing the vectorized tangent image. |
Value
A symmetric matrix, representing the tangent vector.
Vectorize at Identity Matrix (Euclidean)
Description
Converts a symmetric matrix into a vector representation.
Usage
euclidean_vec(sigma, v)
Arguments
sigma |
A symmetric matrix. |
v |
A vector. |
Value
A numeric vector, representing the vectorized tangent image.
Half-underscore operation for use in the log-Cholesky metric
Description
Half-underscore operation for use in the log-Cholesky metric
Usage
half_underscore(x)
Arguments
x |
A symmetric matrix (object of class dsyMatrix) |
Value
The strictly lower triangular part of the matrix, plus half its diagonal part
Create an Identity Matrix
Description
Create an Identity Matrix
Usage
id_matr(sigma)
Arguments
sigma |
A matrix. |
Value
An identity matrix of the same dimensions as sigma
.
Compute the Log-Cholesky Exponential
Description
This function computes the Riemannian exponential map using the Log-Cholesky metric for symmetric positive-definite matrices. The map operates by transforming the tangent vector via Cholesky decomposition of the reference point.
Usage
log_cholesky_exp(sigma, v)
Arguments
sigma |
A symmetric positive-definite matrix of class |
v |
A symmetric matrix of class |
Value
A symmetric positive-definite matrix of class dppMatrix
, representing the point on the manifold.
Compute the Log-Cholesky Logarithm
Description
This function computes the Riemannian logarithmic map using the Log-Cholesky metric for symmetric positive-definite matrices. The Log-Cholesky metric operates by transforming matrices via their Cholesky decomposition.
Usage
log_cholesky_log(sigma, lambda)
Arguments
sigma |
A symmetric positive-definite matrix of class |
lambda |
A symmetric positive-definite matrix of class |
Value
A symmetric matrix of class dspMatrix
, representing the tangent space image of lambda
at sigma
.
Compute the Log-Cholesky Inverse Vectorization
Description
Compute the Log-Cholesky Inverse Vectorization
Usage
log_cholesky_unvec(sigma, w)
Arguments
sigma |
A symmetric positive-definite matrix of class |
w |
A numeric vector representing the vectorized tangent image |
Value
A symmetric matrix of class dspMatrix
Compute the Log-Cholesky Vectorization
Description
Compute the Log-Cholesky Vectorization
Usage
log_cholesky_vec(sigma, v)
Arguments
sigma |
A symmetric positive-definite matrix of class |
v |
A symmetric matrix of class |
Value
A numeric vector representing the vectorized tangent image
Compute the Log-Euclidean Exponential
Description
This function computes the Euclidean exponential map.
Usage
log_euclidean_exp(ref_pt, v)
Arguments
ref_pt |
A reference point. |
v |
A tangent vector to be mapped back to the manifold at |
Value
The point on the manifold corresponding to the tangent vector at ref_pt
.
Compute the Log-Euclidean Logarithm
Description
Compute the Log-Euclidean Logarithm
Usage
log_euclidean_log(sigma, lambda)
Arguments
sigma |
A reference point. |
lambda |
A point on the manifold. |
Value
The tangent space image of lambda
at sigma
.
Compute the Inverse Vectorization (Euclidean)
Description
Converts a vector back into a tangent matrix relative to a reference point using Euclidean metric.
Usage
log_euclidean_unvec(sigma, w)
Arguments
sigma |
A symmetric matrix. |
w |
A numeric vector, representing the vectorized tangent image. |
Value
A symmetric matrix, representing the tangent vector.
Vectorize at Identity Matrix (Euclidean)
Description
Converts a symmetric matrix into a vector representation.
Usage
log_euclidean_vec(sigma, v)
Arguments
sigma |
A symmetric matrix. |
v |
A vector. |
Value
A numeric vector, representing the vectorized tangent image.
Metric Object Constructor
Description
Constructs a metric object that contains the necessary functions for Riemannian operations.
Usage
metric(log, exp, vec, unvec)
Arguments
log |
A function representing the Riemannian logarithmic map. This function should accept a |
exp |
A function representing the Riemannian exponential map. This function should accept a |
vec |
A function representing the vectorization operation for tangent spaces. This function should accept a |
unvec |
A function representing the inverse of the vectorization operation. This function should accept a |
Value
An object of class rmetric
containing the specified functions.
Pre-configured Riemannian metrics for SPD matrices
Description
Ready-to-use metric objects for various Riemannian geometries on the manifold of symmetric positive definite matrices.
Usage
airm
log_euclidean
euclidean
log_cholesky
bures_wasserstein
Format
Objects of class rmetric
containing four functions:
- log
Computes the Riemannian logarithm
- exp
Computes the Riemannian exponential
- vec
Performs vectorization
- unvec
Performs inverse vectorization
An object of class rmetric
of length 4.
An object of class rmetric
of length 4.
An object of class rmetric
of length 4.
An object of class rmetric
of length 4.
An object of class rmetric
of length 4.
Relocate Tangent Representations to a New Reference Point
Description
Changes the reference point for tangent space representations on a Riemannian manifold.
Usage
relocate(old_ref, new_ref, images, met)
Arguments
old_ref |
A reference point on the manifold to be replaced. Must be an object of class |
new_ref |
The new reference point on the manifold. Must be an object of class |
images |
A list of tangent representations relative to the old reference point. Each element in the list must be an object of class |
met |
A metric object of class |
Value
A list of tangent representations relative to the new reference point. Each element in the returned list will be an object of class dspMatrix
.
Examples
if (requireNamespace("Matrix", quietly = TRUE)) {
library(Matrix)
data(airm)
old_ref <- diag(2) |>
Matrix::nearPD() |>
_$mat |>
Matrix::pack()
new_ref <- diag(c(2, 3)) |>
Matrix::nearPD() |>
_$mat |>
Matrix::pack()
images <- list(
diag(2) |> Matrix::symmpart() |> Matrix::pack(),
diag(c(1, 0.5)) |> Matrix::symmpart() |> Matrix::pack()
)
relocate(old_ref, new_ref, images, airm)
}
Generate Random Samples from a Riemannian Normal Distribution
Description
Simulates random samples from a Riemannian normal distribution on symmetric positive definite matrices.
Usage
rspdnorm(n, refpt, disp, met)
Arguments
n |
Number of samples to generate. |
refpt |
Reference point on the manifold, represented as a symmetric positive definite matrix. Must be an object of class |
disp |
Dispersion matrix defining the spread of the distribution. Must be an object of class |
met |
A metric object of class |
Value
An object of class CSample
containing the generated samples.
Examples
if (requireNamespace("Matrix", quietly = TRUE)) {
library(Matrix)
data(airm)
refpt <- diag(2) |>
Matrix::nearPD() |>
_$mat |>
Matrix::pack()
disp <- diag(3) |>
Matrix::nearPD() |>
_$mat |>
Matrix::pack()
rspdnorm(10, refpt, disp, airm)
}
Wrapper for the matrix logarithm
Description
Wrapper for the matrix logarithm
Usage
safe_logm(x)
Arguments
x |
A matrix |
Value
Its matrix logarithm
Solve the Lyapunov Equation
Description
Solves the Lyapunov equation L P + P L = V for L.
Usage
solve_lyapunov(p, v)
Arguments
p |
A symmetric positive-definite matrix of class |
v |
A symmetric matrix of class |
Value
A symmetric matrix of class dspMatrix
, representing the solution L.
Reverse isometry from tangent space at identity to tangent space at P
Description
Reverse isometry from tangent space at identity to tangent space at P
Usage
spd_isometry_from_identity(sigma, v)
Arguments
sigma |
A symmetric positive-definite matrix of class |
v |
A symmetric matrix of class |
Value
A symmetric matrix of class dspMatrix
Isometry from tangent space at P to tangent space at identity
Description
Isometry from tangent space at P to tangent space at identity
Usage
spd_isometry_to_identity(sigma, v)
Arguments
sigma |
A symmetric positive-definite matrix of class |
v |
A symmetric matrix of class |
Value
A symmetric matrix of class dspMatrix
Validate Connections
Description
Validates the connections input.
Usage
validate_conns(conns, tan_imgs, vec_imgs, centered)
Arguments
conns |
List of connection matrices. |
tan_imgs |
List of tangent images. |
vec_imgs |
Matrix of vector images. |
centered |
Logical indicating if the data is centered. |
Value
None. Throws an error if the validation fails.
Validate arguments for Riemannian logarithms
Description
Validate arguments for Riemannian logarithms
Usage
validate_exp_args(sigma, v)
Arguments
sigma |
A dppMatrix object |
v |
A dspMatrix object |
Details
Error if sigma and lambda are not of the same dimensions
Value
None
Validate arguments for Riemannian logarithms
Description
Validate arguments for Riemannian logarithms
Usage
validate_log_args(sigma, lambda)
Arguments
sigma |
A dppMatrix object |
lambda |
A dppMatrix object |
Details
Error if sigma and lambda are not of the same dimensions
Value
None
Validate Metric
Description
Validates that the metric is not NULL.
Usage
validate_metric(metric)
Arguments
metric |
The metric to validate. |
Value
None. Throws an error if the metric is NULL.
Validate Tangent Images
Description
Validates the tangent images input.
Usage
validate_tan_imgs(tan_imgs, vec_imgs, centered)
Arguments
tan_imgs |
List of tangent images. |
vec_imgs |
List of vector images. |
centered |
Logical indicating if the data is centered. |
Value
None. Throws an error if the validation fails.
Validate arguments for inverse vectorization
Description
Validate arguments for inverse vectorization
Usage
validate_unvec_args(sigma, w)
Arguments
sigma |
A dppMatrix object |
w |
A numeric vector |
Details
Error if the dimensionalities don't match
Value
None
Validate arguments for vectorization
Description
Validate arguments for vectorization
Usage
validate_vec_args(sigma, v)
Arguments
sigma |
A dppMatrix object |
v |
A dspMatrix object |
Details
Error if sigma and v are not of the same dimensions
Value
None
Validate Vector Images
Description
Validates the vector images input.
Usage
validate_vec_imgs(vec_imgs, centered)
Arguments
vec_imgs |
List of vector images. |
centered |
Logical indicating if the data is centered. |
Value
None. Throws an error if the validation fails.
Vectorize at Identity Matrix
Description
Converts a symmetric matrix into a vector representation specific to operations at the identity matrix.
Usage
vec_at_id(v)
Arguments
v |
A symmetric matrix of class |
Value
A numeric vector, representing the vectorized tangent image.
Examples
if (requireNamespace("Matrix", quietly = TRUE)) {
library(Matrix)
v <- diag(c(1, sqrt(2))) |>
Matrix::symmpart() |>
Matrix::pack()
vec_at_id(v)
}