{"nbformat_minor": 0, "worksheets": [{"cells": [{"source": ["Matching Pursuits\n", "=================\n", "\n*Important:* Please read the [installation page](http://gpeyre.github.io/numerical-tours/installation_matlab/) for details about how to install the toolboxes.\n", "$\\newcommand{\\dotp}[2]{\\langle #1, #2 \\rangle}$\n", "$\\newcommand{\\enscond}[2]{\\lbrace #1, #2 \\rbrace}$\n", "$\\newcommand{\\pd}[2]{ \\frac{ \\partial #1}{\\partial #2} }$\n", "$\\newcommand{\\umin}[1]{\\underset{#1}{\\min}\\;}$\n", "$\\newcommand{\\umax}[1]{\\underset{#1}{\\max}\\;}$\n", "$\\newcommand{\\umin}[1]{\\underset{#1}{\\min}\\;}$\n", "$\\newcommand{\\uargmin}[1]{\\underset{#1}{argmin}\\;}$\n", "$\\newcommand{\\norm}[1]{\\|#1\\|}$\n", "$\\newcommand{\\abs}[1]{\\left|#1\\right|}$\n", "$\\newcommand{\\choice}[1]{ \\left\\{  \\begin{array}{l} #1 \\end{array} \\right. }$\n", "$\\newcommand{\\pa}[1]{\\left(#1\\right)}$\n", "$\\newcommand{\\diag}[1]{{diag}\\left( #1 \\right)}$\n", "$\\newcommand{\\qandq}{\\quad\\text{and}\\quad}$\n", "$\\newcommand{\\qwhereq}{\\quad\\text{where}\\quad}$\n", "$\\newcommand{\\qifq}{ \\quad \\text{if} \\quad }$\n", "$\\newcommand{\\qarrq}{ \\quad \\Longrightarrow \\quad }$\n", "$\\newcommand{\\ZZ}{\\mathbb{Z}}$\n", "$\\newcommand{\\CC}{\\mathbb{C}}$\n", "$\\newcommand{\\RR}{\\mathbb{R}}$\n", "$\\newcommand{\\EE}{\\mathbb{E}}$\n", "$\\newcommand{\\Zz}{\\mathcal{Z}}$\n", "$\\newcommand{\\Ww}{\\mathcal{W}}$\n", "$\\newcommand{\\Vv}{\\mathcal{V}}$\n", "$\\newcommand{\\Nn}{\\mathcal{N}}$\n", "$\\newcommand{\\NN}{\\mathcal{N}}$\n", "$\\newcommand{\\Hh}{\\mathcal{H}}$\n", "$\\newcommand{\\Bb}{\\mathcal{B}}$\n", "$\\newcommand{\\Ee}{\\mathcal{E}}$\n", "$\\newcommand{\\Cc}{\\mathcal{C}}$\n", "$\\newcommand{\\Gg}{\\mathcal{G}}$\n", "$\\newcommand{\\Ss}{\\mathcal{S}}$\n", "$\\newcommand{\\Pp}{\\mathcal{P}}$\n", "$\\newcommand{\\Ff}{\\mathcal{F}}$\n", "$\\newcommand{\\Xx}{\\mathcal{X}}$\n", "$\\newcommand{\\Mm}{\\mathcal{M}}$\n", "$\\newcommand{\\Ii}{\\mathcal{I}}$\n", "$\\newcommand{\\Dd}{\\mathcal{D}}$\n", "$\\newcommand{\\Ll}{\\mathcal{L}}$\n", "$\\newcommand{\\Tt}{\\mathcal{T}}$\n", "$\\newcommand{\\si}{\\sigma}$\n", "$\\newcommand{\\al}{\\alpha}$\n", "$\\newcommand{\\la}{\\lambda}$\n", "$\\newcommand{\\ga}{\\gamma}$\n", "$\\newcommand{\\Ga}{\\Gamma}$\n", "$\\newcommand{\\La}{\\Lambda}$\n", "$\\newcommand{\\si}{\\sigma}$\n", "$\\newcommand{\\Si}{\\Sigma}$\n", "$\\newcommand{\\be}{\\beta}$\n", "$\\newcommand{\\de}{\\delta}$\n", "$\\newcommand{\\De}{\\Delta}$\n", "$\\newcommand{\\phi}{\\varphi}$\n", "$\\newcommand{\\th}{\\theta}$\n", "$\\newcommand{\\om}{\\omega}$\n", "$\\newcommand{\\Om}{\\Omega}$\n"], "metadata": {}, "cell_type": "markdown"}, {"source": ["This numerical tour explores the use of Matching Pursuits to solve\n", "sparse optimization problems.."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 2, "cell_type": "code", "language": "python", "metadata": {}, "input": ["addpath('toolbox_signal')\n", "addpath('toolbox_general')\n", "addpath('solutions/optim_9_matching_pursuit')"]}, {"source": ["Sparse Regularization\n", "---------------------\n", "We consider the following linear inverse problem\n", "$$ y = \\Phi x_0 + w  \\in \\RR^P$$\n", "where $x_0 \\in \\RR^N$ is the (unknown) signal to recover, $w \\in\n", "\\RR^P$ is a noise vector, and $\\Phi \\in \\RR^{P \\times N}$ models the\n", "acquisition device.\n", "\n", "\n", "Size of the problem."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 3, "cell_type": "code", "language": "python", "metadata": {}, "input": ["N = 200;\n", "P = round(N/4);"]}, {"source": ["We consider here a simple compressed sensing scenario where the operator\n", "$\\Phi$ is a realization of a Gaussian random matrix."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 4, "cell_type": "code", "language": "python", "metadata": {}, "input": ["Phi = randn(P,N);"]}, {"source": ["For simplicity of the exposition of the algorithm, we normalize the\n", "columns of the operator."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 5, "cell_type": "code", "language": "python", "metadata": {}, "input": ["Phi = Phi ./ repmat( sqrt(sum(Phi.^2)), [P 1] );"]}, {"source": ["Sparsity of the input data."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 6, "cell_type": "code", "language": "python", "metadata": {}, "input": ["s = round(P/5);"]}, {"source": ["Sparse vector."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 7, "cell_type": "code", "language": "python", "metadata": {}, "input": ["I = randperm(N); I = I(1:s);\n", "x0 = zeros(N,1); x0(I) = sign(randn(s,1));"]}, {"source": ["Noise level."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 8, "cell_type": "code", "language": "python", "metadata": {}, "input": ["sigma = 0.05 * norm(Phi*x0)/sqrt(P);"]}, {"source": ["Measurements $y=\\Phi x_0 + w$."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 9, "cell_type": "code", "language": "python", "metadata": {}, "input": ["y = Phi*x0 + sigma*randn(P,1);"]}, {"source": ["In the following, we define the $\\ell^0$ pseudo-norm the counting\n", "measure\n", "$$ \\norm{x}_0 = \\abs{\\enscond{i \\in \\{1,\\ldots,N\\}}{x_i \\neq 0}}. $$\n", "\n", "\n", "To recover an approximation of the signal $x_0$, we use $\\ell^0$\n", "sparse optimization, which requires to solve either the sparsity\n", "constraint problem\n", "$$ \\umin{\\norm{x}_0 \\leq s} \\norm{y-\\Phi x}, \\qquad \\Ss(s) $$\n", "or the $\\ell^2$-error constraint problem\n", "$$ \\umin{\\norm{y-\\Phi x} \\leq \\epsilon} \\norm{x}_0, \\qquad \\Ll(\\epsilon)  $$\n", "Note that these are non-convex optimization problems.\n", "\n", "\n", "One needs to set-up either the sparsity parameter $s \\in \\NN$ or the\n", "error contraint $\\epsilon \\geq 0$. These parameters should be\n", "tuned in accordance to the noise level, and when there is not noise,\n", "$w=0$, one should set $ \\epsilon=0 $ and solve the problem\n", "$$ \\umin{\\Phi x = y} \\norm{x}_0.  $$\n", "\n", "\n", "In the special case where $ \\Phi \\in \\RR^{N \\times N} $ is an\n", "orthogonal basis, i.e. $N=P$ and $ \\Phi \\Phi^* = \\text{Id}_N $, one\n", "can compute in closed form the solution to problems $(\\Ss(s))$ and $(\\Ll(\\epsilon))$\n", "using the hard thresholding operator\n", "$$ x = \\Phi \\circ S_T \\circ \\Phi^* (y)\n", "\\qwhereq\n", "  S_T(u)_i = \\choice{\n", "      u_i \\qifq \\abs{u_i}>T, \\\\\n", "      0 \\quad \\text{otherwise}.\n", "  } $$\n", "where the threshold $T \\geq 0$ should be selected so that\n", "either $\\norm{x}_0=s$ or $ \\norm{y-\\Phi x}=\\epsilon $.\n", "\n", "\n", "In the general case howeve, both problems $(\\Ss(s))$ and\n", "$(\\Ll(\\epsilon))$ are NP-hard to solve, see for instance:\n", "\n", "\n", "S. Mallat and G. Davis, and M. Avelaneda,\n", "_Adaptive Greedy Approximations_,\n", "Jour. of Constructive Approximation, vol. 13, No. 1, pp. 57-98, 1997.\n", "\n", "\n", "One thus needs to use alternative strategies to compute approximate solutions.\n", "\n", "\n", "A first line of methods consists in replacing the $\\ell^0$ non-convex\n", "pseudo-norm by its natural convexification, the $\\ell^1$ norm. Other\n", "numerical tours are devoted to this class of methods.\n", "\n", "\n", "Another line of method develop heuristic greedy optimization schemes,\n", "that are usually called Matching Pursuit algorithms. We describe here two\n", "of such methods.\n", "\n", "\n", "Matching Pursuit\n", "----------------\n", "Matching pursuit is a greedy procedure that progressively identify the\n", "location of the spikes by looking at atoms that maximaly correlated with\n", "the current residual.\n", "\n", "\n", "This algorithm is introduced in\n", "\n", "\n", "S. Mallat and S. Zhang,\n", "_Matching Pursuits With Time-Frequency Dictionaries_,\n", "IEEE Transactions on Signal Processing, December 1993.\n", "\n", "\n", "The algorithm is initialized with $x^{(0)}=0$."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 10, "cell_type": "code", "language": "python", "metadata": {}, "input": ["x = zeros(N,1);"]}, {"source": ["The vector is updated from step $\\ell$\n", "to step $\\ell+1$ as\n", "$$ x^{(\\ell+1)} = x^{(\\ell)} + \\mu $$\n", "by adding a 1-sparse vector $\\mu$ that minimizes the\n", "error decay\n", "$$ \\umin{\\norm{\\mu}_0 = 1} \\norm{ y - \\Phi (x^{(\\ell)} + \\mu) }.  $$\n", "\n", "\n", "This update is computed in closed form as\n", "$$ \\mu = c_i \\de_i $$\n", "where the correlation vector is\n", "$$ c = \\Phi^* (y-\\Phi x^{(\\ell)}) \\in \\RR^N $$\n", "and where $\\de_i$ is the Dirac at position $i$ that maximizes the\n", "absolute correlation\n", "$$ \\umax{ i } \\abs{c_i}  $$\n", "\n", "\n", "Compute the correlation."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 11, "cell_type": "code", "language": "python", "metadata": {}, "input": ["c = Phi'*(y-Phi*x);"]}, {"source": ["Display the correlation."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [{"metadata": {}, "png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGwCAIAAADOgk3lAAAACXBIWXMAAAsSAAALEgHS3X78AAAA\nIXRFWHRTb2Z0d2FyZQBBcnRpZmV4IEdob3N0c2NyaXB0IDguNTRTRzzSAAAQ+UlEQVR4nO3d23aj\nxhqFUdij3/+V2RekFUXCWAeg/lU150WG2+6OEYf6KISkeVmWCQDS/K/1AgDAJwQMgEgCBkAkAQMg\nkoABEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACRBAyASAIG\nQCQBAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYAAtzfM0z60XItOf\n1gtwrtl+AVS3TCUHq2VZWi/CLzoP2NR6G8zzXGEnsBgWo+wyWIy1XLdfXWdttF6E37mECNBSgVql\nEjAAIgkYAJFKXGw9T5GryQA75rnchcSIwdMMDIBIAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACR\nBAyASAIGQCQBAyCSgAEQScAAiCRgAEQSMAAipQZsnufNb66uXx4ALvan9QK87ac+rd9fP4Et4qPY\nAPhGXsBuidr86X3GAOhYXsB23Lp1PwN7SJ22AWyKe/6lq4BtUiyAV9yPlhExS72JA4DB9TADWy8Y\nLstyO2Uw6wLoXmrA7hN1+1q3AFbr+Xzfg6JLiABEEjCADvU991oJGACRBAyASAIGQCQBAyCSgAEQ\nScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZAJAED\nIJKAARBJwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQScAAiCRgMIR5nua59ULAoVIDNv98\nLO78CIBu/Gm9AG/b75N6waZlMQOjN3kzsGVZlmXZ/NE8zz/9CIDO5M3A3vUwJ1M4gE1xV7D6Cdi6\n6m//vYVKsQBecT9aRsSsn4DdVr0LiQAjyHsO7FnEmQIAx0oN2P0c62G+ZfoFMILUgAEwOAEDIJKA\nARBJwACIJGDvccMjQBEC9jYNA6hAwACIJGBv8zIz7vmcLWhFwACIJGDwFTNyaEXAAIgkYABEEjAA\nIgkYAJEEDIBIAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACRBAx4j0+QoQgBAyCSgAFv8yEyVCBg\nAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIiUGrB5661s5r+uXx4ALvan9QK8bb9Py7Ksf2fx\nVgEwtnWoMBJ0LG8GtizLT3ESLYBx5M3AfvUw/XqYsYkcDGJZvGv+e+Kef+kqYOvaf0iUYsEZXKDr\nz86pf015lxD3yRXAIHoI2HqmcPuvGxHhAs4VaS71EuL9TGv92twLYCg9zMAAGJCAARBJwACIJGAA\nRBIwACKl3oUI8LrbK2vcrdwTMzAAIgkY0L914mX61RkBAyCSgHXOm2oBvRKwnqkX0DEBAyCSgPXM\nU9ZAxwQMGNQ8u8yeTcAAiCRg13G6B6W4xp5OwACIJGDXcboHcCABAyCSgAEQycepAFThY1/eYgYG\nUIVuvUXAAGqRsRcJGACRBAy65bXz9E3AoE/SRfcEDPrkeRS6J2DQMxmjYwIGQCQBAyCSgAEQScAA\niCRgAETqMGCz178ADKCrd6OXLoBxdDUDW5Zl8bIXgDF0NQPb9DAtUziATXEXsfoPmGIBvOJ+tIyI\nWVeXEAEYh4ABEKnDS4iuGUKu24UrxzG/MgMDylEvXiFgQCHSxesEDIBIAgZAJAEDIJKAARBJwACI\nJGAARBIwACIJWF3zPCW8nSZAGwJWlHQB7BOworwfAcA+AStNxgB+ImAARBIwACIJGACRBAyASAIG\nQCQBAyDSn9YLAPTs9pJ8rwnhcGZgAEQSsDDeYoos68TL9IszCFiStV4aBjAJWCInswCTgGWRLoAb\nAQMgkoABP/KpqlQmYABEEjDowUlTJU+7HsVc9gwCBkAkAQM4nbnsGQQMemB8ZEACBkCkvHejn/8+\nE7psnXOuP938EQA9CQvYPM+3ON1//fCd5x8B0BmXEAGIFDYD+9XzJcT5vy++MDMD2DSnvVStn4D9\ndHVRsQBesXPqX5NLiHTCOx3AaMJmYMuyPN+FuM63Nn8EQK/CAjZtxcnVQla2PwzFJcQhuLYG3XA4\n3whY/9bd3U4P3XA4r0YP2DjP/Lu8BnRm9ICNYDNd45QbzrAeQa0OIuejq9EDZj8APmYAaWv0gA3L\ngQffOPsIco3kFQIGQCQBK8Q5F3Dz8SRvnJEk74XMwL5XBq/177iSfJJ5tm6vYAZWiD0eOuCVl5cR\nMOjNeib06/mQE6aTvLj+z16AEQgYQCfGefZrJWAARBIwgH6Mc/1wEjDeNdQFCihrtKuFmwSMNzhg\ngDoEDCDPUJcKfyJgvMEx0yUXowglYABEEjAYnYk1oQQMgEgCBkAkAQMgkoABEEnAAIgkYABEErAM\nXmoK8EDAgP9wtkQKActQ7aWmBrhe2bIEETDeto5xRrouradK1U6YYJOAVVf2ek5/Y1zZVQ1sErBv\nraPeUANfkXSNttppyM5Wk4AFKBKM7lnPbJKusvICNv/11o/O4zmDVqxzruEYL+tP6wV4zzzPy9/9\n6P7r9Y/TNK3fefgRAP0JC9iv7jMGQMe6Ctjm5OzhiqK2wc16cDgmWF38/Mv3ugrYJsUCeMXzkzLF\n5d3EUZzbbQGuETYDW5bldl7wcMFw80dU4DrVap5PWQlfrl7bhVxhAZu24nT7jm7xotss+bJd5vb+\nW8f+RtN9RuYS4sE09Jl1Mp32WqJq69YldK4kYLR3/ajX2UtTu3kg8BYB42BOwAfnGimXETCO9Nkn\nrZhAsMkH97BPwDieIHEIOxL7BIwjXTPiuFNgKDLGTwRsg/ERoD4B4yIHnhY4JQcmAfuJIRKgOAHj\nIs4JgGMJGACRBAw+5GYfaEvAoB+aylAEDD7niT1oSMB+4ZQWoCYBe4mGxRlzk5kRMhQB+0Vnn7sx\niDHrBaMRMF4iCUA1Ava50Z4eC3qwZsxUMNoQcT0B+9bHO+j1e7bDqRpbpI51Wxy4Oa7ZsoPvQgL2\nuc+eHrs/Tq7c877/XaY1dOzw3dvT5xcQsDbO27l/CpXDqSCbo5q4LRK3wMcSsKudusP5CHZgHALW\nlcFPx4Ch/Gm9ABwvPWPzXPohmOBCEWZg1FL/KmjluMJQBOwKg9/q+paUPKQsJ3RMwE4nXR+QB+BX\nAnY6968DnEHAAIgkYABEEjAAIgkYAJHyAjb/tfMXrlweAJoIeyeOeZ6Xv/fz3X99/xe2/tU0uQ8Q\noC9hAdu3Js0MDH7i4KAnXQVs09+e/RO250kb1J+j119COhB39t9PwNZVf/vvLVTrF/MsXTAtyz/v\nlZw2Ur1H7z9zP0hGxKyfgO0/N0YfzhuY6u8y9ZcQLhYWsPunuBQL+EkfQ4Kp5L6wgE1bT2I9fEfM\n+mbzAqu814ExmraX4n0UDg05XdsnYMFGGFu7f4DQXO5IImDwC2fBUJOAZet+bO3+AUIFoQeagDFN\nrtQxjNzLZTwTMP45nh3V8IFriujw3CRg/HP14ONrCDsHsKOOgg68XHbNHu4U8ycCxlkcdaT4eC/9\n8uTvg9/FPQHjAJuHluONCPXPtBxKPxEwzuXYI4IdNZGAAUOTrlx574UIMKDKFzlbMQP7D68RgcGV\nHQHMFJ8J2L/K7rh8xulINy7elOuvq7nzyNg9AfvXlXfEcpmHYaj5wNR8AaAbIwbMCDIIZyTduHgj\nLovdJsNwAZOuUEdtuAoD05fL4AwMVsMFzFl5oiIvNVUOa4BShgsYfGnYEfyDBz7suuIaAnaiyvcy\nZSk1b353MQ7fDa5ZD8/L/MEDnzSMMwlYHl1s5ZtyFKnv9959IN08cAoSsBO9Pm94K0jRI8KA9f1s\nxlbBl3vaGTtqkTVDEd5KKtWynHsw3/7n0b3kMt/sjeu//XVPe7gmuf79F/8tXa4oM7DhtJ3tlXo2\n62yfzaUSV46JUXG9biABK6HgmFVwkahzdfHB9+clr/zb+99i/3xLrytNwPhFzREzQpHxomz2DlFk\nJdOEgI3o9WO+pzuhOxjHDdZwz00c7Dn7VhHSlbo1oOG+Wmo9jEPA+N1mxmqGbX8c6Wl8eWv9n/3A\na+4MdM8lRD7UUwz4WKm7A9ouTJGVMBQBayz91NVB20qpcpCig2eC7wnYKV7cS3q6RWIq80CM6TAI\nAWupp6G2SL3ow4BvOXaNnsacKfEmjvnvTr08bYqdH13srZv3Wi8sQKSwGdg8z8tf81Yidn7EqWSY\nA3mGj1eEBWxf84kXQBNjnrTnXUL81TpLu//jNE3T9M+0TORK8fJPetJqfz7qdrC4a1ddBWwzUesf\n51m6APZsnfqX1tUlxMkEC2in7WuoBxz8wmZg9zdo3Fq1XjNcv1/nRkReMex7Lbp2Ct8LC9i0Vab1\nO98Ua9jRZNgHXooPv65snm2Xunq7hPiNMacCU6UHXmdJzmZMjHDle+V44fYHBOxfA44ppS6dd/bG\nWq8rtRV4ZtOUJWDTZAcFtlw5MhiFPiBgnzDZP4OJCBSxviKp/iiXdxMHjKP+CNKQW5AwA/vE5jFj\nWpbIVoNnKZdDBOwYxw6Cb/3fDMEfq7/eUsaRVqyZwQnYMQ4caIa9Ge968gA3iafCAlbOu6Pq4ONv\n4lEHHMJNHPRAw+BdD0dN4vu6mYHxKGtO4zIgQbIOrvoEDOhZwWYUWZ4OTv4EjEfRO3RzBYfLIqyZ\nqYtmlCJgfTJY3LtsbVjnBalFxwSsH0bPTVeuFufX+6yZG6eYhxCwTjy/esxgsWoSFcMTO+wbRxGw\nThwyTBt267AtOmamfpRxA2Z04HUf5KSn4UlNqWmggK0H4f1x6JikrIL9c7xQzUABWxUcF0qxfp5P\ndKYxVsvmA1+55PUN89fzDBSw58Nv/4C02wFUNlDA4BXDzjaueeA754UdnzIOuDtdo/uAff6p2A/7\nXK+HFkCo7gN2DJ/Rdb3m5+M29xnMRThQ/wE75IBx1I3GKcv1HGUn6fgT3vsP2IEcYFdqu7Zta/rQ\n96mYD7TkFH0cMIkf8RfEum1lXfObZ2lZ+7wZGMcLOgCOsvMiKp6t66rhNHecjdX3XbUCxvFeOWa6\nHEG+fy/KoezMA+AVAgYH+H4U7vu5igd9Twvq62bNCxjNdHMUHcLauFKRtd3ldYgrCRgUUmRgZYfq\n1CFgAM04ZfmGgAG8QXLq6Op1YPPfif1iF4N2er29sNfHlaufgM3zfOvW/ddA96RlTP0M9JsBmz3Z\nCvCR+nXoZwa2qf4GAOAzbuIAIJKAARCpn+fAJnchAoykq4ABMI7Ob+Jo4nki2HBq+HxD5sXL8PB7\n266Ktouxsy2uXJ5qi/H8x8sWY+eXtl0VFy/Gr3tC2YtbAnaK21ZflqXhC9Ruu12rZXg+FBuuioaL\ncf9yjudtcdnWeX5VSZMd9Xkx9tfPBcvQakfdXBVNdtSdPeHKxXiXmziOV2QD19nV5nmuszBNLMtS\n4eE/LEarRXpYjCb7xuYWuX5HrbNjtF6ED5mBnWXwIfvm+aJEK7dzyeZLUooddWVHTdwTBOx4D5Pu\nhstw/9+RbT7VQcEdNW4APVarHbXCnvAZATtF813h+VmWtstDTdV21LYLM7LQlS9gB3uY9KzXuJvf\nw9NqGZ5/78iL8dPCNNk6NXfUqdHaKLKHNFmMV/aECjvGJmc9AERyFyIAkQQMgEgCBkAkAQMgkoAB\nEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACRBAyASAIGQCQB\nAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBI\n/we599iz4n8SXgAAAABJRU5ErkJggg==\n", "output_type": "display_data"}], "prompt_number": 12, "cell_type": "code", "language": "python", "metadata": {}, "input": ["clf; \n", "stem( abs(c), 'b.' );\n", "axis tight;"]}, {"source": ["Extract the coefficient with maximal correlation"], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 13, "cell_type": "code", "language": "python", "metadata": {}, "input": ["[~,i] = max(abs(c));"]}, {"source": ["Update the coefficients"], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 14, "cell_type": "code", "language": "python", "metadata": {}, "input": ["x(i) = x(i) + c(i);"]}, {"source": ["__Exercise 1__\n", "\n", "Perform the Matching pursuit. Record the evolution of the solution $x^{(\\ell)}$\n", "as well as the maximum correlation $ \\la^{(\\ell)} = \\abs{c_i} $.\n", "Display the decay of the error $ E(\\ell) = \\norm{y-\\Phi x^{(\\ell)}}. $"], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [{"metadata": {}, "png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGwCAIAAADOgk3lAAAACXBIWXMAAAsSAAALEgHS3X78AAAA\nIXRFWHRTb2Z0d2FyZQBBcnRpZmV4IEdob3N0c2NyaXB0IDguNTRTRzzSAAAOf0lEQVR4nO3d0XKj\nRhBAUUj5/3+ZPGgtY0CWhIDp7jnnIeW4nNSsVtZ1DyM8TtM0AEA2/7VeAADsIWAApCRgAKQkYACk\nJGAApCRgAKQkYACkJGAApCRgAKQkYACkJGAApCRgAKQkYACkJGAApCRgAKQkYACkJGAApCRgcIpx\nHN/6yvG3t/4P0Kev1gsA/pmmafGv4zguPgncCRic6D5F3Ts0n6ueJkrD4A8CBmeZt+f28eIz66+/\nfyxa8JSAQRSiBW9xiAOAlExgcJbbnuH948Vn1mwhwltcH4Y21lfFHn3NlauCRExgcKnNmUylYAff\nOQCk5BAHACkJGAApCRgAKQkYACkVP4X4feDLXb0B3hP/iF/xgH0L/ReR5RS1dR7LOo9lncdK8dt8\nbCECkJKAAZBS+YAlmIIB2CHHbuxu4zgOw78/YOk/KMCRUlyrKz+BAVCTgAGQUv2AhR+CAdijfsDu\nMryrAYBXdRQwACrpImD3XURDGEAZXQQMgHp6CZghDKCYXgIGQDEdBcwQBlBJRwEDoJK+AmYIAyij\nr4ANGgZQRXcBA6CGHgNmCAMooMeAAVBApwEzhAFk12nAAMiu34AZwgBS6zdgAKRWKmDjtxe/3hAG\nkFedgI3jOH3TMIDy6gTscxoGkMhX6wWcbjGNTfeZ6+cz0gWwfLWMr37A1sXa+pp/DRvH4YUvByho\n/mqZIma2EJcy/K0BUGgCm5/deGXqWv3nP+kyhwHEVydgw65u/f7PjV8AadhC/MWpeoAsBOwhDQOI\nTMCW5tuQGgYQloBt0DCA+ARsm1OIAMEJ2EMOdABEJmB/0TCAsATsVRoGEIqAPeFiGEBMAvacjUSA\ngATsPRoGEISAvcQ7wwCiEbBXTZO9RIBABOw9GgYQhIC9TcMAIhCwj2gYQCsCtoczHQDNCdhOGgbQ\nloDtp2EADQnYRzQMoBUB+5SbJQI0IWAHmB+sN4cBXEPAjmEvEeBiAnYYDQO4koAdScMALiNgB9Mw\ngGsI2PE0DOACAnYKDQM4m4CdRcMATiVgJ9IwgPMI2Lk0DOAkAnY6DQM4g4BdYdEwGQP4nIBdZJqM\nYgBHErBLaRjAUQTsahoGcAgBa0DDAD4nYG1oGMCHBKwZDQP4hIC15Hg9wG4C1pjj9QD7CFgIRjGA\ndwlYFPOGDUYxgGcELJD1dqKMATwiYOEYxQBeIWARGcUAnhKwuBxQBPhDwYCNtV7pNQxgU6mAjeNY\nrF43GgawVipg0zRNiyMQVbgkBrDw1XoBp1vMZHkLd1v4/U8zjsvzigCfSLeDlTVgr2cpb7E2TZOG\nAaeYv1qmiFnWgBXL0ls0DGAodg2sHy6JARQMWCfDmRt2AJ0rGLB+eKcz0DMBS0/DgD4JWAUuiQEd\nErAibCcCvRGwUjQM6IeAVWM7EeiEgBVkOxHogYCVpWFAbQJWmYYBhQlYcRoGVCVg9TnWAZQkYF1w\nrAOoR8A6omFAJQLWFw0DyhCw7rgkBtQgYD3yu8SAAgSsU451ANkJWNc0DMhLwHq32E4EyELA+GmY\nMx1AIgLGMNhLBBISMP7RMCAXAeOHhgGJCBi/aBiQhYCxpGFACgLGBg0D4hMwtrllIhCcgPGQ200B\nkQkYT2gYEJOA8ZyGAQEJGC/RMCAaAeNVGgaEImC8QcOAOASM92gYEISA8TYNAyIQMPbQMKA5AWMn\nDQPaEjD20zCgIQHjIxoGtCJgfErDgCYEjANoGHA9AeMYfv0KcDEB4zB+/QpwJQHjYBoGXEPAOJ6G\nARcQME6hYcDZBIyzaBhwKgHjRBoGnKdUwMZvrRfCDw0DTlIqYMMwTNM0TZOGhaJhwBlKBWyav1IS\niYYBh/tqvYDjjeM4L9liGhO5VqbpJ13jOPh7gGjS7V1lDdhmlm6fXCRKseLQMIjsjx/9Y8oasEdZ\nkqvgNAw4Sp1rYLefFxxEjM/1MOAQWSewNbNXIuYw4HN1JjByMYcBHxIwmtEw4BMCRksaBuwmYDTm\nVzkD+wgY7flVzsAOAkYUGga8RcAIRMOA1wkYsWgY8CIBIxzHOoBXCBgRrY91yBiwIGDEtbjFlIYB\ncwJGaE7YA48IGAnYTgTWBIwcjGLAgoCRiYYBdwJGMhoG3AgY+WgYMAgYSWkYIGBkpWHQOQEjMQ2D\nngkYuWkYdEvASE/DoE8CRgUaBh0SMIrQMOiNgFGHWyZCVwSMUtwyEfohYBRkFIMeCBg1rX8ZpoxB\nMQJGWYvtxMGOItQiYBTnqhhUJWB0wVUxqEfA6IVRDIoRMPqiYVCGgNEd24lQg4DRo/V2ooxBOgJG\nvxyyh9QEjK4ZxSAvAQOjGKQkYDAMDtlDQgIGP2wnQiICBr8YxSALAYMNGgbxCRhs0zAITsDgIQ2D\nyAQM/qJhEJaAwRMaBjEJGDzneD0EJGDwEsfrIZpSARu/tV4IZWkYxFEqYMMwTNM0TZOGcR4NgyC+\nWi/gSNPinqzDMAzDImabXwNvmaafdN0+8LSigHQ/+pcK2LD1F6BYnOH2tJpnzBON7OavlililjVg\nj+aq2wcpHnoKMIpBQ1kDtp6rxnE0bHG9ecMGGYMLZQ3Y2vzshpJxpcV24iBjcIk6ARt0i6Y2M+Yp\nCeepdowe2lq/39kFWTiJgMHxFoOXhsEZBAxO4dZTcDYBgxNpGJxHwOBc81HMJTE4kIDBFYxicDgB\ng4s4nQjHEjC4jtOJcCABg0t5oxgcRcCgAaMYfE7AoA1vFIMPCRi0pGGwm4BBYxoG+wgYtOdYB+wg\nYBCCYx3wLgGDKJywh7cIGMTidCK8SMAgIqMYPCVgEJSrYvA3AYO4bCfCHwQMorOdCJsEDBKwnQhr\nAgY5OGQPCwIGmRjF4E7AIBknO+BGwCAl24kgYJDVejtRxuiKgEFii+3EwY4iPREwSM8BRfokYFCE\njNEbAYNSXBijHwIG1bgwRicEDGqyo0h5AgaVGcUoTMCgOKMYVQkYdMENqKhHwKAjRjEqETDoix1F\nyhAw6JHDHRQgYNApoxjZCRh0TcbIS8AAO4qkJGDAMBjFSEjAgB8yRiICBizZUSQFAQM2GMWIT8CA\nhxajGIRSMGCjHxThOPNRzPcWoVQLmHrBGeYN801GEKUCNo7jZMsDzuGSGNF8tV7A6RYzmcLBbtP0\nq1vj6CJZKel2sLKOLOssbYbKTAZnWLzQ+SarJ8WLZ9YJbP3I3j+T4nGH1G7fYfeMGcVootQ1MOBK\nrorRVsGAGb/gMuv3O8NlCgYMuJhRjCYEDDiA2ydyPQEDjuH2iVxMwIAjGcW4jIABB3Oyg2sIGHAK\n24mcTcCAs7gqxqkEDDiXq2KcRMCA0xnFOIOAARcxinEsAQOuYxTjQAIGXG09iskYOwgY0MBiFBvs\nKPI+AQOa8ZZnPiFgQGOuirGPgAHtGcXYQcCAKIxivEXAgECc7OB1AgbE4r1ivEjAgIhcFeMpAQPi\nMorxBwEDQnNVjEcEDIjOVTE2CRiQg6tiLAgYkImGcSdgQDK2E7kRMCAf24kMAgbkZRTrnIABiRnF\neiZgQHpGsT4JGFCBUaxDAgbUYRTrioABpbj1VD8EDKjGrac6IWBATa6KlSdgQGVGscIEDCjOVbGq\nBAyoz3ZiSQIG9MJ2YjECBnTEKFaJgAHdMYrVIGBAj5zsKEDAgE55v3N2AgZ0zSiWl4ABvXOyIykB\nAxgGJzsSEjCAf1wVy0XAAH5xVSyLr9YLONI4e6JNi+cgwMturx/3V5RxXFaNCEoFbNAt4DjT9Kth\nw2o4o61qW4jjOI4GfuAgropFVnMCG8fxPootemZEA941H8WGujuK6X76zxqwzSxtxkmxgM+tr4oN\n5XYU56+WKWKWNWDrLM2nLoAz9JCxRLIGbG2apvuPDEoGnGe9ozjIWAt1AjboFnCVxSg21L0wFlm1\nU4gAl3FGsS0BA/jI+s4dMnYNAQP41GIUG9yA6hICBnAMO4oXEzCAI8nYZQQM4HgujF1AwABO4cLY\n2QQM4ER2FM8jYACnM4qdQcAArmAUO5yAAVxHxg4kYABXs6N4CAEDaGA9ivEuAQNoRsM+IWAALbkk\ntpuAATTmZMc+AgYQgrtPvUvAAKLYvPuUjD0iYACxLBo2yNgDX60XAMDSvWHzbt0+XuetWyYwgLjc\n0v4PAgYQnWOKmwQMIAfnOxYEDCANxxTnBAwgGRm7ETCAlJzvEDCAxHo+3yFgAOn1uaMoYAAVdLij\nKGAAdXS1oyhgANV0sqMoYAAF9bCjKGAAZdXeURQwgOKq7igKGEB9JW/eIWAAvSiWMQFrb0zy9LHO\nY1nnsazzdWXOdwgYQI8KnO8QMIB+rTOWiIAB9G5rFJsefnUY4zQlWOVuEbabAZL4lYP4cfhqvYBz\n1c4zwBmy/ORfPGAAvCvLT/6ugQGQkoABkJKAAZCSa2BtzI9H3k6a3D8T5+DJOP4cUl0sL9RqN9c5\nDMM0TUHWuV5GzMfzj3UOHs/3ZXk8b7J8v88JWDPzp8LiqdP8WbJ4+8H6mR1ktZtvk9hcW/NH9f5C\ncHvZivl4Dr/XOf/M4PHcJcvjOf9Wivx4LthCbGYcx7BvU5umKdTT9JHNdQZ8YFM8mMODdXo8d8vy\neEbL0utMYM2sfzTjEItNjziy/F0v1hnz8Yy2nkfW64z5eCYlYG2keCHLKOYDu9iHCWu9zrBrzpKB\n4NeQhu+13f8ZcIV/sIXYQPzvuqQiP7BZXhcWl2YbruSRmKtaW68z5sqnb0OeZ+ldst6W8fQYVQRZ\nTiUFX+fiZWtzYdb5lgKnEEOt8yb499EmAQMgJVuIAKQkYACkJGAApCRgAKQkYACkJGAApCRgAKQk\nYACkJGAApCRgAKQkYACkJGAApCRgAKQkYACkJGAApCRgAKQkYACkJGAApCRgAKQkYACkJGAApCRg\nAKQkYACkJGAApCRgAKQkYACkJGAApCRgAKQkYACkJGAApCRgAKQkYACkJGAApCRgAKT0P8ioI/a8\nxT1yAAAAAElFTkSuQmCC\n", "output_type": "display_data"}], "prompt_number": 15, "cell_type": "code", "language": "python", "metadata": {}, "input": ["exo1()"]}, {"collapsed": false, "outputs": [], "prompt_number": 16, "cell_type": "code", "language": "python", "metadata": {}, "input": ["%% Insert your code here."]}, {"source": ["One can shows that the error $\\norm{y-\\Phi x^{(\\ell)}}$ converges toward zero\n", "when $\\ell$ increases.\n", "\n", "\n", "\n", "Display the signal to recover and the signal recovered with a large number of iterations\n", "(thus corresponding to $\\epsilon=0$).\n", "Note how the noise impacts the recovered vector."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [{"metadata": {}, "png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGwCAIAAADOgk3lAAAACXBIWXMAAAsSAAALEgHS3X78AAAA\nIXRFWHRTb2Z0d2FyZQBBcnRpZmV4IEdob3N0c2NyaXB0IDguNTRTRzzSAAAPf0lEQVR4nO3d3Xab\nOqOGUbRH7/+WtQ/8lUWB4J9gpFfMebBGVuPYMpZ5ECZtqbVOAJDm/1oPAAA+IWAARBIwACIJGACR\nBAyASAIGQCQBAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAA\nIgkYAJEEDIBIAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACRBgxYKaX1EAD4uj+tB3Am6QK4j6FW\nYLXWWmvrUXSklEnToXNfep/e4e0/1Aps6/ZrsjrZCNC7L71Pf3u3/a8HBg/Y1Po1KKW0HsBUa20+\njL+DMYy+htHDGAzjkZj5oU8cxuPt/+nPBhz4DnUKESBOB+1OJWAARBowYD2ciwDg27o49fw9nZxb\nb6gUJyigd196n/7mbiN2ngOuwAC4AwEDIJKAARBJwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCS\ngAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZA\nJAEDIJKAARBJwACI9Kf1AN5WSnl8UWv96Vu73wVgJGEBK6XMZVp+PdMtgJsY7RRiKWW5DgNgVGEr\nsKceK7DVQm17AwBW4o7+hwrYbpwUC+AVy71lRMzGOYUYsbkBOEvYCqzWur0K8XHCcPdbAIwqLGDT\nXpzmP9EtgPsY5xQiALciYABEEjAAIgkYAJEEDIBIAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACR\nBAyASAIGQCQBAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAA\nIgkYAJEEDIBIf1oP4EyllMcXtda2IwHg28YJWCll7tbyawCG5BQiAJHGWYH9ZD6veFu33wAQ4Evv\n07Hf/uMHzLlEbqKU6XiyH9zg6c/+5qc+u/MvWQ2mq7F1JeLQ3ylEACKNswKrtboKEeA+xgnYpFsA\nd+IUIgCRBAyASAIGQCQBAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgk\nYABEEjAAIgkYAJEEDIBIAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQ\n6U/rAbytlPL4otb607d2vwvASMICVkqZy7T8eqZbADcx2inEUspyHQbAqMJWYE89VmCrhdr2BgCs\nxB39dx2wd9uzewPFAnjFcm8ZEbOuA/ZWe3Y/EgNgVF0HbKvWur0K8ZGu3W8BMKqwgE17cZr/RLcA\n7mO0qxABuAkBAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAA\nIgkYAJEEDIBIAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQScC+q5TS\negjTZBj/MoyuxjAZxr86GUYEAQMg0oABc/wCcAd/Wg/gTNIFcB+l1tp6DCcr5b8nJWncxmPOH0z4\ngxs8/dnPHnF5sw/u/0uWw/7sid9F/3UYPGBwE49DtYOJf3CDpz/72SMub/bB/X/JctifPXH6kbqv\nXy2tls9CwADuIPUzMIkCuLkBr0IE4A6cbQMgUuopxJ7Nn89tL4a8/nBh/kSw1RhWj9t2U7QdxsFr\nceV4ehvG9n8vG8bBg7bdFBcP4+lMaPi2PSZgXzG/6rXW1Zvk4j3m9nGvHMP2rdhwUzQcxvKao+1r\ncdmrs/2tkiYTdTuM4+1zwRhaTdTdTdFkoh7MhCuH8S6fgZ2vkxe4n6lWSulnME3UWnt4+qthtBrS\nahhN5sbuK3L9RO1nYrQewoeswL7l5rvsWT+/VD4fSzYfSVdM1AcTNXEmCNj5VovuhmNY/vfO/OUs\nuzqcqHE70HO1mqg9zITPCNhXNJ8K209Z2o6HPvU2UdsO5s5CN76AnWy16Hmc425+DU+rMWwf987D\n+GkwTV6dPifq1GhrdDJDmgzjlZnQw8TY5agHgEiuQgQgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBI\nAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQScAAiCRgAEQSMAAiCRgA\nkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZAJAEDIJKAARBJwACIJGAARBIw\nACIJGACRBAyASAIGQCQBAyCSgAEQScAAiCRgAC2VMpXSehCZ/rQewHcV8wLoXZ263FnVWlsP4YnB\nAza1fg1KKT1MAsMwjG7HYBjTNJXy356qn63RegjPOYUIQCQBAyCSgAEM6A7XhnRxsvV7OjmbDHCg\nlOn0HdWjXh/fbcTOM/Uijt2Nu/zUsf9ND/A9tY6/AssL2PG1MboFcBN5n4HVWg8qVUqJuPoTgF/K\nW4Ede7RteYJx1TNLNIBdcUf/QwVsN06KBfCK5d4yImZ5pxB/ErG5ATjLCCuwxwnDWuvcMKsugOGl\nBmyZqPlr3QK4j3FOIQJwKwIGQCQBAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoAB\nEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZAJAEDIJKAARBpwICVUloPAYCv+9N6AGeSLoD7GCpgtdZp\nk7HV/z5uA8BK3BpgqIDtUiyAVyz3lhExG/AzMADuQMAAiCRg5CllSji9AXzXgAHzoRfAHQwYMIbn\nEAWYBAyAUAIGQCQBAyCSgAEQScAAiCRgtOf3uoAPCBgAkQSM9vxeF/ABAQMgkoABEEnAAIgkYABE\nEjAAIgkYAJEEjIv4bWXgXAIGQCQB4yJ+Wxk4l4ABEOlP6wG8rfz9IKVuDunL4jOW7XcBGElYwEop\nc5mWX890C+AmRjuFWEop/V3r5gI8gNOFrcCeeqzAVgu17Q2uJF1AhA6P/o8NFbDdODU/qVirhgEB\nlnvLiJiNcwqx883dOqPcS9/vBjhH2Aqs1rq9CvFxwnD3W3BDj/dBKQ6bGFxYwKa9OM1/olsw/T1r\n7d3A8MY5hQjArQgYAJEE7D1+owugEwIGQCQBe5vPxgF6IGAARMq7jL4380diVmYAV7ICAyCSgP3W\nY+HVz/LLdZLATQjYUKQLuA8BG0o/C0Fg63GCxIHmWQRsQDIGPfMOPYuAAVxEus4lYFc7OIFw5bkF\n5zEghXfrTwQMoDvbaH3WsLE/dROwqx2cQ7jy9IJTGZCit9/V6YeAAfTolGKNHT8B4xzzaYpRT1bQ\nm4HPjPEiAXvV2O+Ws56dhgGXuV3APvtIc+w98inPbnmOYtTzFbzoxKO947sy0z4zzOH4UH8bffn7\nmtSz53WtUyn/++/0d4//mwd5/e+w//1jPTU/u7PuauWCp3Cu5VN4a9j+aYKn3poMH+9kP5tyH/+U\nl7uVcVZgpZT6V/l54q8+0vzNkci7P9vqqOfdRecwR2cnrrY/2yahW/L3Yz5rh/7TBQjnbtXfXGj+\njXPmF8yZYYo71Ars2Orwavs7Fqtl1sPu6mR5m+UXqzs/WNmsHmt5m+3h/+6y5qcnsnpS2//dPtbu\nXR0fjW431+qWyxv8tCzb/sjW4zY/LYmWj3Iw7O04Dx5r/no7GX668+nfV+r4eS0HsBrM9ls/3cn2\nBi/afaztn8zP4unj7j675bd23xQHP7X7oNvNtTu3D+5z2pstWwdvq4M73L39Tz81bV6F7Z//dBR1\nMJKDhzuwfe8s1P4Xl+X0s22tPFZgq68PlmIAHOq9D8OvwOrUaL3c/8HLU0+XWb+529//1MEi5vhR\nfv/SHA/m9LvdvefdxdPuUv71A/bV8fjTDfjiljzxvfDuXT29/WXv049nyO6seLqwO/ipF++wlNL/\nHmz4gDVbYvb/2j/V1VPYDmZ7xu8VX1qTf2lbvXiya3cH9MrPMpKnJ0jHM07Altdu9L7uDdT/Fu1/\nhK/4YIXRudWHard1zdO/20YeJ2CTbvGC19dqF/tsVA0PuvvcjIOxkY8NFTB4RZ/HOR+nqM+nM+t8\neEQTMHjP986JZe3rs0Yb6pe/LDG8cX6RGa4x6o77+PeZoENWYAzo23vhURvGGE65WD+CgAHTpMoE\nEjAGZF8Mr8t9v/gMDEbgsysOjPoBp4ABYxpvf/2x3DXWMacQuZS/l+FLbE+ODTlDrMAAiCRgXGr+\nd7PgXe+eEqzVTBucgHE1+xTgFD4DAzI49GHFCgyASAIGcBeD/WqBgAEQScCAfYMdrTMNd2WmgAEQ\nScAAiJR3GX35e16jblbCZXHKY/td4HXeQPQvLGCllLlMy69nugVwE6OdQiylFB89A9xA2ArsqccK\nbLVQ294AgJW4o/+uA/Zue3ZvoFgAr1juLSNi1nXA3mrP7kdiAOPxL+o9dB2wrVrr9irER7p2vwXA\nqMICNu3Faf4T3YJrWAE0ZMvPRrsKEYCbyFuBAW1ZAdAJKzAAIgkYRx6XxSRcTwvcjoBx5HGyyCkj\noEMCxhPqBfRJwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQM\ngEgC9l2lj39KyzCWDKOrMUyG8a9OhhFhwIB5+QHu4E/rAZxJugDuY6gVWK21+ucXAe6hjLfHL+W/\nJ2VNBvzssaOwl9jXfx1STyGuyvTThu7/BQBaeexF7CVypQbMnAN+yV4k3VCfgQFwHwN+BgbAHViB\nARAp9TOwns0XmGwvhrx+vTtfk9lqDKvHbbsp2g7j4LW4cjy9DWP7v5cN4+BB226Ki4fxdCY0fNse\nE7CvmF/1WuvqTXLxHnP7uFeOYftWbLgpGg5jedHs9rW47NXZ/lZJk4m6Hcbx9rlgDK0m6u6maDJR\nD2bClcN4l1OI5+vkBe5nqpVS+hlME538iv1qGK2GtBpGk7mx+4pcP1H7mRith/AhK7Bvufkue9bP\nL5XPx5LNR9IVE/XBRE2cCQJ2vtWiu+EYlv+9M385y64OJ2rcDvRcrSZqDzPhMwL2Fc2nwvZTlrbj\noU+9TdS2g7mz0I0vYCdbLXoe57ibX8PTagzbx73zMH4aTJNXp8+JOjXaGp3MkCbDeGUm9DAxdjnq\nASCSqxABiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZA\nJAEDIJKAARBJwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQM\ngEgCBkAkAQMgkoABEEnAAIgkYABE+n+DfO8fP5G2HQAAAABJRU5ErkJggg==\n", "output_type": "display_data"}], "prompt_number": 17, "cell_type": "code", "language": "python", "metadata": {}, "input": ["clf;\n", "subplot(2,1,1);\n", "x1 = x0; x1(x1==0) = NaN; stem(x1, 'b.'); axis tight;\n", "subplot(2,1,2);\n", "x1 = x; x1(x1==0) = NaN; stem(x1, 'b.'); axis tight;"]}, {"source": ["Display the evolution of each entry of $x^{(\\ell)}$ (different colors) with\n", "$\\norm{y-\\Phi x^{(\\ell)}}$. We display the evolution of the correct components (those indexed\n", "by the support $I$ of $x_0$) in color, and the wrong components (indexed by\n", "$J=I^c$) in black."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [{"metadata": {}, "png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGwCAIAAADOgk3lAAAACXBIWXMAAAsSAAALEgHS3X78AAAA\nIXRFWHRTb2Z0d2FyZQBBcnRpZmV4IEdob3N0c2NyaXB0IDguNTRTRzzSAAATcUlEQVR4nO3dW5bj\nNrIFUPEuD8wemd0ja8+M90NltpwpKUmKBHECe3/VqspS8omDAEBqmuf5BgBp/u/qDQCAPQQYAJEE\nGACRBBgAkQQYAJEEGACRBBgAkQQYAJEEGACRBBgAkQQYAJEEGACRBBgAkQQYAJEEGACRBBgAkQQY\nAJEEGACRBBgAkQQYAJEEGACRBBgAkQQYAJEEGACRBBgAkQQYAJEEGACRBBgAkQQYAJEEGACRBBgA\nkX67egPONU3TYZ/15+12u93+c9jnAfRsnuerN+EHxQPsdtA5mP6aDvy0ndswTf1fTyuV2ZcyO3Kz\nL126cEeO7P2fxhAiAJEEGACRBBgAkQRYjBpj+ndl9qXMjtzsS5fK7MhJikx1vrIsvjjE/FflYwWw\niFgIowIDIFJAxn5imqbbnyongG0iKrCATfxExDkA6E1E42kIEYBIAgyASAIMgEgCDIBIAgyASAIM\ngEgCDIBIAgyASMW/0HK+3W4JX8sGVNP9U8AFqMAATqDrfL7iFdjk+wiAxkRXKyowgEPpNLciwADO\noRQ7mQADIJIAAziaUcQmUhdxPP2umumhYLd2A7jeNAmz8+QF2PR2WFluAQwibwhxnuc3KTVN0/uE\nA2hBZ/p8eRXYe/dsexxg/JJnSjSgqZxRxLjef6kAexpOEgtgjcfWMiLM8oYQX4k43AAcpUKA3aNr\nnufpH6ou4HpLQ6R7fY7UIcTHiFr+LLcAxlGhAgNgQAIM4DRGEc8kwACIJMAAiCTAAM5kcdlpBBhA\nE6bBjibAAIgkwABOZhTxHAIMoBWjiIcSYABEEmAA5zOKeAIBBtCQUcTjCDAAIgkwgCaMIh5NgAG0\nZRTxIAIMgEgCDKAVo4iHEmAAzRlFPIIAAyCSAANoyCjicQQYwBWMIn5MgAEQSYABtGUU8SACDOAi\nRhE/I8AAiCTAAIgkwACaW6bBjCJ+QIABEEmAARBJgAFcwSjixwQYAJEEGACRBBjARbyS4zMCDOBq\npsF2EWAARBJgANcxivgBAQbQAaOI2wkwACIJMIBLGUXcS4AB9MEo4kYCDIBIAgzgakYRdxFgAN0w\niriFAAMgkgAD6IBRxO0EGEBPjCKuJsAAiCTAAPpgFHEjAQbQGaOI6wgwACIJMIBuGEXcQoAB9Mco\n4goCDIBIAgyASAIMoCfLNJhRxJ8IMAAiCTAAIhUMsEndDUQzirjOb1dvwJFEF8A4SlVg8zzPHgME\nGEOpCuypL2WZhAMCzPOv8cNpavZ6jrhBrPoBJrEA1nhsLSPCrNQQIgDjEGAAXTJ69JOCAWbMECgl\nYTTvEgUDDIARCDCAXhlPekuAAXTPKOIzAgzYLGKNNeXVfw7MnQZnmKbJgqkWliea+UYFBpBAjH1T\nvwLTSYTDGdigByowYCcx1ohe+AsCDCCEHsO/CTAAIgkwYDNTy6054M8IMGA/02CtOeAPBBgAkQQY\nQAKjiN8IMGAP02CXMYr4DwEGfMQ0GFcRYAAhVL3/JsCAnYwiXkbVe7vdBBgAoQQY8CnTYFxCgAHk\nWIZtdRoEGPAJ02BcSIABBzCKSHv1v9DSfcXgFEnVzPOv8cNpGnxhvQoMitOHo6r6FZjuJ8NqE13z\nPMtILqECg7Ia997EWDv65bfbTYDBCERLWWOfWQEGQCQBBpW1GUU003wBx1yAwSDajCIaq7zAwMdc\ngAEQSYABZBp+FFGAQXHLBNWp43umwa406iiiAAOOZBqMZgQYQKyxC18BBvUZ36tvyMJXgMFATINR\niQADDmYarKmB+w0CDIagPKpvvH6DAIOxGEWkDAEGEG7UfoMAA45nGuwagx12AQb9OjYGjO9RzG9X\nb8Dp9ASJNk3T4cFzxmcu5nl2011gnkcrv24qMIBSRoqx+hWYYRNCnVHHtCyPTq3z4KYCg26d2vob\n5aMAAQa9iwsbhdc1lsOedsHsJsCAs8RFL1kEGPTrjFJGeUQZAgwCnFHKKI8KGmwUUYABx1Pn0YAA\ng+G0TBd1HucRYNA1i+nZZqTaV4BBBmHDNgNcMHlv4lhu4+8908c73BA8XMtLETlbWAV2fznN3dN7\nY/nX9tsGJ0lfTC/GWhumAQwLsB9N0+RugfXcL5VVP7l5Q4jv3fuVj28RdX8OrlI57vW4nCqutSxV\ngbm3Yb0G94tb8jK7jvz84PAtOkOdCuxV5zTlTHC4e3eyRtVy6oKIBoeoxlmINE2Fp8TCAuzxNn4c\nJPyyrMOtQlWSABZhAXZ7Fk7L37ixYRMr3Sub5/KLOErNgcEjHZr1zosxZ+F6dWNMgFFfjSIjPQlq\nnAW6kjeEuJXbpqT01vwTpsFYq/ooogoMhmYx/RCKxlj9CszNU8ymxfHFFimkL6aHY6nAGEKlGLvF\n7k7oZmcr3SkRYIRRJRzOIR1Cxd6DAAN+sZieLAKMVCtb22JNZ/ruGEW8QPg184YAg0iSgM3KXTMC\njDz7qhAt/hvphR1jEmAQ5tSwMQ1GEAFGR7Z+ofbgRVXo7odudral91Dr4Nd/kNndQrHHmU/iKBFH\nBUZ3Tn3ZxEmfXImjRIr6FZiR9yDrm87By4XQ3Q/d7CKWF/sW+o5mFRgd2dHbGLxBDN390M2mNwKM\nISjE13CUyCLA6JFpsGs5SjWV66AIMFINXi6E7n7oZldTpYMiwOiLabCtjt39ZgEz+FnjEPVXIbpP\nzpDYj7YEbhPfb1nTshaxBBUYnRI2PxIw7Ffi/qpfgbnJD9RbqOwrqurVFll7pBTmKCowNmjTSpoG\nu5ZpsOJy+jo/EmAMJKhM6YGAKS7//Aow9mjTtGlAfySSGZkAG1qBhPDllndZi+nl7sWqHP/6izjq\nNVXH6nP+f8c8f587ku7so+qsXSz8xb4qMLqm//EjAcCw6ldgbu9X9mVD+hro9O1/49hq5uwDVfhE\nZCjxRHP9AHOTvGrUyrQgngY7lQNVXPIooiFEdjo7/DwNNghnjd3qV2Ajdx5rNA0qgB+VKaZpKn8U\nUQVGkZxbr3AcWkzPHrEtgABjsw7bHU+DnersA+VEsI8Aqy+3dTANtl6HvQoChF82AowAw8bSDlmj\niPQi8xYTYEP7sHnKzRXt8lYnnWsngk/UX4WY28iu9KYJGGpx2lA7W4+1puygAqNrZzdqlTJPALDH\nctkE3gsCrGuPTdL8j+XP6z/nx2Z6azvevq1cv4WVMmkf02AMQoB17bElmv7x/e8v2DLGYxqM3tQP\nsPmbq7eolK7ic9PJdSX0pqtraSyxo4j1A+y7/vPs6Vbt2+A+d3AT02DrnXSsClxFlDRigHXoTTi9\nX2R4yK/+/EPaMA223klHwIGlK/WX0f/owhb8VXOwppl4XHZsCTLRPAJxvcwX+9avwKZnrt6oXzZV\nVx9u9uF73WdkmgaL3ql+7s1xRZ2CQSuwDr/LavcvOmpfei7jzu6h97zvu2V9QTPsUL8Ci/CmaVjZ\nalzYuLT81abBLufAlhXYhxNgfVnTZR5wISLluUo7ktNHEWBHsog5TsljG30dqvBYT4B1auVCxNuW\nZuXVZ+5umFq2/ut/V8lM2sdierZJu3cE2MUOfMxrkGblpGmwQY5ez/Q8OhJyOwiwU+j5ci1hwAgE\n2E6HNxBLOB3+BNjTz3y/DVtJ1m6d9Gb6s59qOO/D+UFU10eAHanZjfclkAbpbp80DTbI0YNNIu6K\nUgF2yYs2zmv+duzIIF1X02BV6UywSZ0Au7934O7CR2uPvQOXT3v1sVv3dJAXSnELX0zPlXJOcZ03\n6Dx9ua0+NcA+929/v3or3qlTgQFwpL7T6zbsy3xXerUg8Pso5f0n//7711/+/vv89/T37Xb74/bH\n8q9vysEvP/CmLF7zMyt/8v0PTH9Nt9tt/mte/19W/uSP/9T/qMD64/DJ/0o5Gjxa7vL08xYxfKUC\nO97vv5/+rsID1+Nd9dvPeFoAGEqdCuyxKjqp07r+VUz38uvVD79vqXvocd/Lr69/KWD+bd+ZOvuZ\nPC5UpvxKUSfAbg2b/qeNyzJ+uLiPH17ix+bvqgdRtcvAUUoF2Bub3nh7L5KOjUMN90o9FKDv7TuV\nLoDylF/t1Q+w5ysU3rYm939t/D7D9a+fb+bLCo5ff+nubMjRhjfqB1hXPd/GD1kfuIJjx2ZbwfGd\nQKpK+XUJqxAPtiygX7xvsy5p0X4oQJ+t4Gjwe7OceuIqHSg4T/0K7I3//veYz/m+fOOL7+3Rq1Ls\nkGbx2uavduNrAozvlF9XUYEdbCm/zluCeOEw1Nm/2ggbsN7QFdgfp0TM1w997H0vf+6/S/50Bccx\nnyyl1nGgIii/LqQCa2fNC6La/LqVP7A4dgVHAac+wgysVL8CO6PVWN+gH/jbjwq5Dldw9F+PfuLy\nr7zhPMqva9UPsPPsDqdTF9NbwXGe2nsHceoHWM9RkTWmlLW1cDbl1+XqB1iHXi3laJkQp06A7fu9\nKQFpAgw6YRFHTR82lx8uQRz2RRvf7d5xgdc55VcPVGAdWdPYfdiudRgkHW7SUynbCeNQgZXyYSP7\nZgnimuDUxH/OMYyg/OqECmy/EQZ5mk2A1TbsjsOpBgqw943Ipi8AO/zbwhrrcOM73KSnNm2ncqok\n5Vc/6gfYyhZnU8PUeWtrBcfhWn8TW98XGHTCHNhAOlwRVyDbGIryqysCrI7zVnA0+O2FXfXIHZQn\nwHjJCo4vTIANTvnVGwE2nMZBoh0HTlJ/EcdwDeifPwwGvv/Xw1dwnPe/mhnuEuIZ5VeH6gfYcP5z\n8Od1ni6d23H0HHBYqX6AmUL/4lWN9aUy23Q0fvzhp2dhkANOAcqvPmU/kPujrY8n3/R//7Gky6aX\nSL16COx9gHV+zHdcGDt2LeVojGnAAIt4XUP9Csz34V6owMEssAt8aMD0SmEVIu806IL138vbrfCu\nQQ/qV2BbaXR2c+g+5xj2RvnVs+ECTANxuBFWcHiEGTpUP8Ak1lUKtOPNdqHAsSpJ+dU5c2AARKpf\ngenbcgmlfzrlV/9UYBzjVXv9vh2PaOWbTYBFHA3oR/0KTKMAbKX8ilA/wHpgGPOLiANiBQd0zhAi\nHE/dH035lUIF1o5GrTy1FLSkAjudRu2ViESP2EgOpPwKIsDgifbdDkkJWwkwWitfkoqiXMqvLPXn\nwPppLvvZEs6w7/y6KmA3FRjXiChTIjaSoyi/4tSvwC5vg3zRcxxVEURQgcGR9vVU9G8up/xKJMDO\npS8/CCca2hNg8FyDqkjsdUL5Far+HJg2AqAkFRj8i+9DGY3yK1f9CkybAlCSCgw+9eEwtT7WhZRf\n0QQYPCFUoH8CDP6n5QSY5UWXU36lE2AARBJg8BGFVCjlVwECDL5qOQFmsg12y1tGv3R4v9/5j31h\n7QJbqaXGofyqISzApmlakunxzwu5xSWs4ID2qg0hTtOkaaAZF1si5VcZYRXYj+4d4S+F2vcfgFdc\nIQwrrkPWdYBtzZ6nP6A9AhbKrzceW8uIMOs6wDZlz9MpMVjpknf4umLhE10H2HfzPH9fhXiPrqf/\nBOeJ6KLySPlVTFiA3Z6F0/I3cosIkg8OkRdgW2ksWE8fqDDlVz3VltFDYzIPrlK/AtO+cIbPK3tX\nZkvKr5JUYABEEmDQlEnZ9pRfVdUfQtRecB7DgHAhFRhQmfKrsPoVmD4yHXJZwudUYEBZyq/aBBi0\nY0YWDiTAgOKUX1UJMKAm5W55AgyoTPlVWP1ViGYd6I0liA2470egAgPK0lWorX4FprcLo1F+DUIF\nBkAkAQaU4uHlcQgwACIJMKAO5ddQBBgAkQQYUITyazQCDIBIAgyoQPk1IAEGQCQBBsRTfo1JgAEQ\nSYAB2ZRfwxJgAEQSYEAw5dfIBBgAkQQYkEr5NTgBBkAkAQZEUn4hwACIJMCAPMovbgIMgFACDAij\n/OJOgAEQSYABSZRfLAQYAJEEGBBD+cUjAQZAJAEGZFB+8YUAAyCSAAMCKL/4ToABEEmAAb1TfvGU\nAAMgkgADuqb84hUBBkAkAQb0S/nFGwIMgEgCDOiU8ov3BBgAkQQY0CPlFz8SYABEKhhg09Jzq6XS\nfpXZlzI7cutsXz4sv7ral0+U2ZGT/Hb1BhzJyQYYR6kKbJ7n2Xg5hDP7xUqlKrCnvpRlEg7gqbhB\nrCm0QX8TS9P0v51KOx3Awt17sf7TIbUC6//IAp9wj/Oj1ABbyS0AUFXqECIAgyu1ChGAcQgwACIV\nnwOLtqy0/DLM+7gCM3QE+HGlaKLv259+Ul5dbBHebHz0eYk+KW0IsE79+2GAr81l7gUd96DJF2+2\nP/ek3N23P7Rv8WbjE3dnEX1SGjCEGGmaptAkSH9bypvtzz0ptxKt/Cu55yX6pLShAoukX9ahAiel\n5MZHn5fQ6G1GgOVJvA/LSz8p94YydC/ebHzoHi2W9L16QzplCDGMS7lDNU5KdFv/dOOjz0v0xjcT\nWVYP4ssapGUMpMDapNDxnMXjuShwUqJfeP1044udl8SNbyO7HQFgWIYQAYgkwACIJMAAiCTAAIgk\nwACIJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTA\nAIgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAA\niCTAAIj0/wFs0jtlpL62AAAAAElFTkSuQmCC\n", "output_type": "display_data"}], "prompt_number": 18, "cell_type": "code", "language": "python", "metadata": {}, "input": ["clf;\n", "I = find(x0~=0); J = setdiff(1:N,I);\n", "clf; hold on;\n", "h = plot(E, X(I,:)', '-'); set(h, 'LineWidth', 2);\n", "h = plot(E, X(J,:)', 'k-'); set(h, 'LineWidth', 2);\n", "axis tight; box on;"]}, {"source": ["Orthogonal Matching Pursuit\n", "---------------------------\n", "Orthogonal matching pursuit improves over matching pursuit by\n", "reducing the error using an orthogonal projection.\n", "\n", "\n", "This algorithm is introduced in:\n", "\n", "\n", "Y. Pati, R. Rezaiifar, P. Krishnaprasad,\n", "_Orthogonal Matching Pursuit : recursive function approximation with application to wavelet decomposition_,\n", "Asilomar Conf. on Signals, Systems and Comput., 1993\n", "\n", "\n", "At an iteration $\\ell$ of the algorithm, one computes a standard MP step\n", "$$ \\tilde x = x^{(\\ell)} + \\mu $$\n", "and the next iteration is computed by projecting $\\tilde x$ on the\n", "known support\n", "$$ x^{(\\ell+1)} = \\uargmin{ I(x)=I(\\tilde x) } \\norm{y-\\Phi x}  $$\n", "where $I(x)$ is the support of the solution\n", "$$ I(x) = \\enscond{ i \\in \\{1,\\ldots,N\\} }{ x_i \\neq 0 }. $$\n", "\n", "\n", "This $\\ell^2$ minimization can be solved in closed form as\n", "$$ x^{(\\ell+1)}_I = \\Phi_I^+ y\n", "   \\qwhereq  A^+ = (A^*A)^{-1} A^* \\qandq I = I(\\tilde x), $$\n", "where we used the shorthand notation $x_I = (x_i)_{i \\in I} \\in\n", "\\RR^{\\abs{I}}$ and $\\Phi_I \\in \\RR^{P \\times \\abs{I}}$\n", "is the sub-matrix obtained by extracting from $\\Phi$ the columns\n", "indexed by $I$.\n", "\n", "\n", "Perform one step of matching pursuit."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 19, "cell_type": "code", "language": "python", "metadata": {}, "input": ["c = Phi'*(y-Phi*x);\n", "[~,i] = max(abs(c));\n", "x(i) = x(i) + c(i);"]}, {"source": ["Perform $\\ell^2$ projection."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 20, "cell_type": "code", "language": "python", "metadata": {}, "input": ["I = find(x~=0);\n", "x(I) = pinv(Phi(:,I))*y;"]}, {"source": ["Note that the algorithm stops after $P$ steps, providing a zero error,\n", "$\\Phi x^{(P)} = y$."], "metadata": {}, "cell_type": "markdown"}, {"source": ["__Exercise 2__\n", "\n", "Perform the Orthogonal Matching Pursuit algorithm,\n", "and compare it with the Matching Pursuit.\n", "-- MP --- %\n", "-- OMP --- %\n", "isplay"], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [{"metadata": {}, "png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGwCAIAAADOgk3lAAAACXBIWXMAAAsSAAALEgHS3X78AAAA\nIXRFWHRTb2Z0d2FyZQBBcnRpZmV4IEdob3N0c2NyaXB0IDguNTRTRzzSAAAf6klEQVR4nO3dWZqd\nOLOFYainBvbXyMpnZjUzzgU2xjRCvWKFvvfKztzJptUihIB127YFAAA1f42eAQAAchBgAABJBBgA\nQBIBBgCQRIABACQRYAAASQQYAEASAQYAkESAAQAkEWAAAEkEGABAEgEGAJBEgAEAJBFgAABJBBgA\nQBIBBgCQRIABACT9PXoGgFms67osy/EO9OO/+z92vCEdiEeAAeOdU40MAyLRhQgAkESAAf0cHYaX\nSmv9hfILiEcXIjAeuQVkoAIDAEgiwICu9l5ESi6gHAcSAEASFRgAQBIBBgCQRIABACQRYAAASQQY\nAEASNzLHWn+s3x8CAC+2H9bHqDuvwNZ6Ri8KAOAPzu8Dq3XH6FF+DTwl8XT3q5tlcbMgC8ti0sAF\nkViH3rsQ/6XrDwB8ct6FWJf9HmEAmIdAkVhCogqOxLIY5GZBFpbFJLoQwwRmsYTENgAAayQaT7oQ\nAQCSnA/i2JZlYQQ8gP7Mly8OUIEBQAOcOrenWoE99s+ebzfef7vysnYA/e1t0bpShzWlF2Dhh2IQ\nVwAMIcNa0utC3LYtkFI89gnAeIRWF3oVWNjPnsNTB+MlzyjRAPSwbXIdiXJn/64C7DGcSCwAY6hl\n2Lm1lAgzvS7ENxKrG8CkaKAa8BBge3Rt23a8+oSqC4AJ57aIDKtNtQvxHFHHv8ktAOYcHYmozUMF\nBgCmHefWJFlVBBgAdESG1UOAAUB7XAxrgAADgC64SF8bAQYAvXAxrCoCDAA6IsPqIcAAYBAyrAwB\nBgB9cTGsEgIMALqjI7EGAgwARiDDihFgADAaGZaFAAOAQbi7uQwBBgDjMKCjAAEGAENxMSwXAQYA\nZpBhKQgwABiNi2FZCDAAMICLYekIMACwgYthiQgwALCHDIvgMMBWNjwAUVwMS/H36BmoiegCIG/b\niK5Iriqwbds2LoQCUMfFsDiuKrBHl7KMhAMg4KjD1rXbAEW5Tiz/AUZiAdDWK8POraVEmLnqQgQA\nPzj5/uK/ApM4jwAU0b3R3IiORCEOKzAOKgB+MKDjnf8KjDwDqts7NtZ15fjqijrsTw4rMADd0EXf\nA3c3vyDAACSj8OqNFf6EAAOQ48gwirBOuBh2Q4ABgBoybFkWAgxANoqw3rgY9if/oxA5tDC5pter\ntm1jRGJXPOr3hAoMcI5zOG+4GPaL85MmzgoxuSO9mh4Ifb4FfzjSq806l2g8qcCAKfSpw6j2Bph4\nnRNggGd9TqLtn6o7xDonwAD3+owVZETiANNfDCPAgImQLt7MnWEEGOBf545EYnKM+VY7AQZMoXO6\nkGH9THx3MzcyA9ZVr5/ajZA+7mtGV7Pe3UwFBsyCjkTPprwY5r8CY4AvdFV/SlPnJz9J3Azr0DTv\nvaQCAwRUrGY6VEiE1hjzXQzTC7D1l8Cv6LuAD62ToEOGcTB2Ndmpg1iA7T0Su8cD4/ht/3kDWmiR\nBBwgns10MczbNbD9ID8fn5wATs5TY611MYw3rQxzDEr0fjFMrAL7FCjOMCcHO4NuR2LPr8CzlDUv\ndxXGVQX2eJxz6jctlYMwRotqpsOpHmeTw2TdGSbXd+UnwN4ObInNgFQxjbjXLiw6EhFlgo5EsS7E\n/WDYXS5uP/4KXk14XuKgIxG9eR/Q4bytJ8z8SX35r7OXBbdYnA6ryNlWUHKOrpSVL9F4ilVgADcY\n7URvbZ58qw3g+u5mAgy2VG/gnDWd7a5U7f/wsZbwB/OFVDaBIrEER6Ooz90yqUvKX/8VHYlIdjSG\n0X3v9jcTFRgkJVUMjssLrY7Eg78NocTRyhfI2BISJxE4iz9JZzTH/o+6i9N6LTnbCmJSBnRINJ5+\n7gN7w7ledUZ268lvMGq9+I0mO/lWG8zdey/pQkSypucE7Tqy6EiMQag45+vOMOcnQZzlVXd/XHKj\nr9gxmiOMjkTkiBjQIdF4CsxiCX/n2i0k7QN9mp6MK2GRs+Sv6Ww6IrHuZB+/ws2GUBJxMUwiwOhC\nhEXx3X32j7Fu6EhELC/blwDzb3uXN7X9H3aq28mH1Ive2uxvQ4hxcTHM/yjEmc8lI5sGm30FHYar\n2VzwDKIjErtNH8/0H1c/SwUm9Iq26nSH82V0JCYVYS4JdST63hBiNJtH/xXY+Xh2mWElrcA26/sG\n/d2N1GhTtl5R/jaEmPOdYYJ12CwVmGOBZiu+Lklt+7q1NY2KMJcarQFWrHNqoXVGgKkycrpqs1Gb\nczTHWaOFalTktZs4osgO6Jg0wApH4lUXmJnjh5fP3N9J/Tbl8GeMrIEAhtTHaz0isTUybDypTeC8\n33nU8RC/Vi+3c57/e5/58GfevvTzjtHsW0p73ova6CG/Lm+nVXw8h8sNIebPu5slrkpOWoGZlVRU\nLXFHO10006IjEQnMx9WdqwBbfxk1AxkHYSCBzp2Hj5O9PJYwe8Gzz7N6nqC1HlLvqd2kIxGZ1C6G\n+QmwveDddTgAHg/mku+9d6EcP4mfbOGCZ/955xZHqLAYRXFEov0Oqykc581jZyOOQC9npMughsvl\nIgBAkm2x3q/opwIDANRkO72WGZ7EUeJSnn4OrPjvv3VZlv/9b/tv/W9Zln+Wf86//RzIHjMQK/J1\nXIUjD9cfv377Y4v8k8iJ5/3KmrxZZXikexFvKZEh0X1FgIVEbsJzemVPKj69MFaH9IKcy6HJRu7D\nT4Cdx250biMu0bWXX2ejxjgMKb8OXsuvDKk7gO+14Y+nwkuLnwBbxh3t//xz/UH8317m2UGD5btM\nLIwWB9sXF6TXQAzi6KewbzC16RxbfiV9r0rBkZfNqUunsjawkF6juarAjHts/qrXK2MnOEnj63vp\nEIn0Go4As6V1/1tkwJzLr98/pPxalqXX2A3jawOklwUEmJhwu1YYBkf/Yd4EVaLIMt9XEH0guuwg\nwDIpttGjAibme42vT4bOY0d6meI/wOZpO2qVX0f/YYfyy3HBkbFoBJ5xpJc1/gPMcRP55mOR/33o\nKmzHQaPM0HkspJdJfh7m+2jC9KqusPyac/hGnz9BNxOml8QLLf1XYBjIQaPsYBFQgmdEWUaAXdFO\nHSJPwRjgcEH55caEhZcWAuyKXsezuq/olONgEZCN9LKPAEM1j2kXfnW15WCgspwZ6SXBf4BZaEeo\n6ibR52GJaI30UuE/wIajeVLE0PlpkV5C/AeYnerHzpwY4ayVp8tRHQMO5fgPMNhkub0mVCZE4aWI\nAOuEptA3yi9ppJcoXmgJ/IFQmQ3ppct/BTa8GZJ4IgtKMOJDF+klTS/AYp6t9/jbUezMCT6VPEQj\n+7swBEM2HBALsHM181jZkBbIVpgo7HtCKLx88HYNbF1XTmxRqMMzDLnSNhDp5YZYBfZpbw4uhdr9\nA8AFiTIJ0itA7uzfdIClZk/So9OBQpRfckivsPMOKdFUmg6w1FtqaA5QiF3IMdLLH9MBdrdt2/3s\ndY+ux18tNElog/JLCAMOvRILsOXpsD9+En5zB4AJUXg55m0UItBTdhVF+dUH6eWbXgWWipYCdtAf\n0BPp5Z7/AKPJACZEes3Af4AB1tAr0BRDNubhP8BoLIB5UHhNhUEcAJwgvWZDgAHwgPSaEAEGQB7p\nNScCDIA20mta/gdxAPCKAYeTI8AASKLwAl2IAPSQXlgIMABySC/sCDAASkgvHLgGBkADQzZwQYAB\nEEDhhTu6EAFYR3rhEQEGwDTSC28IMBmeXmzmZlncLMhidVny0svmsmRwsyCNOAwwNjngA7UXwlwN\n4iC6AB8YcIgYrgJsf3clMQZIo/BCpNXfC4vX9fdCkWWAMg7gkeyng2oFdimz7K9oAEm2bVkWjmuE\nqAZYZGKRawDglcNRiACAGTi8BgYAmAEVGABAkuo1sBkcA1UuVfJ5AItoAX0eKaroPv/qG+VtZ5MQ\nmHnp7SK9UfogwIz682aAa3Opu0Or36UXmH/djbI7bqNUXJDAzCsuzkF6o3RAF6KkdV1Fk2DbNulD\nMTD/uhtlcdHKv9HdLtIbpQ8qMEmclxnkYKO4nHnp7SIavd0QYHoUj0P31DfK3lCKLkVg5kWX6MDj\n8cLoQhTDrmyQj40i3dY/zrz0dpGe+W4ky+pJXMYgHX0gDsYmifbnHM7bwsFGkX4w2+PMO9suijPf\nh3Y7AgCYFl2IAABJBBgAQBIBBgCQRIABACQRYAAASQQYAEASAQYAkESAAQAkEWAAAEkEGABAEgEG\nAJBEgAEAJBFgAABJBBgAQBIBBgCQRIABACQRYAAASQQYAEASAQYAkESAAQAkEWAAAEkEGABAEgEG\nAJBEgAEAJBFgAABJBBgAQBIBBgCQRIABACQRYAAASX+PngHAuXVd939s23b54fGT47/Hhy+fB3BH\ngAENret6TqmYTEr9PDAtuhCBVi4JdCmwABQiwIAxjjy75Nz6C+UXEEYXImALuQVEogIDWrn0GVJU\nAXVRgQENnTPsnl77b0k1IA8HDwBAEl2IAABJBBgAQBIBBgCQRIABACQRYAAASQQYAEDSHPeBnR9A\nN+1tA8dKmHYNAPDF+31gL89OneqJqoENPNV6AJDEfjrMUYHd5GyWx225rt0KmpxHNnw9+/z8pPSM\nWcrm5vETbhZkYVlMGrggEm9OcB5g6/0kInurvP3h488tHDyXGTvP0n2eI1eLheUCgGVZ3AfYg70J\njmmv4z95NzbVAtF1/0nSAt4/TKQBGGS+ANvFN7uPBdxbd2JY4G9LBL435rs+PxNeLgbIABjESU/x\nm/Fd4W+t/9tcfZZu8QVTuwUPz4PrPQqYxPjGM8KsFVg3b9lzH/0RSAULoRX+lsui9ZwZALMiwDq6\nXFT77Pr7zC07CfE4w426TLMZH1VlZ0UBIgiw7rYt9qqVXIt2zPC9IBu1LMZD6yy1txmYHgE2wqVe\n8ddC3QuyngsrVLnG5CvBBrwgwMZx3wD17FcUCq2zpLE8kZ+xuaRAAwQYGgvE2FLW2oqGVoyKwaa7\nEoAvBBi6eLw8tsQVZEnXsXy31xnBZm0oDVAPAYa+wgVZyTRnFr6xYXF9wRUTI8AwwltBlvrneFNS\n8gIiXAXY8fhk+zeQ46e3dpYtWEu7a5DAaK4CbPkVXRIPQcEf2F5NUZDBI1cBRmgBHyjI4IirANtd\nyq/La9kIOYCCDI8kXmJ55irA9rV/iSgSC3hFQYaTwKm/Ta4CbCGugAwUZNDkJ8D28wUGIgL5KMgg\nxU+AkVhAHRRkEOEnwABURkEG2wgwAEHhgmwhyTAMAQYgztuLwkkyDEKAAUgReI7l+SeEGdojwABk\nOUcUZRlGIMCWZVnWH78G3//Yzv8N2D8JYFniyjKSDLU5f+jt/lipmEAqRJ4BV4FHOfhudto3OH3Y\nb9a8B9i4Pcn+tgc6cZFkbmIpnv1GbNIuxCobJrxDX35rf1cAWrHdwThhMrnhvQLr/mIwrp/Bq8iG\nPnb3bl+W9Ukmr4ezxFsVBWaxxPBtED6EvO76cKB16/9756/x1PP13/Jp/IFjc3jjGUNgFkuY2gaB\nFiF8tHy+18DOMkJUUlwdu2vdkNv+r+LEcnEo/WKq8XwjMIslzG6DtyP/nGSF7+OxueAwIi+xOnzX\n9astpFoG/aPPbON5JjCLjx5X7rnF339rfxs8H9v1Dtpj8S/3umEqPeMq2+NM1pyZ/m9otN34hNlv\nPBfFAAu88eu+xiW2wXI5dC/p9dS5n9zl+H6FYLY8Y8jZxWw7QAVJQajQ/jySaDwFZvFRuAL7XXa0\n3Ab3gu/zk5EfW5ZQ5Lz5o/vxx/pQxkVMM7U5i18JAxFaB+KqE/3HQhJgDb0F2KXn8F6LVFnewNWp\n8/TfPvY4D+EkyG+C33ojc8NMaERJ9qiZjOl3C4bUPYHEGkYwwy5Ht51j+Y2rALv/tvpJxH0DZw+1\nOJeJjz//npmvtuz3aLHPmXzMs4JLcQP3q9a5BcQSzLADFVhD4ctdtQIs3O7fZ+DzY8PLl8KRjX+4\nZ94t8LrtXYQWjJLNMAKsobe42n9Sfg0sKboe/zDmM/GTbSQtzyIvy3WMMXILAjQzjAAb75Jz+z8+\nF7nRlbMZ/E6UshjLvuZHbsEiwQwjwMbL6LKTGFZnXCDGfhr34J/IAThAZWoZRoCNV3LJx/ea6eA7\nxna5YZYwRMUkdrDpSGUYATbe5wM73vheLZ3FJtms2NkmopNhBNh4EttgTh3KppFD+YtvroBbIhkm\n0XhO+kJLDBd5E8L9kxIi5/m+1FyC9W/bfmfYulrOMPsIMJgwZ2MdvkeQMHOLDKuEAANMeHsyy/2H\nhJkHZFgNBBhgy+ejW+LvaIRp5wxDFgIMsIs+RueODKMIy0KAARoIM5/IsAIEGKCHMHOFDMtFgAFp\nrL0zKX70x2JgbvGMDMtCgA3Q+dFHSc/Fp4FbEjeQnSEVMS/uIc/sIsPSCdxrXcLUzeRaj+yzs94a\nabc5DK66yIU1OOfTObbU6G1hqvF8IzCLqR5fFTZwZgK/vbzGrJaS90Q/Tq3kz2ud8o+N//BsP86b\n5SOLPDPNRoYNbzxjCMxivMcXWob/JOPdYEIe+5RKFjntnV5lK99OR2skO+8pTUWe2WLjYYkE2BiP\nL7FUV1KrtSh6Mt43XZHxnVauILsgz8YzkGEE2Bh9Aiz+WbThKbQYRlG9AY3pCI3/3ux1Jafzhqj7\nXTFfKrpdNIzIMLkxPv4DrE8nYcwAsKX7DtGiFKh+eq5er3wa2A3QYjVyk1k/Q+swKrAxDHYhmh1I\n0qIa8LdH1dKzc7X14H47Nw84Ny7DCLAx8gIsYz2UjxDpqencSuzrdvTJ/qYZQx3Wz6AMkzioBWax\nhFbGdMNqmUHrjCHD+hl0Pcz+ZhWYxRIS22CsyAEakZ8PY1t01i3D2LLNdc8wicZTYBZLSGwDO/pf\nMmTrtMbFMD/6ZphE4ykwiyUktoFN5ZcPO9+1hjdcDPOjY4ZJNJ4Cs1hCYhvMpvp9CwjjYpgrvTJM\novEUmMUSEtsARu52eORj/+FimCtdMkyi8RSYxRIS2wBhluNt0WmyuRjmSvsMk2g8eR8YrOtzFBmP\nyXLHOwpaN0wSDZ+8bfudYRO/P8z5rsaxhEYUL/xwMcyblnWYROP51+gZACRFPv3SlNbzbL+98+a8\nwkV2wrroQvztsRO/yqNm459GKPc06JltVV8c2kfree7WUYmf5u5LdL6TybUvqXxvPgmKgxe4M8yb\nBn2JEqcgVGDaqj/YicckTqVFI3Uu8iQaQQ9mrcOc71778XM5iuq+QzL1WYJ1p2+H7x0pTLoIW7gz\nzI2qdZjEyYfALF4EDoz7MSmxDc4qxtWQKzRaa7sW0U4z7gxzqF6GSTSeArN4Fn7b8uNPkqbfem2M\nTZTqe6T7lynHE22suRjmUKUMI8Dq+wyw/R8ZL7R8k71+Bj7xL/KrwxMvnH+t/aoKxQzjzjCfamSY\nRIB5G8SR1HMYk3OP49r75OK6rjFfdJ9UZOfhY4/r51+F5+R89f5x3mBK6wEXincaeJA7pkNuYwlk\n7Fm4Art/rOSNHnW3Zfx6tr8PBYL/3mBp7WCFFIuwhYthXpXVYVRgXb2t7ow0KnkVVnbl0Si3WoTx\nZ8H6WYq9TcH+ARMmehsvj0n0aYKx9Xr70+OFrksz/VilfZ4GNo2Q8m+M7HUMfN5+bbfoZ5hotcHF\nMLdy6zCJsw2BWSxhs8mOuTBQPnik8MF3Pfs8384/RIkuCHeGuZWVYQTYeO0CrEXXXNJXB9S9mSz+\nw5/ZGT/uUXq31F0QLoa5deyTBJiQvDvDRDvflohYzRs6H7mTxNxjPluGyS0Fd4a5lZhhBNh4fW7d\n/RTOlcKbr+NnIN7j8MKkPwwsFBlmGRfDPEvJMIkA8zMKsY/Udvzxv+GJ5PXaVZna/a9q3SQUfz8Q\ndw6NxZ1hnh3jEr0MShTI2BL3Ud271qeW8dqt/1q3EqeeMge+N74oES1fLnSXgothnsXVYRIVmMAs\nlkgNlaRCIe8rYqZZRcU2IinDwt87VYZJd5dxMcyziAyTCLC/Rs+ALfFptP5y/9Xb/dQxTfbbNMcq\nHJH/KGlVV/nG/uwf/zFaX5TV3b7CjvUvvvIJsPru3ZXH4Xpv0d7KlFo3V93np3Bqhd+bMZ23SUlo\nUan30Xr9+0h3YS4yjABr6LHquv9wOzn/3FqTV6stjp+OjwzT1TpjdNPdCf0MI8AybTepf77/I1Cu\nVZnJx29pPcH47233RBI7pJvpbjOvuHJc0Vz/BNir1Hy6fzL7mJRo8krmLSmWJNZGJOlF4GKYQ+cj\nUXD9C4wzKREzkCb12Rwxo7MiB1nFPLoiPIUY1YeTlS9dxoypD0qUnv9udzcrrhwPnh6WyChEE9Yv\nMZ2BGcPxs//2PgVrqp81J01E9Dxduo7kYphzsnWY/wD7dD9mwkMqUq8DBT4ZOWCv26CJjGkGJvv5\nvXkdiYEJoh0uhjmnmWGuAiz7PqpACxvuPGydYZZVOStPWnz1DNPd1hdcDPNJMMP8BNi5M7DuAXA5\ntJImXnhYVi/C6oqvsWLmPzXDFElnWOuMUd+4Hpw2gcTGELhMF2n98+XL+78VmwkAsMB+OvipwAAA\nU5n3dSpvYzQ+P/kmfqhxeNBwxpDi+93QkX8Y/jpng5trjQVvPaa8hOV58+3SfjhY9xLdV5MGWNOB\n4FUax4zULPzqwjkxjvRCO/7SS4WfADuP3ci7eflx9EfSPbZJX11e8ZSnl7OUekN6oRGiayw/AbbU\ne2JF0gSHNBxVCq95Og8PLtPL7Iy5R3oN5yrAqtubhs8nIe0Kn4cUmRnt+gxT50RF9cUxtVpIr1Ge\nnr6E3giwb3buT6qYXs5S6k2tC9E2VxfpNQSFlx3+Ayzm1tFAPRQp++7m5akZ+pxarc5S352Ho0bW\n9EF6DUHhZYr/AEt9BkT1B3m8fd3jF0V+dfjZiRl/6w/phboovAzyH2Cp8h48X/EBjNlSOxgdl18t\nhh3aQXr1R+FlEwFWqk8bF77l+dFjnrmJqEj+hh3anCvHKLwsmzTAhj9T9e1O6uw3jJw9LtTnkrrJ\ntioLYjMnbM6VYxRexs0YYOcjP6OyKfzG1N8Wft3bsoTvENBFeqEKCi8J/gMs41Cv8iT7Wk1M3u3V\nlw+HYyx74taYeipYXTbnyivSS4X/AGvZSbgtDXbu+L6+hxmKfkJxeOKKTWSVJt5mTticK5eILi3+\nA6yldbnt8WP5vm83wGt6GZwlx0gvOe4DbN8H16cfnv7/6wenPXi9fOz0GUuRlSvy4p9ckpFeyEB0\nifLzRuZH4XHmv/4RHmER8y2fkfY2lcFZGNj6bwv1dvfYwB2pPGUNRoXBWfKK9Hq0rgLpIDCLFzGP\nQTp+m1gtnUcnPkxz6fKcjiqq3CeQsQ77706kF7IRXQESASbWhXhep4/rt2yNnzoQXxrtDul1G4iR\n+uevU2v6gMfjT3ru9HU3h5HDlfTqg/RyQCzAPkV2Z5WP1qsrMMNP1+dCjo/dJ5n0opZs3XoU6w7c\nMBIVpFcHRJcb3gLs6PWq2FhXuS3sbbIpn3/++XuxmHZkZjSXgXXSuiEmvZCH9AoYfiqfynSAlT+d\n9v5DC0/dTZpsREP2OJ09dNsen5FX2t62Y+oaKHkGf2CujEQF6dUBj4YKK7niMITpAEs6jN+qrqab\nIbWh6Rif668/b3K39UXSC9UqXoqrMnDDAtKrNQovlwTGmVzcT5yP6Hr81ecEza4B9REK1kLica+w\nsPWtzY8/FF4ZJEYhCsxiCYlt8CmvgYu/kaupKjHW6O0BFvYN0qspCq9sEo2nwCyWkNgGrb03+p3W\nzee4xNSLnZ/fFWP4jkF0tUbhVUKi8RSYxRIS26CPdV3eH/zxcxWxqrohvZqi8Con0XiaHsSBirZt\nOYYm3pJs/+92r17M78CSSK+mSK95EGDTOSXZY4z9/NTPH5U9BwQ7I9cj3SO6ZkOAzSs4MiLqScd3\ntCA76eGvotj3JkSAITLJ7h/+nvI8bYrQ4BF/5tnNcEGA4beY+/BjHocReLrV6U/S52+o7BH8JFZT\npNfMCDA8S7r1KnyN5/GX/cOsz43VxFU3RBcIMHx4a5Ezgm3bHgY6LsHeyNRWqf/jP0isIUgvLAQY\nsmUE22X4eMZVtP1n33+WggTSQnThQIChsshHKj+9Pvv6kSqz8/kJ7n4TQnrhjABDD+HhIWX9ftc2\nLLJRC3wnkWYQ0YU7Agy95T2Zt3pH32V64Xmh9RyL9Y9HDgNM4hFeGTwt174s97Js4ALevzmpRPPE\n2l6WMWDVzcHiZkEacRVg1l5AhXgGj9KkSPMkbzFbbEAKL4S5CrDLay2Bup4iTekEuemRETnx+LUl\nfds7+nAVYI8qvmsKkBbe95PCODsLM/6QQ7YbubN/pfPHs0AsnY9Due0BOFLetnD8jmQ/HVQrsMg1\na38DADP7eggLxy9CVAMMgAOcYaKEahciAGByf42eAQAAchBgAABJXAOz6xhCeXtd5PUxuHK0bp+6\nu8+/+kZ529kkBGZeertIb5Q+CDCjLjcDXPZg3R1a/caGwPzrbpTd8RwAxQUJzLzi4hykN0oHdCFK\nWtdVNAkuj0CUE5h/3Y2yuGjl3+huF+mN0gcVmCTOywxysFFczrz0dhGN3m4IMD2Kx6F76htl+NsA\nSgRmXnSJDjzfNYwuRDHsygb52CjSbX3ke8CFSM98N5Jl9SQuY5COPhAHY5NE+3MO523hYKNIP/D6\nceadbRfFme9Dux0BAEyLLkQAgCQCDAAgiQADAEgiwAAAkggwAIAkAgwAIIkAAwBIIsAAAJIIMACA\nJAIMACCJAAMASCLAAACSCDAAgCQCDAAgiQADAEgiwAAAkggwAIAkAgwAIIkAAwBIIsAAAJIIMACA\nJAIMACCJAAMASCLAAACSCDAAgCQCDAAgiQADAEgiwAAAkggwAIAkAgwAIIkAAwBIIsAAAJIIMACA\nJAIMACDp/wGwvCO9Qb4lRQAAAABJRU5ErkJggg==\n", "output_type": "display_data"}], "prompt_number": 21, "cell_type": "code", "language": "python", "metadata": {}, "input": ["exo2()"]}, {"collapsed": false, "outputs": [], "prompt_number": 22, "cell_type": "code", "language": "python", "metadata": {}, "input": ["%% Insert your code here."]}, {"source": ["Display the signal recovered by OMP."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [{"metadata": {}, "png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGwCAIAAADOgk3lAAAACXBIWXMAAAsSAAALEgHS3X78AAAA\nIXRFWHRTb2Z0d2FyZQBBcnRpZmV4IEdob3N0c2NyaXB0IDguNTRTRzzSAAAOEklEQVR4nO3dUVez\nuAKG0eYs//9fzrlgDcMUxFoLyZvsfTHL+aw20pSHULSl1voAgDT/az0AAHiHgAEQScAAiCRgAEQS\nMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZAJAEDIJKAARBJwACI\nJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQScAAiCRgAEQaMGCllNZDAOByX60H8EnSBTCPoVZg\ntdZaa+tRdKSUh6ZD5y56ns7w9B9qBbY3/ZqsPmwE6N1Fz9O/ftv+1wODB+zR+jEopbQewKPW2nwY\n/wzGMPoaRg9jMIwlMetdf3AYy9P/3a8NOPAd6hQiQJwO2p1KwACINGDAejgXAcDVujj1fJ1Ozq03\nVIoTFNC7i56nf/m2ETvPAVdgAMxAwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQScAAiCRg\nAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZAJAEDIJKAARBJ\nwACIJGAARBIwACJ9tR7Ar5VSlg9qrd996vCzAIwkLGCllLVM249XugUwidFOIZZStuswAEYVtgL7\n0bICe1qo7W8AwJO4o/+hAnYYJ8UCeMV2bxkRs3FOIUZsbgA+JWwFVmvdX4W4nDA8/BQAowoL2OMo\nTuu/6BbAPMY5hQjAVAQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZAJAEDIJKA\nARBJwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAk\nAQMgkoABEOmr9QA+qZSyfFBrbTsSAK42TsBKKWu3th8DMCSnEAGINM4K7DvrecVpTb8BIMBFz9Ox\nn/7jB8y5RCZRyuN8sp/c4Mev/ctXvffNL/I0mK7G1pWIQ3+nEAGINM4KrNbqKkSAeYwTsIduAczE\nKUQAIgkYAJEEDIBIAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQScAA\niCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZAJAEDINJX\n6wH8Will+aDW+t2nDj8LwEjCAlZKWcu0/XilWwCTGO0UYilluw4DYFRhK7AfLSuwp4Xa/gYAPIk7\n+u86YL9tz+ENFAvgFdu9ZUTMug7Yr9pz+JIYAKPqOmB7tdb9VYhLug4/BcCowgL2OIrT+i+6BTCP\n0a5CBGASAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQScAAiCRgAEQS\nMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZAJAEDIJKAXauU0noI\nj4dh/JdhdDWGh2H8VyfDiCBgAEQaMGCOXwBm8NV6AJ8kXQDzKLXW1mP4sFL+/aEkjWksc/5kwp/c\n4Mevfe8etzd74/tfZDvs937wWfRfh8EDBpNYDtVOJv7JDX782vfucXuzN77/RbbDfu8Hpx+p+/qn\npdX2pxAwgBmkvgYmUQCTG/AqRABm4GwbAJFSTyH2bH19bn8x5P2HC+srgq3G8HS/bTdF22GcPBZ3\njqe3Yez/97ZhnNxp201x8zB+nAkNn7bnBOwS66Nea316kty8x9zf751j2D8VG26KhsPYXnO0fyxu\ne3T2v1XSZKLuh3G+fW4YQ6uJergpmkzUk5lw5zB+y2tgn9fJA9zPVCul9DOYJmqtPfz4T8NoNaSn\nYTSZG4ePyP0TtZ+J0XoIb7ICu8rku+xVP79Uvh5LNh9JV0zUhYmaOBME7POeFt0Nx7D978z8cZZD\nHU7UuB3oZ7WaqD3MhPcI2CWaT4X9qyxtx0OfepuobQczs9CNL2Af9rToWc5xN7+Gp9UY9vc78zC+\nG0yTR6fPifpotDU6mSFNhvHKTOhhYhxy1ANAJFchAhBJwACIJGAARBIwACIJGACRBAyASAIGQCQB\nAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBI\nAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQScAAiCRgAEQSMAAiCRgA\nkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIn21HsC1SimthwDQSn083t8H1lo/OJQrDB6wR+vH\noJTSwyQwDMPodgyGcd0wSnl/Bxhx9O8UIgCRBAyASAIGQKQuzvlep5OT2gD3K+Xx9v4vYudpBQZA\nJAEDIJKAARBJwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQaaj3A1vfwKb/P+EFwB+NtgKr\ntdZaI96KDYC/GCpgFl4A8xjqFOLi6V0AnlZjIgdwKO7c1VABW7b+U6IUC+AVJ4f+fRrqFOJDrgCm\nMc4KbDlecCEiwCTGCZhiAUxltFOIAExCwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQScAA\nWirlkfCHc3skYABEEjCAlvwZ17cJGACRBAyASAIGQCQBAyCSgAEQKfUdmUsp+7dgLptfpvAGzQBj\nywtYOf2VP90CmETeKcRa60mlSinnhQNgDHkrsHNL27YnGJ96ZokGcCju6H+ogB3GSbEAXrHdW0bE\nLO8U4nciNjcAnzLCCmw5YVhrXRtm1QUwvNSAbRO1fqxbAPMY5xQi8/D+ScBDwAAIJWAARBIw8nit\nE3gIGAChBAyASAIGQCQBAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEGnAgBXv\ndQgwga/WA/gk6Qq1PG7eJAX4laFWYLXWai8IMIehVmCHnpZlCgdwKO4k1vgBU6z+1fpIe+LAgLZ7\ny4iYDXUKEYB5CBg3KcUyC/ikAQPmnCHADAYMGH1yXAF8loABEEnAAIgkYABEEjAAIgkYAJEEDIBI\nAgZAJAEDIJKAARBJwACIJGAARBIwACLNGzDv7gEQbd6AARBt6oB5gw+AXFMHDIBcX60H8Gvln1eu\n9u+8XDYvanlfZmAYpThjdCAsYKWUtUzbj1e6BQxmOTLXsL3RTiGWUoqLC69nG8PN1GsvbAX2o2UF\n9rRQ29+Av3A8CHeq9aZDxrij/6ECdhgnxbqI7QqD2e4tI2I2zinEiM09BukCehC2Aqu17q9CXE4Y\nHn4KgFGFBexxFKf1X3QL5rQcu769A/CCbqhxTiECvGG9KIk4AgbE+8v6aflaK7BEAgZAJAEDIJKA\nARBJwACIJGBczptfA1cQsPuk7MdTxglMTsC4nAuUGZgDvoYE7D4p+/GUcUJz0tWWgA3CYSDczy9B\ntyVgQ9EwYB4CNghHgsBsBAyASAIGQCQBo2suTgG+I2CQR9ThkfiOzCfKP09rb808jFrtrJ/ZILAY\nJ2CllLVb24+hZ599M/ulbeY+kxhnR38YsOJgFeAt/ddhnBXYN7p4CF48Ll5v9uJR+Y832+b77e3w\n2SXC/d9/e0ePd7fD09cuY/77cufwO5x825PbL55G+MdhfPcdtp9a731/46uXg398FNYv//FJ1M/8\nfHHMHxpPwNH/8AFr7/Vp8PFJ2UO8B/b3zXv4HU6+7Yu3X6bc67u5+3+QDxpmkg/zg9zJVYiX8zcy\nOlHr+4/CX772fndOubZbJmGRwIXGWYHVWru9CrGz4TA+U44ZjBOwR3/dYki/PUHHRWx/xj+F6CQD\nn+Wc8Dy62nusR06sxg8YPQt9NqrXVDqZpWbd3vABKx514D29rbazLia6wVCvgRHHs5HOmaI9G34F\nBsCYBAz4tU5eFmJyApbNfgSYltfAgN/xshCdsALrzuu/7eGSJGBmAgZAJKcQu2NRBfAKKzBgOq5+\nGoOAARBJwACI5DUwYDrLK81OJKazAgMgUt4K7ORtl8vmgMqbWwLncncS3k91ERawUspapu3HK90C\nmMRopxBLKaWnE9udDOaPw/jUDzHG1viUHobRwxgehvFfPw7Dn+BZha3AfrSswJ4Wavsb8KI+ntHA\nHTpJ+Ou6Dthv23N4A8X6CxsP5rHdW0bErOuA/ao9hy+JATCqvJ3+/irENV2Hn7p9gDCG5UnkGTSv\n/uuQFzDgBsuxn90DPev6FCLQinTRv9EuowdgEgIGQCQBAyCS18A+7+RiyPsvmTm5RPO2AWzvt+2m\naDuMVy6XvWE8vQ1j/7+3DePkTttuipuH8eNMaPi0PSdgl9j+QZAf/37jddZp12oM+6diw03RcBjb\nX+fYPxa3PTr73yppMlH3wzjfPjeModVEPdwUTSbqyUy4cxi/5RTi53XyAPcz1ZY/UNnJYJqotfbw\n4z8No9WQnobRZG4cPiL3T9R+JkbrIbzJCuwqk++yV/38Uvl6LNl8JF0xURcmauJMELDPe1p0NxzD\n9r8zO/njzjPrcKLG7UA/q9VE7WEmvEfALtF8KuxfZWk7HvrU20RtO5iZhW58Afuwp0XPco67+TU8\nrcawv9+Zh/HdYJo8On1O1EejrdHJDGkyjFdmQg8T45CjHgAiuQoRgEgCBkAkAQMgkoABEEnAAIgk\nYABEEjAAIgkYAJEEDIBIAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQ\nScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZApP8D\n1Fger0HITDQAAAAASUVORK5CYII=\n", "output_type": "display_data"}], "prompt_number": 23, "cell_type": "code", "language": "python", "metadata": {}, "input": ["clf;\n", "subplot(2,1,1);\n", "x1 = x0; x1(x1==0) = NaN; stem(x1, 'b.'); axis tight;\n", "subplot(2,1,2);\n", "x1 = x; x1(x1==0) = NaN; stem(x1, 'b.'); axis tight;"]}], "metadata": {}}], "nbformat": 3, "metadata": {"kernelspec": {"name": "matlab_kernel", "language": "matlab", "display_name": "Matlab"}, "language_info": {"mimetype": "text/x-matlab", "name": "matlab", "file_extension": ".m", "help_links": [{"url": "https://github.com/calysto/metakernel/blob/master/metakernel/magics/README.md", "text": "MetaKernel Magics"}]}}}