Add SmoothLassoLS method

.pre-commit-config.yaml

@@ -17,12 +17,12 @@ repos:
     rev: v1.10.0
     hooks:
       - id: blacken-docs
-        additional_dependencies: [black==19.3b0]
+        additional_dependencies: [black]
         language: python
         language_version: python
 
   - repo: https://github.com/pre-commit/pre-commit-hooks
-    rev: v1.2.3
+    rev: v3.4.0
     hooks:
       - id: check-yaml
         language: python
@@ -32,8 +32,6 @@ repos:
         language: python
       - id: check-docstring-first
         language: python
-      - id: flake8
-        language: python
       - id: end-of-file-fixer
         language: python
         exclude: docs/notebooks/

CHANGELOG

@@ -4,6 +4,7 @@ v0.2.0
 What's new!
 '''''''''''
 
+- Added :func:`~mrinversion.linear_model.SmoothLassoLS` function to evaluate optimum alpha and lambda values using least-squares analysis.
 - Added :func:`~mrinversion.utils.to_Haeberlen_grid` function to convert the 3D :math:`\rho(\delta_\text{iso}, x, y)`
   distribution to :math:`\rho(\delta_\text{iso}, \zeta_\sigma, \eta_\sigma)` distribution.
 

docs/api/SmoothLassoLS.rst

@@ -0,0 +1,13 @@
+Smooth Lasso cross-validation with Least-squares
+================================================
+
+.. currentmodule:: mrinversion.linear_model
+
+.. autoclass:: SmoothLassoLS
+   :show-inheritance:
+
+   .. rubric:: Methods Documentation
+
+   .. automethod:: fit
+   .. automethod:: predict
+   .. automethod:: residuals

docs/getting_started.rst

@@ -151,7 +151,7 @@ csdmpy library as shown below:
 Here, the anisotropic dimension is sampled at 625 Hz for 32 points with an offset of
 -10 kHz.
 
-Similarly, we can create the CSDM dimensions needed for the *x*-*y* inversion grid as
+Similarly, we can create the CSDM dimensions needed for the *x-y* inversion grid as
 shown below:
 
 .. plot::
@@ -228,7 +228,7 @@ instance to generate the spinning sideband kernel, as follows,
 
 Here, ``K`` is the :math:`32\times 625` kernel, where the 32 is the number of samples
 (sideband amplitudes), and 625 is the number of features (subspectra) on the
-:math:`25 \times 25` *x*-*y* grid. The argument *supersampling* is the supersampling
+:math:`25 \times 25` *x-y* grid. The argument *supersampling* is the supersampling
 factor. In a supersampling scheme, each grid cell is averaged over a :math:`n\times n`
 sub-grid, where :math:`n` is the supersampling factor. A supersampling factor of 1 is
 equivalent to no sub-grid averaging.
@@ -425,7 +425,7 @@ The arguments of the :py:class:`~mrinversion.linear_model.SmoothLassoCV` is a li
 of the *alpha* and *lambda* values, along with the standard deviation of the
 noise, *sigma*. The value of the argument *folds* is the number of folds used in the
 cross-validation. As before, to solve the problem, use the
-:meth:`~mrinversion.linear_model.SmoothLassoCV.fit` method, whose arguments are
+:py:meth:`~mrinversion.linear_model.SmoothLassoCV.fit` method, whose arguments are
 the kernel and signal.
 
 The optimum hyperparameters

docs/introduction.rst

@@ -8,10 +8,6 @@ In ``mrinversion``, we solve for the distribution of the second-rank traceless
 symmetric tensor principal components, through an inversion of a pure anisotropic
 NMR spectrum.
 
-.. whose frequency
-.. contributions are assumed to arise predominantly from the second-rank traceless
-.. symmetric tensors.
-
 In the case of the shielding tensors, the pure anisotropic frequency spectra corresponds
 the cross-sections of the 2D isotropic *v.s.* anisotropic correlation spectrum, such as
 the 2D One Pulse (TOP) MAS, phase adjusted spinning sidebands (PASS), magic-angle turning
@@ -201,7 +197,7 @@ second-rank traceless tensor parameters, :math:`\zeta`-:math:`\eta`, defined as
 .. math::
     :label: zeta_eta_def
 
-    r_\zeta = | \zeta_ | ~~~~\text{and}~~~~
+    r_\zeta = | \zeta | ~~~~\text{and}~~~~
     \theta = \left\{ \begin{array}{l r}
                 \frac{\pi}{4} \eta      &: \zeta \le 0, \\
                 \frac{\pi}{2} \left(1 - \frac{\eta}{2} \right) &: \zeta > 0.

docs/notebooks/auto_examples/sideband/plot_2D_KMg0.5O4SiO2.ipynb

@@ -400,7 +400,7 @@
       },
       "outputs": [],
       "source": [
-        "# Isotropic chemical shift projection of the 2D MAT dataset.\ndata_iso = data_object_truncated.sum(axis=0)\ndata_iso /= data_iso.max()  # normalize the projection\n\n# Isotropic chemical shift projection of the tensor distribution dataset.\nf_sol_iso = f_sol.sum(axis=(0, 1))\n\n# Isotropic chemical shift projection of the tensor distribution for the Q4 region.\nQ4_region_iso = Q4_region.sum(axis=(0, 1))\n\n# Isotropic chemical shift projection of the tensor distribution for the Q3 region.\nQ3_region_iso = Q3_region.sum(axis=(0, 1))\n\n# Normalize the three projections.\nf_sol_iso_max = f_sol_iso.max()\nf_sol_iso /= f_sol_iso_max\nQ4_region_iso /= f_sol_iso_max\nQ3_region_iso /= f_sol_iso_max\n\n# The plot of the different projections.\nplt.figure(figsize=(5.5, 3.5))\nax = plt.subplot(projection=\"csdm\")\nax.plot(f_sol_iso, \"--k\", label=\"tensor\")\nax.plot(Q4_region_iso, \"r\", label=\"Q4\")\nax.plot(Q3_region_iso, \"b\", label=\"Q3\")\nax.plot(data_iso, \"-k\", label=\"MAF\")\nax.plot(data_iso - f_sol_iso - 0.1, \"gray\", label=\"residuals\")\nax.set_title(\"Isotropic projection\")\nax.invert_xaxis()\nplt.legend()\nplt.tight_layout()\nplt.show()"
+        "# Isotropic chemical shift projection of the 2D MAT dataset.\ndata_iso = data_object_truncated.sum(axis=0)\ndata_iso /= data_iso.max()  # normalize the projection\n\n# Isotropic chemical shift projection of the tensor distribution dataset.\nf_sol_iso = f_sol.sum(axis=(0, 1))\n\n# Isotropic chemical shift projection of the tensor distribution for the Q4 region.\nQ4_region_iso = Q4_region.sum(axis=(0, 1))\n\n# Isotropic chemical shift projection of the tensor distribution for the Q3 region.\nQ3_region_iso = Q3_region.sum(axis=(0, 1))\n\n# Normalize the three projections.\nf_sol_iso_max = f_sol_iso.max()\nf_sol_iso /= f_sol_iso_max\nQ4_region_iso /= f_sol_iso_max\nQ3_region_iso /= f_sol_iso_max\n\n# The plot of the different projections.\nplt.figure(figsize=(5.5, 3.5))\nax = plt.subplot(projection=\"csdm\")\nax.plot(f_sol_iso, \"--k\", label=\"tensor\")\nax.plot(Q4_region_iso, \"r\", label=\"Q4\")\nax.plot(Q3_region_iso, \"b\", label=\"Q3\")\nax.plot(data_iso, \"-k\", label=\"MAT\")\nax.plot(data_iso - f_sol_iso - 0.1, \"gray\", label=\"residuals\")\nax.set_title(\"Isotropic projection\")\nax.invert_xaxis()\nplt.legend()\nplt.tight_layout()\nplt.show()"
       ]
     },
     {

docs/notebooks/auto_examples/synthetic/plot_1D_1_sideband.ipynb

@@ -195,7 +195,7 @@
       },
       "outputs": [],
       "source": [
-        "# guess alpha and lambda values.\ns_lasso = SmoothLasso(alpha=5e-5, lambda1=5e-6, inverse_dimension=inverse_dimension)\ns_lasso.fit(K=compressed_K, s=compressed_s)\nf_sol = s_lasso.f"
+        "# guess alpha and lambda values.\ns_lasso = SmoothLasso(\n    hyperparameters={\"alpha\": 5e-5, \"lambda\": 5e-6},\n    inverse_dimension=inverse_dimension,\n)\ns_lasso.fit(K=compressed_K, s=compressed_s)\nf_sol = s_lasso.f"
       ]
     },
     {

docs/notebooks/auto_examples/synthetic/plot_1D_2_sideband_bimodal.ipynb

@@ -177,7 +177,7 @@
       },
       "outputs": [],
       "source": [
-        "# guess alpha and lambda values.\ns_lasso = SmoothLasso(alpha=5e-5, lambda1=5e-6, inverse_dimension=inverse_dimension)\ns_lasso.fit(K=compressed_K, s=compressed_s)\nf_sol = s_lasso.f"
+        "# guess alpha and lambda values.\ns_lasso = SmoothLasso(\n    hyperparameters={\"alpha\": 5e-5, \"lambda\": 5e-6},\n    inverse_dimension=inverse_dimension,\n)\ns_lasso.fit(K=compressed_K, s=compressed_s)\nf_sol = s_lasso.f"
       ]
     },
     {

docs/notebooks/auto_examples/synthetic/plot_1D_3_MAF.ipynb

@@ -33,7 +33,7 @@
       },
       "outputs": [],
       "source": [
-        "import csdmpy as cp\nimport matplotlib.pyplot as plt\nimport numpy as np\n\nfrom mrinversion.kernel.nmr import ShieldingPALineshape\nfrom mrinversion.linear_model import SmoothLasso, SmoothLassoCV, TSVDCompression\nfrom mrinversion.utils import get_polar_grids\n\n# Setup for the matplotlib figures\n\n\n# function for 2D x-y plot.\ndef plot2D(ax, csdm_object, title=\"\"):\n    # convert the dimension coordinates of the csdm_object from Hz to pmm.\n    _ = [item.to(\"ppm\", \"nmr_frequency_ratio\") for item in csdm_object.dimensions]\n\n    levels = (np.arange(9) + 1) / 10\n    plt.figure(figsize=(4.5, 3.5))\n    ax.contourf(csdm_object, cmap=\"gist_ncar\", levels=levels)\n    ax.grid(None)\n    ax.set_title(title)\n    get_polar_grids(ax)\n    ax.set_aspect(\"equal\")"
+        "import csdmpy as cp\nimport matplotlib.pyplot as plt\nimport numpy as np\n\nfrom mrinversion.kernel.nmr import ShieldingPALineshape\nfrom mrinversion.linear_model import SmoothLasso, SmoothLassoCV, TSVDCompression\nfrom mrinversion.utils import get_polar_grids\n\n# Setup for the matplotlib figures\n\n\n# function for 2D x-y plot.\ndef plot2D(ax, csdm_object, title=\"\"):\n    # convert the dimension coordinates of the csdm_object from Hz to pmm.\n    _ = [item.to(\"ppm\", \"nmr_frequency_ratio\") for item in csdm_object.dimensions]\n\n    levels = (np.arange(9) + 1) / 10\n    plt.figure(figsize=(4.5, 3.5))\n    ax.contourf(csdm_object, cmap=\"gist_ncar\", levels=levels)\n    ax.set_xlim(0, 100)\n    ax.set_ylim(0, 100)\n    ax.set_title(title)\n    get_polar_grids(ax)\n    ax.set_aspect(\"equal\")"
       ]
     },
     {
@@ -195,7 +195,7 @@
       },
       "outputs": [],
       "source": [
-        "# guess alpha and lambda values.\ns_lasso = SmoothLasso(alpha=5e-5, lambda1=5e-6, inverse_dimension=inverse_dimension)\ns_lasso.fit(K=compressed_K, s=compressed_s)\nf_sol = s_lasso.f"
+        "# guess alpha and lambda values.\ns_lasso = SmoothLasso(\n    hyperparameters={\"alpha\": 5e-5, \"lambda\": 5e-6},\n    inverse_dimension=inverse_dimension,\n)\ns_lasso.fit(K=compressed_K, s=compressed_s)\nf_sol = s_lasso.f"
       ]
     },
     {
@@ -249,7 +249,7 @@
       },
       "outputs": [],
       "source": [
-        "lambdas = 10 ** (-5.2 - 1 * (np.arange(6) / 5))\nalphas = 10 ** (-4 - 2 * (np.arange(6) / 5))\n\n# set up cross validation smooth lasso method\ns_lasso_cv = SmoothLassoCV(\n    alphas=alphas,\n    lambdas=lambdas,\n    inverse_dimension=inverse_dimension,\n    sigma=0.005,\n    folds=10,\n)\n# run the fit using the compressed kernel and compressed signal.\ns_lasso_cv.fit(compressed_K, compressed_s)"
+        "lambdas = 10 ** (-4.8 - 1.5 * (np.arange(6) / 5))\nalphas = 10 ** (-4 - 2 * (np.arange(6) / 5))\n\n# set up cross validation smooth lasso method\ns_lasso_cv = SmoothLassoCV(\n    alphas=alphas,\n    lambdas=lambdas,\n    inverse_dimension=inverse_dimension,\n    sigma=0.005,\n    folds=10,\n)\n# run the fit using the compressed kernel and compressed signal.\ns_lasso_cv.fit(compressed_K, compressed_s)"
       ]
     },
     {

docs/notebooks/auto_examples/synthetic/plot_1D_4_MAF_bimodal.ipynb

@@ -177,7 +177,7 @@
       },
       "outputs": [],
       "source": [
-        "# guess alpha and lambda values.\ns_lasso = SmoothLasso(alpha=5e-5, lambda1=5e-6, inverse_dimension=inverse_dimension)\ns_lasso.fit(K=compressed_K, s=compressed_s)\nf_sol = s_lasso.f"
+        "# guess alpha and lambda values.\ns_lasso = SmoothLasso(\n    hyperparameters={\"alpha\": 5e-5, \"lambda\": 5e-6},\n    inverse_dimension=inverse_dimension,\n)\ns_lasso.fit(K=compressed_K, s=compressed_s)\nf_sol = s_lasso.f"
       ]
     },
     {

docs/referenceAPI.rst

@@ -8,5 +8,6 @@ API-Reference
     api/shielding_kernel
     api/SmoothLasso
     api/SmoothLassoCV
+    api/SmoothLassoLS
     api/TSVDCompression
     api/utils

docs/requirement.rst

@@ -12,7 +12,7 @@ The ``mrinversion`` library depends on the following packages:
 - `csdmpy>=0.4 <https://csdmpy.readthedocs.io/en/stable/>`_
 - `mrsimulator>=0.6 <https://mrsimulator.readthedocs.io/en/stable/>`_ (for generating
   the NMR line-shape)
-- `scikit-learn>=0.22.1 <https://scikit-learn.org/stable/>`_ (for statistical leaning methods)
+- `scikit-learn>=0.24.1 <https://scikit-learn.org/stable/>`_ (for statistical leaning methods)
 
 **Other packages**
 

docs/requirements.txt

@@ -3,4 +3,4 @@ sphinx>=2.1.2
 sphinx-gallery
 pillow
 matplotlib
-scikit-learn>=0.22
+scikit-learn>=0.24.1

examples/sideband/plot_2D_KMg0.5O4SiO2.py

@@ -394,7 +394,7 @@ def plot2D(csdm_object, **kwargs):
 ax.plot(f_sol_iso, "--k", label="tensor")
 ax.plot(Q4_region_iso, "r", label="Q4")
 ax.plot(Q3_region_iso, "b", label="Q3")
-ax.plot(data_iso, "-k", label="MAF")
+ax.plot(data_iso, "-k", label="MAT")
 ax.plot(data_iso - f_sol_iso - 0.1, "gray", label="residuals")
 ax.set_title("Isotropic projection")
 ax.invert_xaxis()

examples/synthetic/plot_1D_1_sideband.py

@@ -160,7 +160,10 @@ def plot2D(ax, csdm_object, title=""):
 # accordingly.
 
 # guess alpha and lambda values.
-s_lasso = SmoothLasso(alpha=5e-5, lambda1=5e-6, inverse_dimension=inverse_dimension)
+s_lasso = SmoothLasso(
+    hyperparameters={"alpha": 5e-5, "lambda": 5e-6},
+    inverse_dimension=inverse_dimension,
+)
 s_lasso.fit(K=compressed_K, s=compressed_s)
 f_sol = s_lasso.f
 

examples/synthetic/plot_1D_2_sideband_bimodal.py

@@ -143,7 +143,10 @@ def plot2D(ax, csdm_object, title=""):
 # accordingly.
 
 # guess alpha and lambda values.
-s_lasso = SmoothLasso(alpha=5e-5, lambda1=5e-6, inverse_dimension=inverse_dimension)
+s_lasso = SmoothLasso(
+    hyperparameters={"alpha": 5e-5, "lambda": 5e-6},
+    inverse_dimension=inverse_dimension,
+)
 s_lasso.fit(K=compressed_K, s=compressed_s)
 f_sol = s_lasso.f
 

examples/synthetic/plot_1D_3_MAF.py

@@ -33,7 +33,8 @@ def plot2D(ax, csdm_object, title=""):
     levels = (np.arange(9) + 1) / 10
     plt.figure(figsize=(4.5, 3.5))
     ax.contourf(csdm_object, cmap="gist_ncar", levels=levels)
-    ax.grid(None)
+    ax.set_xlim(0, 100)
+    ax.set_ylim(0, 100)
     ax.set_title(title)
     get_polar_grids(ax)
     ax.set_aspect("equal")
@@ -159,7 +160,10 @@ def plot2D(ax, csdm_object, title=""):
 # accordingly.
 
 # guess alpha and lambda values.
-s_lasso = SmoothLasso(alpha=5e-5, lambda1=5e-6, inverse_dimension=inverse_dimension)
+s_lasso = SmoothLasso(
+    hyperparameters={"alpha": 5e-5, "lambda": 5e-6},
+    inverse_dimension=inverse_dimension,
+)
 s_lasso.fit(K=compressed_K, s=compressed_s)
 f_sol = s_lasso.f
 
@@ -209,7 +213,7 @@ def plot2D(ax, csdm_object, title=""):
 # hyperparameters.
 # The following code generates a range of :math:`\lambda` and :math:`\alpha` values
 # that are uniformly sampled on the log scale.
-lambdas = 10 ** (-5.2 - 1 * (np.arange(6) / 5))
+lambdas = 10 ** (-4.8 - 1.5 * (np.arange(6) / 5))
 alphas = 10 ** (-4 - 2 * (np.arange(6) / 5))
 
 # set up cross validation smooth lasso method

examples/synthetic/plot_1D_4_MAF_bimodal.py

@@ -142,7 +142,10 @@ def plot2D(ax, csdm_object, title=""):
 # accordingly.
 
 # guess alpha and lambda values.
-s_lasso = SmoothLasso(alpha=5e-5, lambda1=5e-6, inverse_dimension=inverse_dimension)
+s_lasso = SmoothLasso(
+    hyperparameters={"alpha": 5e-5, "lambda": 5e-6},
+    inverse_dimension=inverse_dimension,
+)
 s_lasso.fit(K=compressed_K, s=compressed_s)
 f_sol = s_lasso.f
 

mrinversion/__init__.py

@@ -1,3 +1,3 @@
 # -*- coding: utf-8 -*-
 
-__version__ = "0.2.0"
+__version__ = "0.3.0dev0"

mrinversion/kernel/base.py

@@ -15,12 +15,27 @@
 class BaseModel:
     """Base kernel class."""
 
-    def __init__(self, kernel_dimension, inverse_kernel_dimension, n_dir, n_inv):
+    def __init__(
+        self,
+        kernel_dimension: list,
+        inverse_kernel_dimension: list,
+        n_dir: int,
+        n_inv: int,
+    ):
+        """Initialize the class
+
+        Args:
+            kernel_dimension: A list of kernel CSDM dimensions.
+            inverse_kernel_dimension: A list of inverse kernel CSDM dimensions.
+            n_dir: Number of kernel dimensions.
+            n_inv: Number of inverse dimensions.
+        """
         kernel = self.__class__.__name__
         message = (
             f"Exactly {n_inv} inverse dimension(s) is/are required for the "
             f"{kernel} kernel."
         )
+
         if isinstance(inverse_kernel_dimension, list):
             if len(inverse_kernel_dimension) != n_inv:
                 raise ValueError(message)
@@ -158,10 +173,9 @@ def __init__(
     def _get_zeta_eta(self, supersampling):
         """Return zeta and eta coordinates over x-y grid"""
 
-        zeta, eta = _x_y_to_zeta_eta_distribution(
+        return _x_y_to_zeta_eta_distribution(
             self.inverse_kernel_dimension, supersampling
         )
-        return zeta, eta
 
 
 def _check_csdm_dimension(dimensions, dimension_id):

mrinversion/kernel/csa_aniso.py

@@ -178,7 +178,7 @@ def __init__(
             inverse_dimension,
             channel,
             magnetic_flux_density,
-            "54.735 deg",
+            "54.7356 deg",
             None,
             None,
         )