Welcome to Mrinversion documentation!

Deployment

PyPI version PyPI - Python Version

Build Status

GitHub Workflow Status Documentation Status

License

License

Metrics

Language grade: Python https://codecov.io/gh/DeepanshS/mrinversion/branch/master/graph/badge.svg Total alerts

GitHub

GitHub issues

About

The mrinversion python package is based on the statistical learning technique for determining the distribution of the magnetic resonance (NMR) tensor parameters from two-dimensional NMR spectra correlating the isotropic to anisotropic frequencies. The library utilizes the mrsimulator package for generating solid-state NMR spectra and scikit-learn package for statistical learning.


Features

The mrinversion package includes the inversion of a two-dimensional solid-state NMR spectrum of dilute spin-systems to a three-dimensional distribution of tensor parameters. At present, we support the inversion of

  • Magic angle turning (MAT), Phase adjusted spinning sidebands (PASS), and similar spectra correlating the isotropic chemical shift resonances to pure anisotropic spinning sideband resonances into a three-dimensional distribution of nuclear shielding tensor parameters, \(\rho(\delta_\text{iso}, \zeta_\sigma, \eta_\sigma)\), where \(\delta_\text{iso}\) is the isotropic chemical shift, and \(\zeta_\sigma\) and \(\eta_\sigma\), are the shielding anisotropy and asymmetry parameters, respectively, defined using the Haeberlen convention.

  • Magic angle flipping (MAF) spectra correlating the isotropic chemical shift resonances to pure anisotropic resonances into a three-dimensional distribution of nuclear shielding tensor parameters, \(\rho(\delta_\text{iso}, \zeta_\sigma, \eta_\sigma)\), where \(\delta_\text{iso}\) is the isotropic chemical shift, and \(\zeta_\sigma\) and \(\eta_\sigma\), are the shielding anisotropy and asymmetry parameters, respectively, defined using the Haeberlen convention.


View our example gallery

https://img.shields.io/badge/View-Example%20Gallery-Purple?s=small

Getting Started

Installation

Requirements

mrinversion has the following strict requirements:

See Package dependencies for a full list of requirements.

Make sure you have the required version of python by typing the following in the terminal,

Tip

You may also click the copy-button located at the top-right corner of the code cell area in the HTML docs, to copy the code lines without the prompts and then paste it as usual. Thanks to Sphinx-copybutton)

$ python --version

For Mac users, python version 3 may be installed under the name python3. You may replace python for python3 in the above command and all subsequent python statements.

Installing mrinversion

On Google Colab Notebook

Colaboratory is a Google research project. It is a Jupyter notebook environment that runs entirely in the cloud. Launch a new notebook on Colab. To install the mrinversion package, type

!pip install mrinversion

in the first cell, and execute. All done! You may now proceed to the next section and start using the library.

On Local machine (Using pip)

The mrinversion package utilizes the mrsimulator package for generating the NMR line-shapes.

For Linux and Mac users, type the following in the terminal to install the package.

$ pip install mrinversion

For Windows users, first, install the mrsimulator package and then install the mrinversion package using the above command.

If you get a PermissionError, it usually means that you do not have the required administrative access to install new packages to your Python installation. In this case, you may consider adding the --user option, at the end of the statement, to install the package into your home directory. You can read more about how to do this in the pip documentation.

$ pip install mrinversion --user
Upgrading to a newer version

To upgrade, type the following in the terminal/Prompt,

$ pip install mrinversion -U

Package dependencies

The mrinversion library depends on the following packages:

Required packages

Other packages

  • pytest>=4.5.0 for unit tests.

  • pre-commit for code formatting

  • sphinx>=2.0 for generating the documentation

  • sphinxjp.themes.basicstrap for documentation.

  • sphinx-copybutton

Introduction

Objective

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.

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 (MAT), extended chemical shift modulation (XCS), magic-angle hopping (MAH), magic-angle flipping (MAF), and Variable Angle Correlation Spectroscopy (VACSY). A key feature of all these 2D isotropic/anisotropic correlation spectra—–either as acquired or after a shear transformation—–is that the anisotropic cross-section can be modeled as a linear combination of subspectra,

()\[s(\nu| \delta_\text{iso}) = \int_{\bf R} \mathcal{K}(\nu, {\bf R}) f({\bf R} | \delta_\text{iso}) d{\bf R},\]

where \(s(\nu| \delta_\text{iso})\) is the observed anisotropic cross-section at a given isotropic shift, \(\delta_\text{iso}\), \(\mathcal{K}(\nu, {\bf R})\) represents a simulated subspectrum of a nuclear spin system with a given set of parameters, \({\bf R}\), and \(f({\bf R} | \delta_\text{iso})\) is the probability of the respective set of parameters. In Eq. (1), \({\bf R}\) represents the anisotropic and asymmetry parameters of the shielding tensor.

Note, Eq. (1) is a Fredholm integral of the first kind.

Generic Linear problem

Linear inverse problems on Fredholm integral of the first kind are frequently encountered in the scientific community and have the following generic form

()\[ {\bf s} = {\bf K \cdot f},\]

where \({\bf K} \in \mathbb{R}^{m\times n}\) is the transforming kernel (matrix), \({\bf f} \in \mathbb{R}^n\) is the unknown and desired solution, and \({\bf s} \in \mathbb{R}^m\) is the known signal, which includes the measurement noise. When the matrix \({\bf K}\) is non-singular and \(m=n\), the solution to the problem in Eq. (2) has a simple closed-form solution,

()\[{\bf f} = {\bf K}^{-1} \cdot {\bf s}.\]

The deciding factor whether the solution \({\bf f}\) exists in Eq. (3) is whether or not the kernel \({\bf K}\) is invertible. Often, most scientific problems with practical applications suffer from singular, near-singular, or ill-conditioned kernels, where \({\bf K}^{-1}\) doesn’t exist. Such types of problems are termed as ill-posed. The inversion of a purely anisotropic NMR spectrum to the distribution of the tensorial parameters is one such ill-posed problem.

Regularized linear problem

A common approach in solving ill-posed problems is to employ the regularization methods of form

()\[{\bf f^\dagger} = \underset{{\bf f} > 0}{\text{argmin}} \left( \|{\bf K \cdot f} - {\bf s}\|^2_2 + g({\bf f}) \right),\]

where \(\|{\bf z}\|_2\) is the l-2 norm of \({\bf z}\), \(g({\bf f})\) is the regularization term, and \({\bf f}^\dagger\) is the regularized solution. The choice of the regularization term, \(g({\bf f})\), is often based on prior knowledge of the system for which the linear problem is defined. For anisotropic NMR spectrum inversion, we choose the smooth-LASSO regularization.

