{"id":2027,"date":"2016-10-17T12:00:08","date_gmt":"2016-10-17T17:00:08","guid":{"rendered":"https:\/\/blogs.mathworks.com\/cleve\/?p=2027"},"modified":"2016-10-04T15:46:57","modified_gmt":"2016-10-04T20:46:57","slug":"compare-gram-schmidt-and-householder-orthogonalization-algorithms-2","status":"publish","type":"post","link":"https:\/\/blogs.mathworks.com\/cleve\/2016\/10\/17\/compare-gram-schmidt-and-householder-orthogonalization-algorithms-2\/","title":{"rendered":"Compare Gram-Schmidt and Householder Orthogonalization Algorithms"},"content":{"rendered":"\r\n<div class=\"content\"><!--introduction--><p>This is a follow-up to <a href=\"https:\/\/blogs.mathworks.com\/cleve\/2016\/10\/03\/householder-reflections-and-the-qr-decomposition\/\">my previous post.<\/a> Classical Gram-Schmidt and Modified Gram-Schmidt are two algorithms for orthogonalizing a set of vectors.  Householder elementary reflectors can be used for the same task.  The three algorithms have very different roundoff error properties.<\/p><!--\/introduction--><h3>Contents<\/h3><div><ul><li><a href=\"#a4fae5cd-143c-4c56-8fa1-99f8e1e35948\">Pete Stewart<\/a><\/li><li><a href=\"#8fd9c75b-6725-45e0-88a6-0491930faf0e\">Classic Gram-Schmidt<\/a><\/li><li><a href=\"#fedd5fd2-c6fa-46ef-925d-2a9d81d10032\">Modified Gram-Schmidt<\/a><\/li><li><a href=\"#3dffecd9-0a4f-4d62-a596-917f9cb8e5fa\">Householder Reflections<\/a><\/li><li><a href=\"#7521cdb1-27e2-49e0-a448-a8fa59fe41a4\">Comparison<\/a><\/li><li><a href=\"#8443e232-bbaf-404c-aee5-31913b28ef9f\">Well conditioned.<\/a><\/li><li><a href=\"#0e6f3a20-0e09-4dcb-8ef1-469ff049a6b6\">Poorly conditioned.<\/a><\/li><li><a href=\"#bfc9eb46-0f52-4e86-9f17-e6addfbbd866\">Singular.<\/a><\/li><li><a href=\"#8bb95823-f903-4ef5-a0cf-827be493ce3e\">Conclusion<\/a><\/li><li><a href=\"#a957cbb2-575f-4b36-ac76-bcab5c8e2a96\">Reference<\/a><\/li><\/ul><\/div><h4>Pete Stewart<a name=\"a4fae5cd-143c-4c56-8fa1-99f8e1e35948\"><\/a><\/h4><p>As I did in my previous post, I am using Pete Stewart's book <i>Matrix Algorithms, Volume I: Basic Decompositions<\/i>. His pseudocode is MATLAB ready.<\/p><h4>Classic Gram-Schmidt<a name=\"8fd9c75b-6725-45e0-88a6-0491930faf0e\"><\/a><\/h4><p>The classic Gram-Schmidt algorithm is the first thing you might think of for producing an orthogonal set of vectors.  For each vector in your data set, remove its projection onto the data set, normalize what is left, and include it in the orthogonal set.  Here is the code. <tt>X<\/tt> is the original set of vectors, <tt>Q<\/tt> is the resulting set of orthogonal vectors, and <tt>R<\/tt> is the set of coefficients, organized into an upper triangular matrix.<\/p><pre class=\"codeinput\">   type <span class=\"string\">gs<\/span>\r\n<\/pre><pre class=\"codeoutput\">\r\nfunction [Q,R] = gs(X)\r\n    % Classical Gram-Schmidt.  [Q,R] = gs(X);\r\n    % G. W. Stewart, \"Matrix Algorithms, Volume 1\", SIAM, 1998.\r\n    [n,p] = size(X);\r\n    Q = zeros(n,p);\r\n    R = zeros(p,p);\r\n    for k = 1:p\r\n        Q(:,k) = X(:,k);\r\n        if k ~= 1\r\n            R(1:k-1,k) = Q(:,k-1)'*Q(:,k);\r\n            Q(:,k) = Q(:,k) - Q(:,1:k-1)*R(1:k-1,k);\r\n        end\r\n        R(k,k) = norm(Q(:,k));\r\n        Q(:,k) = Q(:,k)\/R(k,k);\r\n     end\r\nend\r\n<\/pre><h4>Modified Gram-Schmidt<a name=\"fedd5fd2-c6fa-46ef-925d-2a9d81d10032\"><\/a><\/h4><p>This is a rather different algorithm, not just a simple modification of classical Gram-Schmidt.  The idea is to orthogonalize against the emerging set of vectors instead of against the original set.  There are two variants, a column-oriented one and a row-oriented one.  They produce the same results, in different order.  Here is the column version,<\/p><pre class=\"codeinput\">   type <span class=\"string\">mgs<\/span>\r\n<\/pre><pre class=\"codeoutput\">\r\nfunction [Q,R] =  mgs(X)\r\n    % Modified Gram-Schmidt.  [Q,R] = mgs(X);\r\n    % G. W. Stewart, \"Matrix Algorithms, Volume 1\", SIAM, 1998.\r\n    [n,p] = size(X);\r\n    Q = zeros(n,p);\r\n    R = zeros(p,p);\r\n    for k = 1:p\r\n        Q(:,k) = X(:,k);\r\n        for i = 1:k-1\r\n            R(i,k) = Q(:,i)'*Q(:,k);\r\n            Q(:,k) = Q(:,k) - R(i,k)*Q(:,i);\r\n        end\r\n        R(k,k) = norm(Q(:,k))';\r\n        Q(:,k) = Q(:,k)\/R(k,k);\r\n    end\r\nend\r\n<\/pre><h4>Householder Reflections<a name=\"3dffecd9-0a4f-4d62-a596-917f9cb8e5fa\"><\/a><\/h4><p>Compute <tt>R<\/tt> by applying Householder reflections to <tt>X<\/tt> a column at a time.  Do not actually compute <tt>Q<\/tt>, just save the vectors that generate the reflections.  See the description and codes from <a href=\"https:\/\/blogs.mathworks.com\/cleve\/2016\/10\/03\/householder-reflections-and-the-qr-decomposition\/\">my previous post<\/a>.<\/p><pre class=\"codeinput\">   type <span class=\"string\">house_qr<\/span>\r\n<\/pre><pre class=\"codeoutput\">\r\nfunction [U,R] = house_qr(A)\r\n    % Householder reflections for QR decomposition.\r\n    % [U,R] = house_qr(A) returns\r\n    % U, the reflector generators for use by house_apply.  \r\n    % R, the upper triangular factor.\r\n    H = @(u,x) x - u*(u'*x);\r\n    [m,n] = size(A);\r\n    U = zeros(m,n);\r\n    R = A;\r\n    for j = 1:min(m,n)\r\n        u = house_gen(R(j:m,j));\r\n        U(j:m,j) = u;\r\n        R(j:m,j:n) = H(u,R(j:m,j:n));\r\n        R(j+1:m,j) = 0;\r\n    end\r\nend\r\n<\/pre><h4>Comparison<a name=\"7521cdb1-27e2-49e0-a448-a8fa59fe41a4\"><\/a><\/h4><p>For various matrices <tt>X<\/tt>, let's check how the three algorithms perform on two tasks, accuracy and orthogonality.  How close is <tt>Q*R<\/tt> to <tt>X<\/tt>? And, how close is <tt>Q'*Q<\/tt> to <tt>I<\/tt>.<\/p><pre class=\"codeinput\">   type <span class=\"string\">compare.m<\/span>\r\n<\/pre><pre class=\"codeoutput\">\r\nfunction compare(X);\r\n% compare(X). Compare three QR decompositions.\r\n\r\nI = eye(size(X));\r\nqrerr = zeros(1,3);\r\northerr = zeros(1,3);\r\n\r\n%% Classic Gram Schmidt\r\n[Q,R] = gs(X);\r\nqrerr(1) = norm(Q*R-X,inf)\/norm(X,inf);\r\northerr(1) = norm(Q'*Q-I,inf);\r\n\r\n%% Modified Gram Schmidt\r\n[Q,R] = mgs(X);\r\nqrerr(2) = norm(Q*R-X,inf)\/norm(X,inf);\r\northerr(2) = norm(Q'*Q-I,inf);\r\n\r\n%% Householder QR Decomposition\r\n[U,R] = house_qr(X);\r\nQR = house_apply(U,R);\r\nQQ = house_apply_transpose(U,house_apply(U,I));\r\nqrerr(3) = norm(QR-X,inf)\/norm(X,inf);\r\northerr(3) = norm(QQ-I,inf);\r\n\r\n%% Report results \r\nfprintf('\\n                 Classic   Modified  Householder\\n')\r\nfprintf('QR error      %10.2e %10.2e %10.2e\\n',qrerr)\r\nfprintf('Orthogonality %10.2e %10.2e %10.2e\\n',ortherr)\r\n<\/pre><h4>Well conditioned.<a name=\"8443e232-bbaf-404c-aee5-31913b28ef9f\"><\/a><\/h4><p>First, try a well conditioned matrix, a magic square of odd order.<\/p><pre class=\"codeinput\">   n = 7;\r\n   X = magic(n)\r\n<\/pre><pre class=\"codeoutput\">X =\r\n    30    39    48     1    10    19    28\r\n    38    47     7     9    18    27    29\r\n    46     6     8    17    26    35    37\r\n     5    14    16    25    34    36    45\r\n    13    15    24    33    42    44     4\r\n    21    23    32    41    43     3    12\r\n    22    31    40    49     2    11    20\r\n<\/pre><p>Check its condition.<\/p><pre class=\"codeinput\">   kappa = condest(X)\r\n<\/pre><pre class=\"codeoutput\">kappa =\r\n    8.5681\r\n<\/pre><p>Do the comparison.<\/p><pre class=\"codeinput\">   compare(X);\r\n<\/pre><pre class=\"codeoutput\">\r\n                 Classic   Modified  Householder\r\nQR error        1.73e-16   6.09e-17   5.68e-16\r\nOrthogonality   3.20e+00   1.53e-15   1.96e-15\r\n<\/pre><p>All three algorithms do well with accuracy, but classic Gram-Schmidt fails with orthogonality.<\/p><h4>Poorly conditioned.<a name=\"0e6f3a20-0e09-4dcb-8ef1-469ff049a6b6\"><\/a><\/h4><p>Next, try a Hilbert matrix that is poorly conditioned, but not exactly singular.<\/p><pre class=\"codeinput\">   n = 7;\r\n   X = hilb(n)\r\n<\/pre><pre class=\"codeoutput\">X =\r\n    1.0000    0.5000    0.3333    0.2500    0.2000    0.1667    0.1429\r\n    0.5000    0.3333    0.2500    0.2000    0.1667    0.1429    0.1250\r\n    0.3333    0.2500    0.2000    0.1667    0.1429    0.1250    0.1111\r\n    0.2500    0.2000    0.1667    0.1429    0.1250    0.1111    0.1000\r\n    0.2000    0.1667    0.1429    0.1250    0.1111    0.1000    0.0909\r\n    0.1667    0.1429    0.1250    0.1111    0.1000    0.0909    0.0833\r\n    0.1429    0.1250    0.1111    0.1000    0.0909    0.0833    0.0769\r\n<\/pre><p>Check its condition.<\/p><pre class=\"codeinput\">   kappa = condest(X)\r\n<\/pre><pre class=\"codeoutput\">kappa =\r\n   9.8519e+08\r\n<\/pre><p>Do the comparison.<\/p><pre class=\"codeinput\">   compare(X)\r\n<\/pre><pre class=\"codeoutput\">\r\n                 Classic   Modified  Householder\r\nQR error        5.35e-17   5.35e-17   8.03e-16\r\nOrthogonality   5.21e+00   1.22e-08   1.67e-15\r\n<\/pre><p>All three algorithms do well with accuracy. Classic Gram-Schmidt fails completely with orthogonality. The orthogonality of MGS depends upon <tt>kappa<\/tt>. Householder does well with orthogonality.<\/p><h4>Singular.<a name=\"bfc9eb46-0f52-4e86-9f17-e6addfbbd866\"><\/a><\/h4><p>Finally, an exactly singular matrix, a magic square of even order.<\/p><pre class=\"codeinput\">   n = 8;\r\n   X = magic(n)\r\n<\/pre><pre class=\"codeoutput\">X =\r\n    64     2     3    61    60     6     7    57\r\n     9    55    54    12    13    51    50    16\r\n    17    47    46    20    21    43    42    24\r\n    40    26    27    37    36    30    31    33\r\n    32    34    35    29    28    38    39    25\r\n    41    23    22    44    45    19    18    48\r\n    49    15    14    52    53    11    10    56\r\n     8    58    59     5     4    62    63     1\r\n<\/pre><p>Check its rank.<\/p><pre class=\"codeinput\">   rankX = rank(X)\r\n<\/pre><pre class=\"codeoutput\">rankX =\r\n     3\r\n<\/pre><p>Do the comparison.<\/p><pre class=\"codeinput\">   compare(X)\r\n<\/pre><pre class=\"codeoutput\">\r\n                 Classic   Modified  Householder\r\nQR error        1.43e-16   8.54e-17   4.85e-16\r\nOrthogonality   5.41e+00   2.16e+00   1.30e-15\r\n<\/pre><p>Again, all three algorithms do well with accuracy. Both Gram-Schmidts fail completely with orthogonality. Householder still does well with orthogonality.<\/p><h4>Conclusion<a name=\"8bb95823-f903-4ef5-a0cf-827be493ce3e\"><\/a><\/h4><p>All three of these algorithms provide <tt>Q<\/tt> and <tt>R<\/tt> that do a good job of reproducing the data <tt>X<\/tt>, that is<\/p><div><ul><li><tt>Q<\/tt> * <tt>R<\/tt> is always close to <tt>X<\/tt> for all three algorithms.<\/li><\/ul><\/div><p>On the other hand, their behavior is very different when it comes to producing orthogonality.<\/p><div><ul><li>Classic Gram-Schmidt.  <tt>Q'*Q<\/tt> is almost never close to <tt>I<\/tt>.<\/li><li>Modified Gram-Schmidt.  <tt>Q'*Q<\/tt> depends upon condition of <tt>X<\/tt>   and fails completely when <tt>X<\/tt> is singular.<\/li><li>Householder triangularization.  <tt>Q'*Q<\/tt> is always close to <tt>I<\/tt><\/li><\/ul><\/div><h4>Reference<a name=\"a957cbb2-575f-4b36-ac76-bcab5c8e2a96\"><\/a><\/h4><p>G. W. Stewart, <i>Matrix Algorithms: Volume 1: Basic Decompositions<\/i>, SIAM, xix+458, 1998. <a href=\"http:\/\/epubs.siam.org\/doi\/book\/10.1137\/1.9781611971408\">&lt;http:\/\/epubs.siam.org\/doi\/book\/10.1137\/1.9781611971408<\/a>&gt;<\/p><script language=\"JavaScript\"> <!-- \r\n    function grabCode_0e62a7e371c643b1a14335ce8266b39b() {\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='0e62a7e371c643b1a14335ce8266b39b ' + '##### ' + 'SOURCE BEGIN' + ' #####';\r\n        t2='##### ' + 'SOURCE END' + ' #####' + ' 0e62a7e371c643b1a14335ce8266b39b';\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 2016 The MathWorks, Inc.';\r\n\r\n        w = window.open();\r\n        d = w.document;\r\n        d.write('<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>\\n');\r\n\r\n        d.title = title + ' (MATLAB code)';\r\n        d.close();\r\n    }   \r\n     --> <\/script><p style=\"text-align: right; font-size: xx-small; font-weight:lighter;   font-style: italic; color: gray\"><br><a href=\"javascript:grabCode_0e62a7e371c643b1a14335ce8266b39b()\"><span style=\"font-size: x-small;        font-style: italic;\">Get \r\n      the MATLAB code <noscript>(requires JavaScript)<\/noscript><\/span><\/a><br><br>\r\n      Published with MATLAB&reg; R2016a<br><\/p><\/div><!--\r\n0e62a7e371c643b1a14335ce8266b39b ##### SOURCE BEGIN #####\r\n%% Compare Gram-Schmidt and Householder Orthogonalization Algorithms\r\n% This is a follow-up to\r\n% <https:\/\/blogs.mathworks.com\/cleve\/2016\/10\/03\/householder-reflections-and-the-qr-decomposition\/\r\n% my previous post.>\r\n% Classical Gram-Schmidt and Modified Gram-Schmidt are two algorithms for\r\n% orthogonalizing a set of vectors.  Householder elementary reflectors can\r\n% be used for the same task.  The three algorithms have very\r\n% different roundoff error properties.\r\n\r\n%% Pete Stewart\r\n% As I did in my previous post, I am using Pete Stewart's book\r\n% _Matrix Algorithms, Volume I: Basic Decompositions_.\r\n% His pseudocode is MATLAB ready.\r\n\r\n%% Classic Gram-Schmidt\r\n% The classic Gram-Schmidt algorithm is the first thing you might think\r\n% of for producing an orthogonal set of vectors.  For each vector in\r\n% your data set, remove its projection onto the data set, normalize what\r\n% is left, and include it in the orthogonal set.  Here is the code.\r\n% |X| is the original set of vectors, |Q| is the resulting set of\r\n% orthogonal vectors, and |R| is the set of coefficients, organized into\r\n% an upper triangular matrix.\r\n\r\n   type gs \r\n\r\n%% Modified Gram-Schmidt\r\n% This is a rather different algorithm, not just a simple modification\r\n% of classical Gram-Schmidt.  The idea is to orthogonalize against the\r\n% emerging set of vectors instead of against the original set.  There\r\n% are two variants, a column-oriented one and a row-oriented one.  They\r\n% produce the same results, in different order.  Here is the column\r\n% version,\r\n\r\n   type mgs\r\n   \r\n%% Householder Reflections\r\n% Compute |R| by applying Householder reflections to |X| a column at a\r\n% time.  Do not actually compute |Q|, just save the vectors that generate\r\n% the reflections.  See the description and codes from\r\n% <https:\/\/blogs.mathworks.com\/cleve\/2016\/10\/03\/householder-reflections-and-the-qr-decomposition\/\r\n% my previous post>.\r\n\r\n   type house_qr\r\n\r\n%% Comparison\r\n% For various matrices |X|, let's check how the three algorithms perform\r\n% on two tasks, accuracy and orthogonality.  How close is |Q*R| to |X|?\r\n% And, how close is |Q'*Q| to |I|.\r\n\r\n   type compare.m\r\n   \r\n%% Well conditioned.\r\n% First, try a well conditioned matrix, a magic square of odd order.\r\n\r\n   n = 7;\r\n   X = magic(n)\r\n\r\n%%\r\n% Check its condition.\r\n\r\n   kappa = condest(X)\r\n   \r\n%%\r\n% Do the comparison.\r\n\r\n   compare(X);\r\n   \r\n%%\r\n% All three algorithms do well with accuracy, but classic Gram-Schmidt\r\n% fails with orthogonality.\r\n   \r\n%% Poorly conditioned.\r\n% Next, try a Hilbert matrix that is poorly conditioned, but not exactly\r\n% singular.\r\n\r\n   n = 7;\r\n   X = hilb(n)\r\n  \r\n%%\r\n% Check its condition.\r\n\r\n   kappa = condest(X)\r\n   \r\n%%\r\n% Do the comparison.\r\n\r\n   compare(X)\r\n\r\n%%\r\n% All three algorithms do well with accuracy. \r\n% Classic Gram-Schmidt fails completely with orthogonality.\r\n% The orthogonality of MGS depends upon |kappa|.\r\n% Householder does well with orthogonality.\r\n\r\n%% Singular.\r\n% Finally, an exactly singular matrix, a magic square of even order.\r\n\r\n   n = 8;\r\n   X = magic(n)\r\n\r\n%%\r\n% Check its rank.\r\n\r\n   rankX = rank(X)\r\n   \r\n%%\r\n% Do the comparison.\r\n\r\n   compare(X)\r\n   \r\n%%\r\n% Again, all three algorithms do well with accuracy. \r\n% Both Gram-Schmidts fail completely with orthogonality.\r\n% Householder still does well with orthogonality.\r\n   \r\n%% Conclusion\r\n% All three of these algorithms provide |Q| and |R| that do a good job\r\n% of reproducing the data |X|, that is\r\n%\r\n% * |Q| * |R| is always close to |X| for all three algorithms.\r\n\r\n%%\r\n% On the other hand, their behavior is very different when it comes to\r\n% producing orthogonality.\r\n%\r\n% * Classic Gram-Schmidt.  |Q'*Q| is almost never close to |I|.\r\n% * Modified Gram-Schmidt.  |Q'*Q| depends upon condition of |X|\r\n%   and fails completely when |X| is singular.\r\n% * Householder triangularization.  |Q'*Q| is always close to |I|\r\n   \r\n%% Reference\r\n% G. W. Stewart,\r\n% _Matrix Algorithms: Volume 1: Basic Decompositions_,\r\n% SIAM, xix+458, 1998.\r\n% <http:\/\/epubs.siam.org\/doi\/book\/10.1137\/1.9781611971408\r\n% http:\/\/epubs.siam.org\/doi\/book\/10.1137\/1.9781611971408>\r\n\r\n   \r\n\r\n\r\n##### SOURCE END ##### 0e62a7e371c643b1a14335ce8266b39b\r\n-->","protected":false},"excerpt":{"rendered":"<!--introduction--><p>This is a follow-up to <a href=\"https:\/\/blogs.mathworks.com\/cleve\/2016\/10\/03\/householder-reflections-and-the-qr-decomposition\/\">my previous post.<\/a> Classical Gram-Schmidt and Modified Gram-Schmidt are two algorithms for orthogonalizing a set of vectors.  Householder elementary reflectors can be used for the same task.  The three algorithms have very different roundoff error properties.... <a class=\"read-more\" href=\"https:\/\/blogs.mathworks.com\/cleve\/2016\/10\/17\/compare-gram-schmidt-and-householder-orthogonalization-algorithms-2\/\">read more >><\/a><\/p>","protected":false},"author":78,"featured_media":0,"comment_status":"open","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":[],"categories":[11,9,6,16],"tags":[],"_links":{"self":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/2027"}],"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=2027"}],"version-history":[{"count":1,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/2027\/revisions"}],"predecessor-version":[{"id":2028,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/2027\/revisions\/2028"}],"wp:attachment":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/media?parent=2027"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/categories?post=2027"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/tags?post=2027"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}