Usage of multinmrfit package
1. General information
This notebook showcases usage of the Spectrum object of multinmrfit.
Content
2. Prepare environment
Download and install Anaconda (>= 3.7) on your computer: https://www.anaconda.com/distribution/
Install
multinmrfit:run
Anaconda promptfrom the start menuinstall multinmrfit with:
python -m pip install git+https://github.com/NMRTeamTBI/MultiNMRFit.git@branch_name(wherebranch_nameis the git branch you want to install)
Start Jupyter from the start menu
Open the jupyter notbook
Import required packages.
[1]:
import logging
import sys
import pandas as pd
Initialize logger.
[2]:
logging.basicConfig(format='%(asctime)s | %(levelname)s : %(message)s', level=logging.DEBUG, stream=sys.stdout)
Load multinmrfit.
[3]:
import multinmrfit.base.spectrum as spectrum
import multinmrfit.base.io as io
3. Load NMR data
Data can be loaded from a TSV file containing columns ‘ppm’ and ‘intensity’.
[4]:
test_synthetic_dataset = pd.read_table("./data/data_sim_nmrfit.csv", sep="\t")
[5]:
test_synthetic_dataset.columns
[5]:
Index(['ppm', 'intensity'], dtype='object')
Data can also be loaded from TopSpin files by providing all required information in a dictionnary, assuming one-dimensional spectrum if rowno is None or two-dimensional spectrum if rowno is provided.
[6]:
test_topspin_dataset = {"data_path":"C:/Bruker/TopSpin4.0.7/data",
"dataset":"CFE_test",
"expno":"991",
"procno":"1",
"rowno":"3"}
The window of interest can be provided as a tuple containing lower and upper boundaries.
[7]:
window = (-0.2, 0.2)
We can then load the data in a Spectrum object.
[8]:
sp = spectrum.Spectrum(data=test_synthetic_dataset, window=window)
2024-12-20 09:44:20,814 | DEBUG : create Spectrum object
We can then work directly with this object. For instance, to view the experimental spectrum, use the plot() method with exp=True.
[9]:
fig = sp.plot(exp=True)
2024-12-20 09:44:25,632 | DEBUG : create plot
This method returns a plotly.graph_objects.Figure object that can be updated (e.g. to change the layout as shown below fot the size) and plotted.
[10]:
type(fig)
[10]:
plotly.graph_objs._figure.Figure
[11]:
fig.update_layout(autosize=False, width=900, height=400)
fig.show()
4. Peak picking
Use the peak_picking() method with argument threshold.
[12]:
peak_table = sp.peak_picking(1e6)
2024-12-20 09:44:33,517 | DEBUG : peak peaking
2024-12-20 09:44:33,543 | DEBUG : peak table
ppm intensity cID
0 0.0 1.250091e+06
To visualize the spectrum and the identified peaks, use the plot() method with the peak_table provided as argument pp.
[13]:
fig = sp.plot(pp=peak_table)
2024-12-20 09:44:36,253 | DEBUG : create plot
[14]:
fig.update_layout(autosize=False, width=900, height=400)
fig.show()
5. Spectrum simulation and fitting
5.1. Load models of NMR signals
Load models of signals implemented in multinmrfit.
[15]:
available_models = io.IoHandler.get_models()
2024-12-20 09:44:43,764 | DEBUG : add model from file 'model_acetate.py'
2024-12-20 09:44:43,785 | DEBUG : model name: 13C-acetate (CH3)
2024-12-20 09:44:43,787 | DEBUG : add model from file 'model_doublet.py'
2024-12-20 09:44:43,791 | DEBUG : model name: doublet
2024-12-20 09:44:43,792 | DEBUG : add model from file 'model_doublet_of_doublet.py'
2024-12-20 09:44:43,797 | DEBUG : model name: doublet of doublet
2024-12-20 09:44:43,798 | DEBUG : add model from file 'model_quartet.py'
2024-12-20 09:44:43,802 | DEBUG : model name: quartet
2024-12-20 09:44:43,803 | DEBUG : add model from file 'model_singlet.py'
2024-12-20 09:44:43,807 | DEBUG : model name: singlet
2024-12-20 09:44:43,808 | DEBUG : add model from file 'model_triplet.py'
2024-12-20 09:44:43,813 | DEBUG : model name: triplet
io.IoHandler.get_models() returns a dict with model_name-model_object as key-value pairs.
[16]:
available_models.keys()
[16]:
dict_keys(['13C-acetate (CH3)', 'doublet', 'doublet of doublet', 'quartet', 'singlet', 'triplet'])
[17]:
available_models["doublet"]
[17]:
multinmrfit.models.model_doublet.SignalModel
5.2. Model construction
To simulate or fit a spectrum, we need to provide a list of signals containing the type of signal (e.g. singlet or doublet) and the corresponding parameters (chemical shift, coupling constant, linewidth, intensity, etc). Signals must be provided as a dictionary.
[18]:
signals = {"singlet_TSP": {"model":"singlet", "par": {"x0": {"ini":0.0, "lb":-0.05, "ub":0.05}}}}
#signals = {"singlet_TSP": {"model":"singlet", "par": {"x0": {"ini":0.0, "lb":-0.05, "ub":0.05}}},
# "doublet_TSP": {"model":"doublet", "par": {"x0": {"ini":-0.01, "lb":-0.01, "ub":0.01}, "J": {"ini":0.147, "lb":0.14, "ub":0.15}, "lw": {"ini":0.001}}}}
Then we can build a model of the spectrum with the build_model() method.
[19]:
sp.build_model(signals=signals, available_models=available_models)
2024-12-20 09:45:01,755 | DEBUG : build Model for signal 'singlet_TSP'
2024-12-20 09:45:01,762 | DEBUG : parameters
signal_id model par ini lb ub
0 singlet_TSP singlet x0 1.000 0.0000 1.000000e+01
1 singlet_TSP singlet intensity 1000000.000 1.0000 1.000000e+15
2 singlet_TSP singlet lw 0.001 0.0001 3.000000e-02
3 singlet_TSP singlet gl 0.500 0.0000 1.000000e+00
Parameters can be accessed via the params attibute.
[20]:
sp.params
[20]:
| signal_id | model | par | ini | lb | ub | |
|---|---|---|---|---|---|---|
| 0 | singlet_TSP | singlet | x0 | 0.000 | -0.0500 | 5.000000e-02 |
| 1 | singlet_TSP | singlet | intensity | 1000000.000 | 1.0000 | 1.000000e+15 |
| 2 | singlet_TSP | singlet | lw | 0.001 | 0.0001 | 3.000000e-02 |
| 3 | singlet_TSP | singlet | gl | 0.500 | 0.0000 | 1.000000e+00 |
We can update parameters and offser using the update_params() method.
[21]:
sp.update_params({"singlet_TSP": {"par": {"intensity": {"ini":1e6, "ub":1e12}}}})
Similarly, we can update the offset with the update_offset() method. If offset=None, the offset is removed. To set an offset, provide a dictionary (if empty, offset is initialized to default values).
[22]:
sp.update_offset(offset={})
print(sp.params)
signal_id model par ini lb ub
0 singlet_TSP singlet x0 0.000 -0.050000 5.000000e-02
1 singlet_TSP singlet intensity 1000000.000 1.000000 1.000000e+12
2 singlet_TSP singlet lw 0.001 0.000100 3.000000e-02
3 singlet_TSP singlet gl 0.500 0.000000 1.000000e+00
4 full_spectrum None offset 0.000 -250018.225728 2.500182e+05
[23]:
sp.update_offset(offset=None)
print(sp.params)
signal_id model par ini lb ub
0 singlet_TSP singlet x0 0.000 -0.0500 5.000000e-02
1 singlet_TSP singlet intensity 1000000.000 1.0000 1.000000e+12
2 singlet_TSP singlet lw 0.001 0.0001 3.000000e-02
3 singlet_TSP singlet gl 0.500 0.0000 1.000000e+00
5.3. Simulation
Simulate spectrum using the simulate() method which returns a list of simulated intensities. If params is not given as argument, initial parameters values are used.
[24]:
sim_spectrum = sp.simulate()
[25]:
display(sim_spectrum)
0 12.499688
1 12.906144
2 13.332754
3 13.780870
4 14.251964
...
122 14.251964
123 13.780870
124 13.332754
125 12.906144
126 12.499688
Name: ppm, Length: 127, dtype: float64
To view the spectrum simulated from initial values, just use the plot() method.
[26]:
fig = sp.plot(ini=True)
fig.update_layout(autosize=False, width=900, height=700)
fig.show()
2024-12-20 09:45:16,641 | DEBUG : create plot
5.4. Fitting
To fit experimental spectrum, use the fit() method.
[27]:
sp.fit()
2024-12-20 09:45:20,402 | DEBUG : fit spectrum
2024-12-20 09:45:21,197 | DEBUG : parameters
signal_id model par ini lb ub \
0 singlet_TSP singlet x0 0.000 -0.0500 5.000000e-02
1 singlet_TSP singlet intensity 1000000.000 1.0000 1.000000e+12
2 singlet_TSP singlet lw 0.001 0.0001 3.000000e-02
3 singlet_TSP singlet gl 0.500 0.0000 1.000000e+00
opt opt_sd integral
0 8.048877e-10 0.000004 38655.427821
1 1.250012e+06 62.154093 38655.427821
2 1.000263e-02 0.000002 38655.427821
3 1.000000e+00 0.000177 38655.427821
Estimated parameters, standard deviations and integrals are now in the params attributes (columns opt, opt_sd and integral, respectively).
[28]:
sp.params
[28]:
| signal_id | model | par | ini | lb | ub | opt | opt_sd | integral | |
|---|---|---|---|---|---|---|---|---|---|
| 0 | singlet_TSP | singlet | x0 | 0.000 | -0.0500 | 5.000000e-02 | 8.048877e-10 | 0.000004 | 38655.427821 |
| 1 | singlet_TSP | singlet | intensity | 1000000.000 | 1.0000 | 1.000000e+12 | 1.250012e+06 | 62.154093 | 38655.427821 |
| 2 | singlet_TSP | singlet | lw | 0.001 | 0.0001 | 3.000000e-02 | 1.000263e-02 | 0.000002 | 38655.427821 |
| 3 | singlet_TSP | singlet | gl | 0.500 | 0.0000 | 1.000000e+00 | 1.000000e+00 | 0.000177 | 38655.427821 |
Fitting results can be viewed using the plot method with fit=True.
[29]:
fig = sp.plot(ini=True, fit=True)
fig.update_layout(autosize=False, width=900, height=700)
fig.show()
2024-12-20 09:45:24,827 | DEBUG : create plot
The fit will be better with an offset, we thus add it and fit again the data.
[30]:
sp.update_offset(offset={})
sp.fit()
2024-12-20 09:45:29,255 | DEBUG : fit spectrum
2024-12-20 09:45:30,136 | DEBUG : parameters
signal_id model par ini lb \
0 singlet_TSP singlet x0 0.000 -0.050000
1 singlet_TSP singlet intensity 1000000.000 1.000000
2 singlet_TSP singlet lw 0.001 0.000100
3 singlet_TSP singlet gl 0.500 0.000000
4 full_spectrum None offset 0.000 -250018.225728
ub opt opt_sd integral
0 5.000000e-02 3.940879e-09 7.231758e-07 38644.900682
1 1.000000e+12 1.249999e+06 3.013674e+02 38644.900682
2 3.000000e-02 9.999972e-03 1.117001e-06 38644.900682
3 1.000000e+00 9.999999e-01 9.555780e-05 38644.900682
4 2.500182e+05 9.180373e+01 8.011525e+01 NaN
[31]:
fig = sp.plot(ini=True, fit=True)
fig.update_layout(autosize=False, width=900, height=700)
fig.show()
2024-12-20 09:45:33,141 | DEBUG : create plot