.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "fitting/1D_fitting/plot_2_13C_glycine_multi_spectra_fit.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note Click :ref:`here ` to download the full example code .. rst-class:: sphx-glr-example-title .. _sphx_glr_fitting_1D_fitting_plot_2_13C_glycine_multi_spectra_fit.py: ¹³C MAS NMR of Glycine (CSA) multi-spectra fit ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ .. GENERATED FROM PYTHON SOURCE LINES 7-11 The following is a multi-dataset least-squares fitting example of :math:`^{13}\text{C}` MAS NMR spectrum of Glycine spinning at 5 kHz, 1.94 kHz, and 960 Hz. Before trying multi-dataset fitting, we recommend that you first try individual fits. The experimental datasets are part of DMFIT [#f1]_ examples. .. GENERATED FROM PYTHON SOURCE LINES 11-23 .. code-block:: default import csdmpy as cp import matplotlib.pyplot as plt from lmfit import Minimizer from mrsimulator import Simulator, SpinSystem, Site from mrsimulator.method.lib import BlochDecaySpectrum 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 .. GENERATED FROM PYTHON SOURCE LINES 25-30 Import the datasets ------------------- Import the datasets and assign the standard deviation of noise for each dataset. Here, ``sigma1``, ``sigma2``, and ``sigma3`` are the noise standard deviation for the dataset acquired at 5 kHz, 1.94 kHz, and 960 Hz spinning speeds, respectively. .. GENERATED FROM PYTHON SOURCE LINES 30-58 .. code-block:: default host = "https://nmr.cemhti.cnrs-orleans.fr/Dmfit/Help/csdm/" filename1 = "13C MAS 5000Hz - Glycine.csdf" filename2 = "13C MAS 1940Hz - Glycine.csdf" filename3 = "13C MAS 960Hz - Glycine.csdf" experiment1 = cp.load(host + filename1).real experiment2 = cp.load(host + filename2).real experiment3 = cp.load(host + filename3).real experiments = [experiment1, experiment2, experiment3] # standard deviation of noise from the dataset sigma1 = 1.97637 sigma2 = 3.822249 sigma3 = 3.982936 sigmas = [sigma1, sigma2, sigma3] # Convert the coordinates along each dimension from Hz to ppm for each dataset fig, ax = plt.subplots(1, 3, figsize=(12, 3), subplot_kw={"projection": "csdm"}) for i, experiment in enumerate(experiments): _ = [item.to("ppm", "nmr_frequency_ratio") for item in experiment.dimensions] # plot of the dataset. ax[i].plot(experiment, color="black", linewidth=0.5, label="Experiment") ax[i].set_xlim(280, -10) ax[i].grid() plt.tight_layout() plt.show() .. image-sg:: /fitting/1D_fitting/images/sphx_glr_plot_2_13C_glycine_multi_spectra_fit_001.png :alt: plot 2 13C glycine multi spectra fit :srcset: /fitting/1D_fitting/images/sphx_glr_plot_2_13C_glycine_multi_spectra_fit_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 59-64 Create a fitting model ---------------------- **Spin System**: The objective of a multi-dataset fitting is to optimize the spin system parameters using multiple datasets. In this example, we create two single-site spin systems, which are then shared by three method objects. .. GENERATED FROM PYTHON SOURCE LINES 64-77 .. code-block:: default C1 = Site( isotope="13C", isotropic_chemical_shift=176.0, # in ppm shielding_symmetric=SymmetricTensor(zeta=60, eta=0.6), # zeta in Hz ) C2 = Site( isotope="13C", isotropic_chemical_shift=43.0, # in ppm shielding_symmetric=SymmetricTensor(zeta=30, eta=0.5), # zeta in Hz ) spin_systems = [SpinSystem(sites=[C1], name="C1"), SpinSystem(sites=[C2], name="C2")] .. GENERATED FROM PYTHON SOURCE LINES 78-79 **Method**: Create the three MAS method objects with respective MAS spinning speeds. .. GENERATED FROM PYTHON SOURCE LINES 79-113 .. code-block:: default # Get the spectral dimension parameters from the respective experiment and setup the # corresponding method. # Method for dataset 1 spectral_dims1 = get_spectral_dimensions(experiment1) MAS1 = BlochDecaySpectrum( channels=["13C"], magnetic_flux_density=7.05, # in T rotor_frequency=5000, # in Hz spectral_dimensions=spectral_dims1, experiment=experiment1, # add experimental dataset 1 ) # Method for dataset 2 spectral_dims2 = get_spectral_dimensions(experiment2) MAS2 = BlochDecaySpectrum( channels=["13C"], magnetic_flux_density=7.05, # in T rotor_frequency=1940, # in Hz spectral_dimensions=spectral_dims2, experiment=experiment2, # add experimental dataset 2 ) # Method for dataset 3 spectral_dims3 = get_spectral_dimensions(experiment3) MAS3 = BlochDecaySpectrum( channels=["13C"], magnetic_flux_density=7.05, # in T rotor_frequency=960, # in Hz spectral_dimensions=spectral_dims3, experiment=experiment3, # add experimental dataset 3 ) .. GENERATED FROM PYTHON SOURCE LINES 114-115 **Guess Model Spectrum** .. GENERATED FROM PYTHON SOURCE LINES 115-183 .. code-block:: default # Simulation # ---------- # Add the spin systems and the three methods to the simulator object. sim = Simulator(spin_systems=spin_systems, methods=[MAS1, MAS2, MAS3]) sim.config.decompose_spectrum = "spin_system" sim.run() # Post Simulation Processing # -------------------------- # Add signal processing to simulation dataset from the three methods. # Processor for dataset 1 processor1 = sp.SignalProcessor( operations=[ sp.IFFT(), sp.apodization.Exponential(FWHM="20 Hz", dv_index=0), # spin system 0 sp.apodization.Exponential(FWHM="200 Hz", dv_index=1), # spin system 1 sp.FFT(), sp.Scale(factor=10), # dataset is scaled independently using scale factor. ] ) # Processor for dataset 2 processor2 = sp.SignalProcessor( operations=[ sp.IFFT(), sp.apodization.Exponential(FWHM="30 Hz", dv_index=0), # spin system 0 sp.apodization.Exponential(FWHM="300 Hz", dv_index=1), # spin system 1 sp.FFT(), sp.Scale(factor=100), # dataset is scaled independently using scale factor. ] ) # Processor for dataset 3 processor3 = sp.SignalProcessor( operations=[ sp.IFFT(), sp.apodization.Exponential(FWHM="10 Hz", dv_index=0), # spin system 0 sp.apodization.Exponential(FWHM="150 Hz", dv_index=1), # spin system 1 sp.FFT(), sp.Scale(factor=50), # dataset is scaled independently using scale factor. ] ) processors = [processor1, processor2, processor3] processed_dataset = [] for i, proc in enumerate(processors): processed_dataset.append( proc.apply_operations(dataset=sim.methods[i].simulation).real ) # Plot of the guess Spectrum # -------------------------- fig, ax = plt.subplots(1, 3, figsize=(12, 3), subplot_kw={"projection": "csdm"}) for i, exp_dataset in enumerate(experiments): ax[i].plot(exp_dataset, color="black", linewidth=0.5, label="Experiment") ax[i].plot(processed_dataset[i], linewidth=2, alpha=0.6) ax[i].set_xlim(280, -10) ax[i].grid() plt.legend() plt.tight_layout() plt.show() .. image-sg:: /fitting/1D_fitting/images/sphx_glr_plot_2_13C_glycine_multi_spectra_fit_002.png :alt: plot 2 13C glycine multi spectra fit :srcset: /fitting/1D_fitting/images/sphx_glr_plot_2_13C_glycine_multi_spectra_fit_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 184-194 Least-squares minimization with LMFIT ------------------------------------- Use the :func:`~mrsimulator.utils.spectral_fitting.make_LMFIT_params` for a quick setup of the fitting parameters. Note, the first two arguments of this function is the simulator object and a list of SignalProcessor objects, ``processors``. The fitting parameters corresponding to the signal processor objects are generated using ``SP_i_operation_j_FunctionName_FunctionArg``, where *i* is the *ith* signal processor within the list, *j* is the operation index of the *ith* processor, and *FunctionName* and *FunctionArg* are the operation function name and function argument, respectively. .. GENERATED FROM PYTHON SOURCE LINES 194-197 .. code-block:: default params = sf.make_LMFIT_params(sim, processors, include={"rotor_frequency"}) print(params.pretty_print(columns=["value", "min", "max", "vary", "expr"])) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none Name Value Min Max Vary Expr SP_0_operation_1_Exponential_FWHM 20 -inf inf True None SP_0_operation_2_Exponential_FWHM 200 -inf inf True None SP_0_operation_4_Scale_factor 10 -inf inf True None SP_1_operation_1_Exponential_FWHM 30 -inf inf True None SP_1_operation_2_Exponential_FWHM 300 -inf inf True None SP_1_operation_4_Scale_factor 100 -inf inf True None SP_2_operation_1_Exponential_FWHM 10 -inf inf True None SP_2_operation_2_Exponential_FWHM 150 -inf inf True None SP_2_operation_4_Scale_factor 50 -inf inf True None mth_0_rotor_frequency 5000 4900 5100 True None mth_1_rotor_frequency 1940 1840 2040 True None mth_2_rotor_frequency 960 860 1060 True None sys_0_abundance 50 0 100 True None sys_0_site_0_isotropic_chemical_shift 176 -inf inf True None sys_0_site_0_shielding_symmetric_eta 0.6 0 1 True None sys_0_site_0_shielding_symmetric_zeta 60 -inf inf True None sys_1_abundance 50 0 100 False 100-sys_0_abundance sys_1_site_0_isotropic_chemical_shift 43 -inf inf True None sys_1_site_0_shielding_symmetric_eta 0.5 0 1 True None sys_1_site_0_shielding_symmetric_zeta 30 -inf inf True None None .. GENERATED FROM PYTHON SOURCE LINES 198-199 **Solve the minimizer using LMFIT** .. GENERATED FROM PYTHON SOURCE LINES 199-203 .. code-block:: default minner = Minimizer(sf.LMFIT_min_function, params, fcn_args=(sim, processors, sigmas)) result = minner.minimize() result .. raw:: html

