Factorizing a dataset

Classes:

`ADMMVars`(auxes, duals)
`DiagnosticMetrics`(rec_errors, ...)

Functions:

`cmf_aoadmm`(matrices, rank[, init, ...])	Fit a regularized coupled matrix factorization model with AO-ADMM
`compute_feasibility_gaps`(cmf, regs, ...)	Compute all feasibility gaps.
`parafac2_aoadmm`(matrices, rank[, init, ...])	Alias for cmf_aoadmm with PARAFAC2 constraint (constant cross-product) on mode 1 (B mode)

class matcouply.decomposition.ADMMVars(auxes, duals)[source]

Attributes:

`auxes`	Length three tuple containing a list of auxiliary factor matrices for each mode
`duals`	Length three tuple containing a list of dual variables for each mode

auxes: tuple: Length three tuple containing a list of auxiliary factor matrices for each mode

duals: tuple: Length three tuple containing a list of dual variables for each mode

class matcouply.decomposition.DiagnosticMetrics(rec_errors, feasibility_gaps, regularized_loss, satisfied_stopping_condition, satisfied_feasibility_condition, n_iter, message)[source]

Attributes:

`feasibility_gaps`	List of length equal to the number of iterations plus one containing the feasibility gaps
`message`	Convergence message
`n_iter`	Number of iterations ran
`rec_errors`	List of length equal to the number of iterations plus one containing the reconstruction errors
`regularized_loss`	List of length equal to the number of iterations plus one containing the regularized loss
`satisfied_feasibility_condition`	Boolean specifying whether the feasibility conditions were satisfied, None if no tolerance is set
`satisfied_stopping_condition`	Boolean specifying whether the stopping conditions were satisfied, None if no tolerance is set

feasibility_gaps: list: List of length equal to the number of iterations plus one containing the feasibility gaps

message: str: Convergence message

n_iter: int: Number of iterations ran

rec_errors: list: List of length equal to the number of iterations plus one containing the reconstruction errors

regularized_loss: list: List of length equal to the number of iterations plus one containing the regularized loss

satisfied_feasibility_condition: bool | None: Boolean specifying whether the feasibility conditions were satisfied, None if no tolerance is set

satisfied_stopping_condition: bool | None: Boolean specifying whether the stopping conditions were satisfied, None if no tolerance is set

matcouply.decomposition.cmf_aoadmm(matrices, rank, init='random', n_iter_max=1000, l2_penalty=None, tv_penalty=None, l1_penalty=None, non_negative=None, unimodal=None, generalized_l2_penalty=None, l2_norm_bound=None, lower_bound=None, upper_bound=None, parafac2=None, regs=None, feasibility_penalty_scale=1, constant_feasibility_penalty=False, aux_init='random_uniform', dual_init='random_uniform', svd='truncated_svd', init_params=None, random_state=None, tol=1e-08, absolute_tol=1e-10, feasibility_tol=0.0001, inner_tol=None, inner_n_iter_max=5, update_A=True, update_B_is=True, update_C=True, return_admm_vars=False, return_errors=False, verbose=False)[source]

Fit a regularized coupled matrix factorization model with AO-ADMM

Regularization parameters can be:

None (no regularization)

A single regularization parameter (same regularization in all modes)

A list with specific regularization parameters (one for each mode)

A dictionary with mode-index (0, 1 or 2) as keys and regularization parameter as values (regularization only in specific modes)

\(\mathbf{A}\) is represented by mode-index 0, \(\{\mathbf{B}^{(i)}\}_{i=1}^I\) is represented by mode-index 1 and \(\mathbf{C}\) is represented by mode-index 2.

Parameters:

