add utils and requirements

Author: DKuppens

requirements.txt

@@ -0,0 +1,5 @@
+numpy
+torchio
+torch
+os
+logging

utils/data_loading/load_ivim_subject.py

@@ -0,0 +1,43 @@
+import os
+import logging
+import torch
+import numpy as np
+import torchio as tio
+
+
+def load_ivim_subject(study_subject_path):
+    """
+    loads torchio subject with IVIM data (signals and bvalues)
+    Args:
+        study_subject_path: path in which subject data is located
+
+    Returns:
+
+    """
+    subject_dict = {}
+
+    # find all files that match study path and subject id
+    for file in os.listdir(study_subject_path):
+        file_path = os.path.join(study_subject_path, file)
+        logging.info(f'start loading data from {file_path}')
+
+        # Check file extension for image file
+        if file_path[-2:] == "gz" or file_path[-2:] == "ii":
+
+            # load nifti image
+            image = tio.Image(file_path)
+            image.set_data(image.data.to(dtype=torch.float32))
+            subject_dict['signals'] = image
+
+        # Check if file contains bvalues
+        elif file_path[-2:] == "al":
+            text_file = np.genfromtxt(file_path)
+            bvals = np.array(text_file)
+            subject_dict["xvals"] = tio.Image(tensor=torch.Tensor(np.reshape(bvals, (bvals.shape[0], 1, 1, 1))))
+
+        else:
+            logging.info(f'skipping loading of file {file_path}, no appropriate file extension. ')
+
+        # Create subject
+    if 'xvals' in subject_dict.keys() and 'signals' in subject_dict.keys():
+        return tio.Subject(subject_dict)

utils/data_processing/processors/AverageSignalsOfEqualXvals.py

@@ -0,0 +1,48 @@
+import numpy as np
+import torch
+import torchio
+from torchio.transforms import Transform
+
+
+class AverageSignalsOfEqualXvals(Transform):
+
+    def __init__(self,  **kwargs):
+        super().__init__(**kwargs)
+
+    def apply_transform(self, subject):
+        """
+        normalize signals
+        Args:
+            signals: signals array to normalize
+            xvals: xval array
+
+        Returns:
+            normalized_signals: normalized signals array
+
+        """
+        images_dict = subject.get_images_dict()
+        signals = images_dict['signals'].numpy()
+        xvals = np.squeeze(images_dict['xvals'].numpy())
+        signals, xvals = self.average_signal_of_equal_xvals(signals, xvals)
+        subject.add_image(torchio.Image(tensor=torch.Tensor(signals)), 'signals')
+        subject.add_image(torchio.Image(tensor=torch.Tensor(np.reshape(xvals, (xvals.shape[0], 1, 1, 1)))), 'xvals')
+        return subject
+
+    @staticmethod
+    def average_signal_of_equal_xvals(signals, xvals):
+        """
+        average the signal of equal xvals
+        Args:
+            signals: signal matrix [signals X xvals]
+            xvals: array of xvals
+
+        Returns:
+            averaged_signal_matrix: averaged signal matrix [signals X unique_xvals]
+            unique_xval_arrays: unique xvals in averaged signal matrix
+
+        """
+        unique_xvals = np.unique(xvals)
+        averaged_signals = np.zeros((unique_xvals.shape[0], *signals.shape[1:]))
+        for xval_idx, unique_xval in enumerate(unique_xvals):
+            averaged_signals[xval_idx, ...] = np.squeeze(np.mean(signals[np.where(xvals == unique_xval), ...], axis=1))
+        return averaged_signals, unique_xvals

utils/data_processing/processors/FlattenImageData.py

