Analysis

Utilities to aid in performing and evaluating image registration.

This module provides functions to compute displacements of image coordinates under a transformation, useful for assessing the accuracy of image registration processes.

References

[Power2012] (1,2,3)

Power, JD. et al. (2012). “Spurious but systematic correlations in functional connectivity MRI networks arise from subject motion.” NeuroImage, 59(3):2142-2154. doi:10.1016/j.neuroimage.2011.10.018.

nitransforms.analysis.utils.DEFAULT_FD_RADIUS = 50.0: Default radius (in mm) of a sphere where framewise displacements are calculated. The choice was proposed by [Power2012], and it represents approximately the mean distance from the cerebral cortex to the center of the head.

nitransforms.analysis.utils.compute_fd_from_motion(motion_parameters: ndarray, radius: float = 50.0) → ndarray[source]

Compute framewise displacement (FD) from motion parameters.

The framewise displacement is the sum of the magnitudes of the translational and rotational motion, computed from the frame-to-frame differences along the three spatial axes [Power2012].

Each row in the motion parameters represents one frame, and columns represent each coordinate axis x, y`, and z. Translation parameters are followed by rotation parameters column-wise.

Parameters:

motion_parameters (ndarray) – Motion parameters.
radius (float, optional) – Radius (in mm) of a sphere mimicking the size of a typical human brain.

Returns:

The framewise displacement (FD) at each timepoint as the L1 norm of frame-to-frame displacement across translations and rotation-derived displacements.

Return type:

ndarray

:raises exc:ValueError: If motion_parameters is not a 2D array with shape (T, 6).

nitransforms.analysis.utils.compute_fd_from_transform(img: SpatialImage, xfm: TransformBase, xfm_prev: TransformBase | None = None, radius: float = 50.0, n_vertices: int = 8) → float[source]

Compute the framewise displacement (FD) for a given transformation.

This implementation varies with respect to the original formulation by [Power2012] in that the FD is computed as the average across a number of vertices sampled over the sphere. See sample_unit_sphere() for details about the vertex sampling method.

For n_vertices == 1, FD is computed from rigid-body parameter increments (translation L1 + radius-scaled rotation L1) instead of averaging displacements over sampled sphere points. See compute_fd_from_motion() for direct comparability.

Parameters:

img (SpatialImage) – The reference image. Used to extract the center coordinates.
xfm (TransformBase) – The transformation to test. Applied to coordinates around the image center.
xfm_prev (TransformBase, optional) – A previous transformation to compare with. If None, the identity transformation is assumed (no transformation).
radius (float, optional) – The radius (in mm) of the spherical neighborhood around the center of the image.
n_vertices (int, optional) – The number of vertices to sample on the sphere.

Returns:

The average framewise displacement (FD) for the test transformation.

Return type:

float

:raises exc:ValueError: If n_vertices < 1.

nitransforms.analysis.utils.displacements_within_mask(mask_img: SpatialImage, xfm: TransformBase, xfm_prev: TransformBase | None = None) → ndarray[source]

Compute the distance between voxel coordinates mapped through two transforms.

Parameters:

mask_img (SpatialImage) – A mask image that defines the region of interest. Voxel coordinates within the mask are transformed.
xfm (TransformBase) – The transformation to test. This transformation is applied to the voxel coordinates.
xfm_prev (TransformBase, optional) – A previous (reference) transformation to compare with. If None, the identity transformation is assumed (no transformation).

Returns:

An array of displacements (in mm) for each voxel within the mask.

Return type:

ndarray

nitransforms.analysis.utils.euler_from_matrix(affine: ndarray, degrees: bool = True) → ndarray[source]

Extract XYZ Euler angles from affine or rotation matrices using SciPy.

affinendarray
Array with shape (…, 4, 4) or (…, 3, 3).

degreesbool, optional
If True, return degrees; otherwise radians.

ndarray

Array of shape (…, 3), Euler angles in ‘xyz’ convention.

ValueError: If affine does not end with shape (3, 3) or (4, 4).

nitransforms.analysis.utils.extract_motion_parameters(affine: ndarray) → Tuple[ndarray, ndarray][source]

Extract translation (mm) and rotation (degrees) parameters from an affine matrix.

Parameters:: affine (ndarray) – The affine transformation matrix.
Returns:: Extracted translation and rotation parameters.
Return type:: tuple

nitransforms.analysis.utils.sample_unit_sphere(n_points: int = 8) → ndarray[source]

Returns \(N\) evenly distributed points on a unit-radius sphere.

This function returns a deterministic, quasi-uniform point set on the surface of the unit sphere \(S^2 \subset \mathbb{R}^3\).

Notes

There is no unique notion of “evenly distributed” for arbitrary \(N\) on a sphere. This function uses:
- Platonic solids for certain small \(N\) (high symmetry; e.g. \(N = 6\) gives the \(\pm\) axis points).
- A Fibonacci / golden-angle spiral otherwise (fast, simple, good coverage).

Parameters:

n_points (int) – Number of points on the sphere.

Returns:

An array of shape (n_points, 3) whose rows have unit norm.

Return type:

ndarray

Raises:

TypeError – If n_points is a boolean or not an integer type.
ValueError – If n_points < 1.

Examples

Basic shape + unit norm:

>>> import numpy as np
>>> values = (1, 2, 8, 10, 12, 20)
>>> for n_pts in values:
...     X = sample_unit_sphere(n_pts)
...     X.shape, bool(np.allclose(np.linalg.norm(X, axis=1), 1.0))
((1, 3), True)
((2, 3), True)
((8, 3), True)
((10, 3), True)
((12, 3), True)
((20, 3), True)

Visualization of sampled points for each case:

import numpy as np
import matplotlib.pyplot as plt
import mpl_toolkits.mplot3d

from nitransforms.analysis.utils import sample_unit_sphere

values = (1, 2, 8, 10, 12, 20)

fig = plt.figure(figsize=(10, 6))
for i, n_pts in enumerate(values, start=1):
    X = sample_unit_sphere(n_pts)
    ax = fig.add_subplot(2, 3, i, projection="3d")
    ax.scatter(X[:, 0], X[:, 1], X[:, 2], s=30)
    ax.set_title(f"n={n_pts}")
    ax.set_xlabel("x")
    ax.set_ylabel("y")
    ax.set_zlabel("z")
    ax.set_box_aspect((1, 1, 1))
fig.tight_layout()

(Source code, png, hires.png, pdf)

For \(N = 6\), return the \(\pm\) axis points (octahedron vertices):

>>> X = sample_unit_sphere(6)
>>> # Each row has exactly one coordinate with magnitude 1, others 0
>>> bool(np.all((np.abs(X) == 1.0).sum(axis=1) == 1))
True
>>> bool(np.all((np.abs(X) == 1.0).sum(axis=0) == 2))  # each axis appears twice (±)
True

For \(N = 4\), the tetrahedron has constant pairwise dot product -1/3 off-diagonal:

>>> X = sample_unit_sphere(4)
>>> D = X @ X.T
>>> off = D[~np.eye(4, dtype=bool)]
>>> bool(np.allclose(off, -1/3))
True

For a quasi-uniform set, the second moment matrix is close to I/3, and the minimum angular separation is non-trivial:

>>> X = sample_unit_sphere(200)
>>> M = (X.T @ X) / len(X)
>>> bool(np.allclose(M, np.eye(3) / 3, atol=1e-3))
True
>>> def min_angle_rad(Y):
...     dots = np.clip(Y @ Y.T, -1.0, 1.0)
...     np.fill_diagonal(dots, 1.0)
...     ang = np.arccos(dots)
...     np.fill_diagonal(ang, np.inf)
...     return float(ang.min())
>>> min_angle_rad(X) > 0.18
True

Improper inputs:

>>> sample_unit_sphere(True)
Traceback (most recent call last):
...
TypeError: n_points must be a positive integer

>>> sample_unit_sphere(0)
Traceback (most recent call last):
...
ValueError: n_points must be 1 or greater