¹¹B MAS NMR of Lithium orthoborate crystal

The following is a quadrupolar lineshape fitting example for the 11B MAS NMR of lithium orthoborate crystal. The dataset was shared by Dr. Nathan Barrow.

import csdmpy as cp
import matplotlib.pyplot as plt
from lmfit import Minimizer

from mrsimulator import Simulator, Site, SpinSystem
from mrsimulator.method.lib import BlochDecayCTSpectrum
from mrsimulator import signal_processor as sp
from mrsimulator.utils import spectral_fitting as sf
from mrsimulator.utils import get_spectral_dimensions
from mrsimulator.spin_system.tensors import SymmetricTensor

Import the dataset

host = "https://ssnmr.org/sites/default/files/mrsimulator/"
filename = "11B_lithum_orthoborate.csdf"
experiment = cp.load(host + filename)

# standard deviation of noise from the dataset
sigma = 0.08078374

# For spectral fitting, we only focus on the real part of the complex dataset
experiment = experiment.real

# Convert the coordinates along each dimension from Hz to ppm.
_ = [item.to("ppm", "nmr_frequency_ratio") for item in experiment.dimensions]

# plot of the dataset.
plt.figure(figsize=(4.25, 3.0))
ax = plt.subplot(projection="csdm")
ax.plot(experiment, "k", alpha=0.5)
ax.set_xlim(100, -100)
plt.grid()
plt.tight_layout()
plt.show()
plot 4 11B Quad NMR

Create a fitting model

Spin System

B11 = Site(
    isotope="11B",
    isotropic_chemical_shift=20.0,  # in ppm
    quadrupolar=SymmetricTensor(Cq=2.3e6, eta=0.03),  # Cq in Hz
)
spin_systems = [SpinSystem(sites=[B11])]

Method

# Get the spectral dimension parameters from the experiment.
spectral_dims = get_spectral_dimensions(experiment)

MAS_CT = BlochDecayCTSpectrum(
    channels=["11B"],
    magnetic_flux_density=14.1,  # in T
    rotor_frequency=12500,  # in Hz
    spectral_dimensions=spectral_dims,
    experiment=experiment,  # add the measurement to the method.
)

# Optimize the script by pre-setting the transition pathways for each spin system from
# the method.
for sys in spin_systems:
    sys.transition_pathways = MAS_CT.get_transition_pathways(sys)

Guess Model Spectrum

# Simulation
# ----------
sim = Simulator(spin_systems=spin_systems, methods=[MAS_CT])
sim.run()

# Post Simulation Processing
# --------------------------
processor = sp.SignalProcessor(
    operations=[
        sp.IFFT(),
        sp.apodization.Exponential(FWHM="100 Hz"),
        sp.FFT(),
        sp.Scale(factor=200),
    ]
)
processed_dataset = processor.apply_operations(dataset=sim.methods[0].simulation).real

# Plot of the guess Spectrum
# --------------------------
plt.figure(figsize=(4.25, 3.0))
ax = plt.subplot(projection="csdm")
ax.plot(experiment, "k", linewidth=1, label="Experiment")
ax.plot(processed_dataset, "r", alpha=0.75, linewidth=1, label="guess spectrum")
ax.set_xlim(100, -100)
plt.grid()
plt.legend()
plt.tight_layout()
plt.show()
plot 4 11B Quad NMR

Least-squares minimization with LMFIT

Use the make_LMFIT_params() for a quick setup of the fitting parameters.

params = sf.make_LMFIT_params(sim, processor)
params.pop("sys_0_abundance")
print(params.pretty_print(columns=["value", "min", "max", "vary", "expr"]))

Out:

Name                                      Value      Min      Max     Vary     Expr
SP_0_operation_1_Exponential_FWHM           100     -inf      inf     True     None
SP_0_operation_3_Scale_factor               200     -inf      inf     True     None
sys_0_site_0_isotropic_chemical_shift        20     -inf      inf     True     None
sys_0_site_0_quadrupolar_Cq             2.3e+06     -inf      inf     True     None
sys_0_site_0_quadrupolar_eta               0.03        0        1     True     None
None

Solve the minimizer using LMFIT

Fit Statistics

fitting methodleastsq
# function evals119
# data points4096
# variables5
chi-square 203470.421
reduced chi-square 49.7361087
Akaike info crit. 16006.9680
Bayesian info crit. 16038.5568

Variables

name value standard error relative error initial value min max vary
sys_0_site_0_isotropic_chemical_shift 20.0533817 0.00147191 (0.01%) 20.0 -inf inf True
sys_0_site_0_quadrupolar_Cq 2731196.61 417.270375 (0.02%) 2300000.0 -inf inf True
sys_0_site_0_quadrupolar_eta 0.04712757 7.3037e-04 (1.55%) 0.03 0.00000000 1.00000000 True
SP_0_operation_1_Exponential_FWHM 70.9385144 0.52813486 (0.74%) 100.0 -inf inf True
SP_0_operation_3_Scale_factor 199.362369 0.20570628 (0.10%) 200.0 -inf inf True

Correlations (unreported correlations are < 0.100)

sys_0_site_0_isotropic_chemical_shiftsys_0_site_0_quadrupolar_Cq0.7932
sys_0_site_0_quadrupolar_etaSP_0_operation_1_Exponential_FWHM-0.6610
SP_0_operation_1_Exponential_FWHMSP_0_operation_3_Scale_factor0.5298
sys_0_site_0_quadrupolar_etaSP_0_operation_3_Scale_factor-0.2422
sys_0_site_0_quadrupolar_Cqsys_0_site_0_quadrupolar_eta-0.1331
sys_0_site_0_quadrupolar_CqSP_0_operation_3_Scale_factor0.1264


The best fit solution

best_fit = sf.bestfit(sim, processor)[0].real
residuals = sf.residuals(sim, processor)[0].real

# Plot the spectrum
plt.figure(figsize=(4.25, 3.0))
ax = plt.subplot(projection="csdm")
ax.plot(experiment, "k", linewidth=1, label="Experiment")
ax.plot(best_fit, "r", alpha=0.75, linewidth=1, label="Best Fit")
ax.plot(residuals, alpha=0.75, linewidth=1, label="Residuals")
ax.set_xlim(100, -100)
plt.grid()
plt.legend()
plt.tight_layout()
plt.show()
plot 4 11B Quad NMR

Total running time of the script: ( 0 minutes 2.450 seconds)

Gallery generated by Sphinx-Gallery