@@ -0,0 +1,46 @@
+import numpy as np
+import torch
+import torchio
+
+from torchio.transforms import Transform
+
+
+class FlattenImageData(Transform):
+
+    def __init__(self, **kwargs):
+        super().__init__(**kwargs)
+
+    def apply_transform(self, subject):
+        """
+        flattens image data of signals image of subject
+        Args:
+            signals: signals array to normalize
+            xvals: xval array
+
+        Returns:
+            normalized_signals: normalized signals array
+
+        """
+        images_dict = subject.get_images_dict(include=self.include, exclude=self.exclude)
+        for image_key, image in images_dict.items():
+            flattened_array = self.flatten_image_data(image.numpy())
+            subject.add_image(torchio.Image(tensor=torch.Tensor(np.reshape(flattened_array,
+                                                                           (flattened_array.shape[0],
+                                                                            flattened_array.shape[1], 1, 1)))),
+                              image_key)
+        return subject
+
+    @staticmethod
+    def flatten_image_data(signals):
+        """
+        Flattens 4D array into 2D array
+        Args:
+            signals: signals array to normalize
+
+        Returns:
+            normalized_signals: normalized signals array
+        """
+        bvals, x, y, z = signals.shape
+        signals_array = np.reshape(signals, (bvals, x * y * z))
+        return signals_array
+

utils/data_processing/processors/NormalizeMaxSignal.py

@@ -0,0 +1,43 @@
+import numpy as np
+import torch
+import torchio
+from torchio.transforms import Transform
+
+
+class NormalizeMaxSignal(Transform):
+
+    def __init__(self, **kwargs):
+        super().__init__(**kwargs)
+
+    def apply_transform(self, subject):
+        """
+        normalize signals Image of subject
+        Args:
+            subject: Subject
+
+        Returns:
+            subject: Subject
+        """
+        images_dict = subject.get_images_dict()
+        signals = images_dict['signals'].numpy()
+        signals = self.normalize_signals(signals)
+        subject.add_image(torchio.Image(tensor=torch.Tensor(signals)), 'signals')
+        return subject
+
+    @staticmethod
+    def normalize_signals(signals):
+        """
+         normalize signals
+         Args:
+             signals: signals array to normalize
+
+         Returns:
+             normalized_signals: normalized signals array
+
+         """
+        maxsignal = np.nanmax(signals, axis=0)
+        signals /= maxsignal
+
+        return signals
+
+

utils/data_processing/processors/NormalizeSignals.py

@@ -0,0 +1,52 @@
+import numpy as np
+import torch
+import torchio
+from torchio.transforms import Transform
+
+
+class NormalizeSignals(Transform):
+
+    def __init__(self, xval_threshold, **kwargs):
+        self.xval_threshold = xval_threshold
+        super().__init__(**kwargs)
+
+    def apply_transform(self, subject):
+        """
+        normalize xvals Image of subject
+        Args:
+            subject: Subject
+
+        Returns:
+            subject: Subject
+        """
+        images_dict = subject.get_images_dict()
+        signals = images_dict['signals'].numpy()
+        xvals = np.squeeze(images_dict['xvals'].numpy())
+        signals = self.normalize_signals(signals, xvals, self.xval_threshold)
+        subject.add_image(torchio.Image(tensor=torch.Tensor(signals)), 'signals')
+        return subject
+
+    @staticmethod
+    def normalize_signals(signals, xvals, xval_threshold):
+        """
+         normalize signals
+         Args:
+             signals: signals array to normalize
+             xvals: xval array
+             xval_threshold: threshold below which bvals are considered b0
+
+         Returns:
+             normalized_signals: normalized signals array
+
+         """
+        # get average b0 signals and set signals with S0 of 0 to nan
+        mean_S0 = np.nanmean(signals[xvals <= xval_threshold, :, :, :], axis=0)
+        signals[:, mean_S0 == 0] = np.nan
+
+        # normalize signals to S0 intensity
+        normalized_signals = signals / mean_S0
+        normalized_signals[np.isnan(normalized_signals)] = 0
+
+        return normalized_signals
+
+

utils/data_processing/processors/NormalizeXvals.py

@@ -0,0 +1,41 @@
+import torch
+import torchio
+import numpy as np
+
+from torchio.transforms import Transform
+
+
+class NormalizeXvals(Transform):
+
+    def __init__(self, normalization_factor, **kwargs):
+        self.normalization_factor = normalization_factor
+        super().__init__(**kwargs)
+
+    def apply_transform(self, subject):
+        """
+        normalize xvals Image of subject
+        Args:
+            subject: Subject
+
+        Returns:
+            subject: Subject
+        """
+        images_dict = subject.get_images_dict()
+        xvals = np.squeeze(images_dict['xvals'].numpy())
+        xvals = self.normalize_xvals(xvals, self.normalization_factor)
+        subject.add_image(torchio.Image(tensor=torch.Tensor(np.reshape(xvals, (xvals.shape[0], 1, 1, 1)))), 'xvals')
+        return subject
+
+    @staticmethod
+    def normalize_xvals(xvals, normalization_factor):
+        """
+        normalize signals
+        Args:
+            xvals: xvalue array
+            normalization_factor: factor to multiply xvals with
+
+        Returns:
+            normalized_xvals: normalized signals array
+        """
+        normalized_xvals = xvals * normalization_factor
+        return normalized_xvals

