Home > manopt > manifolds > fixedrank > fixedrankfactory_2factors.m

fixedrankfactory_2factors

PURPOSE

Manifold of m-by-n matrices of rank k with balanced quotient geometry.

SYNOPSIS

function M = fixedrankfactory_2factors(m, n, k)

DESCRIPTION

``` Manifold of m-by-n matrices of rank k with balanced quotient geometry.

function M = fixedrankfactory_2factors(m, n, k)

The first-order geometry follows the balanced quotient geometry described
in the paper,
"Linear regression under fixed-rank constraints: a Riemannian approach",
G. Meyer, S. Bonnabel and R. Sepulchre, ICML 2011.

The second-order geometry follows from the paper
"Fixed-rank matrix factorizations and Riemannian low-rank optimization",
B. Mishra, R. Meyer, S. Bonnabel and R. Sepulchre,
Computational Statistics, 29(3 - 4), pp. 591 - 621, 2014.

A point X on the manifold is represented as a structure with two
fields: L and R. The matrices L (mxk) and R (nxk) are full column-rank
matrices such that X = L*R'.

Tangent vectors are represented as a structure with two fields: L, R.

For first-order geometry, please cite the Manopt paper as well as the research paper:
@InProceedings{meyer2011linear,
Title        = {Linear regression under fixed-rank constraints: a {R}iemannian approach},
Author       = {Meyer, G. and Bonnabel, S. and Sepulchre, R.},
Booktitle    = {{28th International Conference on Machine Learning}},
Year         = {2011},
Organization = {{ICML}}
}

For second-order geometry, please cite the Manopt paper as well as the research paper:
@Article{mishra2014fixedrank,
Title   = {Fixed-rank matrix factorizations and {Riemannian} low-rank optimization},
Author  = {Mishra, B. and Meyer, G. and Bonnabel, S. and Sepulchre, R.},
Journal = {Computational Statistics},
Year    = {2014},
Number  = {3-4},
Pages   = {591--621},
Volume  = {29},
Doi     = {10.1007/s00180-013-0464-z}
}

CROSS-REFERENCE INFORMATION

This function calls:
• hashmd5 Computes the MD5 hash of input data.
• lincomb Computes a linear combination of tangent vectors in the Manopt framework.
• sylvester_nochecks Solve Sylvester equation without input checks.
This function is called by:

SOURCE CODE

```0001 function M = fixedrankfactory_2factors(m, n, k)
0002 % Manifold of m-by-n matrices of rank k with balanced quotient geometry.
0003 %
0004 % function M = fixedrankfactory_2factors(m, n, k)
0005 %
0006 % The first-order geometry follows the balanced quotient geometry described
0007 % in the paper,
0008 % "Linear regression under fixed-rank constraints: a Riemannian approach",
0009 % G. Meyer, S. Bonnabel and R. Sepulchre, ICML 2011.
0010 %
0011 % Paper link: http://www.icml-2011.org/papers/350_icmlpaper.pdf.
0012 %
0013 % The second-order geometry follows from the paper
0014 % "Fixed-rank matrix factorizations and Riemannian low-rank optimization",
0015 % B. Mishra, R. Meyer, S. Bonnabel and R. Sepulchre,
0016 % Computational Statistics, 29(3 - 4), pp. 591 - 621, 2014.
0017 %
0018 % A point X on the manifold is represented as a structure with two
0019 % fields: L and R. The matrices L (mxk) and R (nxk) are full column-rank
0020 % matrices such that X = L*R'.
0021 %
0022 % Tangent vectors are represented as a structure with two fields: L, R.
0023 %
0024 % For first-order geometry, please cite the Manopt paper as well as the research paper:
0025 %     @InProceedings{meyer2011linear,
0026 %       Title        = {Linear regression under fixed-rank constraints: a {R}iemannian approach},
0027 %       Author       = {Meyer, G. and Bonnabel, S. and Sepulchre, R.},
0028 %       Booktitle    = {{28th International Conference on Machine Learning}},
0029 %       Year         = {2011},
0030 %       Organization = {{ICML}}
0031 %     }
0032 %
0033 % For second-order geometry, please cite the Manopt paper as well as the research paper:
0034 %     @Article{mishra2014fixedrank,
0035 %       Title   = {Fixed-rank matrix factorizations and {Riemannian} low-rank optimization},
0036 %       Author  = {Mishra, B. and Meyer, G. and Bonnabel, S. and Sepulchre, R.},
0037 %       Journal = {Computational Statistics},
0038 %       Year    = {2014},
0039 %       Number  = {3-4},
0040 %       Pages   = {591--621},
0041 %       Volume  = {29},
0042 %       Doi     = {10.1007/s00180-013-0464-z}
0043 %     }
0044 %
0045 %
0046 % See also fixedrankembeddedfactory fixedrankfactory_3factors fixedrankfactory_2factors_preconditioned
0047
0048 % This file is part of Manopt: www.manopt.org.
0049 % Original author: Bamdev Mishra, Dec. 30, 2012.
0050 % Contributors:
0051 % Change log:
0052 %
0053 %   July 10, 2013 (NB):
0054 %       Added vec, mat, tangent, tangent2ambient.
0055 %
0056 %    July 03, 2015 (BM):
0057 %      Cosmetic changes including avoiding storing the inverse of a
0058 %       k-by-k matrix.
0059 %
0060 %   Apr. 17, 2018 (NB):
0061 %      Using built-in sylvester instead of lyap, which requires a toolbox.
0062 %
0063 %   Sep.  6, 2018 (NB):
0064 %       Removed M.exp() as it was not implemented.
0065
0066
0067     M.name = @() sprintf('LR'' quotient manifold of %dx%d matrices of rank %d', m, n, k);
0068
0069     M.dim = @() (m+n-k)*k;
0070
0071     % Some precomputations at the point X to be used in the inner product
0072     % (and pretty much everywhere else).
0073     function X = prepare(X)
0074         if ~all(isfield(X,{'LtL','RtR'}))
0075             L = X.L;
0076             R = X.R;
0077             X.LtL = L'*L;
0078             X.RtR = R'*R;
0079         end
0080     end
0081
0082     % Choice of the metric is motivated by the symmetry present in the
0083     % space. The metric is the natural Grassmannian metric on L and R.
0084     M.inner = @iproduct;
0085     function ip = iproduct(X, eta, zeta)
0086         X = prepare(X);
0087         ip = trace(X.LtL\(eta.L'*zeta.L)) + trace( X.RtR\(eta.R'*zeta.R));
0088     end
0089
0090     M.norm = @(X, eta) sqrt(M.inner(X, eta, eta));
0091
0092     M.dist = @(x, y) error('fixedrankfactory_2factors.dist not implemented yet.');
0093
0094     M.typicaldist = @() 10*k;
0095
0096     symm = @(M) .5*(M+M');
0097
0100         X = prepare(X);
0103     end
0104
0105     M.ehess2rhess = @ehess2rhess;
0106     function Hess = ehess2rhess(X, egrad, ehess, eta)
0107         X = prepare(X);
0108
0109         % Riemannian gradient computation.
0111
0112         % Directional derivative of the Riemannian gradient.
0113         Hess.L = ehess.L*X.LtL + 2*egrad.L*symm(eta.L'*X.L);
0114         Hess.R = ehess.R*X.RtR + 2*egrad.R*symm(eta.R'*X.R);
0115
0116         % We need a correction term for the non-constant metric.
0117         Hess.L = Hess.L - rgrad.L*(X.LtL\(symm(X.L'*eta.L))) - eta.L*(X.LtL\(symm(X.L'*rgrad.L))) + X.L*(X.LtL\(symm(eta.L'*rgrad.L)));
0118         Hess.R = Hess.R - rgrad.R*(X.RtR\(symm(X.R'*eta.R))) - eta.R*(X.RtR\(symm(X.R'*rgrad.R))) + X.R*(X.RtR\(symm(eta.R'*rgrad.R)));
0119
0120         % Projection onto the horizontal space.
0121         Hess = M.proj(X, Hess);
0122     end
0123
0124     M.proj = @projection;
0125     % Projection of the vector eta in the ambient space onto the horizontal space.
0126     function etaproj = projection(X, eta)
0127         X = prepare(X);
0128
0129         SS = (X.LtL)*(X.RtR);
0130         AS = (X.LtL)*(X.R'*eta.R) - (eta.L'*X.L)*(X.RtR);
0131         Omega = sylvester_nochecks(SS, SS, AS); % could be sped up since SS appears twice
0132         etaproj.L = eta.L + X.L*Omega';
0133         etaproj.R = eta.R - X.R*Omega;
0134     end
0135
0136     M.tangent = M.proj;
0137     M.tangent2ambient = @(X, eta) eta;
0138
0139     M.retr = @retraction;
0140     function Y = retraction(X, eta, t)
0141         if nargin < 3
0142             t = 1.0;
0143         end
0144
0145         Y.L = X.L + t*eta.L;
0146         Y.R = X.R + t*eta.R;
0147
0148         % Numerical conditioning step: A simpler version.
0149         % We need to ensure that L and R do not have very relative
0150         % skewed norms.
0151
0152         scaling = norm(X.L, 'fro')/norm(X.R, 'fro');
0153         scaling = sqrt(scaling);
0154         Y.L = Y.L / scaling;
0155         Y.R = Y.R * scaling;
0156
0157         % These are reused in the computation of the gradient and Hessian.
0158         Y = prepare(Y);
0159     end
0160
0161
0162     M.hash = @(X) ['z' hashmd5([X.L(:) ; X.R(:)])];
0163
0164     M.rand = @random;
0165     function X = random()
0166         % A random point on the total space.
0167         X.L = randn(m, k);
0168         X.R = randn(n, k);
0169         X = prepare(X);
0170     end
0171
0172     M.randvec = @randomvec;
0173     function eta = randomvec(X)
0174         % A random vector in the horizontal space.
0175         eta.L = randn(m, k);
0176         eta.R = randn(n, k);
0177         eta = projection(X, eta);
0178         nrm = M.norm(X, eta);
0179         eta.L = eta.L / nrm;
0180         eta.R = eta.R / nrm;
0181     end
0182
0183     M.lincomb = @lincomb;
0184
0185     M.zerovec = @(X) struct('L', zeros(m, k),'R', zeros(n, k));
0186
0187     M.transp = @(x1, x2, d) projection(x2, d);
0188
0189     % vec and mat are not isometries, because of the unusual inner metric.
0190     M.vec = @(X, U) [U.L(:) ; U.R(:)];
0191     M.mat = @(X, u) struct('L', reshape(u(1:(m*k)), m, k), ...
0192         'R', reshape(u((m*k+1):end), n, k));
0193     M.vecmatareisometries = @() false;
0194
0195 end
0196
0197 % Linear combination of tangent vectors.
0198 function d = lincomb(x, a1, d1, a2, d2) %#ok<INUSL>
0199
0200     if nargin == 3
0201         d.L = a1*d1.L;
0202         d.R = a1*d1.R;
0203     elseif nargin == 5
0204         d.L = a1*d1.L + a2*d2.L;
0205         d.R = a1*d1.R + a2*d2.R;
0206     else
0207         error('Bad use of fixedrankfactory_2factors.lincomb.');
0208     end
0209
0210 end
0211
0212
0213
0214
0215```

Generated on Mon 10-Sep-2018 11:48:06 by m2html © 2005