Smooth-LASSO regularization

Our prior assumption for the distribution of the tensorial parameters is that it should be smooth and continuous for disordered and sparse and discrete for crystalline materials. Therefore, we employ the smooth-lasso method, which is a linear model that is trained with the combined l1 and l2 priors as the regularizer. The method minimizes the objective function,

()\[\| {\bf K \cdot f - s} \|^2_2 + \alpha \sum_{i=1}^{d} \| {\bf J}_i \cdot {\bf f} \|_2^2 + \lambda \| {\bf f} \|_1 ,\]

where \(\alpha\) and \(\lambda\) are the hyperparameters controlling the smoothness and sparsity of the solution \({\bf f}\). The matrix \({\bf J}_i\) typically reflects some underlying geometry or the structure in the true solution. Here, \({\bf J}_i\) is defined to promote smoothness along the \(\text{i}^\text{th}\) dimension of the solution \({\bf f}\) and is given as

()\[{\bf J}_i = {\bf I}_{n_1} \otimes \cdots \otimes {\bf A}_{n_i} \otimes \cdots \otimes {\bf I}_{n_{d}},\]

where \({\bf I}_{n_i} \in \mathbb{R}^{n_i \times n_i}\) is the identity matrix, and \({\bf A}_{n_i}\) is the first difference matrix given as

()\[\begin{split}{\bf A}_{n_i} = \left(\begin{array}{ccccc} 1 & -1 & 0 & \cdots & \vdots \\ 0 & 1 & -1 & \cdots & \vdots \\ \vdots & \vdots & \vdots & \vdots & 0 \\ 0 & \cdots & 0 & 1 & -1 \end{array}\right) \in \mathbb{R}^{(n_i-1)\times n_i}.\end{split}\]

The symbol \(\otimes\) is the Kronecker product. The terms, \(\left(n_1, n_2, \cdots, n_d\right)\), are the number of points along the respective dimensions, with the constraint that \(\prod_{i=1}^{d}n_i = n\), where \(d\) is the total number of dimensions in the solution \({\bf f}\), and \(n\) is the total number of features in the kernel, \({\bf K}\).

Understanding the x-y plot

A second-rank symmetric tensor, \({\bf S}\), in a three-dimensional space, is described by three principal components, \(s_{xx}\), \(s_{yy}\), and \(s_{zz}\), in the principal axis system (PAS). Often, depending on the context of the problem, the three principal components are expressed with three new parameters following a convention. One such convention is the Haeberlen convention, which defines \(\delta_\text{iso}\), \(\zeta\), and \(\eta\), as the isotropic shift, anisotropy, and asymmetry parameters, respectively. Here, the parameters \(\zeta\) and \(\eta\) contribute to the purely anisotropic frequencies, and determining the distribution of these two parameters is the focus of this library.

Defining the inverse grid

When solving any linear inverse problem, one needs to define an inverse grid before solving the problem. A familiar example is the inverse Fourier transform, where the inverse grid is defined following the Nyquist–Shannon sampling theorem. Unlike inverse Fourier transform, however, there is no well-defined sampling grid for the second-rank traceless symmetric tensor parameters. One obvious choice is to define a two-dimensional \(\zeta\)-\(\eta\) Cartesian grid.

As far as the inversion problem is concerned, \(\zeta\) and \(\eta\) are just labels for the subspectra. In simplistic terms, the inversion problem solves for the probability of each subspectrum, from a given pre-defined basis of subspectra, that describes the observed spectrum. If the subspectra basis is defined over a \(\zeta\)-\(\eta\) Cartesian grid, multiple \((\zeta, \eta)\) coordinates points to the same subspectra. For example, the subspectra from coordinates \((\zeta, \eta=1)\) and \((-\zeta, \eta=1)\) are identical, therefore, distinguishing these coordinates from the subspectra becomes impossible.

The issue of multiple coordinates pointing to the same object is not new. It is a common problem when representing polar coordinates in the Cartesian basis. Try describing the coordinates of the south pole using latitudes and longitudes. You can define the latitude, but describing longitude becomes problematic. A similar situation arises in the context of second-rank traceless tensor parameters when the anisotropy goes to zero. You can specify the anisotropy as zero, but defining asymmetry becomes problematic.

Introducing the \(x\)-\(y\) grid

A simple fix to this issue is to define the \((\zeta, \eta)\) coordinates in a polar basis. We, therefore, introduce a piece-wise polar grid representation of the second-rank traceless tensor parameters, \(\zeta\)-\(\eta\), defined as

()\[\begin{split}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. \end{array} \right.\end{split}\]

Because Cartesian grids are more manageable in computation, we re-express the above polar piece-wise grid as the x-y Cartesian grid following,

()\[x = r_\zeta \cos\theta ~~~~\text{and}~~~~ y = r_\zeta \sin\theta.\]

In the x-y grid system, the basis subspectra are relatively distinguishable. The mrinversion library provides a utility function to render the piece-wise polar grid for your matplotlib figures. Copy-paste the following code in your script.

>>> import matplotlib.pyplot as plt 
>>> from mrinversion.utils import get_polar_grids 
...
>>> plt.figure(figsize=(4, 3.5)) 
>>> ax=plt.gca() 
>>> # add your plots/contours here.
>>> get_polar_grids(ax) 
>>> ax.set_xlabel('x / ppm') 
>>> ax.set_ylabel('y / ppm') 
>>> plt.tight_layout() 
>>> plt.show() 

(png, hires.png, pdf)

_images/introduction-1.png
_images/null.png

The figure depicts the piece-wise polar \(\zeta\)-\(\eta\) grid represented on an x-y grid. The radial and angular grid lines represent the magnitude of \(\zeta\) and \(\eta\), respectively. The blue and red shading represents the positive and negative values of \(\zeta\), respectively. The radian grid lines are drawn at every 0.2 ppm increments of \(\zeta\), and the angular grid lines are drawn at every 0.2 increments of \(\eta\). The x and y-axis are \(\eta=0\), and the diagonal \(x=y\) is \(\eta=1\).

If you are familiar with the matplotlib library, you may notice that most code lines are the basic matplotlib statements, except for the line that says get_polar_grids(ax). The get_polar_grids() is a utility function that generates the piece-wise polar grid for your figures.

Here, the shielding anisotropy parameter, \(\zeta\), is the radial dimension, and the asymmetry parameter, \(\eta\), is the angular dimension, defined using Eqs. (8) and (9). The region in blue and red corresponds to the positive and negative values of \(\zeta\), where the magnitude of the anisotropy increases radially. The x and the y-axis are \(\eta=0\) for the negative and positive \(\zeta\), respectively. When moving towards the diagonal from x or y-axes, the asymmetry parameter, \(\eta\), uniformly increase, where the diagonal is \(\eta=1\).

