Linear manifold structure for optimization over the ShapeFit search space function M = shapefitfactory(VJt) Input: VJt is a matrix of size dxn, such that VJt * ones(n, 1) = 0. Returns M, a structure describing the Euclidean space of d-by-n matrices equipped with the standard Frobenius distance and associated trace inner product, as a manifold for Manopt. Matrices on M, denoted by T, have size dxn and obey T*ones(n, 1) = 0 (centered columns) and <VJt, T> = 1, where <A, B> = Trace(A' * B). See this paper: http://arxiv.org/abs/1506.01437 ShapeFit: Exact location recovery from corrupted pairwise directions, 2015 Paul Hand, Choongbum Lee, Vladislav Voroninski See also: shapefit_smoothed
0001 function M = shapefitfactory(VJt) 0002 % Linear manifold structure for optimization over the ShapeFit search space 0003 % 0004 % function M = shapefitfactory(VJt) 0005 % 0006 % Input: VJt is a matrix of size dxn, such that VJt * ones(n, 1) = 0. 0007 % 0008 % Returns M, a structure describing the Euclidean space of d-by-n matrices 0009 % equipped with the standard Frobenius distance and associated trace inner 0010 % product, as a manifold for Manopt. Matrices on M, denoted by T, have size 0011 % dxn and obey T*ones(n, 1) = 0 (centered columns) and <VJt, T> = 1, where 0012 % <A, B> = Trace(A' * B). 0013 % 0014 % See this paper: http://arxiv.org/abs/1506.01437 0015 % ShapeFit: Exact location recovery from corrupted pairwise directions, 2015 0016 % Paul Hand, Choongbum Lee, Vladislav Voroninski 0017 % 0018 % See also: shapefit_smoothed 0019 0020 % This file is part of Manopt: www.manopt.org. 0021 % Original author: Nicolas Boumal, June 18, 2015. 0022 % Contributors: 0023 % Change log: 0024 % 0025 % Jan. 25, 2017 (NB): 0026 % M.tangent = M.proj now, instead of being identity. This is notably 0027 % necessary so that checkgradient will pick up on gradients that do 0028 % not lie in the appropriate tangent space. 0029 % 0030 % Jan. 4, 2021 (NB): 0031 % Changes for compatibility with Octave 6.1.0. 0032 0033 [d, n] = size(VJt); 0034 0035 M.name = @() sprintf('ShapeFit space of size %d x %d', d, n); 0036 0037 M.dim = @() d*n - d - 1; 0038 0039 M.inner = @(x, d1, d2) d1(:).'*d2(:); 0040 0041 M.norm = @(x, d) norm(d, 'fro'); 0042 0043 M.dist = @(x, y) norm(x-y, 'fro'); 0044 0045 M.typicaldist = @() sqrt(d*n); 0046 0047 VJt_normed = VJt / norm(VJt, 'fro'); 0048 M.proj = @(T, U) projection(U, VJt_normed); 0049 function PU = projection(U, VJt_normed) 0050 % Center the columns 0051 PU = bsxfun(@minus, U, mean(U, 2)); 0052 % Remove component along VJt 0053 % Note: these two actions can be executed separately, without 0054 % interference, owing to VJt having centered columns itself. 0055 PU = PU - (VJt_normed(:)'*U(:))*VJt_normed; 0056 end 0057 0058 M.egrad2rgrad = M.proj; 0059 0060 M.ehess2rhess = @(x, eg, eh, d) projection(eh, VJt_normed); 0061 0062 M.tangent = M.proj; 0063 0064 M.exp = @exp; 0065 function y = exp(x, d, t) 0066 if nargin == 3 0067 y = x + t*d; 0068 else 0069 y = x + d; 0070 end 0071 end 0072 0073 M.retr = M.exp; 0074 0075 M.log = @(x, y) y-x; 0076 0077 M.hash = @(x) ['z' hashmd5(x(:))]; 0078 0079 M.randvec = @(x) randvec(VJt, VJt_normed); 0080 function u = randvec(VJt, VJt_normed) 0081 u = projection(randn(size(VJt)), VJt_normed); 0082 u = u / norm(u, 'fro'); 0083 end 0084 0085 % We exploit the fact that VJt_normed belongs to the manifold 0086 M.rand = @() VJt_normed + randn(1) * randvec(VJt, VJt_normed); 0087 0088 M.lincomb = @matrixlincomb; 0089 0090 M.zerovec = @(x) zeros(d, n); 0091 0092 M.transp = @(x1, x2, d) d; 0093 0094 M.pairmean = @(x1, x2) .5*(x1+x2); 0095 0096 M.vec = @(x, u_mat) u_mat(:); 0097 M.mat = @(x, u_vec) reshape(u_vec, [d, n]); 0098 M.vecmatareisometries = @() true; 0099 0100 end