.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "gallery/0200_cs/matching_pursuit.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code .. rst-class:: sphx-glr-example-title .. _sphx_glr_gallery_0200_cs_matching_pursuit.py: .. _gallery:cs:mp:1: Matching Pursuit =================== .. contents:: :depth: 2 :local: This is a very simple example of using the matching pursuit algorithm. .. GENERATED FROM PYTHON SOURCE LINES 14-19 .. code-block:: default # Configure JAX to work with 64-bit floating point precision. from jax.config import config config.update("jax_enable_x64", True) .. GENERATED FROM PYTHON SOURCE LINES 20-21 Let's import necessary libraries .. GENERATED FROM PYTHON SOURCE LINES 21-37 .. code-block:: default # random number generator from jax import random # numpy import numpy as np import jax.numpy as jnp # utilities import cr.nimble as crn # sample data import cr.sparse.data as crdata # linear operators import cr.sparse.lop as crlop # matching pursuit algorithm import cr.sparse.pursuit.mp as mp import matplotlib.pyplot as plt .. GENERATED FROM PYTHON SOURCE LINES 38-39 Some random number generation keys .. GENERATED FROM PYTHON SOURCE LINES 39-43 .. code-block:: default key = random.PRNGKey(3) keys = random.split(key, 5) .. GENERATED FROM PYTHON SOURCE LINES 44-46 Problem setup ---------------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 46-54 .. code-block:: default # Ambient dimension n = 400 # Number of non-zero entries in the sparse model k = 20 # Number of compressive measurements m = 200 .. GENERATED FROM PYTHON SOURCE LINES 55-57 Spikes as sample data -------------------------- .. GENERATED FROM PYTHON SOURCE LINES 57-62 .. code-block:: default x, omega = crdata.sparse_spikes(keys[0], n, k) plt.figure(figsize=(8,6), dpi= 100, facecolor='w', edgecolor='k') plt.plot(x) .. image-sg:: /gallery/0200_cs/images/sphx_glr_matching_pursuit_001.png :alt: matching pursuit :srcset: /gallery/0200_cs/images/sphx_glr_matching_pursuit_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none [] .. GENERATED FROM PYTHON SOURCE LINES 63-65 Gaussian sensing matrix linear operator ---------------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 65-71 .. code-block:: default Phi = crlop.gaussian_dict(keys[1], m, n) # Make sure that the linear operator is JIT compiled for efficiency. Phi = crlop.jit(Phi) .. GENERATED FROM PYTHON SOURCE LINES 72-75 Compressive sensing/measurements ---------------------------------------- Clean measurements .. GENERATED FROM PYTHON SOURCE LINES 75-85 .. code-block:: default y0 = Phi.times(x) # Noise sigma = 0.01 noise = sigma * random.normal(keys[2], (m,)) # Noisy measurements y = y0 + noise plt.figure(figsize=(8,6), dpi= 100, facecolor='w', edgecolor='k') plt.plot(y) print(f'Measurement noise: {crn.signal_noise_ratio(y0, y):.2f} dB') .. image-sg:: /gallery/0200_cs/images/sphx_glr_matching_pursuit_002.png :alt: matching pursuit :srcset: /gallery/0200_cs/images/sphx_glr_matching_pursuit_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none Measurement noise: 30.07 dB .. GENERATED FROM PYTHON SOURCE LINES 86-88 Reconstruction using matching pursuit ---------------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 88-93 .. code-block:: default sol = mp.solve(Phi, y, max_iters=k*2) print(sol) # solution vector x_hat = sol.x .. rst-class:: sphx-glr-script-out .. code-block:: none iterations=40 m=200, n=400, k=22 r_norm=4.318368e+00 x_norm=4.387438e+00 .. GENERATED FROM PYTHON SOURCE LINES 94-96 Solution ------------------------------------------------ .. GENERATED FROM PYTHON SOURCE LINES 96-102 .. code-block:: default plt.figure(figsize=(8,6), dpi= 100, facecolor='w', edgecolor='k') plt.subplot(211) plt.stem(x) plt.subplot(212) plt.stem(x_hat) .. image-sg:: /gallery/0200_cs/images/sphx_glr_matching_pursuit_003.png :alt: matching pursuit :srcset: /gallery/0200_cs/images/sphx_glr_matching_pursuit_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 103-105 Metrics ------------------------------------------------ .. GENERATED FROM PYTHON SOURCE LINES 105-110 .. code-block:: default snr = crn.signal_noise_ratio(x, x_hat) prd = crn.percent_rms_diff(x, x_hat) n_rmse = crn.normalized_root_mse(x, x_hat) print(f'SNR: {snr:.2f} dB, PRD: {prd:.2f} %, N-RMSE: {n_rmse:.2e}') .. rst-class:: sphx-glr-script-out .. code-block:: none SNR: 27.54 dB, PRD: 4.20 %, N-RMSE: 4.20e-02 .. GENERATED FROM PYTHON SOURCE LINES 111-113 Verifying the support recovery ''''''''''''''''''''''''''''''''''''' .. GENERATED FROM PYTHON SOURCE LINES 113-121 .. code-block:: default print('Support of original signal: ', omega) print('Support of reconstructed signal: ', sol.I) # check if every index in the original support is # also there in the reconstruction support print(np.all(np.in1d(omega, sol.I))) .. rst-class:: sphx-glr-script-out .. code-block:: none Support of original signal: [ 19 77 154 155 192 223 235 236 261 274 277 314 323 342 347 351 369 376 377 380] Support of reconstructed signal: [ 19 77 85 154 155 192 216 223 235 236 261 274 277 314 323 342 347 351 369 376 377 380] True .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 1.792 seconds) .. _sphx_glr_download_gallery_0200_cs_matching_pursuit.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: matching_pursuit.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: matching_pursuit.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_