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.

 Paper link: http://www.icml-2011.org/papers/350_icmlpaper.pdf.

 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}
     }


 See also fixedrankembeddedfactory fixedrankfactory_3factors fixedrankfactory_2factors_preconditioned

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

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     
0061     M.name = @() sprintf('LR'' quotient manifold of %dx%d matrices of rank %d', m, n, k);
0062     
0063     M.dim = @() (m+n-k)*k;
0064     
0065     % Some precomputations at the point X to be used in the inner product
0066     % (and pretty much everywhere else).
0067     function X = prepare(X)
0068         if ~all(isfield(X,{'LtL','RtR'}))
0069             L = X.L;
0070             R = X.R;
0071             X.LtL = L'*L;
0072             X.RtR = R'*R;
0073         end
0074     end
0075     
0076     % Choice of the metric is motivated by the symmetry present in the
0077     % space. The metric is the natural Grassmannian metric on L and R.
0078     M.inner = @iproduct;
0079     function ip = iproduct(X, eta, zeta)
0080         X = prepare(X);
0081         ip = trace(X.LtL\(eta.L'*zeta.L)) + trace( X.RtR\(eta.R'*zeta.R));
0082     end
0083     
0084     M.norm = @(X, eta) sqrt(M.inner(X, eta, eta));
0085     
0086     M.dist = @(x, y) error('fixedrankfactory_2factors.dist not implemented yet.');
0087     
0088     M.typicaldist = @() 10*k;
0089     
0090     symm = @(M) .5*(M+M');
0091     
0092     M.egrad2rgrad = @egrad2rgrad;
0093     function rgrad = egrad2rgrad(X, egrad)
0094         X = prepare(X);
0095         rgrad.L = egrad.L*X.LtL;
0096         rgrad.R = egrad.R*X.RtR;
0097     end
0098     
0099     M.ehess2rhess = @ehess2rhess;
0100     function Hess = ehess2rhess(X, egrad, ehess, eta)
0101         X = prepare(X);
0102         
0103         % Riemannian gradient computation.
0104         rgrad = egrad2rgrad(X, egrad);
0105         
0106         % Directional derivative of the Riemannian gradient.
0107         Hess.L = ehess.L*X.LtL + 2*egrad.L*symm(eta.L'*X.L);
0108         Hess.R = ehess.R*X.RtR + 2*egrad.R*symm(eta.R'*X.R);
0109         
0110         % We need a correction term for the non-constant metric.
0111         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)));
0112         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)));
0113         
0114         % Projection onto the horizontal space.
0115         Hess = M.proj(X, Hess);
0116     end
0117     
0118     M.proj = @projection;
0119     % Projection of the vector eta in the ambient space onto the horizontal space.
0120     function etaproj = projection(X, eta)
0121         X = prepare(X);
0122         
0123         SS = (X.LtL)*(X.RtR);
0124         AS = (X.LtL)*(X.R'*eta.R) - (eta.L'*X.L)*(X.RtR);
0125         Omega = lyap(SS, SS,-AS);
0126         etaproj.L = eta.L + X.L*Omega';
0127         etaproj.R = eta.R - X.R*Omega;
0128     end
0129     
0130     M.tangent = M.proj;
0131     M.tangent2ambient = @(X, eta) eta;
0132     
0133     M.retr = @retraction;
0134     function Y = retraction(X, eta, t)
0135         if nargin < 3
0136             t = 1.0;
0137         end
0138         
0139         Y.L = X.L + t*eta.L;
0140         Y.R = X.R + t*eta.R;
0141         
0142         % Numerical conditioning step: A simpler version.
0143         % We need to ensure that L and R do not have very relative
0144         % skewed norms.
0145         
0146         scaling = norm(X.L, 'fro')/norm(X.R, 'fro');
0147         scaling = sqrt(scaling);
0148         Y.L = Y.L / scaling;
0149         Y.R = Y.R * scaling;
0150         
0151         % These are reused in the computation of the gradient and Hessian.
0152         Y = prepare(Y);
0153     end
0154     
0155     M.exp = @exponential;
0156     function Y = exponential(X, eta, t)
0157         if nargin < 3
0158             t = 1.0;
0159         end
0160         
0161         Y = retraction(X, eta, t);
0162         warning('manopt:fixedrankfactory_2factors:exp', ...
0163             ['Exponential for fixed rank ' ...
0164             'manifold not implemented yet. Used retraction instead.']);
0165     end
0166     
0167     M.hash = @(X) ['z' hashmd5([X.L(:) ; X.R(:)])];
0168     
0169     M.rand = @random;
0170     function X = random()
0171         % A random point on the total space.
0172         X.L = randn(m, k);
0173         X.R = randn(n, k);
0174         X = prepare(X);
0175     end
0176     
0177     M.randvec = @randomvec;
0178     function eta = randomvec(X)
0179         % A random vector in the horizontal space.
0180         eta.L = randn(m, k);
0181         eta.R = randn(n, k);
0182         eta = projection(X, eta);
0183         nrm = M.norm(X, eta);
0184         eta.L = eta.L / nrm;
0185         eta.R = eta.R / nrm;
0186     end
0187     
0188     M.lincomb = @lincomb;
0189     
0190     M.zerovec = @(X) struct('L', zeros(m, k),'R', zeros(n, k));
0191     
0192     M.transp = @(x1, x2, d) projection(x2, d);
0193     
0194     % vec and mat are not isometries, because of the unusual inner metric.
0195     M.vec = @(X, U) [U.L(:) ; U.R(:)];
0196     M.mat = @(X, u) struct('L', reshape(u(1:(m*k)), m, k), ...
0197         'R', reshape(u((m*k+1):end), n, k));
0198     M.vecmatareisometries = @() false;
0199     
0200 end
0201 
0202 % Linear combination of tangent vectors.
0203 function d = lincomb(x, a1, d1, a2, d2) %#ok<INUSL>
0204     
0205     if nargin == 3
0206         d.L = a1*d1.L;
0207         d.R = a1*d1.R;
0208     elseif nargin == 5
0209         d.L = a1*d1.L + a2*d2.L;
0210         d.R = a1*d1.R + a2*d2.R;
0211     else
0212         error('Bad use of fixedrankfactory_2factors.lincomb.');
0213     end
0214     
0215 end
0216 
0217 
0218 
0219 
0220

Generated on Sat 12-Nov-2016 14:11:22 by m2html © 2005