Before getting started

Prepping the 2D dataset for inversion

The following is a list of some requirements and recommendations to help prepare the 2D dataset for inversion.

Common recommendations/requirements

Dataset shear

The inversion method assumes that the 2D dataset is sheared, such that one of the dimensions corresponds to a pure anisotropic spectrum. The anisotropic cross-sections are centered at 0 Hz.

Required: Apply a shear transformation before proceeding.

Calculate the noise standard deviation

Use the noise region of your spectrum to calculate the standard deviation of the noise. You will require this value when implementing cross-validation.

Spinning Sideband correlation dataset specific recommendations

Data-repeat operation

A data-repeat operation on the time-domain signal corresponding to the sideband dimension makes the spinning sidebands look like a stick spectrum after a Fourier transformation, a visual, which most NMR spectroscopists are familiar from the 1D magic-angle spinning spectrum. Like a zero-fill operation, a spinning sideband data-repeat operation is purely cosmetic and adds no information. In terms of computation, however, a data-repeated spinning-sideband spectrum will take longer to solve.

Strongly recommended: Avoid data-repeat operation.

Magic angle flipping dataset specific recommendations

Isotropic shift correction along the anisotropic dimension

Ordinarily, after shear, a MAF spectrum is a 2D isotropic v.s. pure anisotropic frequency correlation spectrum. In certain conditions, this is not true. In a MAF experiment, the sample holder (rotor) physically swaps between two angles (\(90^\circ \leftrightarrow 54.735^\circ\)). It is possible to have a slightly different external magnetic fields at the two angles, in which case, there is an isotropic component along the anisotropic dimension, which is not removed by the shear transformation.

Required: Correct for the isotropic offset along the anisotropic dimension by adding an appropriate coordinates-offset.

Zero-fill operation

Zero filling the time domain dataset is purely cosmetic. It makes the spectrum look visually appealing, but adds no information, that is, a zero-filled dataset contains the same information as a non-zero filled dataset. In terms of computation, however, a zero-filled spectrum will take longer to solve.

Recommendation: If zero-filled, try to keep the total number of points along the anisotropic dimension in the range of 120 - 150 points.

Sinc wiggles artifacts

Kernel correction for spectrum with sinc wiggle artifacts is coming soon.

Getting started with mrinversion

We have put together a set of guidelines for using the mrinversion package. We encourage our users to follow these guidelines for consistency.

Let’s examine the inversion of a purely anisotropic MAS sideband spectrum into a 2D distribution of nuclear shielding anisotropy parameters. For illustrative purposes, we use a synthetic one-dimensional purely anisotropic spectrum. Think of this as a cross-section of your 2D MAT/PASS dataset.

Import relevant modules

>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> from matplotlib import rcParams
>>> from mrinversion.utils import get_polar_grids
...
>>> rcParams['pdf.fonttype'] = 42 # for exporting figures as illustrator editable pdf.
...
>>> # a function to plot the 2D tensor parameter distribution
>>> def plot2D(ax, csdm_object, title=''):
...     # convert the dimension from `Hz` to `ppm`.
...     _ = [item.to('ppm', 'nmr_frequency_ratio') for item in csdm_object.dimensions]
...
...     levels = (np.arange(9)+1)/10
...     ax.contourf(csdm_object, cmap='gist_ncar', levels=levels)
...     ax.grid(None)
...     ax.set_title(title)
...     ax.set_aspect("equal")
...
...     # The get_polar_grids method place a polar zeta-eta grid on the background.
...     get_polar_grids(ax)

Import the dataset

The first step is getting the sideband spectrum. In this example, we get the data from a CSDM 1 compliant file-format. Import the csdmpy module and load the dataset as follows,

Note

The CSDM file-format is a new open-source universal file format for multi-dimensional datasets. It is supported by NMR programs such as SIMPSON 2, DMFIT 3, and RMN 4. A python package supporting CSDM file-format, csdmpy, is also available.

>>> import csdmpy as cp
...
>>> filename = "https://osu.box.com/shared/static/xnlhecn8ifzcwx09f83gsh27rhc5i5l6.csdf"
>>> data_object = cp.load(filename) # load the CSDM file with the csdmpy module

Here, the variable data_object is a CSDM object. The NMR spectroscopic dimension is a frequency dimension. NMR spectroscopists, however, prefer to view the spectrum on a dimensionless scale. If the dataset dimension within the CSDM object is in frequency, you may convert it into ppm as follows,

>>> # convert the dimension coordinates from `Hz` to `ppm`.
>>> data_object.dimensions[0].to('ppm', 'nmr_frequency_ratio')

In the above code, we convert the dimension at index 0 from Hz to ppm. For multi-dimensional datasets, use the appropriate indexing to convert individual dimensions to ppm.

For comparison, let’s also include the true probability distribution from which the synthetic spinning sideband dataset is derived.

>>> datafile = "https://osu.box.com/shared/static/lufeus68orw1izrg8juthcqvp7w0cpzk.csdf"
>>> true_data_object = cp.load(datafile) # the true solution for comparison

The following is the plot of the spinning sideband spectrum as well as the corresponding true probability distribution.

>>> _, ax = plt.subplots(1, 2, figsize=(9, 3.5), subplot_kw={'projection': 'csdm'}) 
>>> ax[0].plot(data_object) 
>>> ax[0].set_xlabel('frequency / ppm') 
>>> ax[0].invert_xaxis() 
>>> ax[0].set_title('Pure anisotropic MAS spectrum') 
...
>>> plot2D(ax[1], true_data_object, title='True distribution') 
>>> plt.tight_layout() 
>>> plt.savefig('filename.pdf') # to save figure as editable pdf 
>>> plt.show() 

(png, hires.png, pdf)

_images/getting_started-5.png
_images/null.png

The figure on the left is the pure anisotropic MAS sideband amplitude spectrum corresponding to the nuclear shielding tensor distribution shown on the right.

Dimension Setup

For the inversion, we need to define (1) the coordinates associated with the pure anisotropic dimension, and (2) the two-dimensional x-y coordinates associated with the anisotropic tensor parameters, i.e., the inversion solution grid.

In mrinversion, the anisotropic spectrum dimension is initialized with a Dimension object from the csdmpy package. This object holds the frequency coordinates of the pure anisotropic spectrum. Because the example NMR dataset is imported as a CSDM object, the anisotropic spectrum dimension is already available as a CSDM Dimension object in the CSDM object and can be copied from there. Alternatively, we can create and initialize a anisotropic spectrum dimension using the csdmpy library as shown below:

