Add plot, update changelog etc. Move back to numpy 1.19 as 1.20 only availible...

Add plot, update changelog etc. Move back to numpy 1.19 as 1.20 only availible from condaforge. Update scripts to pass tests.

CHANGELOG.rst

 This document contains the FSL-MRS release history in reverse chronological order.
 
+WIP
+---
+- SNR measurements should cope with negative peak amplitudes correctly
+- Quantification now takes into account T1 relaxation
+- Quantification module now fits the water reference FID to deal with corruption of first FID points.
+- Restructure of internal quantification code.
+
 1.1.0 (Thursday 18th February 2021)
 -----------------------------------
 - Support for NIfTI-MRS format.

fsl_mrs/scripts/fsl_mrs

@@ -101,6 +101,15 @@ def main():
                           nargs='+',
                           help='Metabolite(s) used as an internal reference.'
                                ' Defaults to tCr (Cr+PCr).')
+    optional.add_argument('--wref_metabolite', type=str, default=None,
+                          help='Metabolite(s) used as an the reference for water scaling.'
+                               ' Uses internal defaults otherwise.')
+    optional.add_argument('--ref_protons', type=int, default=None,
+                          help='Number of protons that reference metabolite is equivalent to.'
+                               ' No effect without setting --wref_metabolite.')
+    optional.add_argument('--ref_int_limits', type=int, default=None, nargs=2,
+                          help='Reference spectrum integration limits (low, high).'
+                               ' No effect without setting --wref_metabolite.')
     optional.add_argument('--h2o_scale', type=float, default=1.0,
                           help='Additional scaling modifier for external water referencing.')
     optional.add_argument('--report', action="store_true",
@@ -307,7 +316,10 @@ def main():
         q_info = quantify.QuantificationInfo(echotime,
                                              repetition_time,
                                              mrs.names,
-                                             mrs.centralFrequency / 1E6)
+                                             mrs.centralFrequency / 1E6,
+                                             water_ref_metab=args.wref_metabolite,
+                                             water_ref_metab_protons=args.ref_protons,
+                                             water_ref_metab_limits=args.ref_int_limits)
 
         if args.tissue_frac:
             q_info.set_fractions(args.tissue_frac)

fsl_mrs/scripts/fsl_mrsi

@@ -94,6 +94,15 @@ def main():
                           nargs='+',
                           help='Metabolite(s) used as an internal reference.'
                                ' Defaults to tCr (Cr+PCr).')
+    optional.add_argument('--wref_metabolite', type=str, default=None,
+                          help='Metabolite(s) used as an the reference for water scaling.'
+                               ' Uses internal defaults otherwise.')
+    optional.add_argument('--ref_protons', type=int, default=None,
+                          help='Number of protons that reference metabolite is equivalent to.'
+                               ' No effect without setting --wref_metabolite.')
+    optional.add_argument('--ref_int_limits', type=int, default=None, nargs=2,
+                          help='Reference spectrum integration limits (low, high).'
+                               ' No effect without setting --wref_metabolite.')
     optional.add_argument('--report', action="store_true",
                           help='output html report')
     optional.add_argument('--verbose', action="store_true",
@@ -440,7 +449,10 @@ def runvoxel(mrs_in, args, Fitargs, echotime, repetition_time):
         q_info = quantify.QuantificationInfo(echotime,
                                              repetition_time,
                                              mrs.names,
-                                             mrs.centralFrequency / 1E6)
+                                             mrs.centralFrequency / 1E6,
+                                             water_ref_metab=args.wref_metabolite,
+                                             water_ref_metab_protons=args.ref_protons,
+                                             water_ref_metab_limits=args.ref_int_limits)
 
         if tissue_seg:
             q_info.set_fractions(tissue_seg)

fsl_mrs/utils/plotting.py

@@ -896,6 +896,81 @@ def plot_indiv(mrs, res, ppmlim=(.2, 4.2)):
     return fig
 
 