Fit Statistics

fitting methodleastsq
# function evals182
# data points12288
# variables19
chi-square 15689.6963
reduced chi-square 1.27880808
Akaike info crit. 3040.95414
Bayesian info crit. 3181.86533

Variables

name value standard error relative error initial value min max vary expression
sys_0_site_0_isotropic_chemical_shift 176.116887 5.4044e-04 (0.00%) 176.0 -inf inf True
sys_0_site_0_shielding_symmetric_zeta 72.3189698 0.15407965 (0.21%) 60.0 -inf inf True
sys_0_site_0_shielding_symmetric_eta 0.91811712 0.00396939 (0.43%) 0.6 0.00000000 1.00000000 True
sys_0_abundance 56.9212657 0.13659678 (0.24%) 50.0 0.00000000 100.000000 True
sys_1_site_0_isotropic_chemical_shift 43.3550919 0.00384271 (0.01%) 43.0 -inf inf True
sys_1_site_0_shielding_symmetric_zeta 22.7833633 0.20784235 (0.91%) 30.0 -inf inf True
sys_1_site_0_shielding_symmetric_eta 0.44484360 0.05491572 (12.34%) 0.5 0.00000000 1.00000000 True
sys_1_abundance 43.0787343 0.13659678 (0.32%) 50.0 0.00000000 100.000000 False 100-sys_0_abundance
mth_0_rotor_frequency 5033.17026 0.41764162 (0.01%) 5000.0 4900.00000 5100.00000 True
mth_1_rotor_frequency 1940.76008 0.05775086 (0.00%) 1940.0 1840.00000 2040.00000 True
mth_2_rotor_frequency 962.164634 0.04765175 (0.00%) 960.0 860.000000 1060.00000 True
SP_0_operation_1_Exponential_FWHM 27.2770184 0.19323483 (0.71%) 20.0 -inf inf True
SP_0_operation_2_Exponential_FWHM 109.089988 1.31094156 (1.20%) 200.0 -inf inf True
SP_0_operation_4_Scale_factor 38.6045004 0.17629942 (0.46%) 10.0 -inf inf True
SP_1_operation_1_Exponential_FWHM 29.3880214 0.16504812 (0.56%) 30.0 -inf inf True
SP_1_operation_2_Exponential_FWHM 136.816220 0.99234027 (0.73%) 300.0 -inf inf True
SP_1_operation_4_Scale_factor 170.686810 0.57431206 (0.34%) 100.0 -inf inf True
SP_2_operation_1_Exponential_FWHM 34.4272669 0.30018680 (0.87%) 10.0 -inf inf True
SP_2_operation_2_Exponential_FWHM 170.832502 1.84941368 (1.08%) 150.0 -inf inf True
SP_2_operation_4_Scale_factor 170.583254 0.91131011 (0.53%) 50.0 -inf inf True