>>> anisotropic_dimension = cp.LinearDimension(count=32, increment='625Hz', coordinates_offset='-10kHz')
>>> print(anisotropic_dimension)
LinearDimension([-10000.  -9375.  -8750.  -8125.  -7500.  -6875.  -6250.  -5625.  -5000.
  -4375.  -3750.  -3125.  -2500.  -1875.  -1250.   -625.      0.    625.
   1250.   1875.   2500.   3125.   3750.   4375.   5000.   5625.   6250.
   6875.   7500.   8125.   8750.   9375.] Hz)

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 shown below:

>>> inverse_dimension = [
...     cp.LinearDimension(count=25, increment='370 Hz', label='x'),  # the x-coordinates
...     cp.LinearDimension(count=25, increment='370 Hz', label='y')   # the y-coordinates
... ]

Both dimensions are sampled at every 370 Hz for 25 points. The inverse dimension at index 0 and 1 are the x and y dimensions, respectively.

Generating the kernel

Import the ShieldingPALineshape class and generate the kernel as follows,

>>> from mrinversion.kernel.nmr import ShieldingPALineshape
>>> lineshapes = ShieldingPALineshape(
...     anisotropic_dimension=anisotropic_dimension,
...     inverse_dimension=inverse_dimension,
...     channel='29Si',
...     magnetic_flux_density='9.4 T',
...     rotor_angle='54.735°',
...     rotor_frequency='625 Hz',
...     number_of_sidebands=32
... )

In the above code, the variable lineshapes is an instance of the ShieldingPALineshape class. The three required arguments of this class are the anisotropic_dimension, inverse_dimension, and channel. We have already defined the first two arguments in the previous subsection. The value of the channel attribute is the observed nucleus. The remaining optional arguments are the metadata that describes the environment under which the spectrum is acquired. In this example, these arguments describe a \(^{29}\text{Si}\) pure anisotropic spinning-sideband spectrum acquired at 9.4 T magnetic flux density and spinning at the magic angle (\(54.735^\circ\)) at 625 Hz. The value of the rotor_frequency argument is the effective anisotropic modulation frequency. For measurements like PASS 5, the anisotropic modulation frequency is the physical rotor frequency. For measurements like the extended chemical shift modulation sequences (XCS) 6, or its variants, where the effective anisotropic modulation frequency is lower than the physical rotor frequency, then it should be set accordingly.

The argument number_of_sidebands is the maximum number of sidebands that will be computed per line-shape within the kernel. For most two-dimensional isotropic vs. pure anisotropic spinning-sideband correlation spectra, the sampling along the sideband dimension is the rotor or the effective anisotropic modulation frequency. Therefore, the number_of_sidebands argument is usually the number of points along the sideband dimension. In this example, this value is 32.

Once the ShieldingPALineshape instance is created, use the kernel() method of the instance to generate the spinning sideband kernel, as follows,

>>> K = lineshapes.kernel(supersampling=1)
>>> print(K.shape)
(32, 625)

Here, K is the \(32\times 625\) kernel, where the 32 is the number of samples (sideband amplitudes), and 625 is the number of features (subspectra) on the \(25 \times 25\) x-y grid. The argument supersampling is the supersampling factor. In a supersampling scheme, each grid cell is averaged over a \(n\times n\) sub-grid, where \(n\) is the supersampling factor. A supersampling factor of 1 is equivalent to no sub-grid averaging.

Data compression (optional)

Often when the kernel, K, is ill-conditioned, the solution becomes unstable in the presence of the measurement noise. An ill-conditioned kernel is the one whose singular values quickly decay to zero. In such cases, we employ the truncated singular value decomposition method to approximately represent the kernel K onto a smaller sub-space, called the range space, where the sub-space kernel is relatively well-defined. We refer to this sub-space kernel as the compressed kernel. Similarly, the measurement data on the sub-space is referred to as the compressed signal. The compression also reduces the time for further computation. To compress the kernel and the data, import the TSVDCompression class and follow,

>>> from mrinversion.linear_model import TSVDCompression
>>> new_system = TSVDCompression(K=K, s=data_object)
compression factor = 1.032258064516129
>>> compressed_K = new_system.compressed_K
>>> compressed_s = new_system.compressed_s

Here, the variable new_system is an instance of the TSVDCompression class. If no truncation index is provided as the argument, when initializing the TSVDCompression class, an optimum truncation index is chosen using the maximum entropy method 7, which is the default behavior. The attributes compressed_K and compressed_s holds the compressed kernel and signal, respectively. The shape of the original signal v.s. the compressed signal is

>>> print(data_object.shape, compressed_s.shape)
(32,) (31,)

Setting up the inverse problem

When setting up the inversion, we solved the smooth LASSO 8 problem. Read the Smooth-LASSO regularization section for further details.

Import the SmoothLasso class and follow,

>>> from mrinversion.linear_model import SmoothLasso
>>> s_lasso = SmoothLasso(alpha=0.01, lambda1=1e-04, inverse_dimension=inverse_dimension)

Here, the variable s_lasso is an instance of the SmoothLasso class. The required arguments of this class are alpha and lambda1, corresponding to the hyperparameters \(\alpha\) and \(\lambda\), respectively, in the Eq. (5). At the moment, we don’t know the optimum value of the alpha and lambda1 parameters. We start with a guess value.

To solve the smooth lasso problem, use the fit() method of the s_lasso instance as follows,

>>> s_lasso.fit(K=compressed_K, s=compressed_s)

The two arguments of the fit() method are the kernel, K, and the signal, s. In the above example, we set the value of K as compressed_K, and correspondingly the value of s as compressed_s. You may also use the uncompressed values of the kernel and signal in this method, if desired.

The solution to the smooth lasso is accessed using the f attribute of the respective object.

>>> f_sol = s_lasso.f

The plot of the solution is

>>> _, ax = plt.subplots(1, 2, figsize=(9, 3.5), subplot_kw={'projection': 'csdm'}) 
>>> plot2D(ax[0], f_sol/f_sol.max(), title='Guess distribution') 
>>> plot2D(ax[1], true_data_object, title='True distribution') 
>>> plt.tight_layout() 
>>> plt.show() 

(png, hires.png, pdf)

_images/getting_started-15.png
_images/null.png

The figure on the left is the guess solution of the nuclear shielding tensor distribution derived from the inversion of the spinning sideband dataset. The figure on the right is the true nuclear shielding tensor distribution.

You may also evaluate the residuals corresponding to the solution using the residuals() method of the object as follows,

>>> residuals = s_lasso.residuals(K=K, s=data_object)
>>> # the plot of the residuals
>>> plt.figure(figsize=(5, 3.5)) 
>>> ax = plt.subplot(projection='csdm') 
>>> ax.plot(residuals, color='black') 
>>> plt.tight_layout() 
>>> plt.show() 

