{"id":940,"date":"2014-04-01T00:01:17","date_gmt":"2014-04-01T05:01:17","guid":{"rendered":"https:\/\/blogs.mathworks.com\/cleve\/?p=940"},"modified":"2016-12-05T13:58:41","modified_gmt":"2016-12-05T18:58:41","slug":"reverse-singular-value-decomposition","status":"publish","type":"post","link":"https:\/\/blogs.mathworks.com\/cleve\/2014\/04\/01\/reverse-singular-value-decomposition\/","title":{"rendered":"Reverse Singular Value Decomposition"},"content":{"rendered":"<div class=\"content\"><!--introduction-->Employing a factorization based on the least significant singular values provides a matrix approximation with many surprisingly useful properties. This <i>Reverse Singular Value Decomposition<\/i>, <i>RSVD<\/i>, is also referred to as <i>Subordinate Component Analysis<\/i>, <i>SCA<\/i>, to distinguish it from Principal Component Analysis.\r\n\r\n<!--\/introduction-->\r\n<h3>Contents<\/h3>\r\n<div>\r\n<ul>\r\n \t<li><a href=\"#ed766650-5ff8-4177-b225-ee09891fbd9f\">RSVD<\/a><\/li>\r\n \t<li><a href=\"#9b743269-ee58-4ecd-a26e-5ae0a283b6a0\">Roundoff Error<\/a><\/li>\r\n \t<li><a href=\"#a28e6044-4b13-4c9a-8ebe-5ce3e1d6be15\">Text Processing<\/a><\/li>\r\n \t<li><a href=\"#08236e6d-dc65-432d-b537-ce7acfa238f5\">Image Processing<\/a><\/li>\r\n<\/ul>\r\n<\/div>\r\n<h4>RSVD<a name=\"ed766650-5ff8-4177-b225-ee09891fbd9f\"><\/a><\/h4>\r\nThe Singular Value Decomposition of a matrix $A$ is computed by\r\n<pre>    [U,S,V] = svd(A);<\/pre>\r\nThis generates two orthogonal matrices $U$ and $V$ and a diagonal matrix $S$ with diagonal elements $\\sigma_k$ that provide the expansion\r\n\r\n$$ A = \\sigma_1 E_1 + \\sigma_2 E_2 + ... + \\sigma_n E_n $$\r\n\r\nwhere $E_k$ is the rank one outer product\r\n\r\n$$ E_k = U(:,k) V(:,k)' $$\r\n\r\nTraditionally, the singular values are arranged in descending order. In contrast, the Reverse Singular Value Decomposition Approximation of rank $r$ is obtained by arranging the singular values in ascending order,\r\n\r\n$$ 0 \\le \\sigma_1 \\le \\sigma_2 \\le ... $$\r\n\r\nand then using the first $r$ terms from the expansion.\r\n\r\nHere is a function that computes the RSVD for square or rectangular matrices with at least as many rows as columns.\r\n<pre class=\"codeinput\">   type <span class=\"string\">rsvd<\/span>\r\n<\/pre>\r\n<pre class=\"codeoutput\">function X = rsvd(A,r)\r\n% RSVD  Approximation by the Reverse Singular Value Decomposition\r\n%   rsvd(A,r) approximates A by a matrix of rank r obtained from\r\n%   the r least significant singular values in ascending order.\r\n\r\n   [m,n] = size(A);\r\n   [U,S,V] = svd(A,'econ');\r\n   k = n:-1:n-r+1;\r\n   X = U(:,k)*S(k,k)*V(:,k)';\r\n<\/pre>\r\n<h4>Roundoff Error<a name=\"9b743269-ee58-4ecd-a26e-5ae0a283b6a0\"><\/a><\/h4>\r\nIn certain situations, the RSVD can reduce or even eliminate roundoff error. For example, according to its <tt>help<\/tt> entry the <tt>elmat<\/tt> function <tt>hilb<\/tt> attempts to compute\r\n<pre>    hilb(N) is the N by N matrix with elements 1\/(i+j-1)<\/pre>\r\nBut the function can only succeed to within roundoff error. Its results are binary floating point numbers approximating the reciprocals of integers described in the <tt>help<\/tt> entry. Here is the printed output with <tt>n = 5<\/tt> and the default <tt>format short<\/tt>.\r\n<pre class=\"codeinput\">   format <span class=\"string\">short<\/span>\r\n   H = hilb(5)\r\n<\/pre>\r\n<pre class=\"codeoutput\">H =\r\n\r\n    1.0000    0.5000    0.3333    0.2500    0.2000\r\n    0.5000    0.3333    0.2500    0.2000    0.1667\r\n    0.3333    0.2500    0.2000    0.1667    0.1429\r\n    0.2500    0.2000    0.1667    0.1429    0.1250\r\n    0.2000    0.1667    0.1429    0.1250    0.1111\r\n\r\n<\/pre>\r\nWe are seeing the effects of the output conversion as well as the underlying binary approximation. Perhaps surprisingly, the reverse singular value decomposition can eliminate both sources of error and produce the exact rational entries of the theoretical Hilbert matrix.\r\n<pre class=\"codeinput\">   format <span class=\"string\">rat<\/span>\r\n   H = rsvd(H,5)\r\n<\/pre>\r\n<pre class=\"codeoutput\">H =\r\n\r\n       1              1\/2            1\/3            1\/4            1\/5     \r\n       1\/2            1\/3            1\/4            1\/5            1\/6     \r\n       1\/3            1\/4            1\/5            1\/6            1\/7     \r\n       1\/4            1\/5            1\/6            1\/7            1\/8     \r\n       1\/5            1\/6            1\/7            1\/8            1\/9     \r\n\r\n<\/pre>\r\n<h4>Text Processing<a name=\"a28e6044-4b13-4c9a-8ebe-5ce3e1d6be15\"><\/a><\/h4>\r\nThe RSVD is capable of uncovering spelling and typographical errors in text. My web site for <i>Numerical Computing with MATLAB<\/i> has a file with the text of Lincoln's Gettysburg Address.\r\n\r\n<a href=\"https:\/\/www.mathworks.com\/content\/dam\/mathworks\/mathworks-dot-com\/moler\/ncm\/gettysburg.txt\">gettysburg.txt<\/a>\r\n\r\nI've used this file for years whenever I want to experiment with text processing. You can download the file and then use the MATLAB data import wizard with column delimeters set to spaces, commas and periods to produce a cell array, <tt>gettysburg<\/tt>, of individual words. There are 278 words. The longest word has 11 characters. So\r\n<pre>    A = cell2mat(gettysburg)<\/pre>\r\nconverts the cell array of strings to a 278-by-11 matrix of doubles.\r\n\r\nIt turns out that an RSVD approximation of rank nine finds three spelling errors in the original text.\r\n<pre>    B = rsvd(A,9);\r\n    k = find(sum(A,2)~=sum(B,2))\r\n    disp([char(A(k,:))  char(B(k,:))])<\/pre>\r\n<pre>    weather      whether\r\n    consicrated  consecrated\r\n    govenment    government<\/pre>\r\nIn all the years that I have been using this data, I have never noticed these errors.\r\n<h4>Image Processing<a name=\"08236e6d-dc65-432d-b537-ce7acfa238f5\"><\/a><\/h4>\r\nWe have also found that the RSVD is capable of softening the appearance of aggression in photographs. A matrix is obtained from the JPEG image by stacking the RGB components vertically. Roughly 90% of the small singular values then produce a pleasant result.\r\n<pre>    RGB = imread('creature.jpg');\r\n    A = [RGB(:,:,1); RGB(:,:,2); RGB(:,:,3)];\r\n    [m,n] = size(A);\r\n    B = rsvd(A,ceil(0.90*n));\r\n    m = m\/3;\r\n    C = cat(3,B(1:m,:),B(m+1:2*m,:),B(2*m+1:3*m,:))\r\n    image(C)<\/pre>\r\n<img decoding=\"async\" src=\"https:\/\/blogs.mathworks.com\/images\/cleve\/creature.jpg\" alt=\"\" hspace=\"5\" vspace=\"5\" \/>\r\n\r\n<img decoding=\"async\" src=\"https:\/\/blogs.mathworks.com\/images\/cleve\/celebrity.jpg\" alt=\"\" hspace=\"5\" vspace=\"5\" \/>\r\n\r\n<script>\/\/ <![CDATA[\r\nfunction grabCode_ec92434ac6f7485e9c4447d400ba9c6e() {\r\n        \/\/ Remember the title so we can use it in the new page\r\n        title = document.title;\r\n\r\n        \/\/ Break up these strings so that their presence\r\n        \/\/ in the Javascript doesn't mess up the search for\r\n        \/\/ the MATLAB code.\r\n        t1='ec92434ac6f7485e9c4447d400ba9c6e ' + '##### ' + 'SOURCE BEGIN' + ' #####';\r\n        t2='##### ' + 'SOURCE END' + ' #####' + ' ec92434ac6f7485e9c4447d400ba9c6e';\r\n    \r\n        b=document.getElementsByTagName('body')[0];\r\n        i1=b.innerHTML.indexOf(t1)+t1.length;\r\n        i2=b.innerHTML.indexOf(t2);\r\n \r\n        code_string = b.innerHTML.substring(i1, i2);\r\n        code_string = code_string.replace(\/REPLACE_WITH_DASH_DASH\/g,'--');\r\n\r\n        \/\/ Use \/x3C\/g instead of the less-than character to avoid errors \r\n        \/\/ in the XML parser.\r\n        \/\/ Use '\\x26#60;' instead of '<' so that the XML parser\r\n        \/\/ doesn't go ahead and substitute the less-than character. \r\n        code_string = code_string.replace(\/\\x3C\/g, '\\x26#60;');\r\n\r\n        copyright = 'Copyright 2014 The MathWorks, Inc.';\r\n\r\n        w = window.open();\r\n        d = w.document;\r\n        d.write('\r\n\r\n<pre>\\n');\r\n        d.write(code_string);\r\n\r\n        \/\/ Add copyright line at the bottom if specified.\r\n        if (copyright.length > 0) {\r\n            d.writeln('');\r\n            d.writeln('%%');\r\n            if (copyright.length > 0) {\r\n                d.writeln('% _' + copyright + '_');\r\n            }\r\n        }\r\n\r\n        d.write('<\/pre>\r\n\r\n\r\n\\n');\r\n\r\n        d.title = title + ' (MATLAB code)';\r\n        d.close();\r\n    }\r\n\/\/ ]]><\/script>\r\n<p style=\"text-align: right; font-size: xx-small; font-weight: lighter; font-style: italic; color: gray;\">\r\n<a><span style=\"font-size: x-small; font-style: italic;\">Get\r\nthe MATLAB code<noscript>(requires JavaScript)<\/noscript><\/span><\/a>\r\n\r\nPublished with MATLAB\u00ae R2014a<\/p>\r\n\r\n<\/div>\r\n<!--\r\nec92434ac6f7485e9c4447d400ba9c6e ##### SOURCE BEGIN #####\r\n%% Reverse Singular Value Decomposition\r\n% Employing a factorization based on the least significant singular values\r\n% provides a matrix approximation with many surprisingly useful properties.\r\n% This _Reverse Singular Value Decomposition_, _RSVD_,\r\n% is also referred to as _Subordinate Component Analysis_, _SCA_, to\r\n% distinguish it from Principal Component Analysis.\r\n\r\n%% RSVD\r\n% The Singular Value Decomposition of a matrix $A$ is computed by\r\n%\r\n%      [U,S,V] = svd(A);\r\n%\r\n% This generates two orthogonal matrices $U$ and $V$ and a diagonal matrix $S$\r\n% with diagonal elements $\\sigma_k$ that provide the expansion\r\n%\r\n% $$ A = \\sigma_1 E_1 +  \\sigma_2 E_2 + ... +  \\sigma_n E_n $$\r\n%\r\n% where $E_k$ is the rank one outer product\r\n%\r\n% $$ E_k = U(:,k) V(:,k)' $$\r\n\r\n%%\r\n% Traditionally, the singular values are arranged in descending order.\r\n% In contrast, the Reverse Singular Value Decomposition Approximation of\r\n% rank $r$ is obtained by arranging the singular values in ascending order,\r\n%\r\n% $$ 0 \\le \\sigma_1 \\le \\sigma_2 \\le ... $$\r\n%\r\n% and then using the first $r$ terms from the expansion.\r\n\r\n%%\r\n% Here is a function that computes the RSVD for square or rectangular matrices\r\n% with at least as many rows as columns.\r\n%\r\ntype rsvd\r\n\r\n%% Roundoff Error\r\n% In certain situations, the RSVD can reduce or even eliminate roundoff error.\r\n% For example, according to its |help| entry the |elmat| function |hilb|\r\n% attempts to compute\r\n%\r\n%      hilb(N) is the N by N matrix with elements 1\/(i+j-1)\r\n%\r\n% But the function can only succeed to within roundoff error.  Its results\r\n% are binary floating point numbers approximating the reciprocals of integers\r\n% described in the |help| entry.  Here is the printed output with |n = 5|\r\n% and the default |format short|.\r\n\r\nformat short\r\nH = hilb(5)\r\n\r\n%%\r\n% We are seeing the effects of the output conversion as well as the underlying\r\n% binary approximation.  Perhaps surprisingly, the reverse singular value\r\n% decomposition can eliminate both sources of error and produce the\r\n% exact rational entries of the theoretical Hilbert matrix.\r\n\r\nformat rat\r\nH = rsvd(H,5)\r\n\r\n%% Text Processing\r\n% The RSVD is capable of uncovering spelling and typographical errors in text.\r\n% My web site for _Numerical Computing with MATLAB_ has a file with the text\r\n% of Lincoln's Gettysburg Address.\r\n%\r\n% <https:\/\/www.mathworks.com\/moler.htmlncm\/gettysburg.txt gettysburg.txt>\r\n%\r\n% I've used this file for years whenever I want to experiment with text\r\n% processing.  You can download the file and then use the MATLAB data import\r\n% wizard with column delimeters set to spaces, commas and periods to produce\r\n% a cell array, |gettysburg|, of individual words.  There are 278 words.\r\n% The longest word has 11 characters.  So\r\n%\r\n%      A = cell2mat(gettysburg)\r\n%\r\n% converts the cell array of strings to a 278-by-11 matrix of doubles.\r\n\r\n%%\r\n% It turns out that an RSVD approximation of rank nine finds three spelling\r\n% errors in the original text.\r\n%\r\n%      B = rsvd(A,9);\r\n%      k = find(sum(A,2)~=sum(B,2))\r\n%      disp([char(A(k,:))  char(B(k,:))])\r\n%\r\n%      weather      whether\r\n%      consicrated  consecrated\r\n%      govenment    government\r\n%\r\n% In all the years that I have been using this data, I have never noticed\r\n% these errors.\r\n\r\n%% Image Processing\r\n% We have also found that the RSVD is capable of softening the appearance\r\n% of aggression in photographs.  A matrix is obtained from the JPEG image\r\n% by stacking the RGB components vertically.  Roughly 90% of the small\r\n% singular values then produce a pleasant result.\r\n%\r\n%      RGB = imread('creature.jpg');\r\n%      A = [RGB(:,:,1); RGB(:,:,2); RGB(:,:,3)];\r\n%      [m,n] = size(A);\r\n%      B = rsvd(A,ceil(0.90*n));\r\n%      m = m\/3;\r\n%      C = cat(3,B(1:m,:),B(m+1:2*m,:),B(2*m+1:3*m,:))\r\n%      image(C)\r\n%\r\n% <<creature.jpg>>\r\n%\r\n%\r\n% <<celebrity.jpg>>\r\n\r\n##### SOURCE END ##### ec92434ac6f7485e9c4447d400ba9c6e\r\n-->","protected":false},"excerpt":{"rendered":"<div class=\"overview-image\"><img decoding=\"async\"  class=\"img-responsive\" src=\"https:\/\/blogs.mathworks.com\/images\/cleve\/creature.jpg\" onError=\"this.style.display ='none';\" \/><\/div><!--introduction-->Employing a factorization based on the least significant singular values provides a matrix approximation with many surprisingly useful properties. This <i>Reverse Singular Value Decomposition<\/i>, <i>RSVD<\/i>, is also referred to as <i>Subordinate Component Analysis<\/i>, <i>SCA<\/i>, to distinguish it from Principal Component Analysis.\r\n\r\n<!--\/introduction-->... <a class=\"read-more\" href=\"https:\/\/blogs.mathworks.com\/cleve\/2014\/04\/01\/reverse-singular-value-decomposition\/\">read more >><\/a><\/p>","protected":false},"author":78,"featured_media":0,"comment_status":"open","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":[],"categories":[13,5,6],"tags":[],"_links":{"self":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/940"}],"collection":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/users\/78"}],"replies":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/comments?post=940"}],"version-history":[{"count":22,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/940\/revisions"}],"predecessor-version":[{"id":2183,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/940\/revisions\/2183"}],"wp:attachment":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/media?parent=940"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/categories?post=940"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/tags?post=940"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}