Commit 30cbd3e9 authored by William Clarke's avatar William Clarke
Browse files

Fix quantification section in interactive notebook.

parent 75ca898f
......@@ -408,6 +408,8 @@
"outputs": [],
"source": [
"import pandas as pd\n",
"from fsl_mrs.utils import quantify\n",
"\n",
"combinationList = [['NAA','NAAG'],\n",
" ['Glu','Gln'],\n",
" ['GPC','PCh'],\n",
......@@ -415,10 +417,18 @@
" ['Glc','Tau']]\n",
"\n",
"res.combine(combinationList)\n",
"res.calculateConcScaling(mrs,waterRefFID=mrs.H2O,\n",
" tissueFractions={'WM':0.4,'GM':0.4,'CSF':0.1},\n",
" TE=11e-3,\n",
" T2='Default')\n",
"\n",
"te = final_data.hdr_ext['EchoTime']\n",
"tr = final_data.hdr_ext['RepetitionTime']\n",
"q_info = quantify.QuantificationInfo(te,\n",
" tr,\n",
" mrs.names,\n",
" mrs.centralFrequency / 1E6)\n",
"q_info.set_fractions({'WM':0.45,'GM':0.45,'CSF':0.1})\n",
" \n",
"res.calculateConcScaling(mrs,\n",
" quant_info=q_info,\n",
" internal_reference=['Cr', 'PCr'])\n",
"\n",
"internal = res.getConc(scaling='internal',function=None).mean().multiply(8)\n",
"molarity = res.getConc(scaling='molarity',function=None).mean()\n",
......@@ -495,9 +505,9 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.7.7"
"version": "3.8.8-final"
}
},
"nbformat": 4,
"nbformat_minor": 4
}
}
\ No newline at end of file
%% Cell type:markdown id: tags:
# Example of SVS processing - interactive notebook
%% Cell type:markdown id: tags:
This notebook demos the process of fitting a single voxel scan in an interactive notebook using the underlying python libraries in FSL-MRS.
To view the plots in this notebook in Jupyter please consult the [plotly getting-started guide](https://plotly.com/python/getting-started/#jupyterlab-support-python-35).
### Contents:
1. [File conversion using spec2nii](#1.-File-conversion-using-spec2nii)
2. [Interactive preprocessing](#2.-Interactive-preprocessing)
3. [Fitting of the resultant spectrum](#3.-Fitting)
4. [Display of fitting results in a notebook](#4.-Display)
Will Clarke
June 2020
University of Oxford
%% Cell type:markdown id: tags:
## 1. File conversion using spec2nii
__THIS IS DONE ON THE COMMAND LINE__
Run spec2nii twix -v to establish the file contents, then using this information run spec2nii twix -e with the appropriate flags to extract the required scans. The -q flag suppresses text output.
This dataset uses a modified versions of the CMRR spectro package sequences on a Siemens scanner. It has three sets of water reference scans. The first is tagged as a phase correction, ans is collected at the start of the main suppressed water scan and will be used for eddy current correction.
The second is collected in a separate scan with only the RF portion of the water suppression disabled. This is used for coil combination (and could be used for eddy current correction). The third is collected with the OVS and all aspects of the water suppression disabled. It therefore experiences eddy currents unlike all the other scans. It will be used for final concentration scaling.
%% Cell type:code id: tags:
``` python
%%bash
spec2nii twix -v -q example_data/meas_MID310_STEAM_metab_FID115673.dat
```
%% Cell type:markdown id: tags:
From the final lines of the output of the first cell we can see that this "twix" file contains two groups of data.
The first is tagged as "image" and contains 64 repetitions (sets) of 4096 point FIDs collected on 32 channels.
The second has a single FID on each channel.
We now call spec2nii again specifying the group we want to extract each time. Each call to spec2nii will generate a NIfTI MRS file with a size of 1x1x1x4096x32xNdyn, where Ndyn is the number of dynamics (repeats).
We repeat this for the water reference scans, extracting just the image data.
%% Cell type:code id: tags:
``` python
%%bash
spec2nii twix -e image -f steam_metab_raw -o data -j -q example_data/meas_MID310_STEAM_metab_FID115673.dat
spec2nii twix -e phasecor -f steam_ecc_raw -o data -j -q example_data/meas_MID310_STEAM_metab_FID115673.dat
spec2nii twix -e image -f steam_wref_comb_raw -o data -j -q example_data/meas_MID311_STEAM_wref1_FID115674.dat
spec2nii twix -e image -f steam_wref_quant_raw -o data -j -q example_data/meas_MID312_STEAM_wref3_FID115675.dat
```
%% Cell type:markdown id: tags:
## 2. Interactive preprocessing
In this section we will preprocess the data using functions in the preproc package in fsl_mrs. This example could be used as a template to construct your own preprocessing script in python.
#### Description of steps
0. Load the data
1. Take averages of water references used for combination across files
2. Coil combine the metab data, the ecc data and the quantification data using the "comb" data as the reference.
3. Phase and frequency align the data where there are multiple transients.
4. Combine data across those transients by taking the mean.
5. Run eddy current correction using the appropriate reference.
6. In this data an additional FID point is collected before the echo centre. Remove this.
7. Run HLSVD on the data to remove the residual water in the water suppressed data.
6. Phase the data by a single peak as a crude zero-order phase correction.
%% Cell type:code id: tags:
``` python
import fsl_mrs.utils.mrs_io as mrs_io
```
%% Cell type:markdown id: tags:
#### 0. Load data
Load all the data into lists of data using the mrs_io.read_FID function
%% Cell type:code id: tags:
``` python
# Load the raw metabolite data
supp_data = mrs_io.read_FID('data/steam_metab_raw.nii.gz')
print(f'Loaded water supressed data with shape {supp_data.shape} and dimensions {supp_data.dim_tags}.')
# Load water ref with eddy currents (for coil combination)
ref_data = mrs_io.read_FID('data/steam_wref_comb_raw.nii.gz')
print(f'Loaded unsupressed data with shape {ref_data.shape} and dimensions {ref_data.dim_tags}.')
# Load water ref without eddy currents (for quantification)
quant_data = mrs_io.read_FID('data/steam_wref_quant_raw.nii.gz')
print(f'Loaded unsupressed data with shape {quant_data.shape} and dimensions {quant_data.dim_tags}.')
# Load phasecor scan (for Eddy)
ecc_data = mrs_io.read_FID('data/steam_ecc_raw.nii.gz')
print(f'Loaded unsupressed data with shape {ecc_data.shape} and dimensions {ecc_data.dim_tags}.')
```
%% Cell type:markdown id: tags:
#### 1. Take averages of reference data for coil combination
Each water reference scan cointained two averages. Calculate the average for use as a coil combination reference.
%% Cell type:code id: tags:
``` python
from fsl_mrs.utils.preproc import nifti_mrs_proc as proc
avg_ref_data = proc.average(ref_data, 'DIM_DYN', figure=True)
```
%% Cell type:markdown id: tags:
#### 2. Coil combination
Coil combine the metab data, the ecc data and the quantification data using the "comb" data as the reference.
%% Cell type:code id: tags:
``` python
supp_data = proc.coilcombine(supp_data, reference=avg_ref_data, figure=True)
quant_data = proc.coilcombine(quant_data, reference=avg_ref_data)
ecc_data = proc.coilcombine(ecc_data, reference=avg_ref_data)
```
%% Cell type:markdown id: tags:
#### Additional step to give resonable display phase for metabolites
Phase using single peak (Cr at 3.03 ppm)
%% Cell type:code id: tags:
``` python
supp_data = proc.apply_fixed_phase(supp_data, 180.0, figure=True)
quant_data = proc.apply_fixed_phase(quant_data, 180.0)
ecc_data = proc.apply_fixed_phase(ecc_data, 180.0)
```
%% Cell type:markdown id: tags:
#### 3. Phase and freq alignment
Phase and frequency align the data where there are multiple transients.
%% Cell type:code id: tags:
``` python
supp_data = proc.align(supp_data, 'DIM_DYN', ppmlim=(0, 4.2), figure=True)
# Alignment for water scans
quant_data = proc.align(quant_data, 'DIM_DYN', ppmlim=(0, 8))
```
%% Cell type:markdown id: tags:
#### 4. Combine scans
Combine data across transients by taking the mean.
%% Cell type:code id: tags:
``` python
supp_data = proc.average(supp_data, 'DIM_DYN', figure=True)
quant_data = proc.average(quant_data, 'DIM_DYN')
```
%% Cell type:markdown id: tags:
#### 5. ECC
Run eddy current correction using the appropriate reference.
%% Cell type:code id: tags:
``` python
supp_data = proc.ecc(supp_data, ecc_data, figure=True)
quant_data = proc.ecc(quant_data, quant_data)
```
%% Cell type:markdown id: tags:
#### 6. Truncation
In this data an additional FID point is collected before the echo centre. Remove this point.
%% Cell type:code id: tags:
``` python
supp_data = proc.truncate_or_pad(supp_data, -1, 'first', figure=True)
quant_data = proc.truncate_or_pad(quant_data, -1, 'first')
```
%% Cell type:markdown id: tags:
#### 7. Remove residual water
Run HLSVD on the data to remove the residual water in the water suppressed data.
%% Cell type:code id: tags:
``` python
limits = [-0.15,0.15]
limunits = 'ppm'
supp_data = proc.remove_peaks(supp_data, limits, limit_units=limunits, figure=True)
```
%% Cell type:markdown id: tags:
#### 8. Shift to reference
Ensure peaks appear at correct frequencies after alignment and ecc.
%% Cell type:code id: tags:
``` python
supp_data = proc.shift_to_reference(supp_data, 3.027, (2.9, 3.1), figure=True)
```
%% Cell type:markdown id: tags:
#### 9. Phasing
Phase the data by a single peak as a basic zero-order phase correction.
%% Cell type:code id: tags:
``` python
final_data = proc.phase_correct(supp_data, (2.9, 3.1), figure=True)
final_wref = proc.phase_correct(quant_data, (4.55, 4.7), hlsvd=False, figure=True)
```
%% Cell type:markdown id: tags:
## 3. Fitting
### Load into MRS objects
- Read pre-baked basis file (this one was generated with `fsl_mrs_sim`).
- Create main MRS object.
- Prepare the data for fitting (this does additional checks such as whether the data needs to be conjugated, and scales the data to improve fitting robustness)
%% Cell type:code id: tags:
``` python
import matplotlib.pyplot as plt
# Create main MRS Object
mrs = final_data.mrs(basis_file='example_data/steam_11ms',
ref_data=final_wref)
mrs.processForFitting()
# Quick plots of the Metab and Water spectra
mrs.plot()
plt.show()
mrs.plot_ref()
plt.show()
plt.figure(figsize=(10,10))
mrs.plot_basis()
plt.show()
```
%% Cell type:markdown id: tags:
### Fitting
Here we show a typical model fitting and some of the parameters that can be user-set.
%% Cell type:code id: tags:
``` python
from fsl_mrs.utils import fitting, misc, plotting
# Separate macromolecule from the rest (it will have its own lineshape parameters)
metab_groups = misc.parse_metab_groups(mrs,'Mac')
# Fit with Newton algorithm
Fitargs = {'ppmlim':[0.2,4.2],
'method':'Newton','baseline_order':4,
'metab_groups':metab_groups,
'model':'voigt'}
res = fitting.fit_FSLModel(mrs,**Fitargs)
# Quick sanity-plot of the fit (see further down for interactive plotting)
_ = plotting.plot_fit(mrs,pred=res.pred,baseline=res.baseline)
```
%% Cell type:markdown id: tags:
### Quantification
Internal and water referencing.
Output is a pandas series.
%% Cell type:code id: tags:
``` python
import pandas as pd
from fsl_mrs.utils import quantify
combinationList = [['NAA','NAAG'],
['Glu','Gln'],
['GPC','PCh'],
['Cr','PCr'],
['Glc','Tau']]
res.combine(combinationList)
res.calculateConcScaling(mrs,waterRefFID=mrs.H2O,
tissueFractions={'WM':0.4,'GM':0.4,'CSF':0.1},
TE=11e-3,
T2='Default')
te = final_data.hdr_ext['EchoTime']
tr = final_data.hdr_ext['RepetitionTime']
q_info = quantify.QuantificationInfo(te,
tr,
mrs.names,
mrs.centralFrequency / 1E6)
q_info.set_fractions({'WM':0.45,'GM':0.45,'CSF':0.1})
res.calculateConcScaling(mrs,
quant_info=q_info,
internal_reference=['Cr', 'PCr'])
internal = res.getConc(scaling='internal',function=None).mean().multiply(8)
molarity = res.getConc(scaling='molarity',function=None).mean()
print(pd.concat([internal.rename('/Cr+PCr',inplace=True), molarity.rename('molarity (mM)',inplace=True)], axis=1))
```
%% Cell type:markdown id: tags:
## 4. Display
Results can be displayed with various plotting functions or rendered into an interactive HTML.
### In notebook
%% Cell type:code id: tags:
``` python
fig = plotting.plotly_fit(mrs,res)
fig.show()
```
%% Cell type:markdown id: tags:
### HTML report
%% Cell type:code id: tags:
``` python
from fsl_mrs.utils import report
import os
import datetime
output = '.'
report.create_report(mrs,res,
filename=os.path.join(output,'report.html'),
fidfile='meas_MID310_STEAM_metab_FID115673.dat',
basisfile='example_data/steam_11ms',
h2ofile='meas_MID311_STEAM_wref1_FID115674.dat',
date=datetime.datetime.now().strftime("%Y-%m-%d %H:%M"))
import webbrowser
current_path = os.getcwd()
# generate a URL
url = os.path.join('file:///'+current_path,'report.html')
webbrowser.open(url)
```
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment