import numpy as np from scipy.stats import entropy from ..util.dtype import dtype_range from .._shared.utils import _supported_float_type, check_shape_equality, warn __all__ = [ 'mean_squared_error', 'normalized_root_mse', 'peak_signal_noise_ratio', 'normalized_mutual_information', ] def _as_floats(image0, image1): """ Promote im1, im2 to nearest appropriate floating point precision. """ float_type = _supported_float_type((image0.dtype, image1.dtype)) image0 = np.asarray(image0, dtype=float_type) image1 = np.asarray(image1, dtype=float_type) return image0, image1 def mean_squared_error(image0, image1): """ Compute the mean-squared error between two images. Parameters ---------- image0, image1 : ndarray Images. Any dimensionality, must have same shape. Returns ------- mse : float The mean-squared error (MSE) metric. Notes ----- .. versionchanged:: 0.16 This function was renamed from ``skimage.measure.compare_mse`` to ``skimage.metrics.mean_squared_error``. """ check_shape_equality(image0, image1) image0, image1 = _as_floats(image0, image1) return np.mean((image0 - image1) ** 2, dtype=np.float64) def normalized_root_mse(image_true, image_test, *, normalization='euclidean'): """ Compute the normalized root mean-squared error (NRMSE) between two images. Parameters ---------- image_true : ndarray Ground-truth image, same shape as im_test. image_test : ndarray Test image. normalization : {'euclidean', 'min-max', 'mean'}, optional Controls the normalization method to use in the denominator of the NRMSE. There is no standard method of normalization across the literature [1]_. The methods available here are as follows: - 'euclidean' : normalize by the averaged Euclidean norm of ``im_true``:: NRMSE = RMSE * sqrt(N) / || im_true || where || . || denotes the Frobenius norm and ``N = im_true.size``. This result is equivalent to:: NRMSE = || im_true - im_test || / || im_true ||. - 'min-max' : normalize by the intensity range of ``im_true``. - 'mean' : normalize by the mean of ``im_true`` Returns ------- nrmse : float The NRMSE metric. Notes ----- .. versionchanged:: 0.16 This function was renamed from ``skimage.measure.compare_nrmse`` to ``skimage.metrics.normalized_root_mse``. References ---------- .. [1] https://en.wikipedia.org/wiki/Root-mean-square_deviation """ check_shape_equality(image_true, image_test) image_true, image_test = _as_floats(image_true, image_test) # Ensure that both 'Euclidean' and 'euclidean' match normalization = normalization.lower() if normalization == 'euclidean': denom = np.sqrt(np.mean((image_true * image_true), dtype=np.float64)) elif normalization == 'min-max': denom = image_true.max() - image_true.min() elif normalization == 'mean': denom = image_true.mean() else: raise ValueError("Unsupported norm_type") return np.sqrt(mean_squared_error(image_true, image_test)) / denom def peak_signal_noise_ratio(image_true, image_test, *, data_range=None): """ Compute the peak signal to noise ratio (PSNR) for an image. Parameters ---------- image_true : ndarray Ground-truth image, same shape as im_test. image_test : ndarray Test image. data_range : int, optional The data range of the input image (distance between minimum and maximum possible values). By default, this is estimated from the image data-type. Returns ------- psnr : float The PSNR metric. Notes ----- .. versionchanged:: 0.16 This function was renamed from ``skimage.measure.compare_psnr`` to ``skimage.metrics.peak_signal_noise_ratio``. References ---------- .. [1] https://en.wikipedia.org/wiki/Peak_signal-to-noise_ratio """ check_shape_equality(image_true, image_test) if data_range is None: if image_true.dtype != image_test.dtype: warn( "Inputs have mismatched dtype. Setting data_range based on " "image_true." ) dmin, dmax = dtype_range[image_true.dtype.type] true_min, true_max = np.min(image_true), np.max(image_true) if true_max > dmax or true_min < dmin: raise ValueError( "image_true has intensity values outside the range expected " "for its data type. Please manually specify the data_range." ) if true_min >= 0: # most common case (255 for uint8, 1 for float) data_range = dmax else: data_range = dmax - dmin image_true, image_test = _as_floats(image_true, image_test) err = mean_squared_error(image_true, image_test) data_range = float(data_range) # prevent overflow for small integer types return 10 * np.log10((data_range**2) / err) def _pad_to(arr, shape): """Pad an array with trailing zeros to a given target shape. Parameters ---------- arr : ndarray The input array. shape : tuple The target shape. Returns ------- padded : ndarray The padded array. Examples -------- >>> _pad_to(np.ones((1, 1), dtype=int), (1, 3)) array([[1, 0, 0]]) """ if not all(s >= i for s, i in zip(shape, arr.shape)): raise ValueError( f'Target shape {shape} cannot be smaller than input' f'shape {arr.shape} along any axis.' ) padding = [(0, s - i) for s, i in zip(shape, arr.shape)] return np.pad(arr, pad_width=padding, mode='constant', constant_values=0) def normalized_mutual_information(image0, image1, *, bins=100): r"""Compute the normalized mutual information (NMI). The normalized mutual information of :math:`A` and :math:`B` is given by:: .. math:: Y(A, B) = \frac{H(A) + H(B)}{H(A, B)} where :math:`H(X) := - \sum_{x \in X}{x \log x}` is the entropy. It was proposed to be useful in registering images by Colin Studholme and colleagues [1]_. It ranges from 1 (perfectly uncorrelated image values) to 2 (perfectly correlated image values, whether positively or negatively). Parameters ---------- image0, image1 : ndarray Images to be compared. The two input images must have the same number of dimensions. bins : int or sequence of int, optional The number of bins along each axis of the joint histogram. Returns ------- nmi : float The normalized mutual information between the two arrays, computed at the granularity given by ``bins``. Higher NMI implies more similar input images. Raises ------ ValueError If the images don't have the same number of dimensions. Notes ----- If the two input images are not the same shape, the smaller image is padded with zeros. References ---------- .. [1] C. Studholme, D.L.G. Hill, & D.J. Hawkes (1999). An overlap invariant entropy measure of 3D medical image alignment. Pattern Recognition 32(1):71-86 :DOI:`10.1016/S0031-3203(98)00091-0` """ if image0.ndim != image1.ndim: raise ValueError( f'NMI requires images of same number of dimensions. ' f'Got {image0.ndim}D for `image0` and ' f'{image1.ndim}D for `image1`.' ) if image0.shape != image1.shape: max_shape = np.maximum(image0.shape, image1.shape) padded0 = _pad_to(image0, max_shape) padded1 = _pad_to(image1, max_shape) else: padded0, padded1 = image0, image1 hist, bin_edges = np.histogramdd( [np.reshape(padded0, -1), np.reshape(padded1, -1)], bins=bins, density=True, ) H0 = entropy(np.sum(hist, axis=0)) H1 = entropy(np.sum(hist, axis=1)) H01 = entropy(np.reshape(hist, -1)) return (H0 + H1) / H01