(png, hires.png, pdf)

_images/getting_started-16.png
_images/null.png

The residuals between the 1D MAS sideband spectrum and the predicted spectrum from the guess shielding tensor parameter distribution.

The argument of the residuals method is the kernel and the signal data. We provide the original kernel, K, and signal, s, because we desire the residuals corresponding to the original data and not the compressed data.

Statistical learning of tensor parameters

The solution from a linear model trained with the combined l1 and l2 priors, such as the smooth LASSO estimator used here, depends on the choice of the hyperparameters. The solution shown in the above figure is when \(\alpha=0.01\) and \(\lambda=1\times 10^{-4}\). Although it’s a solution, it is unlikely that this is the best solution. For this, we employ the statistical learning-based model, such as the n-fold cross-validation.

The SmoothLassoCV class is designed to solve the smooth-lasso problem for a range of \(\alpha\) and \(\lambda\) values and determine the best solution using the n-fold cross-validation. Here, we search the best model on a \(10 \times 10\) pre-defined \(\alpha\)-\(\lambda\) grid, using a 10-fold cross-validation statistical learning method. The \(\lambda\) and \(\alpha\) values are sampled uniformly on a logarithmic scale as,

>>> lambdas = 10 ** (-4 - 2 * (np.arange(10) / 9))
>>> alphas = 10 ** (-3 - 2 * (np.arange(10) / 9))
Smooth-LASSO CV Setup

Setup the smooth lasso cross-validation as follows

>>> from mrinversion.linear_model import SmoothLassoCV
>>> s_lasso_cv = SmoothLassoCV(
...     alphas=alphas,
...     lambdas=lambdas,
...     inverse_dimension=inverse_dimension,
...     sigma=0.005,
...     folds=10
... )
>>> s_lasso_cv.fit(K=compressed_K, s=compressed_s)

The arguments of the SmoothLassoCV is a list 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 fit() method, whose arguments are the kernel and signal.

The optimum hyperparameters

The optimized hyperparameters may be accessed using the hyperparameters attribute of the class instance,

>>> alpha = s_lasso_cv.hyperparameters['alpha']
>>> lambda_1 = s_lasso_cv.hyperparameters['lambda']
The cross-validation surface

The cross-validation error metric is the mean square error metric. You may access this data using the cross_validation_curve attribute.

>>> plt.figure(figsize=(5, 3.5)) 
>>> ax = plt.subplot(projection='csdm') 
>>> ax.contour(np.log10(s_lasso_cv.cross_validation_curve), levels=25) 
>>> ax.scatter(-np.log10(s_lasso_cv.hyperparameters['alpha']),
...         -np.log10(s_lasso_cv.hyperparameters['lambda']),
...         marker='x', color='k') 
>>> plt.tight_layout() 
>>> plt.show() 

(png, hires.png, pdf)

_images/getting_started-20.png
_images/null.png

The ten-folds cross-validation prediction error surface as a function of the hyperparameters \(\alpha\) and \(\beta\).

The optimum solution

The best model selection from the cross-validation method may be accessed using the f attribute.

>>> f_sol_cv = s_lasso_cv.f  # best model selected using the 10-fold cross-validation

The plot of the selected tensor parameter distribution is shown below.

>>> _, ax = plt.subplots(1, 2, figsize=(9, 3.5), subplot_kw={'projection': 'csdm'}) 
>>> plot2D(ax[0], f_sol_cv/f_sol_cv.max(), title='Optimum distribution') 
>>> plot2D(ax[1], true_data_object, title='True distribution') 
>>> plt.tight_layout() 
>>> plt.show() 

(png, hires.png, pdf)

_images/getting_started-22.png
_images/null.png

The figure on the left is the optimum solution selected by the 10-folds cross-validation method. The figure on the right is the true model of the nuclear shielding tensor distribution.

1

Srivastava, D. J., Vosegaard, T., Massiot, D., Grandinetti, P. J., Core Scientific Dataset Model: A lightweight and portable model and file format for multi-dimensional scientific data. PLOS ONE, 15, 1-38, (2020). DOI:10.1371/journal.pone.0225953

2

Bak M., Rasmussen J. T., Nielsen N.C., SIMPSON: A General Simulation Program for Solid-State NMR Spectroscopy. J Magn Reson. 147, 296–330, (2000). DOI:10.1006/jmre.2000.2179

3

Massiot D., Fayon F., Capron M., King I., Le Calvé S., Alonso B., et al. Modelling one- and two-dimensional solid-state NMR spectra. Magn Reson Chem. 40, 70–76, (2002) DOI:10.1002/mrc.984

4

PhySy Ltd. RMN 2.0; 2019. Available from: https://www.physyapps.com/rmn.

5

Dixon, W. T., Spinning sideband free and spinning sideband only NMR spectra in spinning samples. J. Chem. Phys, 77, 1800, (1982). DOI:10.1063/1.444076

6

Gullion, T., Extended chemical shift modulation. J. Mag. Res., 85, 3, (1989). DOI:10.1016/0022-2364(89)90253-9

7

Varshavsky R., Gottlieb A., Linial M., Horn D., Novel unsupervised feature filtering of biological data. Bioinformatics, 22, e507–e513, (2006). DOI:10.1093/bioinformatics/btl214.

8

Hebiri M, Sara A. Van De Geer, The Smooth-Lasso and other l1+l2-penalized methods, arXiv, (2010). arXiv:1003.4885v2

API-Reference

Pure anisotropic Nuclear Shielding Kernel

Generalized Class
class mrinversion.kernel.nmr.ShieldingPALineshape(anisotropic_dimension, inverse_dimension, channel, magnetic_flux_density='9.4 T', rotor_angle='54.735 deg', rotor_frequency='14 kHz', number_of_sidebands=1)[source]

Bases: mrinversion.kernel.base.LineShape

A generalized class for simulating the pure anisotropic NMR nuclear shielding line-shape kernel.

Parameters
  • anisotropic_dimension – A Dimension object, or an equivalent dictionary object. This dimension must represent the pure anisotropic dimension.

  • inverse_dimension – A list of two Dimension objects, or equivalent dictionary objects representing the x-y coordinate grid.

  • channel – The channel is an isotope symbol of the nuclei given as the atomic number followed by the atomic symbol, for example, 1H, 13C, and 29Si. This nucleus must correspond to the recorded frequency resonances.

  • magnetic_flux_density – The magnetic flux density of the external static magnetic field. The default value is 9.4 T.

  • rotor_angle – The angle of the sample holder (rotor) relative to the direction of the external magnetic field. The default value is 54.735 deg (magic angle).

  • rotor_frequency – The effective sample spin rate. Depending on the NMR sequence, this value may be less than the physical sample rotation frequency. The default is 14 kHz.

  • number_of_sidebands – The number of sidebands to simulate along the anisotropic dimension. The default value is 1.

