This practical covers processing unedited and edited single-voxel data. 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.
We will continue using the data we collected yesterday.
In case of issues with this data we have some pre-packaged Siemens data. If required, then click this button to go to the practical using the backup data.
Backup practical
This dataset contains a two sets of single voxel (SVS) data:
press_acc.
mega_acc.
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 any individual coils and 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 marbel_data/press_acc
Sometimes you need to deal with data collected on a phased array receive coil that is uncombined. In this data this isn't the case, the scanner has already performed this step. If the scanner hadn't done it, then we would need to use the high SNR water-reference data to calculate the complex (phase and amplitude) weights needed to appropriately sum each individual coil together.
To see an example of this, look at the (very similar) practical we have using data from a Siemens scanner, available here.
With this data we have a single dimension to deal with in our metabolite and water-reference data:
DIM_DYN. As we saw in the spectrum visualisation practical, the independent acquisitions that make up
this dimension are slightly 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.
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
align subcommand, and we want to run it on the DIM_DYN dimension of the data.
fsl_mrs_proc unlike can exclude elements
that fall outside a similarity metric threshold. We will see this command in the next step.
Alignment is done using an algorithm called spectral registration.
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 press_nifti/metab_raw.nii.gz --ppm 0.2 4.2 --output fsl_mrs_proc -r --filename metab_align fsl_mrs_proc align --file press_nifti/wref_raw.nii.gz --ppm 0 9 --output fsl_mrs_proc --filename wref_align
Notice that we have also specified an output file name --filename
and output location --output.
In the first call above we also 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 alignment 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
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/metab_align.nii.gz --ppm 0.2 4.2 --output fsl_mrs_proc -r --filename 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)
If you once again run mrs_tools info on the output file you will see that the
DIM_DYN dimension has been reduced and the shape is now (1, 1, 1, 1024, 126).
*Reveal command
The previous steps 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 to achieve this.
fsl_mrs_proc average --file fsl_mrs_proc/metab_unlike.nii.gz --dim DIM_DYN --output fsl_mrs_proc -r --filename metab_avg fsl_mrs_proc average --file fsl_mrs_proc/wref_align.nii.gz --dim DIM_DYN --output fsl_mrs_proc --filename 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:
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).
Data exported after online processing, e.g., as DICOM, might not need step #4 running. The online processing will have already applied the appropriate time shift (by removing points). Let's figure out if we need to remove any points in this data. Try removing a single point first and seeing the effect on the spectrum.
fsl_mrs_proc truncate --file fsl_mrs_proc/metab_avg.nii.gz --points -1 --pos first --filename metab_trunc --output fsl_mrs_proc -r
With this command we have removed 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.
Now load the reports, and find the most recent to examine the effect on the phase of the peaks.
firefox fsl_mrs_proc/*.html
Note: the visual improvement might be hard to see until phase correction (step #8) has been run.
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/metab_avg.nii.gz --output fsl_mrs_proc -r --filename 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.
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/metab_hlsvd.nii.gz --reference fsl_mrs_proc/wref_avg.nii.gz --output fsl_mrs_proc -r --filename metab_ecc fsl_mrs_proc ecc --file fsl_mrs_proc/wref_avg.nii.gz --reference fsl_mrs_proc/wref_avg.nii.gz --output fsl_mrs_proc -r --filename 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 might 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 in step #8.
Our penultimate 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/metab_ecc.nii.gz --output fsl_mrs_proc -r --filename metab_phs fsl_mrs_proc phase --ppm 4.6 4.7 --file fsl_mrs_proc/wref_ecc.nii.gz --output fsl_mrs_proc -r --filename wref_phs
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. We will have to do this for water, setting a peak position of
4.65 ppm.
fsl_mrs_proc fshift --shiftRef --file fsl_mrs_proc/metab_phs.nii.gz --output fsl_mrs_proc -r --filename metab fsl_mrs_proc fshift --shiftRef --ppm 4 5 --target 4.65 --file fsl_mrs_proc/wref_phs.nii.gz --output fsl_mrs_proc --filename wref
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.
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!
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 press_nifti/metab_raw.nii.gz\
--reference press_nifti/wref_raw.nii.gz\
--remove-water\
--t1 ../ANAT_defaced_2_1.nii \
--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 ../mega_acc
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 consistent 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 nifti_edited/metab_reshaped.nii.gz\
--reference nifti_edited/wref_low.nii.gz\
--align_ppm_edit 2.5 4\
--remove-water\
--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 alignment 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 phase --file fsl_mrs_proc_edit/diff.nii.gz --output fsl_mrs_proc_edit --ppm 1.8 2.2 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.