Correlations (unreported correlations are < 0.100)

sys_1_site_0_shielding_symmetric_zetasys_1_site_0_shielding_symmetric_eta-0.6387
SP_0_operation_1_Exponential_FWHMSP_0_operation_4_Scale_factor0.5932
SP_2_operation_1_Exponential_FWHMSP_2_operation_4_Scale_factor0.5890
SP_1_operation_1_Exponential_FWHMSP_1_operation_4_Scale_factor0.5196
SP_1_operation_2_Exponential_FWHMSP_1_operation_4_Scale_factor0.4924
SP_0_operation_2_Exponential_FWHMSP_0_operation_4_Scale_factor0.4602
sys_0_abundanceSP_1_operation_2_Exponential_FWHM-0.4373
sys_0_site_0_shielding_symmetric_zetasys_0_site_0_shielding_symmetric_eta-0.4317
SP_2_operation_2_Exponential_FWHMSP_2_operation_4_Scale_factor0.3709
sys_0_abundanceSP_0_operation_2_Exponential_FWHM-0.3607
sys_0_abundanceSP_1_operation_1_Exponential_FWHM0.3026
sys_0_abundanceSP_0_operation_4_Scale_factor-0.2782
sys_0_abundancesys_1_site_0_shielding_symmetric_zeta-0.2019
SP_0_operation_1_Exponential_FWHMSP_0_operation_2_Exponential_FWHM0.1930
sys_0_abundanceSP_2_operation_2_Exponential_FWHM-0.1880
sys_0_abundanceSP_1_operation_4_Scale_factor-0.1870
SP_2_operation_1_Exponential_FWHMSP_2_operation_2_Exponential_FWHM0.1854
sys_1_site_0_shielding_symmetric_zetaSP_2_operation_4_Scale_factor0.1681
sys_0_abundanceSP_2_operation_1_Exponential_FWHM0.1672
SP_0_operation_2_Exponential_FWHMSP_1_operation_2_Exponential_FWHM0.1660
sys_0_abundanceSP_2_operation_4_Scale_factor-0.1613
sys_0_site_0_shielding_symmetric_zetasys_0_abundance0.1549
sys_0_abundanceSP_0_operation_1_Exponential_FWHM0.1530
sys_0_site_0_shielding_symmetric_zetaSP_1_operation_4_Scale_factor0.1457
SP_0_operation_4_Scale_factorSP_1_operation_2_Exponential_FWHM0.1315
SP_1_operation_2_Exponential_FWHMSP_2_operation_2_Exponential_FWHM0.1306
sys_1_site_0_shielding_symmetric_zetaSP_1_operation_4_Scale_factor0.1293
mth_1_rotor_frequencySP_1_operation_1_Exponential_FWHM0.1196
SP_1_operation_1_Exponential_FWHMSP_1_operation_2_Exponential_FWHM0.1142
SP_0_operation_2_Exponential_FWHMSP_1_operation_1_Exponential_FWHM-0.1128
sys_0_abundancesys_1_site_0_shielding_symmetric_eta-0.1090
sys_1_site_0_shielding_symmetric_etaSP_2_operation_4_Scale_factor0.1074
sys_1_site_0_shielding_symmetric_zetaSP_2_operation_2_Exponential_FWHM-0.1031
SP_1_operation_2_Exponential_FWHMSP_2_operation_1_Exponential_FWHM-0.1018