kernel(supersampling=1)[source]

Return the NMR nuclear shielding anisotropic line-shape kernel.

Parameters

supersampling – An integer. Each cell is supersampled by the factor supersampling along every dimension.

Returns

A numpy array containing the line-shape kernel.

Specialized Classes
Magic Angle Flipping
class mrinversion.kernel.nmr.MAF(anisotropic_dimension, inverse_dimension, channel, magnetic_flux_density='9.4 T')[source]

Bases: mrinversion.kernel.csa_aniso.ShieldingPALineshape

A specialized class for simulating the pure anisotropic NMR nuclear shielding line-shape kernel resulting from the 2D MAF spectra.

Parameters
  • anisotropic_dimension – A Dimension object, or an equivalent dictionary object. This dimension must represent the pure anisotropic dimension.

  • inverse_dimension – A list of two Dimension objects, or equivalent dictionary objects representing the x-y coordinate grid.

  • channel – The isotope symbol of the nuclei given as the atomic number followed by the atomic symbol, for example, 1H, 13C, and 29Si. This nucleus must correspond to the recorded frequency resonances.

  • magnetic_flux_density – The magnetic flux density of the external static magnetic field. The default value is 9.4 T.

Assumptions: The simulated line-shapes correspond to an infinite speed spectrum spinning at \(90^\circ\).

kernel(supersampling=1)

Return the NMR nuclear shielding anisotropic line-shape kernel.

Parameters

supersampling – An integer. Each cell is supersampled by the factor supersampling along every dimension.

Returns

A numpy array containing the line-shape kernel.

Spinning Sidebands
class mrinversion.kernel.nmr.SpinningSidebands(anisotropic_dimension, inverse_dimension, channel, magnetic_flux_density='9.4 T')[source]

Bases: mrinversion.kernel.csa_aniso.ShieldingPALineshape

A specialized class for simulating the pure anisotropic spinning sideband amplitudes of the nuclear shielding resonances resulting from a 2D sideband separation spectra.

Parameters
  • anisotropic_dimension – A Dimension object, or an equivalent dictionary object. This dimension must represent the pure anisotropic dimension.

  • inverse_dimension – A list of two Dimension objects, or equivalent dictionary objects representing the x-y coordinate grid.

  • channel – The isotope symbol of the nuclei given as the atomic number followed by the atomic symbol, for example, 1H, 13C, and 29Si. This nucleus must correspond to the recorded frequency resonances.

  • magnetic_flux_density – The magnetic flux density of the external static magnetic field. The default value is 9.4 T.

Assumption: The simulated line-shapes correspond to a finite speed spectrum spinning at the magic angle, \(54.735^\circ\), where the spin rate is the increment along the anisotropic dimension.

kernel(supersampling=1)

Return the NMR nuclear shielding anisotropic line-shape kernel.

Parameters

supersampling – An integer. Each cell is supersampled by the factor supersampling along every dimension.

Returns

A numpy array containing the line-shape kernel.

Smooth Lasso

class mrinversion.linear_model.SmoothLasso(alpha, lambda1, inverse_dimension, max_iterations=10000, tolerance=1e-05, positive=True, method='gradient_decent')[source]

Bases: mrinversion.linear_model._base_l1l2.GeneralL2Lasso

The linear model trained with the combined l1 and l2 priors as the regularizer. The method minimizes the objective function,

()\[\| {\bf Kf - s} \|^2_2 + \alpha \sum_{i=1}^{d} \| {\bf J}_i {\bf f} \|_2^2 + \lambda \| {\bf f} \|_1 ,\]

where \({\bf K} \in \mathbb{R}^{m \times n}\) is the kernel, \({\bf s} \in \mathbb{R}^{m \times m_\text{count}}\) is the known (measured) signal, and \({\bf f} \in \mathbb{R}^{n \times m_\text{count}}\) is the desired solution. The parameters, \(\alpha\) and \(\lambda\), are the hyperparameters controlling the smoothness and sparsity of the solution \({\bf f}\). The matrix \({\bf J}_i\) is given as

()\[{\bf J}_i = {\bf I}_{n_1} \otimes \cdots \otimes {\bf A}_{n_i} \otimes \cdots \otimes {\bf I}_{n_{d}},\]

where \({\bf I}_{n_i} \in \mathbb{R}^{n_i \times n_i}\) is the identity matrix,

()\[\begin{split}{\bf A}_{n_i} = \left(\begin{array}{ccccc} 1 & -1 & 0 & \cdots & \vdots \\ 0 & 1 & -1 & \cdots & \vdots \\ \vdots & \vdots & \vdots & \vdots & 0 \\ 0 & \cdots & 0 & 1 & -1 \end{array}\right) \in \mathbb{R}^{(n_i-1)\times n_i},\end{split}\]

and the symbol \(\otimes\) is the Kronecker product. The terms, \(\left(n_1, n_2, \cdots, n_d\right)\), are the number of points along the respective dimensions, with the constraint that \(\prod_{i=1}^{d}n_i = n\), where \(d\) is the total number of dimensions.

Parameters
  • alpha (float) – The hyperparameter, \(\alpha\).

  • lambda1 (float) – The hyperparameter, \(\lambda\).

  • inverse_dimension (list) – A list of csdmpy Dimension objects representing the inverse space.

  • max_iterations (int) – The maximum number of iterations allowed when solving the problem. The default value is 10000.

  • tolerance (float) – The tolerance at which the solution is considered converged. The default value is 1e-5.

  • positive (bool) – If True, the amplitudes in the solution, \({\bf f}\), is contrained to only positive values, else the solution may contain positive and negative amplitudes. The default is True.

f

A ndarray of shape (m_count, nd, …, n1, n0) representing the solution, \({\bf f} \in \mathbb{R}^{m_\text{count} \times n_d \times \cdots n_1 \times n_0}\).

Type

ndarray or CSDM object.

n_iter

The number of iterations required to reach the specified tolerance.

Type

int

Methods Documentation

fit(K, s)

Fit the model using the coordinate descent method from scikit-learn.

Parameters
  • K (ndarray) – The \(m \times n\) kernel matrix, \({\bf K}\). A numpy array of shape (m, n).

  • s (ndarray or CSDM object.) – A csdm object or an equivalent numpy array holding the signal, \({\bf s}\), as a \(m \times m_\text{count}\) matrix.

predict(K)

Predict the signal using the linear model.

Parameters

K (ndarray) – A \(m \times n\) kernel matrix, \({\bf K}\). A numpy array of shape (m, n).