utils/data_processing/processors/SignalCuration.py

@@ -0,0 +1,61 @@
+import numpy as np
+import torch
+import torchio
+from torchio.transforms import Transform
+
+
+class SignalCuration(Transform):
+
+    def __init__(self, qmri_application, **kwargs):
+        self.qmri_application = qmri_application
+        super().__init__(**kwargs)
+
+    def apply_transform(self, subject):
+        """
+        curates signals Image of subject
+        Args:
+            subject: Subject
+
+        Returns:
+            subject: Subject
+        """
+        images_dict = subject.get_images_dict()
+        if self.qmri_application == 'IVIM' or 'ivim':
+            signals = images_dict['signals'].numpy()
+            xvals = np.squeeze(images_dict['xvals'].numpy())
+            valid_mask = self.ivim_selection(signals, xvals)
+            subject.add_image(torchio.Image(tensor=torch.Tensor(np.expand_dims(valid_mask, 0))), 'valid_mask')
+            return subject
+        else:
+            raise NotImplementedError
+
+    @staticmethod
+    def ivim_selection(signals, xvals):
+        """
+        returns only those signals exhibiting ivim decay
+        Args:
+            signals: signals for corresponding xvals
+            xvals: xvals
+
+        Returns:
+            normalized_valid_signals: normalized_valid_signals that exhibit ivim-like decay
+            masked_signals: normalized_signals where signals not exhibiting ivim-like decay are set to 0
+        """
+
+        # get average b0 signals and set signals with S0 of 0 to nan
+        mean_S0 = np.nanmean(signals[xvals <= 0.0001, ...], axis=0)
+
+        # select only those voxels with average S0 larger than half of median S0 of voxels with S0 larger than 0
+        valid_idcs_median_value = mean_S0 > (0.5 * np.nanmedian(mean_S0[mean_S0 > 0]))
+
+        # check if signal is ivim like
+        signals[signals > 1.5] = 1.5
+        valid_idcs_ivim_curve1 = np.percentile(signals[xvals * 100 < 50, ...], 95,
+                                               axis=0) < 1.3
+        valid_idcs_ivim_curve2 = np.percentile(signals[xvals * 100 > 50, ...], 95,
+                                               axis=0) < 1.2
+        valid_idcs_ivim_curve3 = np.percentile(signals[xvals * 100 > 150, ...], 95,
+                                               axis=0) < 1.0
+        mask_signals = valid_idcs_median_value & valid_idcs_ivim_curve1 & valid_idcs_ivim_curve2 & valid_idcs_ivim_curve3
+
+        return mask_signals

utils/data_processing/processors/SignalMask.py

@@ -0,0 +1,38 @@
+import numpy as np
+import torch
+import torchio
+from torchio.transforms import Transform
+
+
+class SignalMask(Transform):
+
+    def __init__(self, **kwargs):
+        super().__init__(**kwargs)
+
+    def apply_transform(self, subject):
+        """
+        curates signals Image of subject
+        Args:
+            subject: Subject
+
+        Returns:
+            subject: Subject
+        """
+        images_dict = subject.get_images_dict()
+        signals = images_dict['signals'].numpy()
+        signal_mask = self.signal_mask(signals)
+        subject.add_image(torchio.Image(tensor=torch.Tensor(np.expand_dims(signal_mask, 0))), 'signal_mask')
+        return subject
+
+    @staticmethod
+    def signal_mask(signals):
+        """
+        returns mask with nonzero element for signal vectors with nonzero entries
+        Args:
+            signals: signals
+
+        Returns:
+            masked_signals: signals containing nonzero elements
+        """
+
+        return np.sum(signals, axis=0) > 0