Home > manopt > manifolds > fixedrank > fixedrankfactory_3factors.m

fixedrankfactory_3factors

PURPOSE ^

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

SYNOPSIS ^

function M = fixedrankfactory_3factors(m, n, k)

DESCRIPTION ^

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

 function M = fixedrankfactory_3factors(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 three
 fields: L, S and R. The matrices L (mxk) and R (nxk) are orthonormal,
 while the matrix S (kxk) is a symmetric positive definite full rank
 matrix.

 Tangent vectors are represented as a structure with three fields: L, S
 and 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_2factors fixedrankfactory_3factors_preconditioned

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function M = fixedrankfactory_3factors(m, n, k)
0002 % Manifold of m-by-n matrices of rank k with polar quotient geometry.
0003 %
0004 % function M = fixedrankfactory_3factors(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 three
0019 % fields: L, S and R. The matrices L (mxk) and R (nxk) are orthonormal,
0020 % while the matrix S (kxk) is a symmetric positive definite full rank
0021 % matrix.
0022 %
0023 % Tangent vectors are represented as a structure with three fields: L, S
0024 % and R.
0025 %
0026 %
0027 % For first-order geometry, please cite the Manopt paper as well as the research paper:
0028 %     @InProceedings{meyer2011linear,
0029 %       Title        = {Linear regression under fixed-rank constraints: a {R}iemannian approach},
0030 %       Author       = {Meyer, G. and Bonnabel, S. and Sepulchre, R.},
0031 %       Booktitle    = {{28th International Conference on Machine Learning}},
0032 %       Year         = {2011},
0033 %       Organization = {{ICML}}
0034 %     }
0035 % For second-order geometry, please cite the Manopt paper as well as the research paper:
0036 %     @Article{mishra2014fixedrank,
0037 %       Title   = {Fixed-rank matrix factorizations and {Riemannian} low-rank optimization},
0038 %       Author  = {Mishra, B. and Meyer, G. and Bonnabel, S. and Sepulchre, R.},
0039 %       Journal = {Computational Statistics},
0040 %       Year    = {2014},
0041 %       Number  = {3-4},
0042 %       Pages   = {591--621},
0043 %       Volume  = {29},
0044 %       Doi     = {10.1007/s00180-013-0464-z}
0045 %     }
0046 %
0047 %
0048 % See also fixedrankembeddedfactory fixedrankfactory_2factors fixedrankfactory_3factors_preconditioned
0049 
0050 % This file is part of Manopt: www.manopt.org.
0051 % Original author: Bamdev Mishra, Dec. 30, 2012.
0052 % Contributors:
0053 % Change log:
0054     
0055     M.name = @() sprintf('LSR'' quotient manifold of %dx%d matrices of rank %d', m, n, k);
0056     
0057     M.dim = @() (m+n-k)*k;
0058     
0059     % Choice of the metric on the orthnormal space is motivated by the symmetry present in the
0060     % space. The metric on the positive definite space is its natural metric.
0061     M.inner = @(X, eta, zeta) eta.L(:).'*zeta.L(:) + eta.R(:).'*zeta.R(:) ...
0062         + trace( (X.S\eta.S) * (X.S\zeta.S) );
0063     
0064     M.norm = @(X, eta) sqrt(M.inner(X, eta, eta));
0065     
0066     M.dist = @(x, y) error('fixedrankfactory_3factors.dist not implemented yet.');
0067     
0068     M.typicaldist = @() 10*k;
0069     
0070     skew = @(X) .5*(X-X');
0071     symm = @(X) .5*(X+X');
0072     stiefel_proj = @(L, H) H - L*symm(L'*H);
0073     
0074     M.egrad2rgrad = @egrad2rgrad;
0075     function rgrad = egrad2rgrad(X, egrad)
0076         rgrad.L = stiefel_proj(X.L, egrad.L);
0077         rgrad.S = X.S*symm(egrad.S)*X.S;
0078         rgrad.R = stiefel_proj(X.R, egrad.R);
0079     end
0080     
0081     
0082     M.ehess2rhess = @ehess2rhess;
0083     function Hess = ehess2rhess(X, egrad, ehess, eta)
0084         
0085         % Riemannian gradient for the factor S.
0086         rgrad.S = X.S*symm(egrad.S)*X.S;
0087         
0088         % Directional derivatives of the Riemannian gradient.
0089         Hess.L = ehess.L - eta.L*symm(X.L'*egrad.L);
0090         Hess.L = stiefel_proj(X.L, Hess.L);
0091         
0092         Hess.R = ehess.R - eta.R*symm(X.R'*egrad.R);
0093         Hess.R = stiefel_proj(X.R, Hess.R);
0094         
0095         Hess.S = X.S*symm(ehess.S)*X.S +  2*symm(eta.S*symm(egrad.S)*X.S);
0096         
0097         % Correction factor for the non-constant metric on the factor S.
0098         Hess.S = Hess.S - symm(eta.S*(X.S\rgrad.S));
0099         
0100         % Projection onto the horizontal space.
0101         Hess = M.proj(X, Hess);
0102     end
0103     
0104     
0105     M.proj = @projection;
0106     function etaproj = projection(X, eta)
0107         % First, projection onto the tangent space of the total space.
0108         eta.L = stiefel_proj(X.L, eta.L);
0109         eta.R = stiefel_proj(X.R, eta.R);
0110         eta.S = symm(eta.S);
0111         
0112         % Then, projection onto the horizontal space.
0113         SS = X.S*X.S;
0114         AS = X.S*(skew(X.L'*eta.L) + skew(X.R'*eta.R) - 2*skew(X.S\eta.S))*X.S;
0115         omega = lyap(SS, -AS);
0116         
0117         etaproj.L = eta.L - X.L*omega;
0118         etaproj.S = eta.S - (X.S*omega - omega*X.S);
0119         etaproj.R = eta.R - X.R*omega;
0120     end
0121     
0122     M.tangent = M.proj;
0123     M.tangent2ambient = @(X, eta) eta;
0124     
0125     M.retr = @retraction;
0126     function Y = retraction(X, eta, t)
0127         if nargin < 3
0128             t = 1.0;
0129         end
0130         
0131         L = chol(X.S);
0132         Y.S = L'*expm(L'\(t*eta.S)/L)*L;
0133         Y.L = uf(X.L + t*eta.L);
0134         Y.R = uf(X.R + t*eta.R);
0135     end
0136     
0137     M.exp = @exponential;
0138     function Y = exponential(X, eta, t)
0139         if nargin < 3
0140             t = 1.0;
0141         end
0142         Y = retraction(X, eta, t);
0143         warning('manopt:fixedrankfactory_3factors:exp', ...
0144             ['Exponential for fixed rank ' ...
0145             'manifold not implemented yet. Lsed retraction instead.']);
0146     end
0147     
0148     M.hash = @(X) ['z' hashmd5([X.L(:) ; X.S(:) ; X.R(:)])];
0149     
0150     M.rand = @random;
0151     % Factors L and R are on Stiefel manifolds, hence we reuse
0152     % their random generators.
0153     stiefelm = stiefelfactory(m, k);
0154     stiefeln = stiefelfactory(n, k);
0155     function X = random()
0156         X.L = stiefelm.rand();
0157         X.R = stiefeln.rand();
0158         X.S = diag(1+rand(k, 1));
0159     end
0160     
0161     M.randvec = @randomvec;
0162     function eta = randomvec(X)
0163         % A random vector on the horizontal space.
0164         eta.L = randn(m, k);
0165         eta.R = randn(n, k);
0166         eta.S = randn(k, k);
0167         eta = projection(X, eta);
0168         nrm = M.norm(X, eta);
0169         eta.L = eta.L / nrm;
0170         eta.R = eta.R / nrm;
0171         eta.S = eta.S / nrm;
0172     end
0173     
0174     M.lincomb = @lincomb;
0175     
0176     M.zerovec = @(X) struct('L', zeros(m, k), 'S', zeros(k, k), ...
0177         'R', zeros(n, k));
0178     
0179     M.transp = @(x1, x2, d) projection(x2, d);
0180     
0181     % vec and mat are not isometries, because of the scaled inner metric.
0182     M.vec = @(X, U) [U.L(:) ; U.S(:); U.R(:)];
0183     M.mat = @(X, u) struct('L', reshape(u(1:(m*k)), m, k), ...
0184         'S', reshape(u((m*k+1): m*k + k*k), k, k), ...
0185         'R', reshape(u((m*k+ k*k + 1):end), n, k));
0186     M.vecmatareisometries = @() false;
0187     
0188 end
0189 
0190 % Linear combination of tangent vectors.
0191 function d = lincomb(x, a1, d1, a2, d2) %#ok<INLSL>
0192     
0193     if nargin == 3
0194         d.L = a1*d1.L;
0195         d.R = a1*d1.R;
0196         d.S = a1*d1.S;
0197     elseif nargin == 5
0198         d.L = a1*d1.L + a2*d2.L;
0199         d.R = a1*d1.R + a2*d2.R;
0200         d.S = a1*d1.S + a2*d2.S;
0201     else
0202         error('Bad use of fixedrankfactory_3factors.lincomb.');
0203     end
0204     
0205 end
0206 
0207 function A = uf(A)
0208     [L, unused, R] = svd(A, 0); %#ok
0209     A = L*R';
0210 end

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