tensorly.decomposition.parafac2

parafac2(tensor_slices, rank: int, n_iter_max: int = 2000, init='random', svd: Literal['truncated_svd', 'symeig_svd', 'randomized_svd'] = 'truncated_svd', normalize_factors: bool = False, tol: float = 1e-08, nn_modes: Sequence[int] | Literal['all'] | None = None, random_state=None, verbose: bool | int = False, return_errors: bool = False, n_iter_parafac: int = 5, linesearch: bool = True)[source]

PARAFAC2 decomposition [1] of a third order tensor via alternating least squares (ALS)

Computes a rank-rank PARAFAC2 decomposition of the third-order tensor defined by tensor_slices. The decomposition is on the form \((A [B_i] C)\) such that the i-th frontal slice, \(X_i\), of \(X\) is given by

\[X_i = B_i diag(a_i) C^T,\]

where \(diag(a_i)\) is the diagonal matrix whose nonzero entries are equal to the \(i\)-th row of the \(I \times R\) factor matrix \(A\), \(B_i\) is a \(J_i \times R\) factor matrix such that the cross product matrix \(B_{i_1}^T B_{i_1}\) is constant for all \(i\), and \(C\) is a \(K \times R\) factor matrix. To compute this decomposition, we reformulate the expression for \(B_i\) such that

\[B_i = P_i B,\]

where \(P_i\) is a \(J_i \times R\) orthogonal matrix and \(B\) is a \(R \times R\) matrix.

An alternative formulation of the PARAFAC2 decomposition is that the tensor element \(X_{ijk}\) is given by

\[X_{ijk} = \sum_{r=1}^R A_{ir} B_{ijr} C_{kr},\]

with the same constraints hold for \(B_i\) as above.

Parameters:
tensor_slicesndarray or list of ndarrays

Either a third order tensor or a list of second order tensors that may have different number of rows. Note that the second mode factor matrices are allowed to change over the first mode, not the third mode as some other implementations use (see note below).

rankint

Number of components.

n_iter_maxint, optional

(Default: 2000) Maximum number of iteration

Changed in version 0.6.1: Previously, the default maximum number of iterations was 100.

init{‘svd’, ‘random’, CPTensor, Parafac2Tensor}

Type of factor matrix initialization. See initialize_factors.

svdstr, default is ‘truncated_svd’

function to use to compute the SVD, acceptable values in tensorly.SVD_FUNS

normalize_factorsbool (optional)

If True, aggregate the weights of each factor in a 1D-tensor of shape (rank, ), which will contain the norms of the factors. Note that there may be some inaccuracies in the component weights.

tolfloat, optional

(Default: 1e-8) Relative reconstruction error decrease tolerance. The algorithm is considered to have converged when \(\left|\| X - \hat{X}_{n-1} \|^2 - \| X - \hat{X}_{n} \|^2\right| < \epsilon \| X - \hat{X}_{n-1} \|^2\). That is, when the relative change in sum of squared error is less than the tolerance.

Changed in version 0.6.1: Previously, the stopping condition was \(\left|\| X - \hat{X}_{n-1} \| - \| X - \hat{X}_{n} \|\right| < \epsilon\).

nn_modes: None, ‘all’ or array of integers

(Default: None) Used to specify which modes to impose non-negativity constraints on. We cannot impose non-negativity constraints on the the B-mode (mode 1) with the ALS algorithm, so if this mode is among the constrained modes, then a warning will be shown (see notes for more info).

random_state{None, int, np.random.RandomState}
verboseint, optional

Level of verbosity

return_errorsbool, optional

Activate return of iteration errors

n_iter_parafacint, optional

Number of PARAFAC iterations to perform for each PARAFAC2 iteration

linesearchbool, default is False

Whether to perform line search as proposed by Bro in his PhD dissertation [2] (similar to the PLSToolbox line search described in [3]).

Returns:
Parafac2Tensor(weight, factors, projection_matrices)
  • weights1D array of shape (rank, )

    all ones if normalize_factors is False (default), weights of the (normalized) factors otherwise

  • factorsList of factors of the CP decomposition element i is of shape

    (tensor.shape[i], rank)

  • projection_matricesList of projection matrices used to create evolving

    factors.

errorslist

A list of reconstruction errors at each iteration of the algorithms.

Notes

This formulation of the PARAFAC2 decomposition is slightly different from the one in [1]. The difference lies in that, here, the second mode changes over the first mode, whereas in [1], the second mode changes over the third mode. This change allows the function to accept both lists of matrices and a single nd-array as input without any mode reordering.

Because of the reformulation above, \(B_i = P_i B\), the \(B_i\) matrices cannot be constrained to be non-negative with ALS. If this mode is constrained to be non-negative, then \(B\) will be non-negative, but not the orthogonal P_i matrices. Consequently, the B_i matrices are unlikely to be non-negative.

References

[1] (1,2,3)

Kiers, H.A.L., ten Berge, J.M.F. and Bro, R. (1999), PARAFAC2—Part I. A direct fitting algorithm for the PARAFAC2 model. J. Chemometrics, 13: 275-294.

[2]

R. Bro, “Multi-Way Analysis in the Food Industry: Models, Algorithms, and Applications”, PhD., University of Amsterdam, 1998

[3]

H. Yu, D. Augustijn, R. Bro, “Accelerating PARAFAC2 algorithms for non-negative complex tensor decomposition.” Chemometrics and Intelligent Laboratory Systems 214 (2021): 104312.