Returns

A numpy array of shape (m, m_count) with the predicted values

Return type

ndarray

residuals(K, s)

Return the residual as the difference the data and the prediced data(fit), following

()\[\text{residuals} = {\bf s - Kf^*}\]

where \({\bf f^*}\) is the optimum solution.

Parameters
  • K (ndarray.) – A \(m \times n\) kernel matrix, \({\bf K}\). A numpy array of shape (m, n).

  • s (ndarray ot CSDM object.) – A csdm object or a \(m \times m_\text{count}\) signal matrix, \({\bf s}\).

Returns

If s is a csdm object, returns a csdm object with the residuals. If s is a numpy array, return a \(m \times m_\text{count}\) residue matrix. csdm object

Return type

ndarray or CSDM object.

score(K, s, sample_weights=None)

The coefficient of determination, \(R^2\), of the prediction. For more information, read scikit-learn documentation.

Smooth Lasso cross-validation

class mrinversion.linear_model.SmoothLassoCV(alphas, lambdas, inverse_dimension, folds=10, max_iterations=10000, tolerance=1e-05, positive=True, sigma=0.0, randomize=False, times=2, verbose=False, n_jobs=- 1, method='gradient_decent')[source]

Bases: mrinversion.linear_model._base_l1l2.GeneralL2LassoCV

The linear model trained with the combined l1 and l2 priors as the regularizer. The method minimizes the objective function,

()\[\| {\bf Kf - s} \|^2_2 + \alpha \sum_{i=1}^{d} \| {\bf J}_i {\bf f} \|_2^2 + \lambda \| {\bf f} \|_1 ,\]

where \({\bf K} \in \mathbb{R}^{m \times n}\) is the kernel, \({\bf s} \in \mathbb{R}^{m \times m_\text{count}}\) is the known signal containing noise, and \({\bf f} \in \mathbb{R}^{n \times m_\text{count}}\) is the desired solution. The parameters, \(\alpha\) and \(\lambda\), are the hyperparameters controlling the smoothness and sparsity of the solution \({\bf f}\). The matrix \({\bf J}_i\) is given as

()\[{\bf J}_i = {\bf I}_{n_1} \otimes \cdots \otimes {\bf A}_{n_i} \otimes \cdots \otimes {\bf I}_{n_{d}},\]

where \({\bf I}_{n_i} \in \mathbb{R}^{n_i \times n_i}\) is the identity matrix,

()\[\begin{split}{\bf A}_{n_i} = \left(\begin{array}{ccccc} 1 & -1 & 0 & \cdots & \vdots \\ 0 & 1 & -1 & \cdots & \vdots \\ \vdots & \vdots & \vdots & \vdots & 0 \\ 0 & \cdots & 0 & 1 & -1 \end{array}\right) \in \mathbb{R}^{(n_i-1)\times n_i},\end{split}\]

and the symbol \(\otimes\) is the Kronecker product. The terms, \(\left(n_1, n_2, \cdots, n_d\right)\), are the number of points along the respective dimensions, with the constraint that \(\prod_{i=1}^{d}n_i = n\), where \(d\) is the total number of dimensions.

The cross-validation is carried out using a stratified splitting of the signal.

Parameters
  • alphas (ndarray) – A list of \(\alpha\) hyperparameters.

  • lambdas (ndarray) – A list of \(\lambda\) hyperparameters.

  • inverse_dimension (list) – A list of csdmpy Dimension objects representing the inverse space.

  • folds (int) – The number of folds used in cross-validation.The default is 10.

  • max_iterations (int) – The maximum number of iterations allowed when solving the problem. The default value is 10000.

  • tolerance (float) – The tolerance at which the solution is considered converged. The default value is 1e-5.

  • positive (bool) – If True, the amplitudes in the solution, \({\bf f}\), is contrained to only positive values, else the solution may contain positive and negative amplitudes. The default is True.

  • sigma (float) – The standard deviation of the noise in the signal. The default is 0.0.

  • sigma – The standard deviation of the noise in the signal. The default is 0.0.

  • randomize (bool) – If true, the folds are created by randomly assigning the samples to each fold. If false, a stratified sampled is used to generate folds. The default is False.

  • times (int) – The number of times to randomized n-folds are created. Only applicable when randomize attribute is True.

  • verbose (bool) – If true, prints the process.

  • n_jobs (int) – The number of CPUs used for computation. The default is -1, that is, all available CPUs are used.

f

A ndarray of shape (m_count, nd, …, n1, n0). The solution, \({\bf f} \in \mathbb{R}^{m_\text{count} \times n_d \times \cdots n_1 \times n_0}\) or an equivalent CSDM object.

Type

ndarray or CSDM object.

n_iter

The number of iterations required to reach the specified tolerance.

Type

int.

hyperparameters

