# procrustes.symmetric

Symmetric Procrustes Module.

procrustes.symmetric.symmetric(a: numpy.ndarray, b: numpy.ndarray, pad: bool = True, translate: bool = False, scale: bool = False, unpad_col: bool = False, unpad_row: bool = False, check_finite: bool = True, weight: Optional[numpy.ndarray] = None, lapack_driver: str = 'gesvd') [source]

Perform symmetric Procrustes.

Given a matrix $$\mathbf{A}_{m \times n}$$ and a reference matrix $$\mathbf{B}_{m \times n}$$ with $$m \geqslant n$$, find the symmetrix transformation matrix $$\mathbf{X}_{n \times n}$$ that makes $$\mathbf{AX}$$ as close as possible to $$\mathbf{B}$$. In other words,

$\underbrace{\text{min}}_{\left\{\mathbf{X} \left| \mathbf{X} = \mathbf{X}^\dagger \right. \right\}} \|\mathbf{A} \mathbf{X} - \mathbf{B}\|_{F}^2$

This Procrustes method requires the $$\mathbf{A}$$ and $$\mathbf{B}$$ matrices to have the same shape with $$m \geqslant n$$, which is guaranteed with the default pad argument for any given $$\mathbf{A}$$ and $$\mathbf{B}$$ matrices. In preparing the $$\mathbf{A}$$ and $$\mathbf{B}$$ matrices, the (optional) order of operations is: 1) unpad zero rows/columns, 2) translate the matrices to the origin, 3) weight entries of $$\mathbf{A}$$, 4) scale the matrices to have unit norm, 5) pad matrices with zero rows/columns so they have the same shape.

Parameters
• a (ndarray) -- The 2D-array $$\mathbf{A}$$ which is going to be transformed.

• b (ndarray) -- The 2D-array $$\mathbf{B}$$ representing the reference matrix.

• pad (bool, optional) -- Add zero rows (at the bottom) and/or columns (to the right-hand side) of matrices $$\mathbf{A}$$ and $$\mathbf{B}$$ so that they have the same shape.

• translate (bool, optional) -- If True, both arrays are centered at origin (columns of the arrays will have mean zero).

• scale (bool, optional) -- If True, both arrays are normalized with respect to the Frobenius norm, i.e., $$\text{Tr}\left[\mathbf{A}^\dagger\mathbf{A}\right] = 1$$ and $$\text{Tr}\left[\mathbf{B}^\dagger\mathbf{B}\right] = 1$$.

• unpad_col (bool, optional) -- If True, zero columns (with values less than 1.0e-8) on the right-hand side of the intial $$\mathbf{A}$$ and $$\mathbf{B}$$ matrices are removed.

• unpad_row (bool, optional) -- If True, zero rows (with values less than 1.0e-8) at the bottom of the intial $$\mathbf{A}$$ and $$\mathbf{B}$$ matrices are removed.

• check_finite (bool, optional) -- If True, convert the input to an array, checking for NaNs or Infs.

• weight (ndarray, optional) -- The 1D-array representing the weights of each row of $$\mathbf{A}$$. This defines the elements of the diagonal matrix $$\mathbf{W}$$ that is multiplied by $$\mathbf{A}$$ matrix, i.e., $$\mathbf{A} \rightarrow \mathbf{WA}$$.

• lapack_driver ({'gesvd', 'gesdd'}, optional) -- Whether to use the more efficient divide-and-conquer approach ('gesdd') or the more robust general rectangular approach ('gesvd') to compute the singular-value decomposition with scipy.linalg.svd.

Returns

res -- The Procrustes result represented as a class:utils.ProcrustesResult object.

Return type

ProcrustesResult

Notes

The optimal symmetrix matrix is obtained by,

$\mathbf{X}_{\text{opt}} = \arg \underbrace{\text{min}}_{\left\{\mathbf{X} \left| \mathbf{X} = \mathbf{X}^\dagger \right. \right\}} \|\mathbf{A} \mathbf{X} - \mathbf{B}\|_{F}^2 = \underbrace{\text{min}}_{\left\{\mathbf{X} \left| \mathbf{X} = \mathbf{X}^\dagger \right. \right\}} \text{Tr}\left[\left(\mathbf{A}\mathbf{X} - \mathbf{B} \right)^\dagger \left(\mathbf{A}\mathbf{X} - \mathbf{B} \right)\right]$

Considering the singular value decomposition of $$\mathbf{A}$$,

$\mathbf{A}_{m \times n} = \mathbf{U}_{m \times m} \mathbf{\Sigma}_{m \times n} \mathbf{V}_{n \times n}^\dagger$

where $$\mathbf{\Sigma}_{m \times n}$$ is a rectangular diagonal matrix with non-negative singular values $$\sigma_i = [\mathbf{\Sigma}]_{ii}$$ listed in descending order, define

$\mathbf{C}_{m \times n} = \mathbf{U}_{m \times m}^\dagger \mathbf{B}_{m \times n} \mathbf{V}_{n \times n}$

with elements denoted by $$c_{ij}$$. Then we compute the symmetric matrix $$\mathbf{Y}_{n \times n}$$ with

$\begin{split}[\mathbf{Y}]_{ij} = \begin{cases} 0 && i \text{ and } j > \text{rank} \left(\mathbf{A}\right) \\ \frac{\sigma_i c_{ij} + \sigma_j c_{ji}}{\sigma_i^2 + \sigma_j^2} && \text{otherwise} \end{cases}\end{split}$

It is worth noting that the first part of this definition only applies in the unusual case where $$\mathbf{A}$$ has rank less than $$n$$. The $$\mathbf{X}_\text{opt}$$ is given by

$\mathbf{X}_\text{opt} = \mathbf{V Y V}^{\dagger}$

Examples

>>> import numpy as np
>>> a = np.array([[5., 2., 8.],
...               [2., 2., 3.],
...               [1., 5., 6.],
...               [7., 3., 2.]])
>>> b = np.array([[ 52284.5, 209138. , 470560.5],
...               [ 22788.5,  91154. , 205096.5],
...               [ 46139.5, 184558. , 415255.5],
...               [ 22788.5,  91154. , 205096.5]])
>>> res = symmetric(a, b, pad=True, translate=True, scale=True)
>>> res.t   # symmetric transformation array
array([[0.0166352 , 0.06654081, 0.14971682],
[0.06654081, 0.26616324, 0.59886729],
[0.14971682, 0.59886729, 1.34745141]])
>>> res.error   # error
4.483083428047388e-31