matrices (list of ndarray) – Matrices that the coupled matrix factorization should model
rank (int) – The rank of the model
init ({"random", "svd", "threshold_svd", "parafac_als", "parafac_hals", "parafac2_als"} (default="random")) – Initialization method. If equal to "parafac_als", "parafac_hals" or "parafac2_als", then the corresponding methods in TensorLy are used to initialize the model. In that case, you can pass additional keyword-arguments to the decomposition function with the init_params parameter.
n_iter_max (int (default=1000)) – Maximum number of iterations.
l2_penalty (Regularization parameter (default=None)) – Strength of the L2 penalty, imposed as 0.5 * l2_penalty * tl.sum(M**2), where M represents a single factor matrix (note that this differs by a factor \(0.5\) compared to the expression in [RSC+22, RSCA21]).
tv_penalty (Regularization parameter (default=None)) – Strength of the TV penalty. To use this regularizer, you must have the GPL-lisenced library: condat_tv installed.
l1_penalty (Regularization parameter (default=None)) – Strength of the sparsity inducing regularization
non_negative (Regularization parameter (default=None)) – If True, then the corresponding factor matrices are non-negative
unimodal (Regularization parameter (default=None)) – If True, then the corresponding factor matrices have unimodal columns
generalized_l2_penalty (Regularization parameter (must be None, list or dict) (default=None)) – List or dict containing square matrices that parametrize a generalized L2 norm.
l2_norm_bound (Regularization parameter (default=None)) – Maximum L2 norm of the columns
lower_bound (Regularization parameter (default=None)) – Lower value of the columns
upper_bound (Regularization parameter (default=None)) – Upper value of the columns
parafac2 (bool (default=False)) – If True, then the \(\mathbf{B}^{(i)}\)-matrices follow the PARAFAC2 constraint
regs (List of lists of matcouply.penalties.ADMMPenalty (optional, default=None)) – The first element of this list contains a list of matcouply.penalties.ADMMPenalty-instances for \(\mathbf{A}\), the second for \(\{\mathbf{B}^{(i)}\}_{i=1}^I\) and the third for \(\mathbf{C}\).
feasibility_penalty_scale (float (default=1)) – What function to multiply the automatically computed feasibility penalty parameter by (see [RSC+22] for details)
constant_feasibility_penalty (bool, "A" or "B" (default=False)) –
If True, then all rows of \(\mathbf{A}\) will have the same feasibility penalty parameter and all \(\mathbf{B}^{(i)}\)-matrices will have the same feasibility penalty parameter. This makes it possible to use matrix-penalties for \(\mathbf{A}\) and multi-matrix penalties that require the same penalty parameter for all \(\mathbf{B}^{(i)}\)-matrices.

The maximum penalty parameter for all rows of \(\mathbf{A}\) are used as the penalty parameter for \(\mathbf{A}\) and the maximum penalty parameter among all \(\mathbf{B}^{(i)}\)-matrices are used as the penalty parameter for \(\{\mathbf{B}^{(i)}\}_{i=1}^I\).

If "A" or "B", then this option is only enabled for the specified parameters.
aux_init (str (default="random_uniform")) – Initialization scheme for the auxiliary variables. See matcouply.penalties.ADMMPenalty.init_aux() for possible options.
dual_init (str (default="random_uniform")) – Initialization scheme for the dual variables. See matcouply.penalties.ADMMPenalty.init_dual() for possible options.
svd (str (default="truncated_svd")) – String that specifies which SVD algorithm to use. Valid strings are the keys of tensorly.SVD_FUNS.
init_params (dict) – Keyword arguments to use for the initialization.
random_state (None, int or valid tensorly random state) –
tol (float (default=1e-8)) – Relative loss decrease condition. For stopping, the optimization requires that abs(losses[-2] - losses[-1]) / losses[-2] < tol or that losses[-2] < absolute_tol.
absolute_tol (float (default=1e-10)) – Absolute loss value condition. For stopping, the optimization requires that abs(losses[-2] - losses[-1]) / losses[-2] < tol or that losses[-2] < absolute_tol.
feasibility_tol (float (default=1e-4)) – If any feasibility gap is greater than feasibility_tol, then the optimization will not stop (unless the maximum number of iterations are reached).
inner_tol (float (optional, default=None)) – If set, then this specifies the stopping condition for the inner subproblems, solved with ADMM. If it is None, then the ADMM algorithm will always run for inner_n_iter_max-iterations.
inner_n_iter_max (int (default=5)) – Maximum number of iterations for the ADMM subproblems
update_A (bool (default=True)) – If False, then \(\mathbf{A}\) will not be updated.
update_B_is (bool (default=True)) – If False, then the \(\mathbf{B}^{(i)}\)-matrices will not be updated.
update_C (bool (default=True)) – If False, then \(\mathbf{C}\) will not be updated.
return_admm_vars (bool (default=False)) – If True, then the auxiliary and dual variables will also be returned.
return_errors (bool (default=False)) – If True, then the reconstruction error, feasibility gaps and loss per iteration will be returned.
verbose (bool or int (default=False)) – If True, then a message with convergence info will be printed out every iteration. If int > 1, then a message with convergence info will be printed out ever verbose iteration.

Returns:

CoupledMatrixFactorization – The fitted coupled matrix factorization
ADMMVars – (only returned if return_admm_vars=True) NamedTuple containing the auxuliary and dual variables. The feasibility gaps are computed by taking the relative normed differences between the auxiliary variables and the factor matrices of the model.
DiagnosticMetrics – (only returned if return_errors=True) NamedTuple containing lists of the relative norm-error, the feasibility gaps and the regularized loss for each iteration, each with length equal to the number of iterations plus one (the initial values), a stopping message and two boolean values, one signifying whether the stopping conditions (including feasibility gap) were satisfied and one signifying whether the feasibility gaps were sufficiently low (according to feasibility_tol).

Note

