Subspace Power Method for CP Decomposition of Symmetric Tensors

The function cp_spm computes the symmetric CP decomposition of a symmetric tensor using the Subspace Power Method (SPM). The symmetric CP decomposition is described in Symmetric CP Decomposition, while SPM is described in the following reference:

Contents

Create a sample problem

For consistency, we use the same example as in Symmetric CP Decomposition

d = 3; % order
n = 10; % size
r = 2; % true rank

rng(5); % Set random number generator state for consistent results

info = create_problem('Size', n*ones(d,1), 'Num_Factors', r, ...
    'Symmetric', 1:d, 'Factor_Generator', @rand, 'Lambda_Generator', @rand, 'Noise', 0.1);

X = info.Data;
M_true = info.Soln;
S_true = symktensor(M_true); % Convert from ktensor to symktensor

Check that the tensor is symmetric

issymmetric(X)
ans =

  logical

   1

Run CP-SPM

SPM estimates the rank by picking a cut-off of the eigenvalues of the tensor flattening. We plot these eigenvalues as follows:

cp_spm(X, 'rank_sel', 'plot_only');

We observe the first two eigenvalues contain most of the energy, and conclude the tensor is approximately rank 2

rng(5); % Set random number generator state for consistent results

tic
[S, info] = cp_spm(X, 2);
toc

fprintf('\n');
fprintf('Final function value: %.2g\n', fg_explicit(S, X, norm(X)^2));
fprintf('Check similarity score (1=perfect): %.2f\n', score(S, S_true));
fprintf('\n');
Elapsed time is 0.022587 seconds.

Final function value: 0.39
Check similarity score (1=perfect): 0.93

Compare with CP-SYM using L-BFGS from Poblano Toolbox

We compare SPM with cp_sym. Its options are explained in Symmetric CP Decomposition; This is the recommended way to run the method:

optparams = lbfgs('defaults'); % Get the optimization parameters
optparams.RelFuncTol = 1e-10; % Tighten the stopping tolerance
optparams.StopTol = 1e-6; % Tighten the stopping tolerance
rng(5); % Set random number generator state for consistent results

tic
[S,info] = cp_sym(X, 2,'unique',false,'l1param',0,'alg_options',optparams);
toc

fprintf('\n');
fprintf('Final function value: %.2g\n', fg_explicit(S, X, norm(X)^2));
fprintf('Stopping condition: %s\n', info.optout.ExitDescription);
fprintf('Check similarity score (1=perfect): %.2f\n', score(S, S_true));
fprintf('\n');
 Iter  FuncEvals       F(X)          ||G(X)||/N        
------ --------- ---------------- ----------------
     0         1       7.34429496       3.10499258
     1         8       0.40258352       0.08233575
     2        10       0.37608104       0.04310174
     3        12       0.36658247       0.01555076
     4        14       0.36528638       0.00530333
     5        16       0.36488180       0.00311884
     6        18       0.36484022       0.00174885
     7        21       0.36477178       0.00360656
     8        23       0.36472681       0.00208400
     9        26       0.36467068       0.00313366
    10        28       0.36461167       0.00116220
    11        30       0.36460241       0.00025735
    12        32       0.36460199       0.00012652
    13        34       0.36460188       0.00014958
    14        36       0.36460181       0.00004645
    15        38       0.36460180       0.00000968
    16        40       0.36460180       0.00000896
    17        42       0.36460180       0.00000387
    18        44       0.36460180       0.00000212
    19        46       0.36460180       0.00000117
Elapsed time is 0.142301 seconds.

Final function value: 0.36
Stopping condition: Relative change in F < RelFuncTol
Check similarity score (1=perfect): 0.99