Commit a339f232 authored by William Clarke's avatar William Clarke
Browse files

Merge branch 'enh/fsl_mrsi_correlations' into 'master'

Enh/fsl mrsi correlations

Closes #28

See merge request fsl/fsl_mrs!46
parents 0bd0d2f4 67f8d2d9
Pipeline #13860 canceled with stages
in 12 minutes and 34 seconds
This document contains the FSL-MRS release history in reverse chronological order.
1.1.11 (WIP)
-----------------------------
1.1.11 (Monday 4th April 2022)
------------------------------
- Now able to choose the number of workers in fsl_mrs_sim.
- HSLVD routines now use dense (non-sparse) routine for better numerical stability
- Updates and corrections to documentation, references to new FSL Course MRS section added.
- Fixed bugs in LCModel basis set handling.
- Removed divide by zero warnings in quantification of voxels where fitting has failed.
- New outputs from fsl_mrsi script: parameter correlation matrices, group mappings and parameter names
1.1.10 (Thursday 12 January 2022)
---------------------------------
......
......@@ -355,6 +355,27 @@ class MRSI(object):
return data
def list_to_correlation_array(self, data_list, indicies=None, cleanup=True, dtype=float):
'''Convert 5D array of correlation matricies indexed from an mrsi object
to a numpy array with the shape of the first three dimensions matching
that of the MRSI object .'''
if indicies is None:
indicies = self.get_indicies_in_order()
size_m, size_n = data_list[0].shape
if size_m != size_n:
raise ValueError(f'Only symmetric matricies are handled, size is ({size_m},{size_n}).')
data = np.zeros(self.spatial_shape + (size_m, size_n), dtype=dtype)
for d, ind in zip(data_list, indicies):
data[ind] = d
if cleanup:
data[np.isnan(data)] = 0
data[np.isinf(data)] = 0
return data
def check_basis(self, ppmlim=(.2, 4.2)):
"""Check orientation of basis using a single generated mrs object.
......
......@@ -107,6 +107,8 @@ def main():
help='Additional scaling modifier for external water referencing.')
optional.add_argument('--report', action="store_true",
help='output html report')
optional.add_argument('--output_correlations', action="store_true",
help='Output correlation matricies for each fit.')
optional.add_argument('--verbose', action="store_true",
help='spit out verbose info')
optional.add_argument('--overwrite', action="store_true",
......@@ -324,12 +326,14 @@ def main():
nuisance_folder = os.path.join(args.output, 'nuisance')
qc_folder = os.path.join(args.output, 'qc')
fit_folder = os.path.join(args.output, 'fit')
misc_folder = os.path.join(args.output, 'misc')
os.mkdir(concs_folder)
os.mkdir(uncer_folder)
os.mkdir(nuisance_folder)
os.mkdir(qc_folder)
os.mkdir(fit_folder)
os.mkdir(misc_folder)
# Extract concentrations
indicies = [res[1] for res in results]
......@@ -457,7 +461,6 @@ def main():
dtype=float))
# fit
# TODO: check if data has been conjugated, if so conjugate the predictions
mrs_scale = mrsi.get_scalings_in_order()
pred_list = []
for res, scale in zip(results, mrs_scale):
......@@ -492,6 +495,31 @@ def main():
cleanup=False,
dtype=np.complex64))
# Save a parameter mappings of:
# 1) metabolites to groups
results[0][0].metab_in_group_json(
os.path.join(misc_folder, 'metabolite_groups.json'))
# 2) A list of parameters (to go with the correlation matrix)
results[0][0].fit_parameters_json(
os.path.join(misc_folder, 'mrs_fit_parameters.json'))
if args.output_correlations:
# Per voxel correlations of parameters
corr_list = [res[0].corr for res in results]
corr_mats = mrsi.list_to_correlation_array(
corr_list,
indicies=indicies,
cleanup=True)
# Save
file_nm = os.path.join(misc_folder, 'fit_correlations.nii.gz')
corr_img = nib.Nifti1Image(corr_mats, mrsi_data.voxToWorldMat)
corr_img.header.set_intent(
1005, # NIFTI_INTENT_SYMMATRIX
params=(corr_list[0].shape[0], ),
name='MRS fit correlation matrix')
nib.save(corr_img, file_nm)
if args.verbose:
print('\n\n\nDone.')
......
......@@ -34,6 +34,7 @@ def test_fsl_mrsi(tmp_path):
data['seg_wm'],
data['seg_gm'],
data['seg_csf'],
'--output_correlations',
'--overwrite',
'--combine', 'Cr', 'PCr'])
......@@ -56,6 +57,10 @@ def test_fsl_mrsi(tmp_path):
assert (tmp_path / 'fit_out/nuisance/gamma_group0.nii.gz').exists()
assert (tmp_path / 'fit_out/nuisance/sigma_group0.nii.gz').exists()
assert (tmp_path / 'fit_out/misc/metabolite_groups.json').exists()
assert (tmp_path / 'fit_out/misc/mrs_fit_parameters.json').exists()
assert (tmp_path / 'fit_out/misc/fit_correlations.nii.gz').exists()
def test_fsl_mrsi_noh2o(tmp_path):
......@@ -86,3 +91,6 @@ def test_fsl_mrsi_noh2o(tmp_path):
assert (tmp_path / 'fit_out/nuisance/combined_lw_group0.nii.gz').exists()
assert (tmp_path / 'fit_out/nuisance/gamma_group0.nii.gz').exists()
assert (tmp_path / 'fit_out/nuisance/sigma_group0.nii.gz').exists()
assert (tmp_path / 'fit_out/misc/metabolite_groups.json').exists()
assert (tmp_path / 'fit_out/misc/mrs_fit_parameters.json').exists()
......@@ -10,6 +10,7 @@ from fsl_mrs.core import MRS
from fsl_mrs.utils.fitting import fit_FSLModel
from pytest import fixture
import numpy as np
import json
# Set up some synthetic data to use
......@@ -112,3 +113,42 @@ def test_qcOutput(data):
assert np.allclose(FWHM, 10.0, atol=1E0)
assert SNR.size == 3
def test_metabs_in_groups(data):
res = data[0]
met_g = res.metabs_in_groups()
assert met_g == [['Cr', 'PCr', 'NAA']]
def test_metabs_in_group(data):
res = data[0]
met_g = res.metabs_in_group(0)
assert met_g == ['Cr', 'PCr', 'NAA']
def test_metab_in_group_json(data, tmp_path):
res = data[0]
met_g = res.metab_in_group_json()
assert json.loads(met_g) == {'0': ['Cr', 'PCr', 'NAA']}
met_g2 = res.metab_in_group_json(tmp_path / 'test.json')
assert json.loads(met_g2) == {'0': ['Cr', 'PCr', 'NAA']}
assert (tmp_path / 'test.json').is_file()
with open(tmp_path / 'test.json') as fp:
assert json.load(fp) == {'0': ['Cr', 'PCr', 'NAA']}
def test_fit_parameters_json(data, tmp_path):
res = data[0]
res.fit_parameters_json(tmp_path / 'params.json')
with open(tmp_path / 'params.json') as fp:
saved_dict = json.load(fp)
assert saved_dict['parameters'] ==\
['Cr', 'PCr', 'NAA', 'gamma_0', 'sigma_0', 'eps_0', 'Phi0', 'Phi1', 'B_real_0', 'B_imag_0']
assert saved_dict['parameters_inc_comb'] ==\
['Cr', 'PCr', 'NAA', 'gamma_0', 'sigma_0', 'eps_0', 'Phi0', 'Phi1', 'B_real_0', 'B_imag_0', 'Cr+PCr']
assert saved_dict['metabolites'] == ['Cr', 'PCr', 'NAA']
assert saved_dict['metabolites_inc_comb'] == ['Cr', 'PCr', 'NAA', 'Cr+PCr']
......@@ -6,15 +6,17 @@
# Copyright (C) 2020 University of Oxford
# SHBASECOPYRIGHT
from copy import deepcopy
import json
import pandas as pd
import numpy as np
import fsl_mrs.utils.models as models
import fsl_mrs.utils.quantify as quant
import fsl_mrs.utils.qc as qc
from fsl_mrs.utils.misc import FIDToSpec, SpecToFID, calculate_lap_cov
import numpy as np
from copy import deepcopy
class FitRes(object):
"""
......@@ -410,6 +412,65 @@ class FitRes(object):
df.to_csv(filename, index=False, header=True)
def metabs_in_groups(self):
"""Return list of metabolites in each metabolite group
:return: Returns nested list with metabolite strings in each metabolite_group
:rtype: list
"""
group_list = []
for g in range(self.g):
metabs = []
for i, m in enumerate(self.original_metabs):
if self.metab_groups[i] == g:
metabs.append(m)
group_list.append(metabs)
return group_list
def metabs_in_group(self, group):
"""Get list of metabolites in specific group
:param group: Metabolite group index
:type group: int
:return: List of emtabolites in specified group
:rtype: List
"""
if group >= self.g:
raise ValueError(f'Group must be in the range 0 to {self.g - 1}.')
return self.metabs_in_groups()[group]
def metab_in_group_json(self, save_path=None):
"""Generate and (optionally) save a json representation
of the metabolites in each group
:param save_path: Optional path to save json to, defaults to None
:type save_path: str, optional
:return: Returns json formmatted string of metabolite groups
:rtype: str or Pathlib.Path
"""
dict_repr = {idx: ml for idx, ml in enumerate(self.metabs_in_groups())}
if save_path:
with open(save_path, 'w') as jf:
json.dump(dict_repr, jf)
return json.dumps(dict_repr)
else:
return json.dumps(dict_repr)
def fit_parameters_json(self, save_path):
"""Save a list of parameter names
:param save_path: Path to save json file containing parameter names
:type save_path: str or Pathlib.Path
"""
param_dict = {'parameters': self.params_names,
'parameters_inc_comb': self.params_names_inc_comb,
'metabolites': self.original_metabs,
'metabolites_inc_comb': self.metabs}
with open(save_path, 'w') as jf:
json.dump(param_dict, jf)
# Functions to return physically meaningful units from the fitting results
def getConc(self, scaling='raw', metab=None, function='mean'):
if function is None:
......
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