A dictionary with the \(\alpha\) and :math:lambda` hyperparameters.

Type

dict.

cross_validation_curve

The cross-validation error metric determined as the mean square error.

Type

CSDM object.

Methods Documentation

fit(K, s)

Fit the model using the coordinate descent method from scikit-learn for all alpha anf lambda values using the n-folds cross-validation technique. The cross-validation metric is the mean squared error.

Parameters
  • K – A \(m \times n\) kernel matrix, \({\bf K}\). A numpy array of shape (m, n).

  • s – A \(m \times m_\text{count}\) signal matrix, \({\bf s}\) as a csdm object or a numpy array or shape (m, m_count).

predict(K)

Predict the signal using the linear model.

Parameters

K – A \(m \times n\) kernel matrix, \({\bf K}\). A numpy array of shape (m, n).

Returns

A numpy array of shape (m, m_count) with the predicted values.

residuals(K, s)

Return the residual as the difference the data and the prediced data(fit), following

()\[\text{residuals} = {\bf s - Kf^*}\]

where \({\bf f^*}\) is the optimum solution.

Parameters
  • K – A \(m \times n\) kernel matrix, \({\bf K}\). A numpy array of shape (m, n).

  • s – A csdm object or a \(m \times m_\text{count}\) signal matrix, \({\bf s}\).

Returns

If s is a csdm object, returns a csdm object with the residuals. If s is a numpy array, return a \(m \times m_\text{count}\) residue matrix.

TSVDCompression

class mrinversion.linear_model.TSVDCompression(K, s, r=None)[source]

Bases: object

SVD compression.

Parameters
  • K – The kernel.

  • s – The data.

  • r – The number of singular values used in data compression.

truncation_index

The number of singular values retained.

Type

int

compressed_K

The compressed kernel.

Type

ndarray

compressed_s

The compressed data.

Type

ndarray of CSDM object

Utils

mrinversion.kernel.utils.x_y_to_zeta_eta(x, y)[source]

Convert the coordinates \((x,y)\) to \((\zeta, \eta)\) using the following definition,

()\[\begin{split}\left.\begin{array}{rl} \zeta &= \sqrt{x^2 + y^2}, \\ \eta &= \frac{4}{\pi} \tan^{-1} \left| \frac{x}{y} \right| \end{array} {~~~~~~~~} \right\} {~~~~~~~~} |x| \le |y|.\end{split}\]
()\[\begin{split}\left.\begin{array}{rl} \zeta &= -\sqrt{x^2 + y^2}, \\ \eta &= \frac{4}{\pi} \tan^{-1} \left| \frac{y}{x} \right| \end{array} {~~~~~~~~} \right\} {~~~~~~~~} |x| > |y|.\end{split}\]
Parameters
  • x – floats or Quantity object. The coordinate x.

  • y – floats or Quantity object. The coordinate y.

Returns

A list of two ndarrays. The first array is the \(\zeta\) coordinates. The second array is the \(\eta\) coordinates.

mrinversion.kernel.utils.zeta_eta_to_x_y(zeta, eta)[source]

Convert the coordinates \((\zeta,\eta)\) to \((x, y)\) using the following definition,

()\[\begin{split}\left. \begin{array}{rl} x &= |\zeta| \sin\theta, \\ y &= |\zeta| \cos\theta \end{array} {~~~~~~~~} \right\} {~~~~~~~~} \zeta \ge 0\end{split}\]
()\[\begin{split}\left. \begin{array}{rl} x &= |\zeta| \cos\theta, \\ y &= |\zeta| \sin\theta \end{array} {~~~~~~~~} \right\} {~~~~~~~~} \zeta < 0\end{split}\]

where \(\theta = \frac{\pi}{4}\eta\).

Parameters
  • x – ndarray or list of floats. The coordinate x.

  • y – ndarray or list of floats. The coordinate y.

Returns

A list of ndarrays. The first array holds the coordinate \(x\). The second array holds the coordinates \(y\).

mrinversion.utils.get_polar_grids(ax, ticks=None, offset=0)[source]

Generate a piece-wise polar grid of Haeberlen parameters, zeta and eta.

Parameters
  • ax – Matplotlib Axes.

  • ticks – Tick coordinates where radial grids are drawn. The value can be a list or a numpy array. The default value is None.

  • offset – The grid is drawn at an offset away from the origin.

mrinversion.utils.to_Haeberlen_grid(csdm_object, zeta, eta, n=5)[source]

Convert the three-dimensional p(iso, x, y) to p(iso, zeta, eta) tensor distribution.

Parameters
  • csdm_object (CSDM) – A CSDM object containing the 3D p(iso, x, y) distribution.

  • zeta (CSDM.Dimension) – A CSDM dimension object describing the zeta dimension.

  • eta (CSDM.Dimension) – A CSDM dimension object describing the eta dimension.

  • n (int) – An interger used in linear interpolation of the data. The default is 5.

mrinversion.utils.plot_3d(ax, csdm_object, elev=28, azim=-150, x_lim=None, y_lim=None, z_lim=None, max_2d=None, max_1d=None, cmap=<matplotlib.colors.LinearSegmentedColormap object>, box=False, clip_percent=0.0, linewidth=1, alpha=0.15, **kwargs)[source]

Generate a 3D density plot with 2D contour and 1D projections.

Parameters
  • ax – Matplotlib Axes to render the plot.

  • csdm_object – A 3D{1} CSDM object holding the data.

  • elev – (optional) The 3D view angle, elevation angle in the z plane.

  • azim – (optional) The 3D view angle, azimuth angle in the x-y plane.

  • x_lim – (optional) The x limit given as a list, [x_min, x_max].

  • y_lim – (optional) The y limit given as a list, [y_min, y_max].

  • z_lim – (optional) The z limit given as a list, [z_min, z_max].

  • max_2d – (Optional) The normalization factor of the 2D contour projections. The attribute is meaningful when multiple 3D datasets are viewed on the same plot. The value is given as a list, [yz, xz, xy], where ij is the maximum of the projection onto the ij plane, \(i,j \in [x, y, z]\).

  • max_1d – (Optional) The normalization factor of the 1D projections. The attribute is meaningful when multiple 3D datasets are viewed on the same plot. The value is given as a list, [x, y, z], where i is the maximum of the projection onto the i axis, \(i \in [x, y, z]\).

  • cmap – (Optional) The colormap used in rendering the volumetric plot. The same colormap is used for the 2D contour projections. For 1D plots, the first color in the colormap scheme is used for the line color.

  • box – (Optional) If True, draw a box around the 3D data region.

  • clip_percent – (Optional) The amplitudes of the dataset below the given percent is made transparent for the volumetric plot.

  • linewidth – (Optional) The linewidth of the 2D countours, 1D plots and box.

  • alpha – (Optional) The amount of alpha(transparency) applied in rendering the 3D volume.

Examples

Project details

Changelog

v0.2.0

What’s new!
  • Added to_Haeberlen_grid() function to convert the 3D \(\rho(\delta_\text{iso}, x, y)\) distribution to \(\rho(\delta_\text{iso}, \zeta_\sigma, \eta_\sigma)\) distribution.

Changes
  • Update code to comply with updated csdmpy library.

License

Mrinversion License

Mrinversion is licensed under BSD 3-Clause License.

Copyright (c) 2020, Deepansh J. Srivastava,

All rights reserved.

Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:

  • Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.

  • Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.

  • Neither the name of the copyright holder nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.

THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS “AS IS” AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

Acknowledgment

The development of the mrinversion library is supported in part by the US National Science Foundation under Grant No. DIBBS OAC 1640899 and Grant No. CHE 1807922.

How to cite

If you use this work in your publication, please cite the following.

  • Srivastava, D. J.; Grandinetti P. J., Statistical learning of NMR tensors from 2D isotropic/anisotropic correlation nuclear magnetic resonance spectra, J. Chem. Phys. 153, 134201 (2020). https://doi.org/10.1063/5.0023345.

  • Deepansh J. Srivastava, Maxwell Venetos, Philip J. Grandinetti, Shyam Dwaraknath, & Alexis McCarthy. (2021, May 26). mrsimulator: v0.6.0 (Version v0.6.0). Zenodo. http://doi.org/10.5281/zenodo.4814638

Additionally, if you use the CSDM data model, please consider citing

  • Srivastava DJ, Vosegaard T, Massiot D, Grandinetti PJ (2020) Core Scientific Dataset Model: A lightweight and portable model and file format for multi-dimensional scientific data. PLOS ONE 15(1): e0225953. https://doi.org/10.1371/journal.pone.0225953


Indices and tables