.. GENERATED FROM PYTHON SOURCE LINES 204-206 The best fit solution --------------------- .. GENERATED FROM PYTHON SOURCE LINES 206-222 .. code-block:: default all_best_fit = sf.bestfit(sim, processors) # a list of best fit simulations all_residuals = sf.residuals(sim, processors) # a list of residuals # Plot the spectrum fig, ax = plt.subplots(1, 3, figsize=(12, 3), subplot_kw={"projection": "csdm"}) for i, proc in enumerate(processors): ax[i].plot(experiments[i], color="black", linewidth=0.5, label="Experiment") ax[i].plot(all_residuals[i].real, color="gray", linewidth=0.5, label="Residual") ax[i].plot(all_best_fit[i].real, linewidth=2, alpha=0.6) ax[i].set_xlim(280, -10) ax[i].grid() plt.legend() plt.tight_layout() plt.show() .. image-sg:: /fitting/1D_fitting/images/sphx_glr_plot_2_13C_glycine_multi_spectra_fit_003.png :alt: plot 2 13C glycine multi spectra fit :srcset: /fitting/1D_fitting/images/sphx_glr_plot_2_13C_glycine_multi_spectra_fit_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 223-227 .. [#f1] D.Massiot, F.Fayon, M.Capron, I.King, S.Le Calvé, B.Alonso, J.O.Durand, B.Bujoli, Z.Gan, G.Hoatson, 'Modelling one and two-dimensional solid-state NMR spectra.', Magn. Reson. Chem. **40** 70-76 (2002) `DOI: 10.1002/mrc.984 `_ .. rst-class:: sphx-glr-timing **Total running time of the script:** ( 0 minutes 11.623 seconds) .. _sphx_glr_download_fitting_1D_fitting_plot_2_13C_glycine_multi_spectra_fit.py: .. only :: html .. container:: sphx-glr-footer :class: sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_2_13C_glycine_multi_spectra_fit.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_2_13C_glycine_multi_spectra_fit.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_