"""The feature extraction module contains classes for feature extraction."""
import sys
import numpy as np
import pymia.filtering.filter as fltr
import SimpleITK as sitk
[docs]class AtlasCoordinates(fltr.Filter):
"""Represents an atlas coordinates feature extractor."""
[docs] def __init__(self):
"""Initializes a new instance of the AtlasCoordinates class."""
super().__init__()
[docs] def execute(self, image: sitk.Image, params: fltr.FilterParams = None) -> sitk.Image:
"""Executes a atlas coordinates feature extractor on an image.
Args:
image (sitk.Image): The image.
params (fltr.FilterParams): The parameters (unused).
Returns:
sitk.Image: The atlas coordinates image
(a vector image with 3 components, which represent the physical x, y, z coordinates in mm).
Raises:
ValueError: If image is not 3-D.
"""
if image.GetDimension() != 3:
raise ValueError('image needs to be 3-D')
x, y, z = image.GetSize()
# create matrix with homogenous indices in axis 3
coords = np.zeros((x, y, z, 4))
coords[..., 0] = np.arange(x)[:, np.newaxis, np.newaxis]
coords[..., 1] = np.arange(y)[np.newaxis, :, np.newaxis]
coords[..., 2] = np.arange(z)[np.newaxis, np.newaxis, :]
coords[..., 3] = 1
# reshape such that each voxel is one row
lin_coords = np.reshape(coords, [coords.shape[0] * coords.shape[1] * coords.shape[2], 4])
# generate transformation matrix
tmp_mat = image.GetDirection() + image.GetOrigin()
tfm = np.reshape(tmp_mat, [3, 4], order='F')
tfm = np.vstack((tfm, [0, 0, 0, 1]))
atlas_coords = (tfm @ np.transpose(lin_coords))[0:3, :]
atlas_coords = np.reshape(np.transpose(atlas_coords), [z, y, x, 3], 'F')
img_out = sitk.GetImageFromArray(atlas_coords)
img_out.CopyInformation(image)
return img_out
def __str__(self):
"""Gets a printable string representation.
Returns:
str: String representation.
"""
return 'AtlasCoordinates:\n' \
.format(self=self)
[docs]def first_order_texture_features_function(values):
"""Calculates first-order texture features.
Args:
values (np.array): The values to calculate the first-order texture features from.
Returns:
np.array: A vector containing the first-order texture features:
- mean
- variance
- sigma
- skewness
- kurtosis
- entropy
- energy
- snr
- min
- max
- range
- percentile10th
- percentile25th
- percentile50th
- percentile75th
- percentile90th
"""
eps = sys.float_info.epsilon # to avoid division by zero
mean = np.mean(values)
std = np.std(values)
snr = mean / std if std != 0 else 0
min_ = np.min(values)
max_ = np.max(values)
num_values = len(values)
p = values / (np.sum(values) + eps)
return np.array([mean,
np.var(values), # variance
std,
np.sqrt(num_values * (num_values - 1)) / (num_values - 2) * np.sum((values - mean) ** 3) /
(num_values*std**3 + eps), # adjusted Fisher-Pearson coefficient of skewness
np.sum((values - mean) ** 4) / (num_values * std ** 4 + eps), # kurtosis
np.sum(-p * np.log2(p)), # entropy
np.sum(p**2), # energy (intensity histogram uniformity)
snr,
min_,
max_,
max_ - min_,
np.percentile(values, 10),
np.percentile(values, 25),
np.percentile(values, 50),
np.percentile(values, 75),
np.percentile(values, 90)
])
[docs]class RandomizedTrainingMaskGenerator:
"""Represents a training mask generator.
A training mask is an image with intensity values 0 and 1, where 1 represents masked.
Such a mask can be used to sample voxels for training.
"""
[docs] @staticmethod
def get_mask(ground_truth: sitk.Image,
ground_truth_labels: list,
label_percentages: list,
background_mask: sitk.Image = None) -> sitk.Image:
"""Gets a training mask.
Args:
ground_truth (sitk.Image): The ground truth image.
ground_truth_labels (list of int): The ground truth labels,
where 0=background, 1=label1, 2=label2, ..., e.g. [0, 1]
label_percentages (list of float): The percentage of voxels of a corresponding label to extract as mask,
e.g. [0.2, 0.2].
background_mask (sitk.Image): A mask, where intensity 0 indicates voxels to exclude independent of the
label.
Returns:
sitk.Image: The training mask.
"""
# initialize mask
ground_truth_array = sitk.GetArrayFromImage(ground_truth)
mask_array = np.zeros(ground_truth_array.shape, dtype=np.uint8)
# exclude background
if background_mask is not None:
background_mask_array = sitk.GetArrayFromImage(background_mask)
background_mask_array = np.logical_not(background_mask_array)
ground_truth_array = ground_truth_array.astype(float) # convert to float because of np.nan
ground_truth_array[background_mask_array] = np.nan
for label_idx, label in enumerate(ground_truth_labels):
indices = np.transpose(np.where(ground_truth_array == label))
np.random.shuffle(indices)
no_mask_items = int(indices.shape[0] * label_percentages[label_idx])
for no in range(no_mask_items):
x = indices[no][0]
y = indices[no][1]
z = indices[no][2]
mask_array[x, y, z] = 1 # this is a masked item
mask = sitk.GetImageFromArray(mask_array)
mask.SetOrigin(ground_truth.GetOrigin())
mask.SetDirection(ground_truth.GetDirection())
mask.SetSpacing(ground_truth.GetSpacing())
return mask