If the maximum number of iterations is reached, then you should check that the feasibility gaps are sufficiently small (see Fitting coupled matrix factorization models for more details)

Note

If you use norm-dependent regularization (e.g. generalized_l2_penalty) in one mode, then you must use norm-dependent regularization in all modes. You may, for example, use l2_penalty, norm_bound, l1_penalty or lower_bound and upper_bound. See [RSC+22] for more details.

Note

For simplicity, the model used here is a permuted version of that in [RSC+22, RSCA21] (where \(\mathbf{B_k}\)-matrices vary over the rows in \(\mathbf{C}\) instead of the rows in \(\mathbf{A}\)).

Examples

Here is a small example that shows what the diagnostic metrics contain. For more detailed examples, see Gallery of examples

>>> import matcouply
>>> data = matcouply.random.random_coupled_matrices([(10, 5), (12, 5), (11, 5)], rank=2, full=True, random_state=1)
>>> cmf, diagnostics = matcouply.decomposition.cmf_aoadmm(data, 2, n_iter_max=4, non_negative=True, return_errors=True, random_state=10)
>>> diagnostics.message
'MAXIMUM NUMBER OF ITERATIONS REACHED'

>>> diagnostics.satisfied_stopping_condition
False

>>> diagnostics.satisfied_feasibility_condition
True

>>> len(diagnostics.regularized_loss)
5

We see that the length of the regularized loss list is the number of iterations plus one. This is because the initial decomposition is included in this.

>>> len(diagnostics.feasibility_gaps[0])
3

The feasibility gaps list contain tuples (one for each mode) of tuples (one for each penalty for the given mode).

>>> len(diagnostics.feasibility_gaps[0][0])
1

matcouply.decomposition.compute_feasibility_gaps(cmf, regs, A_aux_list, B_aux_list, C_aux_list)[source]

Compute all feasibility gaps.

The feasibility gaps for AO-ADMM are given by

\[\frac{\|\mathbf{x} - \mathbf{z}\|_2}{\|\mathbf{z}\|_2},\]

where \(\mathbf{x}\) is a component vector and \(\mathbf{z}\) is an auxiliary variable that represents a component vector. If a decomposition obtained with AO-ADMM is valid, then all feasibility gaps should be small compared to the components. If any of the feasibility penalties are large, then the decomposition may not satisfy the imposed constraints.

To compute the feasibility gap for the \(\mathbf{A}\) and \(\mathbf{C}\) matrices, we use the frobenius norm, and to compute the feasibility gap for the \(\mathbf{B}^{(i)}\)-matrices, we compute \(\sqrt{\sum_i \|\mathbf{B}^{(i)} - \mathbf{Z}^{(\mathbf{B}^{(i)})}\|^2}\), where \(\mathbf{Z}^{(\mathbf{B}^{(i)})}\|^2\) is the auxiliary variable for \(\mathbf{B}^{(i)}\).

Parameters:

cmf (CoupledMatrixFactorization - (weights, factors)) –
Coupled matrix factorization represented by weights and factors as described in What are coupled matrix factorizations?.
- weights1D array of shape (rank,) or None
  weights of the factors
- factorsList of factors of the coupled matrix decomposition
  List on the form [A, [B_0, B_1, ..., B_i], C], where A represents \(\mathbf{A}\), [B_0, B_1, ..., B_i] represents a list of all \(\mathbf{B}^{(i)}\)-matrices and C represents \(\mathbf{C}\)
regs (list of list of matcouply.penalties.ADMMPenalty) – The regs list should have three elements, the first being a list of penalties imposed on mode 0, the second being a list of penalties imposed on mode 1 and the last being a list of penalties imposed on mode 2.
A_aux_list (list) – A list of all auxiliary variables for the A-matrix
B_aux_list (list) – A list of all auxiliary variables for the B_is-matrices
C_aux_list (list) – A list of all auxiliary variables for the C-matrix

Returns:

list – Feasibility gaps for A
list – Feasibility gaps for B_is
list – Veasibility gaps for C

matcouply.decomposition.parafac2_aoadmm(matrices, rank, init='random', n_iter_max=1000, l2_penalty=0, tv_penalty=None, l1_penalty=None, non_negative=None, unimodal=None, generalized_l2_penalty=None, l2_norm_bound=None, lower_bound=None, upper_bound=None, regs=None, feasibility_penalty_scale=1, constant_feasibility_penalty=False, aux_init='random_uniform', dual_init='random_uniform', svd='truncated_svd', init_params=None, random_state=None, tol=1e-08, absolute_tol=1e-10, feasibility_tol=0.0001, inner_tol=None, inner_n_iter_max=5, update_A=True, update_B_is=True, update_C=True, return_errors=False, return_admm_vars=False, verbose=False)[source]

Alias for cmf_aoadmm with PARAFAC2 constraint (constant cross-product) on mode 1 (B mode)