+def plot_references(mrs, res):
+    """Generate figure showing the areas and metabolites used to estimate the water scaling.
+
+    :param mrs: MRS object with unsuppressed water reference data.
+    :type mrs: fsl_mrs.core.mrs.MRS
+    :param res: FSL-MRS results object.
+    :type res: fsl_mrs.utils.results.FitRes
+    :return: Figure
+    :rtype: go.Figure
+    """
+    fig = make_subplots(rows=1, cols=2, subplot_titles=['Unsuppressed', 'Suppressed'])
+
+    quant_info = res.concScalings['quant_info']
+    min_val = np.min([res.ppmlim[0], quant_info.ref_limits[0], quant_info.h2o_limits[0]])
+    max_val = np.max([res.ppmlim[1], quant_info.ref_limits[1], quant_info.h2o_limits[1]])
+    data_range = (min_val, max_val)
+    axis = mrs.getAxes(ppmlim=data_range)
+    first, last = mrs.ppmlim_to_range(ppmlim=data_range)
+    water_first, water_last = res.concScalings['ref_info']['water_ref'].limits
+    water_axis = mrs.getAxes()[water_first:water_last]
+    metab_first, metab_last = res.concScalings['ref_info']['metab_ref'].limits
+    metab_axis = mrs.getAxes()[metab_first:metab_last]
+
+    y_data = np.real(FIDToSpec(mrs.H2O))[first:last]
+    trace1 = go.Scatter(x=axis, y=y_data,
+                        mode='lines',
+                        name='data',
+                        line=dict(color='rgb(0,0,0)', width=1))
+    y_data = np.real(FIDToSpec(res.concScalings['ref_info']['water_ref'].fid))[water_first:water_last]
+    trace2 = go.Scatter(x=water_axis, y=y_data,
+                        mode='lines',
+                        name='Fitted, integrated water',
+                        fill='tozeroy',
+                        line=dict(color='rgb(255,0,0)', width=1))
+    fig.add_trace(trace1, 1, 1)
+    fig.add_trace(trace2, 1, 1)
+
+    y_data = np.real(FIDToSpec(mrs.FID))[first:last]
+    trace3 = go.Scatter(x=axis, y=y_data,
+                        mode='lines',
+                        name='data',
+                        line=dict(color='rgb(0,0,0)', width=1))
+    y_data = np.real(FIDToSpec(res.concScalings['ref_info']['metab_ref'].original_fid))[first:last]
+    trace4 = go.Scatter(x=axis, y=y_data,
+                        mode='lines',
+                        name='Fitted Reference',
+                        line=dict(color='rgb(0,0,255)', width=1))
+    y_data = np.real(FIDToSpec(res.concScalings['ref_info']['metab_ref'].fid))[metab_first:metab_last]
+    trace5 = go.Scatter(x=metab_axis, y=y_data,
+                        mode='lines',
+                        name='Metabolite integrated',
+                        fill='tozeroy',
+                        line=dict(color='rgb(255,0,0)', width=1))
+    fig.add_trace(trace3, 1, 2)
+    fig.add_trace(trace4, 1, 2)
+    fig.add_trace(trace5, 1, 2)
+
+    fig.layout.xaxis.update(title_text='Chemical shift (ppm)',
+                            tick0=2, dtick=.5,
+                            range=[max_val, min_val])
+    fig.layout.xaxis2.update(title_text='Chemical shift (ppm)',
+                             tick0=2, dtick=.5,
+                             range=[max_val, min_val])
+
+    fig.layout.yaxis2.update(zeroline=True,
+                             zerolinewidth=1,
+                             zerolinecolor='Gray',
+                             showgrid=False)
+    fig.layout.yaxis.update(zeroline=True,
+                            zerolinewidth=1,
+                            zerolinecolor='Gray',
+                            showgrid=False)
+    return fig
+
+
 def create_table(df):
     """
     Generate plotly graphical Table from pandas dataframe

fsl_mrs/utils/quantify.py

@@ -6,7 +6,6 @@
 # Copyright (C) 2020 University of Oxford
 # SHBASECOPYRIGHT
 
-from fsl_mrs.core.nifti_mrs import NIFTI_MRS
 import numpy as np
 import matplotlib.pyplot as plt
 from scipy.optimize import minimize
@@ -94,6 +93,8 @@ class WaterRef(FIDIntegrator):
 class RefIntegral(FIDIntegrator):
     def __init__(self, mrs_obj, res_obj, metab, limits):
         super().__init__(mrs_obj, limits)
+        self.original_fid = res_obj.predictedFID(
+            mrs_obj, mode=metab, noBaseline=False, no_phase=False)
         self.fid = res_obj.predictedFID(
             mrs_obj, mode=metab, noBaseline=True, no_phase=True)
 

fsl_mrs/utils/report.py

@@ -77,13 +77,13 @@ def create_plotly_div(mrs, res):
     divs['indiv'] = to_div(fig)  # plotly.offline.plot(fig, output_type='div',include_plotlyjs='cdn')
 
     # Quantification table
-    if (res.concScalings['molality'] is not None) and hasattr(res.concScalings['quant_info'], 'f_GM'):
+    if (res.concScalings['molality'] is not None) and res.concScalings['quant_info'].f_GM is not None:
         Q = res.concScalings['quant_info']
         quant_df = pd.DataFrame()
         quant_df['Tissue-water densities (g/cm^3)'] = [Q.d_GM, Q.d_WM, Q.d_CSF]
         quant_df['Tissue volume fractions'] = [Q.f_GM, Q.f_WM, Q.f_CSF]
-        quant_df['Water T2 (ms)'] = [Q.t2['H2O_GM'] * 1000, Q.t2['H2O_WM'] * 1000, Q.t2['H2O_CSF'] * 1000]
-        quant_df['Water T1 (s)'] = [Q.t1['H2O_GM'], Q.t1['H2O_WM'], Q.t1['H2O_CSF']]
+        quant_df['Water T2 (ms)'] = np.around([Q.t2['H2O_GM'] * 1000, Q.t2['H2O_WM'] * 1000, Q.t2['H2O_CSF'] * 1000])
+        quant_df['Water T1 (s)'] = np.around([Q.t1['H2O_GM'], Q.t1['H2O_WM'], Q.t1['H2O_CSF']], decimals=2)
         quant_df.index = ['GM', 'WM', 'CSF']
         quant_df.index.name = 'Tissue'
         quant_df.reset_index(inplace=True)
@@ -92,6 +92,10 @@ def create_plotly_div(mrs, res):
         fig.update_layout(height=100, margin=dict(l=0, r=0, t=0, b=0))
         divs['quant_table'] = to_div(fig)
 
+    if res.concScalings['molality'] is not None:
+        fig = plotting.plot_references(mrs, res)
+        divs['refs'] = to_div(fig)
+
     return divs
 
 
@@ -303,7 +307,7 @@ def create_report(mrs, res, filename, fidfile, basisfile, h2ofile, date, locatio
                     <td>{Q.tr} s</td>
                 </tr>
                 <tr>
-                    <td class="titles">T<sub>2</sub> relaxation corrected water concentration:</td>
+                    <td class="titles">Relaxation corrected water concentration:</td>
                     <td>{relax_water_conc:0.0f} mmol/kg</td>
                 </tr>
                 <tr>
@@ -370,6 +374,14 @@ def create_report(mrs, res, filename, fidfile, basisfile, h2ofile, date, locatio
             """
         template += section
 
+    # Real/imag section
+    section = f"""
+    <h2><a name="refs">Water referencing</a></h2>
+    {divs['refs']}
+    <hr>
+    """
+    template += section
+
     # Details of the methods
     # methodsfile="/path/to/file"
     # methods = np.readtxt(methodsfile)

requirements.txt

-numpy>=1.20
-scipy>=1.4
+numpy>=1.19
+scipy>=1.6
 matplotlib
 plotly>=4
 pandas