This practical covers eight example boxes in Chapter 3, section 3: "Magnetic resonance spectroscopy processing: A single voxel spectroscopy example". The boxes appear from page 82 to 93. As with all these mini-practicals, instructions on how to run FSL commands can be found on the Getting Started page.
Please download the dataset for this example here:
Data download
Once it is downloaded, move it to your preferred working directory and unzip it there.
This dataset contains a two sets of single voxel (SVS) data:
7t_unedited.
3t_edited.
Processing SVS data involves multiple steps. In this workflow we will apply the following steps:
Each of these steps are detailed in the following sub-sections. We'll go into detail in each one as to why they are being run. At the end we'll see a way of running all of these steps with a single command.
The overall aim of the pre-processing is to end up with a single water-suppressed 'metabolite' spectrum and a single water-reference spectrum. Both should have been processed in a near identical way to preserve the relative scaling between them. We therefore must apply processing to take our 32 individual coils and 64 dynamics and produce a single spectrum, appropriately phased, and with the minimum linewidth and maximum SNR possible.
Start by moving into the unedited data directory:
cd data_ch3_svs/7t_unedited
We will first tackle the uncombined coil dimension. To do this we will use the high SNR water-reference data to
calculate the complex (phase and amplitude) weights needed to appropriately sum each individual coil together. But
first, we need to form some water-reference data with an appropriate shape. Our raw water reference data also has
a DIM_DYN dimension which needs to be handled for it to be used as a reference. The water signal is
high and fairly constant between each element of this dimension so we simply average to collapse the dimension.
For this, and all the other steps in this section we will use the fsl_mrs_proc command.
fsl_mrs_proc has a series of sub commands that you will see in turn that apply a single
pre-processing operation to a NIfTI-MRS file and return a modified file. In this case we want the
average subcommand, and we want to run it on the DIM_DYN dimension of the raw
water-reference data data/steam_wref_raw.nii.gz.
fsl_mrs_proc average --file data/steam_wref_raw.nii.gz --dim DIM_DYN --output fsl_mrs_proc --filename steam_wref_tmp
DIM_DYN dimension
using mrs_tools split. This would have been a sensible thing to do if we knew one of the elements was
in some way corrupted.
Notice that we have also specified an output file name --filename and output location
--output. If you once again run mrs_tools info on the output file you will see that the
DIM_DYN dimension has been removed and the shape is now (1, 1, 1, 4096, 32).
*Reveal command
Now that we have formed our reference for coil combination we will use it to combine the metabolite
water-suppressed data, and the water-reference data itself. We call the coilcombine function of
fsl_mrs_proc to do this, providing the previously formed data as a reference. The wSVD algorithm will
then be run, combining data and removing the DIM_COIL dimension.
fsl_mrs_proc coilcombine --file data/steam_metab_raw.nii.gz --reference fsl_mrs_proc/steam_wref_tmp.nii.gz --filename steam_metab_comb --output fsl_mrs_proc -r fsl_mrs_proc coilcombine --file data/steam_wref_raw.nii.gz --reference fsl_mrs_proc/steam_wref_tmp.nii.gz --filename steam_wref_comb --output fsl_mrs_proc
In the first call above notice that we provided the -r flag. This requests that
fsl_mrs_proc produces an HTML report on the step. All fsl_mrs_proc commands can produce
their own report, which we can then combine together to form a single complete pre-processing report. We'll see
how to do that later.
Try opening the coil-combination report now in a browser window. This command can be run after each step to
(re)open the reports which generated after each step of pre-processing. The reports are stored in the
fsl_mrs_proc directory and named according to the time they were run.
firefox fsl_mrs_proc/*.html
You will see the results of coil combination on the metabolite water-suppressed data. Notice that the combined
signal SNR is much higher than the individual channels, revealing many peaks in the spectrum. For brevity, by
default only the result of the first index of the remaining DIM_DYN dimension are shown, to change
this you can add the --allreports flag.
At this point we now have a single dimension left to deal with in our metabolite and water-reference data:
DIM_DYN. As we saw in the spectrum visualisation step above the independent acquisitions that make up
this dimension are shifted along the frequency axis, and also have different phase. These shifts mostly arise due
to physiological (respiratory and cardiac) effects. In this step we want to use the similarity of each signal to
re-phase and shift the spectra back into alignment.
fsl_mrs_proc unlike can exclude elements
that fall outside a similarity metric threshold. In this data there isn't any large motion artefacts, and we will
not run this step.
fsl_mrs_proc align. Both the water-reference and metabolite data require alignment, and will have
experienced different shifts as they were acquired at different times. We therefore run the align
command on each separately and direct the algorithm to only consider parts of the spectrum that are relevant to each
case:
fsl_mrs_proc align --file fsl_mrs_proc/steam_metab_comb.nii.gz --ppm 0.2 4.2 --output fsl_mrs_proc -r --filename steam_metab_align fsl_mrs_proc align --file fsl_mrs_proc/steam_wref_comb.nii.gz --ppm 0 9 --output fsl_mrs_proc --filename steam_wref_align
Sometimes, due to subject motion, particular transients are badly corrupted by artefact. Even after frequency and phase alignment the spectrum collected during the corrupted transient remains an outlier. In this case we want to exclude these based on a similarity metric.
In FSL-MRS this can be done using fsl_mrs_proc unlike which calculates a simple
MSE metric, and excludes spectra that fall outside a defined standard deviation from the mean.
fsl_mrs_proc unlike --file fsl_mrs_proc/steam_metab_align.nii.gz --ppm 0.2 4.2 --output fsl_mrs_proc -r --filename steam_metab_unlike
Note: If you have very high quality data, this step may be unnecessary, as it will
exclude reasonable quality data. Consider also running with larger limits by setting the
--SD to a larger value (default = 1.96)
The previous step aligned all the data ready for averaging, but it didn't perform the averaging step. Run
mrs_tools info on the output of the last step to check for yourself. In this step we simply call
fsl_mrs_proc average again to achieve this.
fsl_mrs_proc average --file fsl_mrs_proc/steam_metab_unlike.nii.gz --dim DIM_DYN --output fsl_mrs_proc -r --filename steam_metab_avg fsl_mrs_proc average --file fsl_mrs_proc/steam_wref_align.nii.gz --dim DIM_DYN --output fsl_mrs_proc --filename steam_wref_avg
Now look at (visualise) the metabolite data across the full chemical shift range (0 - 6 ppm) after this step reveal command. Three things are apparent about the data:
Eddy current correction relies on the high SNR water-reference data experiencing the same eddy currents as the low
SNR metabolite data. This is the case if the water-reference acquisition sequence had identical sequence
gradients. If the water-reference data experiences eddy current effects we would expect to see a distinctly
non-Lorentzian peak shape. Take a look at the averaged water-reference data using mrs_tools vis and
answer the following question (you will have to expand the ppm range as above):
We will now use the water reference data to correct the eddy current effects in both the metabolite data and in the water-reference data itself.
fsl_mrs_proc ecc --file fsl_mrs_proc/steam_metab_avg.nii.gz --reference fsl_mrs_proc/steam_wref_avg.nii.gz --output fsl_mrs_proc -r --filename steam_metab_ecc fsl_mrs_proc ecc --file fsl_mrs_proc/steam_wref_avg.nii.gz --reference fsl_mrs_proc/steam_wref_avg.nii.gz --output fsl_mrs_proc --filename steam_wref_ecc
Now compare the water reference data before and after the eddy current correction step and notice the much nicer lineshape. The effect can also be observed on the NAA peak (at 2 ppm) as seen in the most recent report.
firefox fsl_mrs_proc/*.html
This will open several tabs in firefox. To inspect the ECC result, click through the tabs at the top until you find the one titled "ECC".
Note:If the water-reference data is not at the central spectrometer frequency
this step will introduce a global frequency shift. In this data that doesn't happen,
but if it does it can be corrected using fsl_mrs_proc fshift.
Now we will run that one additional step mentioned above. During data acquisition the signal receive hardware has to be activated at the right time to capture the echo. To ensure that the centre of the echo is captured, and to avoid any artefacts caused by the switching of hardware components it is common to start the signal reception a bit too early. We therefore want to remove any data that we captured and happened before the centre. If we don't do this the spectrum will have first-order phase (a phase linearly dependent on the frequency).
In this data a single point must be removed. This comes from explicit knowledge of the sequence implementation.
Other sequences report this correction number in a range of different ways. The lack of consistency in reporting
of important meta-data is a frustrating aspect of MRS data analysis, and improving the situation is one aim of
NIfTI-MRS and spec2nii.
Once again we treat the water-reference and metabolite data identically.
fsl_mrs_proc truncate --file fsl_mrs_proc/steam_metab_ecc.nii.gz --points -1 --pos first --filename steam_metab_trunc --output fsl_mrs_proc -r fsl_mrs_proc truncate --file fsl_mrs_proc/steam_wref_ecc.nii.gz --points -1 --pos first --filename steam_wref_trunc --output fsl_mrs_proc
Here we are removing a point --points -1 from the first time points --pos first. You
can use mrs_tools info to see the effect of this truncation on the data dimensions.
Note: the visual improvement might be hard to see until phase correction (step #8) has been run.
Even after all the processing above, a systematic frequency shift can remain in your spectrum. This will cause poor fitting performance, as the peaks of the metabolites will not align with their corresponding peaks in the basis spectra. These global shifts arise from inexact frequency alignment to water or large frequency drift at acquisition time, or shift induced by processing, such as the eddy current correction.
fsl_mrs_proc fshift command can also be used with the --shiftppm or
--shifthz arguments to shift spectra by a fixed amount.
In FSL-MRS we use fsl_mrs_proc fshift --shiftRef to remove these global shifts.
This command uses a simple peak searching algorithm to find peaks near their expected values,
and then sets them to the exact chemical shift expected. By default, it searches for the 3 ppm
creatine peak, and sets it to 3.027.
fsl_mrs_proc fshift --shiftRef --file fsl_mrs_proc/steam_metab_trunc.nii.gz --output fsl_mrs_proc -r --filename steam_metab_shifted
To choose a different peak and target chemical shift, you can use the --ppm to set
a search range, and --target to set a final position.
Back to our list of remaining steps.
We will now deal with that residual water peak. HLSVD is an algorithm that fits multiple lorentzian peaks to an
MRS spectrum. We can use this in a restricted frequency range and 'delete' any peaks in that range.
fsl_mrs_proc remove --file fsl_mrs_proc/steam_metab_shifted.nii.gz --output fsl_mrs_proc -r --filename steam_metab_hlsvd
Hopefully it is clear why we don't repeat this step to delete the water signal from the water-reference data ;-).
Use mrs_tools vis to look at the metabolite spectra before and after this step.
We need to extend the default x axis limits (in ppm) to cover the water chemical shift (~4.65 ppm).
Use the mrs_tools vis --help option to look at the argument options, and
construct an input
to extend the ppm limits to cover the water region.
Our final preprocessing step is to apply phase correction. This will put the spectrum into a purely positive
absorptive mode. I.e. all the singlet peaks should be positive and have no zero-crossings. We run
fsl_mrs_proc phase on both bits of data to form our final pre-processed files. For the
water-reference data the algorithm is directed just to run phase correction using the water peak between 4.6 and
4.7 ppm (--ppm 4.6 4.7).
fsl_mrs_proc phase --file fsl_mrs_proc/steam_metab_hlsvd.nii.gz --output fsl_mrs_proc -r --filename metab fsl_mrs_proc phase --ppm 4.6 4.7 --file fsl_mrs_proc/steam_wref_trunc.nii.gz --output fsl_mrs_proc -r --filename wref
Now merge the individual HTML reports in the fsl_mrs_proc directory into a single pre-processing
report. Simultaneously we will delete all the individual reports, using --delete.
merge_mrs_reports -d "FSL-MRS pre-processing practical." -o fsl_mrs_proc --delete fsl_mrs_proc/report*.html
If at any point you ran a command twice two reports will have been generated. Both reports will appear in the merged report. (Manual editing can be carried out by deleting the repeated HTML file). Now open the report and look at the last few stages we have run.
firefox fsl_mrs_proc/mergedReports.html
Run through this checklist. Hopefully the answer is yes to all of the questions!
mrs_tools vis).Running all those steps was quite hard work. FSL-MRS pre-processing was designed to be flexible and tunable by
the user. It would be quite easy to take this practical and formulate your own pre-processing bash script.
However, for standard SVS data, FSL-MRS also has a packaged processing pipeline fsl_mrs_preproc. The
below command should run a very similar set of pre-processing steps from a single command-line call. Try it and
compare the reports generated.
fsl_mrs_preproc\
--output fsl_mrs_preproc\
--data data/steam_metab_raw.nii.gz\
--reference data/steam_wref_raw.nii.gz\
--remove-water\
--truncate-fid 1\
--t1 T1.nii.gz \
--overwrite --report
Sometimes multiple types of water reference are acquired. This can be to avoid
magnetisation transfer effects of OVS bands in the water estimation. Water references
acquired for different reasons can be passed to fsl_mrs_preproc using the
following arguments:
--reference--ecc--quant
There are a number of additional options available in fsl_mrs_proc
and fsl_mrs_preproc that can be adjusted to improve performance for
some data.
--help argument. For example fsl_mrs_proc unlike --help.
--align_window {N} and/or the
--align_limits {min} {max} arguments.
The former aligns the data iteratively after windowed averaging, reducing artefact
by including the effect of phase cycling, whilst also boosting SNR.
The latter can be used to exclude artefact containing regions.
Both arguments are exposed for fsl_mrs_proc align
and fsl_mrs_preproc.
fsl_mrs_preproc
with the --noaverage (and optionally --noremoval).
All preprocessing steps will be carried out, as above, except for #5 (averaging),
and #4 (bad transient removal) respectively.
First move to the correct data directory
cd ../3t_edited
Much of the processing required for edited SVS data is identical to that carried out above for unedited data. There are just two additional steps. This is all with the caveat that processing must be applied separately, but consistently for each editing condition (ON and OFF for MEGA editing).
For the reason of consitent processing across the two or more
editing conditions, we recommend using the packaged, one-step
editing command fsl_mrs_preproc_edit.
If you need to customise the processing further, we will also
explain the separated additional steps.
fsl_mrs_preproc_edit\
--output fsl_mrs_preproc_edit\
--data data/mpress_metab_raw.nii.gz\
--reference data/mpress_wref_raw.nii.gz\
--align_ppm_edit 2.5 4\
--remove-water\
--truncate-fid 4\
--overwrite --report
This command works in the same way as fsl_mrs_preproc.
There are a couple of arguments to be aware of:
--align_ppm_dynamic and --align_window_dynamic
control the behaviour of the phase and frequency alignment of
temporal averages.
--align_ppm_edit controls the behaviour of the phase and
frequency alignment of the editing conditions.
If you wish to run the steps separately, the two additional steps required for processing edited data are:
After at least performing basic processing (as described above for unedited data), namely coil combination, and temporal averaging, you can run alignment of the editing dimension. *Reveal basic commands
We once again use the fsl_mrs_proc align command, but this time
we specify that the alignment should be run over the editing dimension (DIM-EDIT)
and that we only want to calculate the metrics for alignment over the region where
the two spectra are most similar (2.5 to 4 ppm).
fsl_mrs_proc align --dim DIM_EDIT --ppm 2.5 4 --file fsl_mrs_proc_edit/basic.nii.gz --filename edit_aligned --output fsl_mrs_proc_edit -r
After alignment of the edit dimension, we can subtract the two spectra to calculate the difference spectrum.
fsl_mrs_proc subtract --dim DIM_EDIT --file fsl_mrs_proc_edit/edit_aligned.nii.gz --filename diff --output fsl_mrs_proc_edit
Note: fsl_mrs_proc subtract can either be used to subtract
two spectra held in a single NIfTI file, or in two separate files.
Depending on the order of data in the file, you may need to apply a 180 degree phase to view the spectrum in a standard convention (NAA peak down).
fsl_mrs_proc fixed_phase --p0 180 --file fsl_mrs_proc_edit/diff.nii.gz --output fsl_mrs_proc_edit mrs_tools vis fsl_mrs_proc_edit/diff.nii.gz
The alignment of edited spectra is a difficult step in processing.
Fundamentally we are trying to align two similar, but different, signals,
without suppressing the true difference, that we wish to measure.
You might have to extensively trial different values of alignment
limits (set using align_ppm_edit), before you don't
see characteristic subtraction artefacts around the 3.0 and 3.2 ppm
creatine and choline singlets.
The alignment of temporal averages in edited data can often be badly affected
by artefactual lipid signal appearing in the 1 to 2 ppm region. Phase cycling
causes these lipids vary in phase every transient. Alignment algorithms, can
end up aligning the lipids, undoing the effect of phase cycling, resulting in
increased artefactual signal and decreased metabolite signal. Use the
--align_ppm_dynamic and --align_window_dynamic arguments
to customise the alignment limits to exclude lipid regions and/or incorporate
phase cycling.
Note that FSL-MRS has a python API that mirrors, almost exactly, all the commands run in the above tutorial. The API also can perform fitting steps, that will be discussed in later practicals.
If you want to try processing data using the API, you can see examples of its use in example python notebooks in the official FSL-MRS repository.