# align

 Bunch(**kwds) floating alias of float32 affine_registration(moving, static[, ...]) Find the affine transformation between two 3D images. read_mapping(disp, domain_img, codomain_img) Read a syn registration mapping from a nifti file. register_dwi_series(data, gtab[, affine, ...]) Register a DWI series to the mean of the B0 images in that series. register_dwi_to_template(dwi, gtab[, ...]) Register DWI data to a template through the B0 volumes. register_series(series, ref[, pipeline, ...]) Register a series to a reference image. resample(moving, static[, moving_affine, ...]) Resample an image (moving) from one space to another (static). streamline_registration(moving, static[, ...]) Register two collections of streamlines ('bundles') to each other. syn_registration(moving, static[, ...]) Register a 2D/3D source image (moving) to a 2D/3D target image (static). write_mapping(mapping, fname) Write out a syn registration mapping to a nifti file.

## Module: align._public

Registration API: simplified API for registration of MRI data and of streamlines.

 AffineMap(affine[, domain_grid_shape, ...]) Methods AffineRegistration([metric, level_iters, ...]) Methods AffineTransform3D Methods CCMetric(dim[, sigma_diff, radius]) Methods DiffeomorphicMap(dim, disp_shape[, ...]) Methods EMMetric(dim[, smooth, inner_iter, ...]) Methods MutualInformationMetric([nbins, ...]) Methods RigidIsoScalingTransform3D Methods RigidScalingTransform3D Methods RigidTransform3D Methods SSDMetric(dim[, smooth, inner_iter, step_type]) Methods StreamlineLinearRegistration([metric, x0, ...]) Methods SymmetricDiffeomorphicRegistration(metric[, ...]) Methods TranslationTransform3D Methods partial partial(func, *args, **keywords) - new function with partial application of the given arguments and keywords. affine_registration(moving, static[, ...]) Find the affine transformation between two 3D images. load_nifti(fname[, return_img, ...]) Load data and other information from a nifti file. load_trk(filename, reference[, to_space, ...]) Load the stateful tractogram of the .trk format read_img_arr_or_path(data[, affine]) Helper function that handles inputs that can be paths, nifti img or arrays read_mapping(disp, domain_img, codomain_img) Read a syn registration mapping from a nifti file. register_dwi_series(data, gtab[, affine, ...]) Register a DWI series to the mean of the B0 images in that series. register_dwi_to_template(dwi, gtab[, ...]) Register DWI data to a template through the B0 volumes. register_series(series, ref[, pipeline, ...]) Register a series to a reference image. resample(moving, static[, moving_affine, ...]) Resample an image (moving) from one space to another (static). save_nifti(fname, data, affine[, hdr, dtype]) Save a data array into a nifti file. set_number_of_points Change the number of points of streamlines streamline_registration(moving, static[, ...]) Register two collections of streamlines ('bundles') to each other. syn_registration(moving, static[, ...]) Register a 2D/3D source image (moving) to a 2D/3D target image (static). transform_centers_of_mass(static, ...) Transformation to align the center of mass of the input images. transform_tracking_output(tracking_output, ...) Apply a linear transformation, given by affine, to streamlines. warn(/, message[, category, stacklevel, source]) Issue a warning, or maybe ignore it or raise an exception. write_mapping(mapping, fname) Write out a syn registration mapping to a nifti file.

## Module: align.bundlemin

 determine_num_threads Determine the effective number of threads to be used for OpenMP calls distance_matrix_mdf Minimum direct flipped distance matrix between two streamline sets

## Module: align.crosscorr

Utility functions used by the Cross Correlation (CC) metric

 compute_cc_backward_step_2d(grad_moving, ...) Gradient of the CC Metric w.r.t. compute_cc_backward_step_3d(grad_moving, ...) Gradient of the CC Metric w.r.t. compute_cc_forward_step_2d(grad_static, ...) Gradient of the CC Metric w.r.t. compute_cc_forward_step_3d(grad_static, ...) Gradient of the CC Metric w.r.t. precompute_cc_factors_2d(static, moving, radius) Precomputations to quickly compute the gradient of the CC Metric precompute_cc_factors_2d_test(static, ...) Precomputations to quickly compute the gradient of the CC Metric precompute_cc_factors_3d(static, moving, radius) Precomputations to quickly compute the gradient of the CC Metric precompute_cc_factors_3d_test(static, ...) Precomputations to quickly compute the gradient of the CC Metric

## Module: align.expectmax

 compute_em_demons_step_2d(delta_field, ...) Demons step for EM metric in 2D compute_em_demons_step_3d(delta_field, ...) Demons step for EM metric in 3D compute_masked_class_stats_2d(mask, v, ...) Computes the mean and std. compute_masked_class_stats_3d(mask, v, ...) Computes the mean and std. quantize_positive_2d(v, num_levels) Quantizes a 2D image to num_levels quantization levels quantize_positive_3d(v, num_levels) Quantizes a 3D volume to num_levels quantization levels

## Module: align.imaffine

Affine image registration module consisting of the following classes:

AffineMap: encapsulates the necessary information to perform affine

transforms between two domains, defined by a static and a moving image. The domain of the transform is the set of points in the static image’s grid, and the codomain is the set of points in the moving image. When we call the transform method, AffineMap maps each point x of the domain (static grid) to the codomain (moving grid) and interpolates the moving image at that point to obtain the intensity value to be placed at x in the resulting grid. The transform_inverse method performs the opposite operation mapping points in the codomain to points in the domain.

ParzenJointHistogram: computes the marginal and joint distributions of

intensities of a pair of images, using Parzen windows [Parzen62] with a cubic spline kernel, as proposed by Mattes et al. [Mattes03]. It also computes the gradient of the joint histogram w.r.t. the parameters of a given transform.

MutualInformationMetric: computes the value and gradient of the mutual

information metric the way Optimizer needs them. That is, given a set of transform parameters, it will use ParzenJointHistogram to compute the value and gradient of the joint intensity histogram evaluated at the given parameters, and evaluate the the value and gradient of the histogram’s mutual information.

AffineRegistration: it runs the multi-resolution registration, putting

all the pieces together. It needs to create the scale space of the images and run the multi-resolution registration by using the Metric and the Optimizer at each level of the Gaussian pyramid. At each level, it will setup the metric to compute value and gradient of the metric with the input images with different levels of smoothing.

### References

[Parzen62] E. Parzen. On the estimation of a probability density

function and the mode. Annals of Mathematical Statistics, 33(3), 1065-1076, 1962.

[Mattes03] Mattes, D., Haynor, D. R., Vesselle, H., Lewellen, T. K.,

& Eubank, W. PET-CT image registration in the chest using free-form deformations. IEEE Transactions on Medical Imaging, 22(1), 120-8, 2003.

 AffineInvalidValuesError AffineInversionError AffineMap(affine[, domain_grid_shape, ...]) Methods AffineRegistration([metric, level_iters, ...]) Methods IsotropicScaleSpace(image, factors, sigmas) Methods MutualInformationMetric([nbins, ...]) Methods Optimizer(fun, x0[, args, method, jac, ...]) Attributes: Methods ScaleSpace(image, num_levels[, ...]) Methods compute_parzen_mi Computes the mutual information and its gradient (if requested) deprecated_params(old_name[, new_name, ...]) Deprecate a renamed or removed function argument. get_direction_and_spacings(affine, dim) Extracts the rotational and spacing components from a matrix interpolate_scalar_2d(image, locations) Bilinear interpolation of a 2D scalar image interpolate_scalar_3d(image, locations) Trilinear interpolation of a 3D scalar image sample_domain_regular Take floor(total_voxels/k) samples from a (2D or 3D) grid transform_centers_of_mass(static, ...) Transformation to align the center of mass of the input images. transform_geometric_centers(static, ...) Transformation to align the geometric center of the input images. transform_origins(static, static_grid2world, ...) Transformation to align the origins of the input images. warn(/, message[, category, stacklevel, source]) Issue a warning, or maybe ignore it or raise an exception.

## Module: align.imwarp

Classes and functions for Symmetric Diffeomorphic Registration

 Bunch(**kwds) DiffeomorphicMap(dim, disp_shape[, ...]) Methods DiffeomorphicRegistration([metric]) Methods ScaleSpace(image, num_levels[, ...]) Methods Streamlines alias of ArraySequence SymmetricDiffeomorphicRegistration(metric[, ...]) Methods floating alias of float32 get_direction_and_spacings(affine, dim) Extracts the rotational and spacing components from a matrix mult_aff(A, B) Returns the matrix product A.dot(B) considering None as the identity

## Module: align.metrics

Metrics for Symmetric Diffeomorphic Registration

 CCMetric(dim[, sigma_diff, radius]) Methods EMMetric(dim[, smooth, inner_iter, ...]) Methods SSDMetric(dim[, smooth, inner_iter, step_type]) Methods Methods floating alias of float32 gradient(f, *varargs[, axis, edge_order]) Return the gradient of an N-dimensional array. v_cycle_2d(n, k, delta_field, ...[, depth]) Multi-resolution Gauss-Seidel solver using V-type cycles v_cycle_3d(n, k, delta_field, ...[, depth]) Multi-resolution Gauss-Seidel solver using V-type cycles

## Module: align.parzenhist

 Methods compute_parzen_mi Computes the mutual information and its gradient (if requested) cubic_spline Evaluates the cubic spline at a set of values cubic_spline_derivative Evaluates the cubic spline derivative at a set of values sample_domain_regular Take floor(total_voxels/k) samples from a (2D or 3D) grid

## Module: align.reslice

 Pool([processes, initializer, initargs, ...]) Returns a process pool object affine_transform(input, matrix[, offset, ...]) Apply an affine transformation. determine_num_processes(num_processes) Determine the effective number of processes for parallelization. reslice(data, affine, zooms, new_zooms[, ...]) Reslice data with new voxel resolution defined by new_zooms.

## Module: align.scalespace

 IsotropicScaleSpace(image, factors, sigmas) Methods ScaleSpace(image, num_levels[, ...]) Methods floating alias of float32 gaussian_filter(input, sigma[, order, ...]) Multidimensional Gaussian filter.

## Module: align.streamlinear

 BundleMinDistanceAsymmetricMetric([num_threads]) Asymmetric Bundle-based Minimum distance. BundleMinDistanceMatrixMetric([num_threads]) Bundle-based Minimum Distance aka BMD BundleMinDistanceMetric([num_threads]) Bundle-based Minimum Distance aka BMD BundleSumDistanceMatrixMetric([num_threads]) Bundle-based Sum Distance aka BMD JointBundleMinDistanceMetric([num_threads]) Bundle-based Minimum Distance for joint optimization. JointStreamlineRegistrationMap(xopt, fopt, ...) Methods Optimizer(fun, x0[, args, method, jac, ...]) Attributes: StreamlineDistanceMetric([num_threads]) Methods StreamlineLinearRegistration([metric, x0, ...]) Methods StreamlineRegistrationMap(matopt, xopt, ...) Methods Streamlines alias of ArraySequence combinations(iterable, r) Return successive r-length combinations of elements in the iterable. bundle_min_distance(t, static, moving) MDF-based pairwise distance optimization function (MIN). MDF-based pairwise distance optimization function (MIN). bundle_min_distance_fast(t, static, moving, ...) MDF-based pairwise distance optimization function (MIN). bundle_sum_distance(t, static, moving[, ...]) MDF distance optimization function (SUM). center_streamlines(streamlines) Move streamlines to the origin compose_matrix([scale, shear, angles, ...]) Return 4x4 transformation matrix from sequence of transformations. compose_matrix44(t[, dtype]) Compose a 4x4 transformation matrix. Compose multiple 4x4 affine transformations in one 4x4 matrix decompose_matrix(matrix) Return sequence of transformations from transformation matrix. decompose_matrix44(mat[, size]) Given a 4x4 homogeneous matrix return the parameter vector. distance_matrix_mdf Minimum direct flipped distance matrix between two streamline sets get_unique_pairs(n_bundle[, pairs]) Make unique pairs from n_bundle bundles. groupwise_slr(bundles[, x0, tol, max_iter, ...]) Function to perform unbiased groupwise bundle registration. length Euclidean length of streamlines progressive_slr(static, moving, metric, x0, ...) Progressive SLR. qbx_and_merge(streamlines, thresholds[, ...]) Run QuickBundlesX and then run again on the centroids of the last layer remove_clusters_by_size(clusters[, min_size]) select_random_set_of_streamlines(...[, rng]) Select a random set of streamlines set_number_of_points Change the number of points of streamlines slr_with_qbx(static, moving[, x0, ...]) Utility function for registering large tractograms. Return the current time in seconds since the Epoch. transform_streamlines(streamlines, mat[, ...]) Apply affine transformation to streamlines unlist_streamlines(streamlines) Return the streamlines not as a list but as an array and an offset whole_brain_slr(static, moving[, x0, ...]) Utility function for registering large tractograms.

## Module: align.sumsqdiff

Utility functions used by the Sum of Squared Differences (SSD) metric

 compute_energy_ssd_2d(delta_field) Sum of squared differences between two 2D images compute_energy_ssd_3d(delta_field) Sum of squared differences between two 3D volumes The residual displacement field to be fit on the next iteration The residual displacement field to be fit on the next iteration compute_ssd_demons_step_2d(delta_field, ...) Demons step for 2D SSD-driven registration compute_ssd_demons_step_3d(delta_field, ...) Demons step for 3D SSD-driven registration One iteration of a large linear system solver for 2D SSD registration One iteration of a large linear system solver for 3D SSD registration solve_2d_symmetric_positive_definite Solves a 2-variable symmetric positive-definite linear system solve_3d_symmetric_positive_definite Solves a 3-variable symmetric positive-definite linear system

## Module: align.transforms

 AffineTransform2D Methods AffineTransform3D Methods RigidIsoScalingTransform2D Methods RigidIsoScalingTransform3D Methods RigidScalingTransform2D Methods RigidScalingTransform3D Methods RigidTransform2D Methods RigidTransform3D Methods RotationTransform2D Methods RotationTransform3D Methods ScalingTransform2D Methods ScalingTransform3D Methods Transform Base class (contract) for all transforms for affine image registration Each transform must define the following (fast, nogil) methods: TranslationTransform2D Methods TranslationTransform3D Methods

## Module: align.vector_fields

 compose_vector_fields_2d(d1, d2, ...) Computes the composition of two 2D displacement fields compose_vector_fields_3d(d1, d2, ...) Computes the composition of two 3D displacement fields create_circle Create a binary 2D image where pixel values are 1 iff their distance to the center of the image is less than or equal to radius. create_harmonic_fields_2d Creates an invertible 2D displacement field create_harmonic_fields_3d Creates an invertible 3D displacement field create_random_displacement_2d Creates a random 2D displacement 'exactly' mapping points of two grids create_random_displacement_3d Creates a random 3D displacement 'exactly' mapping points of two grids create_sphere Create a binary 3D image where voxel values are 1 iff their distance to the center of the image is less than or equal to radius. Down-samples the 2D input vector field by a factor of 2 Down-samples the input vector field by a factor of 2. Down-samples the input 3D vector field by a factor of 2 Down-samples the input 2D image by a factor of 2 Down-samples the input volume by a factor of 2 gradient Gradient of an image in physical space invert_vector_field_fixed_point_2d(d, ...[, ...]) Computes the inverse of a 2D displacement fields invert_vector_field_fixed_point_3d(d, ...[, ...]) Computes the inverse of a 3D displacement fields is_valid_affine reorient_vector_field_2d(d, affine) Linearly transforms all vectors of a 2D displacement field reorient_vector_field_3d(d, affine) Linearly transforms all vectors of a 3D displacement field resample_displacement_field_2d(field, ...) Resamples a 2D vector field to a custom target shape resample_displacement_field_3d(field, ...) Resamples a 3D vector field to a custom target shape simplify_warp_function_2d(d, affine_idx_in, ...) Simplifies a nonlinear warping function combined with an affine transform simplify_warp_function_3d(d, affine_idx_in, ...) Simplifies a nonlinear warping function combined with an affine transform sparse_gradient Gradient of an image in physical space transform_2d_affine(image, ref_shape[, affine]) Transforms a 2D image by an affine transform with bilinear interp. transform_2d_affine_nn(image, ref_shape[, ...]) Transforms a 2D image by an affine transform with NN interpolation transform_3d_affine(volume, ref_shape, affine) Transforms a 3D volume by an affine transform with trilinear interp. transform_3d_affine_nn(volume, ref_shape[, ...]) Transforms a 3D volume by an affine transform with NN interpolation warp_2d(image, d1[, affine_idx_in, ...]) Warps a 2D image using bilinear interpolation warp_2d_nn(image, d1[, affine_idx_in, ...]) Warps a 2D image using nearest neighbor interpolation warp_3d(volume, d1[, affine_idx_in, ...]) Warps a 3D volume using trilinear interpolation warp_3d_nn(volume, d1[, affine_idx_in, ...]) Warps a 3D volume using using nearest-neighbor interpolation warp_coordinates_2d(points, d1, in2world, ...) Parameters: warp_coordinates_3d(points, d1, in2world, ...) Parameters:

### Bunch

class dipy.align.Bunch(**kwds)

Bases: object

__init__(**kwds)

A ‘bunch’ of values (a replacement of Enum)

This is a temporary replacement of Enum, which is not available on all versions of Python 2

### floating

dipy.align.floating

alias of float32

### affine_registration

dipy.align.affine_registration(moving, static, moving_affine=None, static_affine=None, pipeline=None, starting_affine=None, metric='MI', level_iters=None, sigmas=None, factors=None, ret_metric=False, moving_mask=None, static_mask=None, **metric_kwargs)

Find the affine transformation between two 3D images. Alternatively, find the combination of several linear transformations.

Parameters:
movingarray, nifti image or str

Containing the data for the moving object, or full path to a nifti file with the moving data.

staticarray, nifti image or str

Containing the data for the static object, or full path to a nifti file with the moving data.

moving_affine4x4 array, optional

An affine transformation associated with the moving object. Required if data is provided as an array. If provided together with nifti/path, will over-ride the affine that is in the nifti.

static_affine4x4 array, optional

An affine transformation associated with the static object. Required if data is provided as an array. If provided together with nifti/path, will over-ride the affine that is in the nifti.

pipelinelist of str, optional

Sequence of transforms to use in the gradual fitting. Default: gradual fit of the full affine (executed from left to right): ["center_of_mass", "translation", "rigid", "affine"] Alternatively, any other combination of the following registration methods might be used: center_of_mass, translation, rigid, rigid_isoscaling, rigid_scaling and affine.

starting_affine: 4x4 array, optional

Initial guess for the transformation between the spaces. Default: identity.

metricstr, optional.

Currently only supports ‘MI’ for MutualInformationMetric.

AffineRegistration key-word argument: the number of iterations at each scale of the scale space. level_iters[0] corresponds to the coarsest scale, level_iters[-1] the finest, where n is the length of the sequence. By default, a 3-level scale space with iterations sequence equal to [10000, 1000, 100] will be used.

sigmassequence of floats, optional

AffineRegistration key-word argument: custom smoothing parameter to build the scale space (one parameter for each scale). By default, the sequence of sigmas will be [3, 1, 0].

AffineRegistration key-word argument: custom scale factors to build the scale space (one factor for each scale). By default, the sequence of factors will be [4, 2, 1].

ret_metricboolean, optional

Set it to True to return the value of the optimized coefficients and the optimization quality metric.

moving_maskarray, shape (S’, R’, C’) or (R’, C’), optional

moving image mask that defines which pixels in the moving image are used to calculate the mutual information.

static_maskarray, shape (S, R, C) or (R, C), optional

static image mask that defines which pixels in the static image are used to calculate the mutual information.

nbinsint, optional

MutualInformationMetric key-word argument: the number of bins to be used for computing the intensity histograms. The default is 32.

sampling_proportionNone or float in interval (0, 1], optional

MutualInformationMetric key-word argument: There are two types of sampling: dense and sparse. Dense sampling uses all voxels for estimating the (joint and marginal) intensity histograms, while sparse sampling uses a subset of them. If sampling_proportion is None, then dense sampling is used. If sampling_proportion is a floating point value in (0,1] then sparse sampling is used, where sampling_proportion specifies the proportion of voxels to be used. The default is None (dense sampling).

Returns:
transformedarray with moving data resampled to the static space
after computing the affine transformation
affinethe affine 4x4 associated with the transformation.
xoptthe value of the optimized coefficients.
foptthe value of the optimization quality metric.

Notes

Performs a gradual registration between the two inputs, using a pipeline that gradually approximates the final registration. If the final default step (affine) is omitted, the resulting affine may not have all 12 degrees of freedom adjusted.

Read a syn registration mapping from a nifti file.

Parameters:
dispstr or Nifti1Image

A file of image containing the mapping displacement field in each voxel Shape (x, y, z, 3, 2)

domain_imgstr or Nifti1Image
codomain_imgstr or Nifti1Image
Returns:
A DiffeomorphicMap object.

Notes

See write_mapping() for the data format expected.

### register_dwi_series

dipy.align.register_dwi_series(data, gtab, affine=None, b0_ref=0, pipeline=None, static_mask=None)

Register a DWI series to the mean of the B0 images in that series.

all first registered to the first B0 volume

Parameters:
data4D array or nibabel Nifti1Image class instance or str

Diffusion data. Either as a 4D array or as a nifti image object, or as a string containing the full path to a nifti file.

gtaba GradientTable class instance or tuple of strings

If provided as a tuple of strings, these are assumed to be full paths to the bvals and bvecs files (in that order).

affine4x4 array, optional.

Must be provided for data provided as an array. If provided together with Nifti1Image or str data, this input will over-ride the affine that is stored in the data input. Default: use the affine stored in data.

b0_refint, optional.

Which b0 volume to use as reference. Default: 0

pipelinelist of callables, optional.

The transformations to perform in sequence (from left to right): Default: [center_of_mass, translation, rigid, affine]

static_maskarray, shape (S, R, C) or (R, C), optional

static image mask that defines which pixels in the static image are used to calculate the mutual information.

Returns:
xform_img, affine_array: a Nifti1Image containing the registered data and
using the affine of the original data and a list containing the affine
transforms associated with each of the

### register_dwi_to_template

dipy.align.register_dwi_to_template(dwi, gtab, dwi_affine=None, template=None, template_affine=None, reg_method='syn', **reg_kwargs)

Register DWI data to a template through the B0 volumes.

Parameters:
dwi4D array, nifti image or str

Containing the DWI data, or full path to a nifti file with DWI.

The gradients associated with the DWI data, or a sequence with (fbval, fbvec), full paths to bvals and bvecs files.

dwi_affine4x4 array, optional

An affine transformation associated with the DWI. Required if data is provided as an array. If provided together with nifti/path, will over-ride the affine that is in the nifti.

template3D array, nifti image or str

Containing the data for the template, or full path to a nifti file with the template data.

template_affine4x4 array, optional

An affine transformation associated with the template. Required if data is provided as an array. If provided together with nifti/path, will over-ride the affine that is in the nifti.

reg_methodstr,

One of “syn” or “aff”, which designates which registration method is used. Either syn, which uses the syn_registration() function or affine_registration() function. Default: “syn”.

reg_kwargskey-word arguments for syn_registration() or

affine_registration()

Returns:
warped_b0, mapping: The fist is an array with the b0 volume warped to the
template. If reg_method is “syn”, the second is a DiffeomorphicMap class
instance that can be used to transform between the two spaces. Otherwise,
if reg_method is “aff”, this is a 4x4 matrix encoding the affine transform.

Notes

This function assumes that the DWI data is already internally registered. See register_dwi_series().

### register_series

dipy.align.register_series(series, ref, pipeline=None, series_affine=None, ref_affine=None, static_mask=None)

Register a series to a reference image.

Parameters:
series4D array or nib.Nifti1Image class instance or str

The data is 4D with the last dimension separating different 3D volumes

refint or 3D array or nib.Nifti1Image class instance or str

If this is an int, this is the index of the reference image within the series. Otherwise it is an array of data to register with (associated with a ref_affine required) or a nifti img or full path to a file containing one.

pipelinesequence, optional

Sequence of transforms to do for each volume in the series. Default: (executed from left to right): [center_of_mass, translation, rigid, affine]

series_affine, ref_affine4x4 arrays, optional.

The affine. If provided, this input will over-ride the affine provided together with the nifti img or file.

static_maskarray, shape (S, R, C) or (R, C), optional

static image mask that defines which pixels in the static image are used to calculate the mutual information.

Returns:
xformed, affines4D array with transformed data and a (4,4,n) array
with 4x4 matrices associated with each of the volumes of the input moving
data that was used to transform it into register with the static data.

### resample

dipy.align.resample(moving, static, moving_affine=None, static_affine=None, between_affine=None)

Resample an image (moving) from one space to another (static).

Parameters:
movingarray, nifti image or str

Containing the data for the moving object, or full path to a nifti file with the moving data.

moving_affine4x4 array, optional

An affine transformation associated with the moving object. Required if data is provided as an array. If provided together with nifti/path, will over-ride the affine that is in the nifti.

staticarray, nifti image or str

Containing the data for the static object, or full path to a nifti file with the moving data.

static_affine4x4 array, optional

An affine transformation associated with the static object. Required if data is provided as an array. If provided together with nifti/path, will over-ride the affine that is in the nifti.

between_affine: 4x4 array, optional

If an additional affine is needed between the two spaces. Default: identity (no additional registration).

Returns:
A Nifti1Image class instance with the data from the moving object
resampled into the space of the static object.

### streamline_registration

dipy.align.streamline_registration(moving, static, n_points=100, native_resampled=False)

Register two collections of streamlines (‘bundles’) to each other.

Parameters:
moving, staticlists of 3 by n, or str

The two bundles to be registered. Given either as lists of arrays with 3D coordinates, or strings containing full paths to these files.

n_pointsint, optional

How many points to resample to. Default: 100.

native_resampledbool, optional

Whether to return the moving bundle in the original space, but resampled in the static space to n_points.

Returns:
alignedlist

Streamlines from the moving group, moved to be closely matched to the static group.

matrixarray (4, 4)

The affine transformation that takes us from ‘moving’ to ‘static’

### syn_registration

dipy.align.syn_registration(moving, static, moving_affine=None, static_affine=None, step_length=0.25, metric='CC', dim=3, level_iters=None, prealign=None, **metric_kwargs)

Register a 2D/3D source image (moving) to a 2D/3D target image (static).

Parameters:
moving, staticarray or nib.Nifti1Image or str.

Either as a 2D/3D array or as a nifti image object, or as a string containing the full path to a nifti file.

moving_affine, static_affine4x4 array, optional.

Must be provided for data provided as an array. If provided together with Nifti1Image or str data, this input will over-ride the affine that is stored in the data input. Default: use the affine stored in data.

metricstring, optional

The metric to be optimized. One of CC, EM, SSD, Default: ‘CC’ => CCMetric.

dim: int (either 2 or 3), optional

The dimensions of the image domain. Default: 3

level_iterslist of int, optional

the number of iterations at each level of the Gaussian Pyramid (the length of the list defines the number of pyramid levels to be used). Default: [10, 10, 5].

metric_kwargsdict, optional

Parameters for initialization of the metric object. If not provided, uses the default settings of each metric.

Returns:
warped_movingndarray

The data in moving, warped towards the static data.

forwardndarray (…, 3)

The vector field describing the forward warping from the source to the target.

backwardndarray (…, 3)

The vector field describing the backward warping from the target to the source.

### write_mapping

dipy.align.write_mapping(mapping, fname)

Write out a syn registration mapping to a nifti file.

Parameters:
mappinga DiffeomorphicMap object derived from syn_registration()
fnamestr

Full path to the nifti file storing the mapping

Notes

The data in the file is organized with shape (X, Y, Z, 3, 2), such that the forward mapping in each voxel is in data[i, j, k, :, 0] and the backward mapping in each voxel is in data[i, j, k, :, 1].

### AffineMap

class dipy.align._public.AffineMap(affine, domain_grid_shape=None, domain_grid2world=None, codomain_grid_shape=None, codomain_grid2world=None)

Bases: object

Methods

 Return the value of the transformation, not a reference. set_affine(affine) Set the affine transform (operating in physical space). transform(image[, interpolation, ...]) Transform the input image from co-domain to domain space. transform_inverse(image[, interpolation, ...]) Transform the input image from domain to co-domain space.
__init__(affine, domain_grid_shape=None, domain_grid2world=None, codomain_grid_shape=None, codomain_grid2world=None)

AffineMap.

Implements an affine transformation whose domain is given by domain_grid and domain_grid2world, and whose co-domain is given by codomain_grid and codomain_grid2world.

The actual transform is represented by the affine matrix, which operate in world coordinates. Therefore, to transform a moving image towards a static image, we first map each voxel (i,j,k) of the static image to world coordinates (x,y,z) by applying domain_grid2world. Then we apply the affine transform to (x,y,z) obtaining (x’, y’, z’) in moving image’s world coordinates. Finally, (x’, y’, z’) is mapped to voxel coordinates (i’, j’, k’) in the moving image by multiplying (x’, y’, z’) by the inverse of codomain_grid2world. The codomain_grid_shape is used analogously to transform the static image towards the moving image when calling transform_inverse.

If the domain/co-domain information is not provided (None) then the sampling information needs to be specified each time the transform or transform_inverse is called to transform images. Note that such sampling information is not necessary to transform points defined in physical space, such as stream lines.

Parameters:
affinearray, shape (dim + 1, dim + 1)

the matrix defining the affine transform, where dim is the dimension of the space this map operates in (2 for 2D images, 3 for 3D images). If None, then self represents the identity transformation.

domain_grid_shapesequence, shape (dim,), optional

the shape of the default domain sampling grid. When transform is called to transform an image, the resulting image will have this shape, unless a different sampling information is provided. If None, then the sampling grid shape must be specified each time the transform method is called.

domain_grid2worldarray, shape (dim + 1, dim + 1), optional

the grid-to-world transform associated with the domain grid. If None (the default), then the grid-to-world transform is assumed to be the identity.

codomain_grid_shapesequence of integers, shape (dim,)

the shape of the default co-domain sampling grid. When transform_inverse is called to transform an image, the resulting image will have this shape, unless a different sampling information is provided. If None (the default), then the sampling grid shape must be specified each time the transform_inverse method is called.

codomain_grid2worldarray, shape (dim + 1, dim + 1)

the grid-to-world transform associated with the co-domain grid. If None (the default), then the grid-to-world transform is assumed to be the identity.

get_affine()

Return the value of the transformation, not a reference.

Returns:
affinendarray

Copy of the transform, not a reference.

set_affine(affine)

Set the affine transform (operating in physical space).

Also sets self.affine_inv - the inverse of affine, or None if there is no inverse.

Parameters:
affinearray, shape (dim + 1, dim + 1)

the matrix representing the affine transform operating in physical space. The domain and co-domain information remains unchanged. If None, then self represents the identity transformation.

transform(image, interpolation='linear', image_grid2world=None, sampling_grid_shape=None, sampling_grid2world=None, resample_only=False)

Transform the input image from co-domain to domain space.

By default, the transformed image is sampled at a grid defined by self.domain_shape and self.domain_grid2world. If such information was not provided then sampling_grid_shape is mandatory.

Parameters:
image2D or 3D array

the image to be transformed

interpolationstring, either ‘linear’ or ‘nearest’

the type of interpolation to be used, either ‘linear’ (for k-linear interpolation) or ‘nearest’ for nearest neighbor

image_grid2worldarray, shape (dim + 1, dim + 1), optional

the grid-to-world transform associated with image. If None (the default), then the grid-to-world transform is assumed to be the identity.

sampling_grid_shapesequence, shape (dim,), optional

the shape of the grid where the transformed image must be sampled. If None (the default), then self.codomain_shape is used instead (which must have been set at initialization, otherwise an exception will be raised).

sampling_grid2worldarray, shape (dim + 1, dim + 1), optional

the grid-to-world transform associated with the sampling grid (specified by sampling_grid_shape, or by default self.codomain_shape). If None (the default), then the grid-to-world transform is assumed to be the identity.

resample_onlyBoolean, optional

If False (the default) the affine transform is applied normally. If True, then the affine transform is not applied, and the input image is just re-sampled on the domain grid of this transform.

Returns:
transformedarray, shape sampling_grid_shape or

self.codomain_shape

the transformed image, sampled at the requested grid

transform_inverse(image, interpolation='linear', image_grid2world=None, sampling_grid_shape=None, sampling_grid2world=None, resample_only=False)

Transform the input image from domain to co-domain space.

By default, the transformed image is sampled at a grid defined by self.codomain_shape and self.codomain_grid2world. If such information was not provided then sampling_grid_shape is mandatory.

Parameters:
image2D or 3D array

the image to be transformed

interpolationstring, either ‘linear’ or ‘nearest’

the type of interpolation to be used, either ‘linear’ (for k-linear interpolation) or ‘nearest’ for nearest neighbor

image_grid2worldarray, shape (dim + 1, dim + 1), optional

the grid-to-world transform associated with image. If None (the default), then the grid-to-world transform is assumed to be the identity.

sampling_grid_shapesequence, shape (dim,), optional

the shape of the grid where the transformed image must be sampled. If None (the default), then self.codomain_shape is used instead (which must have been set at initialization, otherwise an exception will be raised).

sampling_grid2worldarray, shape (dim + 1, dim + 1), optional

the grid-to-world transform associated with the sampling grid (specified by sampling_grid_shape, or by default self.codomain_shape). If None (the default), then the grid-to-world transform is assumed to be the identity.

resample_onlyBoolean, optional

If False (the default) the affine transform is applied normally. If True, then the affine transform is not applied, and the input image is just re-sampled on the domain grid of this transform.

Returns:
transformedarray, shape sampling_grid_shape or

self.codomain_shape

the transformed image, sampled at the requested grid

### AffineRegistration

class dipy.align._public.AffineRegistration(metric=None, level_iters=None, sigmas=None, factors=None, method='L-BFGS-B', ss_sigma_factor=None, options=None, verbosity=1)

Bases: object

Methods

 optimize(static, moving, transform, params0) Start the optimization process.
__init__(metric=None, level_iters=None, sigmas=None, factors=None, method='L-BFGS-B', ss_sigma_factor=None, options=None, verbosity=1)

Initialize an instance of the AffineRegistration class.

Parameters:
metricNone or object, optional

an instance of a metric. The default is None, implying the Mutual Information metric with default settings.

the number of iterations at each scale of the scale space. level_iters[0] corresponds to the coarsest scale, level_iters[-1] the finest, where n is the length of the sequence. By default, a 3-level scale space with iterations sequence equal to [10000, 1000, 100] will be used.

sigmassequence of floats, optional

custom smoothing parameter to build the scale space (one parameter for each scale). By default, the sequence of sigmas will be [3, 1, 0].

custom scale factors to build the scale space (one factor for each scale). By default, the sequence of factors will be [4, 2, 1].

methodstring, optional

optimization method to be used. If Scipy version < 0.12, then only L-BFGS-B is available. Otherwise, method can be any gradient-based method available in dipy.core.Optimize: CG, BFGS, Newton-CG, dogleg or trust-ncg. The default is ‘L-BFGS-B’.

ss_sigma_factorfloat, optional

If None, this parameter is not used and an isotropic scale space with the given factors and sigmas will be built. If not None, an anisotropic scale space will be used by automatically selecting the smoothing sigmas along each axis according to the voxel dimensions of the given image. The ss_sigma_factor is used to scale the automatically computed sigmas. For example, in the isotropic case, the sigma of the kernel will be $$factor * (2 ^ i)$$ where $$i = 1, 2, ..., n_scales - 1$$ is the scale (the finest resolution image $$i=0$$ is never smoothed). The default is None.

optionsdict, optional

extra optimization options. The default is None, implying no extra options are passed to the optimizer.

verbosity: int (one of {0, 1, 2, 3}), optional

Set the verbosity level of the algorithm: 0 : do not print anything 1 : print information about the current status of the algorithm 2 : print high level information of the components involved in

the registration that can be used to detect a failing component.

3print as much information as possible to isolate the cause

of a bug.

Default: 1

docstring_addendum = 'verbosity: int (one of {0, 1, 2, 3}), optional\n            Set the verbosity level of the algorithm:\n            0 : do not print anything\n            1 : print information about the current status of the algorithm\n            2 : print high level information of the components involved in\n                the registration that can be used to detect a failing\n                component.\n            3 : print as much information as possible to isolate the cause\n                of a bug.\n            Default: 1\n    '

Start the optimization process.

Parameters:
static2D or 3D array

the image to be used as reference during optimization.

moving2D or 3D array

the image to be used as “moving” during optimization. It is necessary to pre-align the moving image to ensure its domain lies inside the domain of the deformation fields. This is assumed to be accomplished by “pre-aligning” the moving image towards the static using an affine transformation given by the ‘starting_affine’ matrix

transforminstance of Transform

the transformation with respect to whose parameters the gradient must be computed

params0array, shape (n,)

parameters from which to start the optimization. If None, the optimization will start at the identity transform. n is the number of parameters of the specified transformation.

static_grid2worldarray, shape (dim+1, dim+1), optional

the voxel-to-space transformation associated with the static image. The default is None, implying the transform is the identity.

moving_grid2worldarray, shape (dim+1, dim+1), optional

the voxel-to-space transformation associated with the moving image. The default is None, implying the transform is the identity.

starting_affinestring, or matrix, or None, optional
If string:

‘mass’: align centers of gravity ‘voxel-origin’: align physical coordinates of voxel (0,0,0) ‘centers’: align physical coordinates of central voxels

If matrix:

array, shape (dim+1, dim+1).

If None:

Start from identity.

The default is None.

ret_metricboolean, optional

if True, it returns the parameters for measuring the similarity between the images (default ‘False’). The metric containing optimal parameters and the distance between the images.

static_maskarray, shape (S, R, C) or (R, C), optional

static image mask that defines which pixels in the static image are used to calculate the mutual information.

moving_maskarray, shape (S’, R’, C’) or (R’, C’), optional

moving image mask that defines which pixels in the moving image are used to calculate the mutual information.

Returns:
affine_mapinstance of AffineMap

the affine resulting affine transformation

xoptoptimal parameters

the optimal parameters (translation, rotation shear etc.)

foptSimilarity metric

the value of the function at the optimal parameters.

### AffineTransform3D

class dipy.align._public.AffineTransform3D

Bases: Transform

Methods

 get_identity_parameters Parameter values corresponding to the identity transform jacobian Jacobian function of this transform param_to_matrix Matrix representation of this transform with the given parameters
 get_dim get_number_of_parameters
__init__()

Affine transform in 3D

### CCMetric

Methods

 Computes one step bringing the static image towards the moving. Computes one step bringing the moving image towards the static. Frees the resources allocated during initialization Numerical value assigned by this metric to the current image pair Prepares the metric to compute one displacement field iteration. set_levels_above(levels) Informs the metric how many pyramid levels are above the current one set_levels_below(levels) Informs the metric how many pyramid levels are below the current one set_moving_image(moving_image, ...) Sets the moving image being compared against the static one. set_static_image(static_image, ...) Sets the static image being compared against the moving one. use_moving_image_dynamics(...) This is called by the optimizer just after setting the moving image use_static_image_dynamics(...) This is called by the optimizer just after setting the static image.

Normalized Cross-Correlation Similarity metric.

Parameters:
dimint (either 2 or 3)

the dimension of the image domain

sigma_diffthe standard deviation of the Gaussian smoothing kernel to

be applied to the update field at each iteration

the radius of the squared (cubic) neighborhood at each voxel to be considered to compute the cross correlation

compute_backward()

Computes one step bringing the static image towards the moving.

Computes the update displacement field to be used for registration of the static image towards the moving image

compute_forward()

Computes one step bringing the moving image towards the static.

Computes the update displacement field to be used for registration of the moving image towards the static image

free_iteration()

Frees the resources allocated during initialization

get_energy()

Numerical value assigned by this metric to the current image pair

Returns the Cross Correlation (data term) energy computed at the largest iteration

initialize_iteration()

Prepares the metric to compute one displacement field iteration.

Pre-computes the cross-correlation factors for efficient computation of the gradient of the Cross Correlation w.r.t. the displacement field. It also pre-computes the image gradients in the physical space by re-orienting the gradients in the voxel space using the corresponding affine transformations.

### DiffeomorphicMap

class dipy.align._public.DiffeomorphicMap(dim, disp_shape, disp_grid2world=None, domain_shape=None, domain_grid2world=None, codomain_shape=None, codomain_grid2world=None, prealign=None)

Bases: object

Methods

 Creates a zero displacement field Inversion error of the displacement fields expand_fields(expand_factors, new_shape) Expands the displacement fields from current shape to new_shape Deformation field to transform an image in the backward direction Deformation field to transform an image in the forward direction Constructs a simplified version of this Diffeomorhic Map Try to interpret obj as a matrix Inverse of this DiffeomorphicMap instance Shallow copy of this DiffeomorphicMap instance transform(image[, interpolation, ...]) Warps an image in the forward direction transform_inverse(image[, interpolation, ...]) Warps an image in the backward direction Composition of this DiffeomorphicMap with a given endomorphism
 transform_points transform_points_inverse
__init__(dim, disp_shape, disp_grid2world=None, domain_shape=None, domain_grid2world=None, codomain_shape=None, codomain_grid2world=None, prealign=None)

DiffeomorphicMap

Implements a diffeomorphic transformation on the physical space. The deformation fields encoding the direct and inverse transformations share the same domain discretization (both the discretization grid shape and voxel-to-space matrix). The input coordinates (physical coordinates) are first aligned using prealign, and then displaced using the corresponding vector field interpolated at the aligned coordinates.

Parameters:
dimint, 2 or 3

the transformation’s dimension

disp_shapearray, shape (dim,)

the number of slices (if 3D), rows and columns of the deformation field’s discretization

disp_grid2worldthe voxel-to-space transform between the def. fields

grid and space

domain_shapearray, shape (dim,)

the number of slices (if 3D), rows and columns of the default discretization of this map’s domain

domain_grid2worldarray, shape (dim+1, dim+1)

the default voxel-to-space transformation between this map’s discretization and physical space

codomain_shapearray, shape (dim,)

the number of slices (if 3D), rows and columns of the images that are ‘normally’ warped using this transformation in the forward direction (this will provide default transformation parameters to warp images under this transformation). By default, we assume that the inverse transformation is ‘normally’ used to warp images with the same discretization and voxel-to-space transformation as the deformation field grid.

codomain_grid2worldarray, shape (dim+1, dim+1)

the voxel-to-space transformation of images that are ‘normally’ warped using this transformation (in the forward direction).

prealignarray, shape (dim+1, dim+1)

the linear transformation to be applied to align input images to the reference space before warping under the deformation field.

allocate()

Creates a zero displacement field

Creates a zero displacement field (the identity transformation).

compute_inversion_error()

Inversion error of the displacement fields

Estimates the inversion error of the displacement fields by computing statistics of the residual vectors obtained after composing the forward and backward displacement fields.

Returns:
residualarray, shape (R, C) or (S, R, C)

the displacement field resulting from composing the forward and backward displacement fields of this transformation (the residual should be zero for a perfect diffeomorphism)

statsarray, shape (3,)

statistics from the norms of the vectors of the residual displacement field: maximum, mean and standard deviation

Notes

Since the forward and backward displacement fields have the same discretization, the final composition is given by

comp[i] = forward[ i + Dinv * backward[i]]

where Dinv is the space-to-grid transformation of the displacement fields

expand_fields(expand_factors, new_shape)

Expands the displacement fields from current shape to new_shape

Up-samples the discretization of the displacement fields to be of new_shape shape.

Parameters:
expand_factorsarray, shape (dim,)

the factors scaling current spacings (voxel sizes) to spacings in the expanded discretization.

new_shapearray, shape (dim,)

the shape of the arrays holding the up-sampled discretization

get_backward_field()

Deformation field to transform an image in the backward direction

Returns the deformation field that must be used to warp an image under this transformation in the backward direction (note the ‘is_inverse’ flag).

get_forward_field()

Deformation field to transform an image in the forward direction

Returns the deformation field that must be used to warp an image under this transformation in the forward direction (note the ‘is_inverse’ flag).

get_simplified_transform()

Constructs a simplified version of this Diffeomorhic Map

The simplified version incorporates the pre-align transform, as well as the domain and codomain affine transforms into the displacement field. The resulting transformation may be regarded as operating on the image spaces given by the domain and codomain discretization. As a result, self.prealign, self.disp_grid2world, self.domain_grid2world and self.codomain affine will be None (denoting Identity) in the resulting diffeomorphic map.

interpret_matrix(obj)

Try to interpret obj as a matrix

Some operations are performed faster if we know in advance if a matrix is the identity (so we can skip the actual matrix-vector multiplication). This function returns None if the given object is None or the ‘identity’ string. It returns the same object if it is a numpy array. It raises an exception otherwise.

Parameters:
objobject

any object

Returns:
objobject

the same object given as argument if obj is None or a numpy array. None if obj is the ‘identity’ string.

inverse()

Inverse of this DiffeomorphicMap instance

Returns a diffeomorphic map object representing the inverse of this transformation. The internal arrays are not copied but just referenced.

Returns:
invDiffeomorphicMap object

the inverse of this diffeomorphic map.

shallow_copy()

Shallow copy of this DiffeomorphicMap instance

Creates a shallow copy of this diffeomorphic map (the arrays are not copied but just referenced)

Returns:
new_mapDiffeomorphicMap object

the shallow copy of this diffeomorphic map

transform(image, interpolation='linear', image_world2grid=None, out_shape=None, out_grid2world=None)

Warps an image in the forward direction

Transforms the input image under this transformation in the forward direction. It uses the “is_inverse” flag to switch between “forward” and “backward” (if is_inverse is False, then transform(…) warps the image forwards, else it warps the image backwards).

Parameters:
imagearray, shape (s, r, c) if dim = 3 or (r, c) if dim = 2

the image to be warped under this transformation in the forward direction

interpolationstring, either ‘linear’ or ‘nearest’

the type of interpolation to be used for warping, either ‘linear’ (for k-linear interpolation) or ‘nearest’ for nearest neighbor

image_world2gridarray, shape (dim+1, dim+1)

the transformation bringing world (space) coordinates to voxel coordinates of the image given as input

out_shapearray, shape (dim,)

the number of slices, rows and columns of the desired warped image

out_grid2worldthe transformation bringing voxel coordinates of the

warped image to physical space

Returns:
warpedarray, shape = out_shape or self.codomain_shape if None

the warped image under this transformation in the forward direction

Notes

See _warp_forward and _warp_backward documentation for further information.

transform_inverse(image, interpolation='linear', image_world2grid=None, out_shape=None, out_grid2world=None)

Warps an image in the backward direction

Transforms the input image under this transformation in the backward direction. It uses the “is_inverse” flag to switch between “forward” and “backward” (if is_inverse is False, then transform_inverse(…) warps the image backwards, else it warps the image forwards)

Parameters:
imagearray, shape (s, r, c) if dim = 3 or (r, c) if dim = 2

the image to be warped under this transformation in the forward direction

interpolationstring, either ‘linear’ or ‘nearest’

the type of interpolation to be used for warping, either ‘linear’ (for k-linear interpolation) or ‘nearest’ for nearest neighbor

image_world2gridarray, shape (dim+1, dim+1)

the transformation bringing world (space) coordinates to voxel coordinates of the image given as input

out_shapearray, shape (dim,)

the number of slices, rows, and columns of the desired warped image

out_grid2worldthe transformation bringing voxel coordinates of the

warped image to physical space

Returns:
warpedarray, shape = out_shape or self.codomain_shape if None

warped image under this transformation in the backward direction

Notes

See _warp_forward and _warp_backward documentation for further information.

transform_points(points, coord2world=None, world2coord=None)
transform_points_inverse(points, coord2world=None, world2coord=None)
warp_endomorphism(phi)

Composition of this DiffeomorphicMap with a given endomorphism

Creates a new DiffeomorphicMap C with the same properties as self and composes its displacement fields with phi’s corresponding fields. The resulting diffeomorphism is of the form C(x) = phi(self(x)) with inverse C^{-1}(y) = self^{-1}(phi^{-1}(y)). We assume that phi is an endomorphism with the same discretization and domain affine as self to ensure that the composition inherits self’s properties (we also assume that the pre-aligning matrix of phi is None or identity).

Parameters:
phiDiffeomorphicMap object

the endomorphism to be warped by this diffeomorphic map

Returns:
compositionthe composition of this diffeomorphic map with the

endomorphism given as input

Notes

The problem with our current representation of a DiffeomorphicMap is that the set of Diffeomorphism that can be represented this way (a pre-aligning matrix followed by a non-linear endomorphism given as a displacement field) is not closed under the composition operation.

Supporting a general DiffeomorphicMap class, closed under composition, may be extremely costly computationally, and the kind of transformations we actually need for Avants’ mid-point algorithm (SyN) are much simpler.

### EMMetric

class dipy.align._public.EMMetric(dim, smooth=1.0, inner_iter=5, q_levels=256, double_gradient=True, step_type='gauss_newton')

Methods

 Computes one step bringing the static image towards the moving. compute_demons_step([forward_step]) Demons step for EM metric Computes one step bringing the reference image towards the static. compute_gauss_newton_step([forward_step]) Computes the Gauss-Newton energy minimization step Frees the resources allocated during initialization The numerical value assigned by this metric to the current image pair Prepares the metric to compute one displacement field iteration. set_levels_above(levels) Informs the metric how many pyramid levels are above the current one set_levels_below(levels) Informs the metric how many pyramid levels are below the current one set_moving_image(moving_image, ...) Sets the moving image being compared against the static one. set_static_image(static_image, ...) Sets the static image being compared against the moving one. This is called by the optimizer just after setting the moving image. This is called by the optimizer just after setting the static image.
__init__(dim, smooth=1.0, inner_iter=5, q_levels=256, double_gradient=True, step_type='gauss_newton')

Expectation-Maximization Metric

Similarity metric based on the Expectation-Maximization algorithm to handle multi-modal images. The transfer function is modeled as a set of hidden random variables that are estimated at each iteration of the algorithm.

Parameters:
dimint (either 2 or 3)

the dimension of the image domain

smoothfloat

smoothness parameter, the larger the value the smoother the deformation field

inner_iterint

number of iterations to be performed at each level of the multi- resolution Gauss-Seidel optimization algorithm (this is not the number of steps per Gaussian Pyramid level, that parameter must be set for the optimizer, not the metric)

q_levelsnumber of quantization levels (equal to the number of hidden

variables in the EM algorithm)

if True, the gradient of the expected static image under the moving modality will be added to the gradient of the moving image, similarly, the gradient of the expected moving image under the static modality will be added to the gradient of the static image.

step_typestring (‘gauss_newton’, ‘demons’)

the optimization schedule to be used in the multi-resolution Gauss-Seidel optimization algorithm (not used if Demons Step is selected)

compute_backward()

Computes one step bringing the static image towards the moving.

Computes the update displacement field to be used for registration of the static image towards the moving image

compute_demons_step(forward_step=True)

Demons step for EM metric

Parameters:
forward_stepboolean

if True, computes the Demons step in the forward direction (warping the moving towards the static image). If False, computes the backward step (warping the static image to the moving image)

Returns:
displacementarray, shape (R, C, 2) or (S, R, C, 3)

the Demons step

compute_forward()

Computes one step bringing the reference image towards the static.

Computes the forward update field to register the moving image towards the static image in a gradient-based optimization algorithm

compute_gauss_newton_step(forward_step=True)

Computes the Gauss-Newton energy minimization step

Computes the Newton step to minimize this energy, i.e., minimizes the linearized energy function with respect to the regularized displacement field (this step does not require post-smoothing, as opposed to the demons step, which does not include regularization). To accelerate convergence we use the multi-grid Gauss-Seidel algorithm proposed by Bruhn and Weickert et al [Bruhn05]

Parameters:
forward_stepboolean

if True, computes the Newton step in the forward direction (warping the moving towards the static image). If False, computes the backward step (warping the static image to the moving image)

Returns:
displacementarray, shape (R, C, 2) or (S, R, C, 3)

the Newton step

References

[Bruhn05] Andres Bruhn and Joachim Weickert, “Towards ultimate motion

estimation: combining highest accuracy with real-time performance”, 10th IEEE International Conference on Computer Vision, 2005. ICCV 2005.

free_iteration()

Frees the resources allocated during initialization

get_energy()

The numerical value assigned by this metric to the current image pair

Returns the EM (data term) energy computed at the largest iteration

initialize_iteration()

Prepares the metric to compute one displacement field iteration.

Pre-computes the transfer functions (hidden random variables) and variances of the estimators. Also pre-computes the gradient of both input images. Note that once the images are transformed to the opposite modality, the gradient of the transformed images can be used with the gradient of the corresponding modality in the same fashion as diff-demons does for mono-modality images. If the flag self.use_double_gradient is True these gradients are averaged.

use_moving_image_dynamics(original_moving_image, transformation)

This is called by the optimizer just after setting the moving image.

EMMetric takes advantage of the image dynamics by computing the current moving image mask from the original_moving_image mask (warped by nearest neighbor interpolation)

Parameters:
original_moving_imagearray, shape (R, C) or (S, R, C)

the original moving image from which the current moving image was generated, the current moving image is the one that was provided via ‘set_moving_image(…)’, which may not be the same as the original moving image but a warped version of it.

transformationDiffeomorphicMap object

the transformation that was applied to the original_moving_image to generate the current moving image

use_static_image_dynamics(original_static_image, transformation)

This is called by the optimizer just after setting the static image.

EMMetric takes advantage of the image dynamics by computing the current static image mask from the originalstaticImage mask (warped by nearest neighbor interpolation)

Parameters:
original_static_imagearray, shape (R, C) or (S, R, C)

the original static image from which the current static image was generated, the current static image is the one that was provided via ‘set_static_image(…)’, which may not be the same as the original static image but a warped version of it (even the static image changes during Symmetric Normalization, not only the moving one).

transformationDiffeomorphicMap object

the transformation that was applied to the original_static_image to generate the current static image

### MutualInformationMetric

class dipy.align._public.MutualInformationMetric(nbins=32, sampling_proportion=None)

Bases: object

Methods

 distance(params) Numeric value of the negative Mutual Information. distance_and_gradient(params) Numeric value of the metric and its gradient at given parameters. gradient(params) Numeric value of the metric's gradient at the given parameters. setup(transform, static, moving[, ...]) Prepare the metric to compute intensity densities and gradients.
__init__(nbins=32, sampling_proportion=None)

Initialize an instance of the Mutual Information metric.

This class implements the methods required by Optimizer to drive the registration process.

Parameters:
nbinsint, optional

the number of bins to be used for computing the intensity histograms. The default is 32.

sampling_proportionNone or float in interval (0, 1], optional

There are two types of sampling: dense and sparse. Dense sampling uses all voxels for estimating the (joint and marginal) intensity histograms, while sparse sampling uses a subset of them. If sampling_proportion is None, then dense sampling is used. If sampling_proportion is a floating point value in (0,1] then sparse sampling is used, where sampling_proportion specifies the proportion of voxels to be used. The default is None.

Notes

Since we use linear interpolation, images are not, in general, differentiable at exact voxel coordinates, but they are differentiable between voxel coordinates. When using sparse sampling, selected voxels are slightly moved by adding a small random displacement within one voxel to prevent sampling points from being located exactly at voxel coordinates. When using dense sampling, this random displacement is not applied.

distance(params)

Numeric value of the negative Mutual Information.

We need to change the sign so we can use standard minimization algorithms.

Parameters:
paramsarray, shape (n,)

the parameter vector of the transform currently used by the metric (the transform name is provided when self.setup is called), n is the number of parameters of the transform

Returns:
neg_mifloat

the negative mutual information of the input images after transforming the moving image by the currently set transform with params parameters

Numeric value of the metric and its gradient at given parameters.

Parameters:
paramsarray, shape (n,)

the parameter vector of the transform currently used by the metric (the transform name is provided when self.setup is called), n is the number of parameters of the transform

Returns:
neg_mifloat

the negative mutual information of the input images after transforming the moving image by the currently set transform with params parameters

the gradient of the negative Mutual Information

Numeric value of the metric’s gradient at the given parameters.

Parameters:
paramsarray, shape (n,)

the parameter vector of the transform currently used by the metric (the transform name is provided when self.setup is called), n is the number of parameters of the transform

Returns:

the gradient of the negative Mutual Information

Prepare the metric to compute intensity densities and gradients.

The histograms will be setup to compute probability densities of intensities within the minimum and maximum values of static and moving

Parameters:
transform: instance of Transform

the transformation with respect to whose parameters the gradient must be computed

staticarray, shape (S, R, C) or (R, C)

static image

movingarray, shape (S’, R’, C’) or (R’, C’)

moving image. The dimensions of the static (S, R, C) and moving (S’, R’, C’) images do not need to be the same.

static_grid2worldarray (dim+1, dim+1), optional

the grid-to-space transform of the static image. The default is None, implying the transform is the identity.

moving_grid2worldarray (dim+1, dim+1)

the grid-to-space transform of the moving image. The default is None, implying the spacing along all axes is 1.

starting_affinearray, shape (dim+1, dim+1), optional

the pre-aligning matrix (an affine transform) that roughly aligns the moving image towards the static image. If None, no pre-alignment is performed. If a pre-alignment matrix is available, it is recommended to provide this matrix as starting_affine instead of manually transforming the moving image to reduce interpolation artifacts. The default is None, implying no pre-alignment is performed.

static_maskarray, shape (S, R, C) or (R, C), optional

static image mask that defines which pixels in the static image are used to calculate the mutual information.

moving_maskarray, shape (S’, R’, C’) or (R’, C’), optional

moving image mask that defines which pixels in the moving image are used to calculate the mutual information.

### RigidIsoScalingTransform3D

class dipy.align._public.RigidIsoScalingTransform3D

Bases: Transform

Methods

 get_identity_parameters Parameter values corresponding to the identity transform jacobian Jacobian function of this transform param_to_matrix Matrix representation of this transform with the given parameters
 get_dim get_number_of_parameters
__init__()

Rigid isoscaling transform in 3D.

(rotation + translation + scaling) The parameter vector theta of length 7 is interpreted as follows: theta[0] : rotation about the x axis theta[1] : rotation about the y axis theta[2] : rotation about the z axis theta[3] : translation along the x axis theta[4] : translation along the y axis theta[5] : translation along the z axis theta[6] : isotropic scaling

### RigidScalingTransform3D

class dipy.align._public.RigidScalingTransform3D

Bases: Transform

Methods

 get_identity_parameters Parameter values corresponding to the identity transform jacobian Jacobian function of this transform param_to_matrix Matrix representation of this transform with the given parameters
 get_dim get_number_of_parameters
__init__()

Rigid Scaling transform in 3D (rotation + translation + scaling).

The parameter vector theta of length 9 is interpreted as follows: theta[0] : rotation about the x axis theta[1] : rotation about the y axis theta[2] : rotation about the z axis theta[3] : translation along the x axis theta[4] : translation along the y axis theta[5] : translation along the z axis theta[6] : scaling in the x axis theta[7] : scaling in the y axis theta[8] : scaling in the z axis

### RigidTransform3D

class dipy.align._public.RigidTransform3D

Bases: Transform

Methods

 get_identity_parameters Parameter values corresponding to the identity transform jacobian Jacobian function of this transform param_to_matrix Matrix representation of this transform with the given parameters
 get_dim get_number_of_parameters
__init__()

Rigid transform in 3D (rotation + translation) The parameter vector theta of length 6 is interpreted as follows: theta[0] : rotation about the x axis theta[1] : rotation about the y axis theta[2] : rotation about the z axis theta[3] : translation along the x axis theta[4] : translation along the y axis theta[5] : translation along the z axis

### SSDMetric

class dipy.align._public.SSDMetric(dim, smooth=4, inner_iter=10, step_type='demons')

Methods

 Computes one step bringing the static image towards the moving. compute_demons_step([forward_step]) Demons step for SSD metric Computes one step bringing the reference image towards the static. compute_gauss_newton_step([forward_step]) Computes the Gauss-Newton energy minimization step Nothing to free for the SSD metric The numerical value assigned by this metric to the current image pair Prepares the metric to compute one displacement field iteration. set_levels_above(levels) Informs the metric how many pyramid levels are above the current one set_levels_below(levels) Informs the metric how many pyramid levels are below the current one set_moving_image(moving_image, ...) Sets the moving image being compared against the static one. set_static_image(static_image, ...) Sets the static image being compared against the moving one. use_moving_image_dynamics(...) This is called by the optimizer just after setting the moving image use_static_image_dynamics(...) This is called by the optimizer just after setting the static image.
__init__(dim, smooth=4, inner_iter=10, step_type='demons')

Sum of Squared Differences (SSD) Metric

Similarity metric for (mono-modal) nonlinear image registration defined by the sum of squared differences (SSD)

Parameters:
dimint (either 2 or 3)

the dimension of the image domain

smoothfloat

smoothness parameter, the larger the value the smoother the deformation field

inner_iterint

number of iterations to be performed at each level of the multi- resolution Gauss-Seidel optimization algorithm (this is not the number of steps per Gaussian Pyramid level, that parameter must be set for the optimizer, not the metric)

step_typestring

the displacement field step to be computed when ‘compute_forward’ and ‘compute_backward’ are called. Either ‘demons’ or ‘gauss_newton’

compute_backward()

Computes one step bringing the static image towards the moving.

Computes the updated displacement field to be used for registration of the static image towards the moving image

compute_demons_step(forward_step=True)

Demons step for SSD metric

Computes the demons step proposed by Vercauteren et al.[Vercauteren09] for the SSD metric.

Parameters:
forward_stepboolean

if True, computes the Demons step in the forward direction (warping the moving towards the static image). If False, computes the backward step (warping the static image to the moving image)

Returns:
displacementarray, shape (R, C, 2) or (S, R, C, 3)

the Demons step

References

[Vercauteren09] Tom Vercauteren, Xavier Pennec, Aymeric Perchant,

Nicholas Ayache, “Diffeomorphic Demons: Efficient Non-parametric Image Registration”, Neuroimage 2009

compute_forward()

Computes one step bringing the reference image towards the static.

Computes the update displacement field to be used for registration of the moving image towards the static image

compute_gauss_newton_step(forward_step=True)

Computes the Gauss-Newton energy minimization step

Minimizes the linearized energy function (Newton step) defined by the sum of squared differences of corresponding pixels of the input images with respect to the displacement field.

Parameters:
forward_stepboolean

if True, computes the Newton step in the forward direction (warping the moving towards the static image). If False, computes the backward step (warping the static image to the moving image)

Returns:
displacementarray, shape = static_image.shape + (3,)

if forward_step==True, the forward SSD Gauss-Newton step, else, the backward step

free_iteration()

Nothing to free for the SSD metric

get_energy()

The numerical value assigned by this metric to the current image pair

Returns the Sum of Squared Differences (data term) energy computed at the largest iteration

initialize_iteration()

Prepares the metric to compute one displacement field iteration.

Pre-computes the gradient of the input images to be used in the computation of the forward and backward steps.

### StreamlineLinearRegistration

class dipy.align._public.StreamlineLinearRegistration(metric=None, x0='rigid', method='L-BFGS-B', bounds=None, verbose=False, options=None, evolution=False, num_threads=None)

Bases: object

Methods

 optimize(static, moving[, mat]) Find the minimum of the provided metric.
__init__(metric=None, x0='rigid', method='L-BFGS-B', bounds=None, verbose=False, options=None, evolution=False, num_threads=None)

Linear registration of 2 sets of streamlines [Garyfallidis15].

Parameters:
metricStreamlineDistanceMetric,

If None and fast is False then the BMD distance is used. If fast is True then a faster implementation of BMD is used. Otherwise, use the given distance metric.

x0array or int or str

Initial parametrization for the optimization.

If 1D array with:

a) 6 elements then only rigid registration is performed with the 3 first elements for translation and 3 for rotation. b) 7 elements also isotropic scaling is performed (similarity). c) 12 elements then translation, rotation (in degrees), scaling and shearing is performed (affine).

Here is an example of x0 with 12 elements: x0=np.array([0, 10, 0, 40, 0, 0, 2., 1.5, 1, 0.1, -0.5, 0])

This has translation (0, 10, 0), rotation (40, 0, 0) in degrees, scaling (2., 1.5, 1) and shearing (0.1, -0.5, 0).

If int:
1. 6

x0 = np.array([0, 0, 0, 0, 0, 0])

2. 7

x0 = np.array([0, 0, 0, 0, 0, 0, 1.])

3. 12

x0 = np.array([0, 0, 0, 0, 0, 0, 1., 1., 1, 0, 0, 0])

If str:
1. “rigid”

x0 = np.array([0, 0, 0, 0, 0, 0])

2. “similarity”

x0 = np.array([0, 0, 0, 0, 0, 0, 1.])

3. “affine”

x0 = np.array([0, 0, 0, 0, 0, 0, 1., 1., 1, 0, 0, 0])

methodstr,

‘L_BFGS_B’ or ‘Powell’ optimizers can be used. Default is ‘L_BFGS_B’.

boundslist of tuples or None,

If method == ‘L_BFGS_B’ then we can use bounded optimization. For example for the six parameters of rigid rotation we can set the bounds = [(-30, 30), (-30, 30), (-30, 30),

(-45, 45), (-45, 45), (-45, 45)]

That means that we have set the bounds for the three translations and three rotation axes (in degrees).

verbosebool, optional.

If True, if True then information about the optimization is shown. Default: False.

optionsNone or dict,

Extra options to be used with the selected method.

evolutionboolean

If True save the transformation for each iteration of the optimizer. Default is False. Supported only with Scipy >= 0.11.

Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus |num_threads + 1| is used (enter -1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.

References

Garyfallidis et al. “Robust and efficient linear registration of white-matter fascicles in the space of streamlines”, NeuroImage, 117, 124–140, 2015

[Garyfallidis14]

Garyfallidis et al., “Direct native-space fiber bundle alignment for group comparisons”, ISMRM, 2014.

[Garyfallidis17]

Garyfallidis et al. Recognition of white matter bundles using local and global streamline-based registration and clustering, Neuroimage, 2017.

optimize(static, moving, mat=None)

Find the minimum of the provided metric.

Parameters:
staticstreamlines

Reference or fixed set of streamlines.

movingstreamlines

Moving set of streamlines.

matarray

Transformation (4, 4) matrix to start the registration. mat is applied to moving. Default value None which means that initial transformation will be generated by shifting the centers of moving and static sets of streamlines to the origin.

Returns:
mapStreamlineRegistrationMap

### SymmetricDiffeomorphicRegistration

class dipy.align._public.SymmetricDiffeomorphicRegistration(metric, level_iters=None, step_length=0.25, ss_sigma_factor=0.2, opt_tol=1e-05, inv_iter=20, inv_tol=0.001, callback=None)

Methods

 Return the resulting diffeomorphic map. optimize(static, moving[, ...]) Starts the optimization set_level_iters(level_iters) Sets the number of iterations at each pyramid level update(current_displacement, ...) Composition of the current displacement field with the given field
__init__(metric, level_iters=None, step_length=0.25, ss_sigma_factor=0.2, opt_tol=1e-05, inv_iter=20, inv_tol=0.001, callback=None)

Symmetric Diffeomorphic Registration (SyN) Algorithm

Performs the multi-resolution optimization algorithm for non-linear registration using a given similarity metric.

Parameters:
metricSimilarityMetric object

the metric to be optimized

level_iterslist of int

the number of iterations at each level of the Gaussian Pyramid (the length of the list defines the number of pyramid levels to be used)

opt_tolfloat

the optimization will stop when the estimated derivative of the energy profile w.r.t. time falls below this threshold

inv_iterint

the number of iterations to be performed by the displacement field inversion algorithm

step_lengthfloat

the length of the maximum displacement vector of the update displacement field at each iteration

ss_sigma_factorfloat

parameter of the scale-space smoothing kernel. For example, the std. dev. of the kernel will be factor*(2^i) in the isotropic case where i = 0, 1, …, n_scales is the scale

inv_tolfloat

the displacement field inversion algorithm will stop iterating when the inversion error falls below this threshold

callbackfunction(SymmetricDiffeomorphicRegistration)

a function receiving a SymmetricDiffeomorphicRegistration object to be called after each iteration (this optimizer will call this function passing self as parameter)

get_map()

Return the resulting diffeomorphic map.

Returns the DiffeomorphicMap registering the moving image towards the static image.

optimize(static, moving, static_grid2world=None, moving_grid2world=None, prealign=None)

Starts the optimization

Parameters:
staticarray, shape (S, R, C) or (R, C)

the image to be used as reference during optimization. The displacement fields will have the same discretization as the static image.

movingarray, shape (S, R, C) or (R, C)

the image to be used as “moving” during optimization. Since the deformation fields’ discretization is the same as the static image, it is necessary to pre-align the moving image to ensure its domain lies inside the domain of the deformation fields. This is assumed to be accomplished by “pre-aligning” the moving image towards the static using an affine transformation given by the ‘prealign’ matrix

static_grid2worldarray, shape (dim+1, dim+1)

the voxel-to-space transformation associated to the static image

moving_grid2worldarray, shape (dim+1, dim+1)

the voxel-to-space transformation associated to the moving image

prealignarray, shape (dim+1, dim+1)

the affine transformation (operating on the physical space) pre-aligning the moving image towards the static

Returns:
static_to_refDiffeomorphicMap object

the diffeomorphic map that brings the moving image towards the static one in the forward direction (i.e. by calling static_to_ref.transform) and the static image towards the moving one in the backward direction (i.e. by calling static_to_ref.transform_inverse).

update(current_displacement, new_displacement, disp_world2grid, time_scaling)

Composition of the current displacement field with the given field

Interpolates new displacement at the locations defined by current_displacement. Equivalently, computes the composition C of the given displacement fields as C(x) = B(A(x)), where A is current_displacement and B is new_displacement. This function is intended to be used with deformation fields of the same sampling (e.g. to be called by a registration algorithm).

Parameters:
current_displacementarray, shape (R’, C’, 2) or (S’, R’, C’, 3)

the displacement field defining where to interpolate new_displacement

new_displacementarray, shape (R, C, 2) or (S, R, C, 3)

the displacement field to be warped by current_displacement

disp_world2gridarray, shape (dim+1, dim+1)

the space-to-grid transform associated with the displacements’ grid (we assume that both displacements are discretized over the same grid)

time_scalingfloat

scaling factor applied to d2. The effect may be interpreted as moving d1 displacements along a factor (time_scaling) of d2.

Returns:
updatedarray, shape (the same as new_displacement)

the warped displacement field

mean_normthe mean norm of all vectors in current_displacement

### TranslationTransform3D

class dipy.align._public.TranslationTransform3D

Bases: Transform

Methods

 get_identity_parameters Parameter values corresponding to the identity transform jacobian Jacobian function of this transform param_to_matrix Matrix representation of this transform with the given parameters
 get_dim get_number_of_parameters
__init__()

Translation transform in 3D

### partial

class dipy.align._public.partial

Bases: object

partial(func, *args, **keywords) - new function with partial application of the given arguments and keywords.

Attributes:
args

tuple of arguments to future partial calls

func

function object to use in future partial calls

keywords

dictionary of keyword arguments to future partial calls

Methods

 __call__(*args, **kwargs) Call self as a function.
__init__(*args, **kwargs)
args

tuple of arguments to future partial calls

func

function object to use in future partial calls

keywords

dictionary of keyword arguments to future partial calls

### affine_registration

dipy.align._public.affine_registration(moving, static, moving_affine=None, static_affine=None, pipeline=None, starting_affine=None, metric='MI', level_iters=None, sigmas=None, factors=None, ret_metric=False, moving_mask=None, static_mask=None, **metric_kwargs)

Find the affine transformation between two 3D images. Alternatively, find the combination of several linear transformations.

Parameters:
movingarray, nifti image or str

Containing the data for the moving object, or full path to a nifti file with the moving data.

staticarray, nifti image or str

Containing the data for the static object, or full path to a nifti file with the moving data.

moving_affine4x4 array, optional

An affine transformation associated with the moving object. Required if data is provided as an array. If provided together with nifti/path, will over-ride the affine that is in the nifti.

static_affine4x4 array, optional

An affine transformation associated with the static object. Required if data is provided as an array. If provided together with nifti/path, will over-ride the affine that is in the nifti.

pipelinelist of str, optional

Sequence of transforms to use in the gradual fitting. Default: gradual fit of the full affine (executed from left to right): ["center_of_mass", "translation", "rigid", "affine"] Alternatively, any other combination of the following registration methods might be used: center_of_mass, translation, rigid, rigid_isoscaling, rigid_scaling and affine.

starting_affine: 4x4 array, optional

Initial guess for the transformation between the spaces. Default: identity.

metricstr, optional.

Currently only supports ‘MI’ for MutualInformationMetric.

AffineRegistration key-word argument: the number of iterations at each scale of the scale space. level_iters[0] corresponds to the coarsest scale, level_iters[-1] the finest, where n is the length of the sequence. By default, a 3-level scale space with iterations sequence equal to [10000, 1000, 100] will be used.

sigmassequence of floats, optional

AffineRegistration key-word argument: custom smoothing parameter to build the scale space (one parameter for each scale). By default, the sequence of sigmas will be [3, 1, 0].

AffineRegistration key-word argument: custom scale factors to build the scale space (one factor for each scale). By default, the sequence of factors will be [4, 2, 1].

ret_metricboolean, optional

Set it to True to return the value of the optimized coefficients and the optimization quality metric.

moving_maskarray, shape (S’, R’, C’) or (R’, C’), optional

moving image mask that defines which pixels in the moving image are used to calculate the mutual information.

static_maskarray, shape (S, R, C) or (R, C), optional

static image mask that defines which pixels in the static image are used to calculate the mutual information.

nbinsint, optional

MutualInformationMetric key-word argument: the number of bins to be used for computing the intensity histograms. The default is 32.

sampling_proportionNone or float in interval (0, 1], optional

MutualInformationMetric key-word argument: There are two types of sampling: dense and sparse. Dense sampling uses all voxels for estimating the (joint and marginal) intensity histograms, while sparse sampling uses a subset of them. If sampling_proportion is None, then dense sampling is used. If sampling_proportion is a floating point value in (0,1] then sparse sampling is used, where sampling_proportion specifies the proportion of voxels to be used. The default is None (dense sampling).

Returns:
transformedarray with moving data resampled to the static space
after computing the affine transformation
affinethe affine 4x4 associated with the transformation.
xoptthe value of the optimized coefficients.
foptthe value of the optimization quality metric.

Notes

Performs a gradual registration between the two inputs, using a pipeline that gradually approximates the final registration. If the final default step (affine) is omitted, the resulting affine may not have all 12 degrees of freedom adjusted.

Load data and other information from a nifti file.

Parameters:
fnamestr

Full path to a nifti file.

return_imgbool, optional

Whether to return the nibabel nifti img object. Default: False

return_voxsize: bool, optional

Whether to return the nifti header zooms. Default: False

return_coordsbool, optional

Whether to return the nifti header aff2axcodes. Default: False

as_ndarray: bool, optional

convert nibabel ArrayProxy to a numpy.ndarray. If you want to save memory and delay this casting, just turn this option to False (default: True)

Returns:
A tuple, with (at the most, if all keyword args are set to True):
(data, img.affine, img, vox_size, nib.aff2axcodes(img.affine))

load_nifti_data

Load the stateful tractogram of the .trk format

Parameters:
filenamestring

Filename with valid extension

referenceNifti or Trk filename, Nifti1Image or TrkFile, Nifti1Header or

trk.header (dict), or ‘same’ if the input is a trk file. Reference that provides the spatial attribute. Typically a nifti-related object from the native diffusion used for streamlines generation

to_spaceEnum (dipy.io.stateful_tractogram.Space)

to_originEnum (dipy.io.stateful_tractogram.Origin)

NIFTI standard, default (center of the voxel) TRACKVIS standard (corner of the voxel)

bbox_valid_checkbool

Verification for negative voxel coordinates or values above the volume dimensions. Default is True, to enforce valid file.

Verification that the reference has the same header as the spatial attributes as the input tractogram when a Trk is loaded

Returns:
outputStatefulTractogram

The tractogram to load (must have been saved properly)

Helper function that handles inputs that can be paths, nifti img or arrays

Parameters:
dataarray or nib.Nifti1Image or str.

Either as a 3D/4D array or as a nifti image object, or as a string containing the full path to a nifti file.

affine4x4 array, optional.

Must be provided for data provided as an array. If provided together with Nifti1Image or str data, this input will over-ride the affine that is stored in the data input. Default: use the affine stored in data.

Returns:
data, affinendarray and 4x4 array

Read a syn registration mapping from a nifti file.

Parameters:
dispstr or Nifti1Image

A file of image containing the mapping displacement field in each voxel Shape (x, y, z, 3, 2)

domain_imgstr or Nifti1Image
codomain_imgstr or Nifti1Image
Returns:
A DiffeomorphicMap object.

Notes

See write_mapping() for the data format expected.

### register_dwi_series

dipy.align._public.register_dwi_series(data, gtab, affine=None, b0_ref=0, pipeline=None, static_mask=None)

Register a DWI series to the mean of the B0 images in that series.

all first registered to the first B0 volume

Parameters:
data4D array or nibabel Nifti1Image class instance or str

Diffusion data. Either as a 4D array or as a nifti image object, or as a string containing the full path to a nifti file.

gtaba GradientTable class instance or tuple of strings

If provided as a tuple of strings, these are assumed to be full paths to the bvals and bvecs files (in that order).

affine4x4 array, optional.

Must be provided for data provided as an array. If provided together with Nifti1Image or str data, this input will over-ride the affine that is stored in the data input. Default: use the affine stored in data.

b0_refint, optional.

Which b0 volume to use as reference. Default: 0

pipelinelist of callables, optional.

The transformations to perform in sequence (from left to right): Default: [center_of_mass, translation, rigid, affine]

static_maskarray, shape (S, R, C) or (R, C), optional

static image mask that defines which pixels in the static image are used to calculate the mutual information.

Returns:
xform_img, affine_array: a Nifti1Image containing the registered data and
using the affine of the original data and a list containing the affine
transforms associated with each of the

### register_dwi_to_template

dipy.align._public.register_dwi_to_template(dwi, gtab, dwi_affine=None, template=None, template_affine=None, reg_method='syn', **reg_kwargs)

Register DWI data to a template through the B0 volumes.

Parameters:
dwi4D array, nifti image or str

Containing the DWI data, or full path to a nifti file with DWI.

The gradients associated with the DWI data, or a sequence with (fbval, fbvec), full paths to bvals and bvecs files.

dwi_affine4x4 array, optional

An affine transformation associated with the DWI. Required if data is provided as an array. If provided together with nifti/path, will over-ride the affine that is in the nifti.

template3D array, nifti image or str

Containing the data for the template, or full path to a nifti file with the template data.

template_affine4x4 array, optional

An affine transformation associated with the template. Required if data is provided as an array. If provided together with nifti/path, will over-ride the affine that is in the nifti.

reg_methodstr,

One of “syn” or “aff”, which designates which registration method is used. Either syn, which uses the syn_registration() function or affine_registration() function. Default: “syn”.

reg_kwargskey-word arguments for syn_registration() or

affine_registration()

Returns:
warped_b0, mapping: The fist is an array with the b0 volume warped to the
template. If reg_method is “syn”, the second is a DiffeomorphicMap class
instance that can be used to transform between the two spaces. Otherwise,
if reg_method is “aff”, this is a 4x4 matrix encoding the affine transform.

Notes

This function assumes that the DWI data is already internally registered. See register_dwi_series().

### register_series

dipy.align._public.register_series(series, ref, pipeline=None, series_affine=None, ref_affine=None, static_mask=None)

Register a series to a reference image.

Parameters:
series4D array or nib.Nifti1Image class instance or str

The data is 4D with the last dimension separating different 3D volumes

refint or 3D array or nib.Nifti1Image class instance or str

If this is an int, this is the index of the reference image within the series. Otherwise it is an array of data to register with (associated with a ref_affine required) or a nifti img or full path to a file containing one.

pipelinesequence, optional

Sequence of transforms to do for each volume in the series. Default: (executed from left to right): [center_of_mass, translation, rigid, affine]

series_affine, ref_affine4x4 arrays, optional.

The affine. If provided, this input will over-ride the affine provided together with the nifti img or file.

static_maskarray, shape (S, R, C) or (R, C), optional

static image mask that defines which pixels in the static image are used to calculate the mutual information.

Returns:
xformed, affines4D array with transformed data and a (4,4,n) array
with 4x4 matrices associated with each of the volumes of the input moving
data that was used to transform it into register with the static data.

### resample

dipy.align._public.resample(moving, static, moving_affine=None, static_affine=None, between_affine=None)

Resample an image (moving) from one space to another (static).

Parameters:
movingarray, nifti image or str

Containing the data for the moving object, or full path to a nifti file with the moving data.

moving_affine4x4 array, optional

An affine transformation associated with the moving object. Required if data is provided as an array. If provided together with nifti/path, will over-ride the affine that is in the nifti.

staticarray, nifti image or str

Containing the data for the static object, or full path to a nifti file with the moving data.

static_affine4x4 array, optional

An affine transformation associated with the static object. Required if data is provided as an array. If provided together with nifti/path, will over-ride the affine that is in the nifti.

between_affine: 4x4 array, optional

If an additional affine is needed between the two spaces. Default: identity (no additional registration).

Returns:
A Nifti1Image class instance with the data from the moving object
resampled into the space of the static object.

### save_nifti

dipy.align._public.save_nifti(fname, data, affine, hdr=None, dtype=None)

Save a data array into a nifti file.

Parameters:
fnamestr

The full path to the file to be saved.

datandarray

The array with the data to save.

affine4x4 array

The affine transform associated with the file.

Returns:
None

### set_number_of_points

dipy.align._public.set_number_of_points()
Change the number of points of streamlines

(either by downsampling or upsampling)

Change the number of points of streamlines in order to obtain nb_points-1 segments of equal length. Points of streamlines will be modified along the curve.

Parameters:
streamlinesndarray or a list or dipy.tracking.Streamlines

If ndarray, must have shape (N,3) where N is the number of points of the streamline. If list, each item must be ndarray shape (Ni,3) where Ni is the number of points of streamline i. If dipy.tracking.Streamlines, its common_shape must be 3.

nb_pointsint

integer representing number of points wanted along the curve.

Returns:
new_streamlinesndarray or a list or dipy.tracking.Streamlines

Results of the downsampling or upsampling process.

Examples

>>> from dipy.tracking.streamline import set_number_of_points
>>> import numpy as np


One streamline, a semi-circle:

>>> theta = np.pi*np.linspace(0, 1, 100)
>>> x = np.cos(theta)
>>> y = np.sin(theta)
>>> z = 0 * x
>>> streamline = np.vstack((x, y, z)).T
>>> modified_streamline = set_number_of_points(streamline, 3)
>>> len(modified_streamline)
3


Multiple streamlines:

>>> streamlines = [streamline, streamline[::2]]
>>> new_streamlines = set_number_of_points(streamlines, 10)
>>> [len(s) for s in streamlines]
[100, 50]
>>> [len(s) for s in new_streamlines]
[10, 10]


### streamline_registration

dipy.align._public.streamline_registration(moving, static, n_points=100, native_resampled=False)

Register two collections of streamlines (‘bundles’) to each other.

Parameters:
moving, staticlists of 3 by n, or str

The two bundles to be registered. Given either as lists of arrays with 3D coordinates, or strings containing full paths to these files.

n_pointsint, optional

How many points to resample to. Default: 100.

native_resampledbool, optional

Whether to return the moving bundle in the original space, but resampled in the static space to n_points.

Returns:
alignedlist

Streamlines from the moving group, moved to be closely matched to the static group.

matrixarray (4, 4)

The affine transformation that takes us from ‘moving’ to ‘static’

### syn_registration

dipy.align._public.syn_registration(moving, static, moving_affine=None, static_affine=None, step_length=0.25, metric='CC', dim=3, level_iters=None, prealign=None, **metric_kwargs)

Register a 2D/3D source image (moving) to a 2D/3D target image (static).

Parameters:
moving, staticarray or nib.Nifti1Image or str.

Either as a 2D/3D array or as a nifti image object, or as a string containing the full path to a nifti file.

moving_affine, static_affine4x4 array, optional.

Must be provided for data provided as an array. If provided together with Nifti1Image or str data, this input will over-ride the affine that is stored in the data input. Default: use the affine stored in data.

metricstring, optional

The metric to be optimized. One of CC, EM, SSD, Default: ‘CC’ => CCMetric.

dim: int (either 2 or 3), optional

The dimensions of the image domain. Default: 3

level_iterslist of int, optional

the number of iterations at each level of the Gaussian Pyramid (the length of the list defines the number of pyramid levels to be used). Default: [10, 10, 5].

metric_kwargsdict, optional

Parameters for initialization of the metric object. If not provided, uses the default settings of each metric.

Returns:
warped_movingndarray

The data in moving, warped towards the static data.

forwardndarray (…, 3)

The vector field describing the forward warping from the source to the target.

backwardndarray (…, 3)

The vector field describing the backward warping from the target to the source.

### transform_centers_of_mass

dipy.align._public.transform_centers_of_mass(static, static_grid2world, moving, moving_grid2world)

Transformation to align the center of mass of the input images.

Parameters:
staticarray, shape (S, R, C)

static image

static_grid2worldarray, shape (dim+1, dim+1)

the voxel-to-space transformation of the static image

movingarray, shape (S, R, C)

moving image

moving_grid2worldarray, shape (dim+1, dim+1)

the voxel-to-space transformation of the moving image

Returns:
affine_mapinstance of AffineMap

the affine transformation (translation only, in this case) aligning the center of mass of the moving image towards the one of the static image

### transform_tracking_output

dipy.align._public.transform_tracking_output(tracking_output, affine, save_seeds=False)

Apply a linear transformation, given by affine, to streamlines.

Parameters:
streamlinesStreamlines generator

Either streamlines (list, ArraySequence) or a tuple with streamlines and seeds together

affinearray (4, 4)

The mapping between voxel indices and the point space for seeds. The voxel_to_rasmm matrix, typically from a NIFTI file.

save_seedsbool, optional

If set, seeds associated to streamlines will be also moved and returned

Returns
——-
streamlinesgenerator

A generator for the sequence of transformed streamlines. If save_seeds is True, also return a generator for the transformed seeds.

### warn

dipy.align._public.warn(/, message, category=None, stacklevel=1, source=None)

Issue a warning, or maybe ignore it or raise an exception.

### write_mapping

dipy.align._public.write_mapping(mapping, fname)

Write out a syn registration mapping to a nifti file.

Parameters:
mappinga DiffeomorphicMap object derived from syn_registration()
fnamestr

Full path to the nifti file storing the mapping

Notes

The data in the file is organized with shape (X, Y, Z, 3, 2), such that the forward mapping in each voxel is in data[i, j, k, :, 0] and the backward mapping in each voxel is in data[i, j, k, :, 1].

Determine the effective number of threads to be used for OpenMP calls

• For num_threads = None, - if the OMP_NUM_THREADS environment variable is set, return that value - otherwise, return the maximum number of cores retrieved by openmp.opm_get_num_procs().

• For num_threads > 0, return this value.

• For num_threads < 0, return the maximal number of threads minus |num_threads + 1|. In particular num_threads = -1 will use as many threads as there are available cores on the machine.

• For num_threads = 0 a ValueError is raised.

Parameters:

Desired number of threads to be used.

### distance_matrix_mdf

dipy.align.bundlemin.distance_matrix_mdf()

Minimum direct flipped distance matrix between two streamline sets

All streamlines need to have the same number of points

Parameters:
streamlines_asequence

of streamlines as arrays, [(N, 3) .. (N, 3)]

streamlines_bsequence

of streamlines as arrays, [(N, 3) .. (N, 3)]

Returns:
DMarray, shape (len(streamlines_a), len(streamlines_b))

distance matrix

### compute_cc_backward_step_2d

Gradient of the CC Metric w.r.t. the backward transformation

Computes the gradient of the Cross Correlation metric for symmetric registration (SyN) [Rb7e0e20c1a76-Avants2008] w.r.t. the displacement associated to the static image (‘forward’ step) as in [Rb7e0e20c1a76-Avants2011]

Parameters:

the gradient of the moving image

factorsarray, shape (R, C, 5)

the precomputed cross correlation terms obtained via precompute_cc_factors_2d

Returns:
outarray, shape (R, C, 2)

the gradient of the cross correlation metric with respect to the displacement associated to the static image

energythe cross correlation energy (data term) at this iteration

References

### compute_cc_backward_step_3d

Gradient of the CC Metric w.r.t. the backward transformation

Computes the gradient of the Cross Correlation metric for symmetric registration (SyN) [Avants08] w.r.t. the displacement associated to the static volume (‘backward’ step) as in [Avants11]

Parameters:
grad_movingarray, shape (S, R, C, 3)

the gradient of the moving volume

factorsarray, shape (S, R, C, 5)

the precomputed cross correlation terms obtained via precompute_cc_factors_3d

the radius of the neighborhood used for the CC metric when computing the factors. The returned vector field will be zero along a boundary of width radius voxels.

Returns:
outarray, shape (S, R, C, 3)

the gradient of the cross correlation metric with respect to the displacement associated to the static volume

energythe cross correlation energy (data term) at this iteration

References

[Avants08] Avants, B. B., Epstein, C. L., Grossman, M., & Gee, J. C. (2008)

Symmetric Diffeomorphic Image Registration with Cross-Correlation: Evaluating Automated Labeling of Elderly and Neurodegenerative Brain, Med Image Anal. 12(1), 26-41.

[Avants11] Avants, B. B., Tustison, N., & Song, G. (2011).

### compute_cc_forward_step_2d

Gradient of the CC Metric w.r.t. the forward transformation

Computes the gradient of the Cross Correlation metric for symmetric registration (SyN) [R1684b0862e45-Avants2008] w.r.t. the displacement associated to the moving image (‘backward’ step) as in [R1684b0862e45-Avants2011]

Parameters:

the gradient of the static image

factorsarray, shape (R, C, 5)

the precomputed cross correlation terms obtained via precompute_cc_factors_2d

Returns:
outarray, shape (R, C, 2)

the gradient of the cross correlation metric with respect to the displacement associated to the moving image

energythe cross correlation energy (data term) at this iteration

Notes

Currently, the gradient of the static image is not being used, but some authors suggest that symmetrizing the gradient by including both, the moving and static gradients may improve the registration quality. We are leaving this parameter as a placeholder for future investigation

References

### compute_cc_forward_step_3d

Gradient of the CC Metric w.r.t. the forward transformation

Computes the gradient of the Cross Correlation metric for symmetric registration (SyN) [R418c58b05724-Avants2008] w.r.t. the displacement associated to the moving volume (‘forward’ step) as in [R418c58b05724-Avants2011]

Parameters:
grad_staticarray, shape (S, R, C, 3)

the gradient of the static volume

factorsarray, shape (S, R, C, 5)

the precomputed cross correlation terms obtained via precompute_cc_factors_3d

the radius of the neighborhood used for the CC metric when computing the factors. The returned vector field will be zero along a boundary of width radius voxels.

Returns:
outarray, shape (S, R, C, 3)

the gradient of the cross correlation metric with respect to the displacement associated to the moving volume

energythe cross correlation energy (data term) at this iteration

References

### precompute_cc_factors_2d

Precomputations to quickly compute the gradient of the CC Metric

Pre-computes the separate terms of the cross correlation metric [Rc27aa752565c-Avants2008] and image norms at each voxel considering a neighborhood of the given radius to efficiently [Rc27aa752565c-Avants2011] compute the gradient of the metric with respect to the deformation field.

Parameters:
staticarray, shape (R, C)

the static volume, which also defines the reference registration domain

movingarray, shape (R, C)

the moving volume (notice that both images must already be in a common reference domain, i.e. the same R, C)

Returns:
factorsarray, shape (R, C, 5)

the precomputed cross correlation terms: factors[:,:,0] : static minus its mean value along the neighborhood factors[:,:,1] : moving minus its mean value along the neighborhood factors[:,:,2] : sum of the pointwise products of static and moving

along the neighborhood

factors[:,:,3] : sum of sq. values of static along the neighborhood factors[:,:,4] : sum of sq. values of moving along the neighborhood

References

### precompute_cc_factors_2d_test

Precomputations to quickly compute the gradient of the CC Metric

This version of precompute_cc_factors_2d is for testing purposes, it directly computes the local cross-correlation without any optimization.

### precompute_cc_factors_3d

Precomputations to quickly compute the gradient of the CC Metric

Pre-computes the separate terms of the cross correlation metric and image norms at each voxel considering a neighborhood of the given radius to efficiently compute the gradient of the metric with respect to the deformation field [Reda22c03a495-Ocegueda2016] [Reda22c03a495-Avants2008] [Reda22c03a495-Avants2011].

Parameters:
staticarray, shape (S, R, C)

the static volume, which also defines the reference registration domain

movingarray, shape (S, R, C)

the moving volume (notice that both images must already be in a common reference domain, i.e. the same S, R, C)

Returns:
factorsarray, shape (S, R, C, 5)

the precomputed cross correlation terms: factors[:,:,:,0] : static minus its mean value along the neighborhood factors[:,:,:,1] : moving minus its mean value along the neighborhood factors[:,:,:,2] : sum of the pointwise products of static and moving

along the neighborhood

factors[:,:,:,3] : sum of sq. values of static along the neighborhood factors[:,:,:,4] : sum of sq. values of moving along the neighborhood

References

### precompute_cc_factors_3d_test

Precomputations to quickly compute the gradient of the CC Metric

This version of precompute_cc_factors_3d is for testing purposes, it directly computes the local cross-correlation factors without any optimization, so it is less error-prone than the accelerated version.

### compute_em_demons_step_2d

Demons step for EM metric in 2D

Computes the demons step [Vercauteren09] for SSD-driven registration ( eq. 4 in [Vercauteren09] ) using the EM algorithm [Arce14] to handle multi-modality images.

In this case, $$\sigma_i$$ in eq. 4 of [Vercauteren] is estimated using the EM algorithm, while in the original version of diffeomorphic demons it is estimated by the difference between the image values at each pixel.

Parameters:
delta_fieldarray, shape (R, C)

contains, at each pixel, the difference between the moving image (warped under the current deformation s(. , .) ) J and the static image I: delta_field[i,j] = J(s(i,j)) - I(i,j). The order is important, changing to delta_field[i,j] = I(i,j) - J(s(i,j)) yields the backward demons step warping the static image towards the moving, which may not be the intended behavior unless the ‘gradient_moving’ passed corresponds to the gradient of the static image

sigma_sq_fieldarray, shape (R, C)

contains, at each pixel (i, j), the estimated variance (not std) of the hidden variable associated to the intensity at static[i,j] (which must have been previously quantized)

the gradient of the moving image

sigma_sq_xfloat

parameter controlling the amount of regularization. It corresponds to $$\sigma_x^2$$ in algorithm 1 of Vercauteren et al.[2]

outarray, shape (R, C, 2)

the resulting demons step will be written to this array

Returns:
demons_steparray, shape (R, C, 2)

the demons step to be applied for updating the current displacement field

energyfloat

the current em energy (before applying the returned demons_step)

References

[Arce14] Arce-santana, E., Campos-delgado, D. U., & Vigueras-g, F. (2014).

Non-rigid Multimodal Image Registration Based on the Expectation-Maximization Algorithm, (168140), 36-47.

[Vercauteren09] Vercauteren, T., Pennec, X., Perchant, A., & Ayache, N.

(2009). Diffeomorphic demons: efficient non-parametric image registration. NeuroImage, 45(1 Suppl), S61-72. doi:10.1016/j.neuroimage.2008.10.040

### compute_em_demons_step_3d

Demons step for EM metric in 3D

Computes the demons step [Vercauteren09] for SSD-driven registration ( eq. 4 in [Vercauteren09] ) using the EM algorithm [Arce14] to handle multi-modality images.

In this case, $$\sigma_i$$ in eq. 4 of [Vercauteren09] is estimated using the EM algorithm, while in the original version of diffeomorphic demons it is estimated by the difference between the image values at each pixel.

Parameters:
delta_fieldarray, shape (S, R, C)

contains, at each pixel, the difference between the moving image (warped under the current deformation s ) J and the static image I: delta_field[k,i,j] = J(s(k,i,j)) - I(k,i,j). The order is important, changing to delta_field[k,i,j] = I(k,i,j) - J(s(k,i,j)) yields the backward demons step warping the static image towards the moving, which may not be the intended behavior unless the ‘gradient_moving’ passed corresponds to the gradient of the static image

sigma_sq_fieldarray, shape (S, R, C)

contains, at each pixel (k, i, j), the estimated variance (not std) of the hidden variable associated to the intensity at static[k,i,j] (which must have been previously quantized)

gradient_movingarray, shape (S, R, C, 2)

the gradient of the moving image

sigma_sq_xfloat

parameter controlling the amount of regularization. It corresponds to $$\sigma_x^2$$ in algorithm 1 of Vercauteren et al.[2].

outarray, shape (S, R, C, 2)

the resulting demons step will be written to this array

Returns:
demons_steparray, shape (S, R, C, 3)

the demons step to be applied for updating the current displacement field

energyfloat

the current em energy (before applying the returned demons_step)

References

[Arce14] Arce-santana, E., Campos-delgado, D. U., & Vigueras-g, F. (2014).

Non-rigid Multimodal Image Registration Based on the Expectation-Maximization Algorithm, (168140), 36-47.

[Vercauteren09] Vercauteren, T., Pennec, X., Perchant, A., & Ayache, N.

(2009). Diffeomorphic demons: efficient non-parametric image registration. NeuroImage, 45(1 Suppl), S61-72. doi:10.1016/j.neuroimage.2008.10.040

Computes the mean and std. for each quantization level.

Computes the mean and standard deviation of the intensities in ‘v’ for each corresponding label in ‘labels’. In other words, for each label L, it computes the mean and standard deviation of the intensities in ‘v’ at pixels whose label in ‘labels’ is L. This is used by the EM metric to compute statistics for each hidden variable represented by the labels.

Parameters:

the mask of pixels that will be taken into account for computing the statistics. All zero pixels in mask will be ignored

varray, shape (R, C)

the image which the statistics will be computed from

num_labelsint

the number of different labels in ‘labels’ (equal to the number of hidden variables in the EM metric)

labelsarray, shape (R, C)

the label assigned to each pixel

Returns:
meansarray, shape (num_labels,)

means[i], 0<=i<num_labels will be the mean intensity in v of all voxels labeled i, or 0 if no voxels are labeled i

variancesarray, shape (num_labels,)

variances[i], 0<=i<num_labels will be the standard deviation of the intensities in v of all voxels labeled i, or infinite if less than 2 voxels are labeled i.

Computes the mean and std. for each quantization level.

Computes the mean and standard deviation of the intensities in ‘v’ for each corresponding label in ‘labels’. In other words, for each label L, it computes the mean and standard deviation of the intensities in ‘v’ at voxels whose label in ‘labels’ is L. This is used by the EM metric to compute statistics for each hidden variable represented by the labels.

Parameters:

the mask of voxels that will be taken into account for computing the statistics. All zero voxels in mask will be ignored

varray, shape (S, R, C)

the volume which the statistics will be computed from

num_labelsint

the number of different labels in ‘labels’ (equal to the number of hidden variables in the EM metric)

labelsarray, shape (S, R, C)

the label assigned to each pixel

Returns:
meansarray, shape (num_labels,)

means[i], 0<=i<num_labels will be the mean intensity in v of all voxels labeled i, or 0 if no voxels are labeled i

variancesarray, shape (num_labels,)

variances[i], 0<=i<num_labels will be the standard deviation of the intensities in v of all voxels labeled i, or infinite if less than 2 voxels are labeled i.

### quantize_positive_2d

dipy.align.expectmax.quantize_positive_2d(v, num_levels)

Quantizes a 2D image to num_levels quantization levels

Quantizes the input image at num_levels intensity levels considering <=0 as a special value. Those input pixels <=0, and only those, will be assigned a quantization level of 0. The positive values are divided into the remaining num_levels-1 uniform quantization levels.

The following are undefined, and raise a ValueError: * Quantizing at zero levels because at least one level must be assigned * Quantizing at one level because positive values should be assigned a

level different from the secial level 0 (at least 2 levels are needed)

Parameters:
varray, shape (R, C)

the image to be quantized

num_levelsint

the number of levels

Returns:
outarray, shape (R, C), same shape as v

the quantized image

levels: array, shape (num_levels,)

the quantization values: levels[0]=0, and levels[i] is the mid-point of the interval of intensities that are assigned to quantization level i, i=1, …, num_levels-1.

hist: array, shape (num_levels,)

histogram: the number of pixels that were assigned to each quantization level

### quantize_positive_3d

dipy.align.expectmax.quantize_positive_3d(v, num_levels)

Quantizes a 3D volume to num_levels quantization levels

Quantizes the input volume at num_levels intensity levels considering <=0 as a special value. Those input voxels <=0, and only those, will be assigned a quantization level of 0. The positive values are divided into the remaining num_levels-1 uniform quantization levels.

The following are undefined, and raise a ValueError: * Quantizing at zero levels because at least one level must be assigned * Quantizing at one level because positive values should be assigned a

level different from the secial level 0 (at least 2 levels are needed)

Parameters:
varray, shape (S, R, C)

the volume to be quantized

num_levelsint

the number of levels

Returns:
outarray, shape (S, R, C), same shape as v

the quantized volume

levels: array, shape (num_levels,)

the quantization values: levels[0]=0, and levels[i] is the mid-point of the interval of intensities that are assigned to quantization level i, i=1, …, num_levels-1.

hist: array, shape (num_levels,)

histogram: the number of voxels that were assigned to each quantization level

### AffineInvalidValuesError

class dipy.align.imaffine.AffineInvalidValuesError

Bases: Exception

Attributes:
args

Methods

 with_traceback Exception.with_traceback(tb) -- set self.__traceback__ to tb and return self.
__init__(*args, **kwargs)

### AffineInversionError

class dipy.align.imaffine.AffineInversionError

Bases: Exception

Attributes:
args

Methods

 with_traceback Exception.with_traceback(tb) -- set self.__traceback__ to tb and return self.
__init__(*args, **kwargs)

### AffineMap

class dipy.align.imaffine.AffineMap(affine, domain_grid_shape=None, domain_grid2world=None, codomain_grid_shape=None, codomain_grid2world=None)

Bases: object

Methods

 Return the value of the transformation, not a reference. set_affine(affine) Set the affine transform (operating in physical space). transform(image[, interpolation, ...]) Transform the input image from co-domain to domain space. transform_inverse(image[, interpolation, ...]) Transform the input image from domain to co-domain space.
__init__(affine, domain_grid_shape=None, domain_grid2world=None, codomain_grid_shape=None, codomain_grid2world=None)

AffineMap.

Implements an affine transformation whose domain is given by domain_grid and domain_grid2world, and whose co-domain is given by codomain_grid and codomain_grid2world.

The actual transform is represented by the affine matrix, which operate in world coordinates. Therefore, to transform a moving image towards a static image, we first map each voxel (i,j,k) of the static image to world coordinates (x,y,z) by applying domain_grid2world. Then we apply the affine transform to (x,y,z) obtaining (x’, y’, z’) in moving image’s world coordinates. Finally, (x’, y’, z’) is mapped to voxel coordinates (i’, j’, k’) in the moving image by multiplying (x’, y’, z’) by the inverse of codomain_grid2world. The codomain_grid_shape is used analogously to transform the static image towards the moving image when calling transform_inverse.

If the domain/co-domain information is not provided (None) then the sampling information needs to be specified each time the transform or transform_inverse is called to transform images. Note that such sampling information is not necessary to transform points defined in physical space, such as stream lines.

Parameters:
affinearray, shape (dim + 1, dim + 1)

the matrix defining the affine transform, where dim is the dimension of the space this map operates in (2 for 2D images, 3 for 3D images). If None, then self represents the identity transformation.

domain_grid_shapesequence, shape (dim,), optional

the shape of the default domain sampling grid. When transform is called to transform an image, the resulting image will have this shape, unless a different sampling information is provided. If None, then the sampling grid shape must be specified each time the transform method is called.

domain_grid2worldarray, shape (dim + 1, dim + 1), optional

the grid-to-world transform associated with the domain grid. If None (the default), then the grid-to-world transform is assumed to be the identity.

codomain_grid_shapesequence of integers, shape (dim,)

the shape of the default co-domain sampling grid. When transform_inverse is called to transform an image, the resulting image will have this shape, unless a different sampling information is provided. If None (the default), then the sampling grid shape must be specified each time the transform_inverse method is called.

codomain_grid2worldarray, shape (dim + 1, dim + 1)

the grid-to-world transform associated with the co-domain grid. If None (the default), then the grid-to-world transform is assumed to be the identity.

get_affine()

Return the value of the transformation, not a reference.

Returns:
affinendarray

Copy of the transform, not a reference.

set_affine(affine)

Set the affine transform (operating in physical space).

Also sets self.affine_inv - the inverse of affine, or None if there is no inverse.

Parameters:
affinearray, shape (dim + 1, dim + 1)

the matrix representing the affine transform operating in physical space. The domain and co-domain information remains unchanged. If None, then self represents the identity transformation.

transform(image, interpolation='linear', image_grid2world=None, sampling_grid_shape=None, sampling_grid2world=None, resample_only=False)

Transform the input image from co-domain to domain space.

By default, the transformed image is sampled at a grid defined by self.domain_shape and self.domain_grid2world. If such information was not provided then sampling_grid_shape is mandatory.

Parameters:
image2D or 3D array

the image to be transformed

interpolationstring, either ‘linear’ or ‘nearest’

the type of interpolation to be used, either ‘linear’ (for k-linear interpolation) or ‘nearest’ for nearest neighbor

image_grid2worldarray, shape (dim + 1, dim + 1), optional

the grid-to-world transform associated with image. If None (the default), then the grid-to-world transform is assumed to be the identity.

sampling_grid_shapesequence, shape (dim,), optional

the shape of the grid where the transformed image must be sampled. If None (the default), then self.codomain_shape is used instead (which must have been set at initialization, otherwise an exception will be raised).

sampling_grid2worldarray, shape (dim + 1, dim + 1), optional

the grid-to-world transform associated with the sampling grid (specified by sampling_grid_shape, or by default self.codomain_shape). If None (the default), then the grid-to-world transform is assumed to be the identity.

resample_onlyBoolean, optional

If False (the default) the affine transform is applied normally. If True, then the affine transform is not applied, and the input image is just re-sampled on the domain grid of this transform.

Returns:
transformedarray, shape sampling_grid_shape or

self.codomain_shape

the transformed image, sampled at the requested grid

transform_inverse(image, interpolation='linear', image_grid2world=None, sampling_grid_shape=None, sampling_grid2world=None, resample_only=False)

Transform the input image from domain to co-domain space.

By default, the transformed image is sampled at a grid defined by self.codomain_shape and self.codomain_grid2world. If such information was not provided then sampling_grid_shape is mandatory.

Parameters:
image2D or 3D array

the image to be transformed

interpolationstring, either ‘linear’ or ‘nearest’

the type of interpolation to be used, either ‘linear’ (for k-linear interpolation) or ‘nearest’ for nearest neighbor

image_grid2worldarray, shape (dim + 1, dim + 1), optional

the grid-to-world transform associated with image. If None (the default), then the grid-to-world transform is assumed to be the identity.

sampling_grid_shapesequence, shape (dim,), optional

the shape of the grid where the transformed image must be sampled. If None (the default), then self.codomain_shape is used instead (which must have been set at initialization, otherwise an exception will be raised).

sampling_grid2worldarray, shape (dim + 1, dim + 1), optional

the grid-to-world transform associated with the sampling grid (specified by sampling_grid_shape, or by default self.codomain_shape). If None (the default), then the grid-to-world transform is assumed to be the identity.

resample_onlyBoolean, optional

If False (the default) the affine transform is applied normally. If True, then the affine transform is not applied, and the input image is just re-sampled on the domain grid of this transform.

Returns:
transformedarray, shape sampling_grid_shape or

self.codomain_shape

the transformed image, sampled at the requested grid

### AffineRegistration

class dipy.align.imaffine.AffineRegistration(metric=None, level_iters=None, sigmas=None, factors=None, method='L-BFGS-B', ss_sigma_factor=None, options=None, verbosity=1)

Bases: object

Methods

 optimize(static, moving, transform, params0) Start the optimization process.
__init__(metric=None, level_iters=None, sigmas=None, factors=None, method='L-BFGS-B', ss_sigma_factor=None, options=None, verbosity=1)

Initialize an instance of the AffineRegistration class.

Parameters:
metricNone or object, optional

an instance of a metric. The default is None, implying the Mutual Information metric with default settings.

the number of iterations at each scale of the scale space. level_iters[0] corresponds to the coarsest scale, level_iters[-1] the finest, where n is the length of the sequence. By default, a 3-level scale space with iterations sequence equal to [10000, 1000, 100] will be used.

sigmassequence of floats, optional

custom smoothing parameter to build the scale space (one parameter for each scale). By default, the sequence of sigmas will be [3, 1, 0].

custom scale factors to build the scale space (one factor for each scale). By default, the sequence of factors will be [4, 2, 1].

methodstring, optional

optimization method to be used. If Scipy version < 0.12, then only L-BFGS-B is available. Otherwise, method can be any gradient-based method available in dipy.core.Optimize: CG, BFGS, Newton-CG, dogleg or trust-ncg. The default is ‘L-BFGS-B’.

ss_sigma_factorfloat, optional

If None, this parameter is not used and an isotropic scale space with the given factors and sigmas will be built. If not None, an anisotropic scale space will be used by automatically selecting the smoothing sigmas along each axis according to the voxel dimensions of the given image. The ss_sigma_factor is used to scale the automatically computed sigmas. For example, in the isotropic case, the sigma of the kernel will be $$factor * (2 ^ i)$$ where $$i = 1, 2, ..., n_scales - 1$$ is the scale (the finest resolution image $$i=0$$ is never smoothed). The default is None.

optionsdict, optional

extra optimization options. The default is None, implying no extra options are passed to the optimizer.

verbosity: int (one of {0, 1, 2, 3}), optional

Set the verbosity level of the algorithm: 0 : do not print anything 1 : print information about the current status of the algorithm 2 : print high level information of the components involved in

the registration that can be used to detect a failing component.

3print as much information as possible to isolate the cause

of a bug.

Default: 1

docstring_addendum = 'verbosity: int (one of {0, 1, 2, 3}), optional\n            Set the verbosity level of the algorithm:\n            0 : do not print anything\n            1 : print information about the current status of the algorithm\n            2 : print high level information of the components involved in\n                the registration that can be used to detect a failing\n                component.\n            3 : print as much information as possible to isolate the cause\n                of a bug.\n            Default: 1\n    '

Start the optimization process.

Parameters:
static2D or 3D array

the image to be used as reference during optimization.

moving2D or 3D array

the image to be used as “moving” during optimization. It is necessary to pre-align the moving image to ensure its domain lies inside the domain of the deformation fields. This is assumed to be accomplished by “pre-aligning” the moving image towards the static using an affine transformation given by the ‘starting_affine’ matrix

transforminstance of Transform

the transformation with respect to whose parameters the gradient must be computed

params0array, shape (n,)

parameters from which to start the optimization. If None, the optimization will start at the identity transform. n is the number of parameters of the specified transformation.

static_grid2worldarray, shape (dim+1, dim+1), optional

the voxel-to-space transformation associated with the static image. The default is None, implying the transform is the identity.

moving_grid2worldarray, shape (dim+1, dim+1), optional

the voxel-to-space transformation associated with the moving image. The default is None, implying the transform is the identity.

starting_affinestring, or matrix, or None, optional
If string:

‘mass’: align centers of gravity ‘voxel-origin’: align physical coordinates of voxel (0,0,0) ‘centers’: align physical coordinates of central voxels

If matrix:

array, shape (dim+1, dim+1).

If None:

Start from identity.

The default is None.

ret_metricboolean, optional

if True, it returns the parameters for measuring the similarity between the images (default ‘False’). The metric containing optimal parameters and the distance between the images.

static_maskarray, shape (S, R, C) or (R, C), optional

static image mask that defines which pixels in the static image are used to calculate the mutual information.

moving_maskarray, shape (S’, R’, C’) or (R’, C’), optional

moving image mask that defines which pixels in the moving image are used to calculate the mutual information.

Returns:
affine_mapinstance of AffineMap

the affine resulting affine transformation

xoptoptimal parameters

the optimal parameters (translation, rotation shear etc.)

foptSimilarity metric

the value of the function at the optimal parameters.

### IsotropicScaleSpace

class dipy.align.imaffine.IsotropicScaleSpace(image, factors, sigmas, image_grid2world=None, input_spacing=None, mask0=False)

Bases: ScaleSpace

Methods

 get_affine(level) Voxel-to-space transformation at a given level. get_affine_inv(level) Space-to-voxel transformation at a given level. get_domain_shape(level) Shape the sub-sampled image must have at a particular level. get_expand_factors(from_level, to_level) Ratio of voxel size from pyramid level from_level to to_level. get_image(level) Smoothed image at a given level. get_scaling(level) Adjustment factor for input-spacing to reflect voxel sizes at level. get_sigmas(level) Smoothing parameters used at a given level. get_spacing(level) Spacings the sub-sampled image must have at a particular level. print_level(level) Prints properties of a pyramid level.
__init__(image, factors, sigmas, image_grid2world=None, input_spacing=None, mask0=False)

IsotropicScaleSpace.

Computes the Scale Space representation of an image using isotropic smoothing kernels for all scales. The scale space is simply a list of images produced by smoothing the input image with a Gaussian kernel with different smoothing parameters.

This specialization of ScaleSpace allows the user to provide custom scale and smoothing factors for all scales.

Parameters:
imagearray, shape (r,c) or (s, r, c) where s is the number of

slices, r is the number of rows and c is the number of columns of the input image.

factorslist of floats

custom scale factors to build the scale space (one factor for each scale).

sigmaslist of floats

custom smoothing parameter to build the scale space (one parameter for each scale).

image_grid2worldarray, shape (dim + 1, dim + 1), optional

the grid-to-space transform of the image grid. The default is the identity matrix.

input_spacingarray, shape (dim,), optional

the spacing (voxel size) between voxels in physical space. The default if 1.0 along all axes.

if True, all smoothed images will be zero at all voxels that are zero in the input image. The default is False.

### MutualInformationMetric

class dipy.align.imaffine.MutualInformationMetric(nbins=32, sampling_proportion=None)

Bases: object

Methods

 distance(params) Numeric value of the negative Mutual Information. distance_and_gradient(params) Numeric value of the metric and its gradient at given parameters. gradient(params) Numeric value of the metric's gradient at the given parameters. setup(transform, static, moving[, ...]) Prepare the metric to compute intensity densities and gradients.
__init__(nbins=32, sampling_proportion=None)

Initialize an instance of the Mutual Information metric.

This class implements the methods required by Optimizer to drive the registration process.

Parameters:
nbinsint, optional

the number of bins to be used for computing the intensity histograms. The default is 32.

sampling_proportionNone or float in interval (0, 1], optional

There are two types of sampling: dense and sparse. Dense sampling uses all voxels for estimating the (joint and marginal) intensity histograms, while sparse sampling uses a subset of them. If sampling_proportion is None, then dense sampling is used. If sampling_proportion is a floating point value in (0,1] then sparse sampling is used, where sampling_proportion specifies the proportion of voxels to be used. The default is None.

Notes

Since we use linear interpolation, images are not, in general, differentiable at exact voxel coordinates, but they are differentiable between voxel coordinates. When using sparse sampling, selected voxels are slightly moved by adding a small random displacement within one voxel to prevent sampling points from being located exactly at voxel coordinates. When using dense sampling, this random displacement is not applied.

distance(params)

Numeric value of the negative Mutual Information.

We need to change the sign so we can use standard minimization algorithms.

Parameters:
paramsarray, shape (n,)

the parameter vector of the transform currently used by the metric (the transform name is provided when self.setup is called), n is the number of parameters of the transform

Returns:
neg_mifloat

the negative mutual information of the input images after transforming the moving image by the currently set transform with params parameters

Numeric value of the metric and its gradient at given parameters.

Parameters:
paramsarray, shape (n,)

the parameter vector of the transform currently used by the metric (the transform name is provided when self.setup is called), n is the number of parameters of the transform

Returns:
neg_mifloat

the negative mutual information of the input images after transforming the moving image by the currently set transform with params parameters

the gradient of the negative Mutual Information

Numeric value of the metric’s gradient at the given parameters.

Parameters:
paramsarray, shape (n,)

the parameter vector of the transform currently used by the metric (the transform name is provided when self.setup is called), n is the number of parameters of the transform

Returns:

the gradient of the negative Mutual Information

Prepare the metric to compute intensity densities and gradients.

The histograms will be setup to compute probability densities of intensities within the minimum and maximum values of static and moving

Parameters:
transform: instance of Transform

the transformation with respect to whose parameters the gradient must be computed

staticarray, shape (S, R, C) or (R, C)

static image

movingarray, shape (S’, R’, C’) or (R’, C’)

moving image. The dimensions of the static (S, R, C) and moving (S’, R’, C’) images do not need to be the same.

static_grid2worldarray (dim+1, dim+1), optional

the grid-to-space transform of the static image. The default is None, implying the transform is the identity.

moving_grid2worldarray (dim+1, dim+1)

the grid-to-space transform of the moving image. The default is None, implying the spacing along all axes is 1.

starting_affinearray, shape (dim+1, dim+1), optional

the pre-aligning matrix (an affine transform) that roughly aligns the moving image towards the static image. If None, no pre-alignment is performed. If a pre-alignment matrix is available, it is recommended to provide this matrix as starting_affine instead of manually transforming the moving image to reduce interpolation artifacts. The default is None, implying no pre-alignment is performed.

static_maskarray, shape (S, R, C) or (R, C), optional

static image mask that defines which pixels in the static image are used to calculate the mutual information.

moving_maskarray, shape (S’, R’, C’) or (R’, C’), optional

moving image mask that defines which pixels in the moving image are used to calculate the mutual information.

### Optimizer

class dipy.align.imaffine.Optimizer(fun, x0, args=(), method='L-BFGS-B', jac=None, hess=None, hessp=None, bounds=None, constraints=(), tol=None, callback=None, options=None, evolution=False)

Bases: object

Attributes:
evolution
fopt
message
nfev
nit
xopt

Methods

 print_summary
__init__(fun, x0, args=(), method='L-BFGS-B', jac=None, hess=None, hessp=None, bounds=None, constraints=(), tol=None, callback=None, options=None, evolution=False)

A class for handling minimization of scalar function of one or more variables.

Parameters:
funcallable

Objective function.

x0ndarray

Initial guess.

argstuple, optional

Extra arguments passed to the objective function and its derivatives (Jacobian, Hessian).

methodstr, optional

Type of solver. Should be one of

• ‘Powell’

• ‘CG’

• ‘BFGS’

• ‘Newton-CG’

• ‘Anneal’

• ‘L-BFGS-B’

• ‘TNC’

• ‘COBYLA’

• ‘SLSQP’

• ‘dogleg’

• ‘trust-ncg’

jacbool or callable, optional

Jacobian of objective function. Only for CG, BFGS, Newton-CG, dogleg, trust-ncg. If jac is a Boolean and is True, fun is assumed to return the value of Jacobian along with the objective function. If False, the Jacobian will be estimated numerically. jac can also be a callable returning the Jacobian of the objective. In this case, it must accept the same arguments as fun.

hess, hesspcallable, optional

Hessian of objective function or Hessian of objective function times an arbitrary vector p. Only for Newton-CG, dogleg, trust-ncg. Only one of hessp or hess needs to be given. If hess is provided, then hessp will be ignored. If neither hess nor hessp is provided, then the hessian product will be approximated using finite differences on jac. hessp must compute the Hessian times an arbitrary vector.

boundssequence, optional

Bounds for variables (only for L-BFGS-B, TNC and SLSQP). (min, max) pairs for each element in x, defining the bounds on that parameter. Use None for one of min or max when there is no bound in that direction.

constraintsdict or sequence of dict, optional

Constraints definition (only for COBYLA and SLSQP). Each constraint is defined in a dictionary with fields:

typestr

Constraint type: ‘eq’ for equality, ‘ineq’ for inequality.

funcallable

The function defining the constraint.

jaccallable, optional

The Jacobian of fun (only for SLSQP).

argssequence, optional

Extra arguments to be passed to the function and Jacobian.

Equality constraint means that the constraint function result is to be zero whereas inequality means that it is to be non-negative. Note that COBYLA only supports inequality constraints.

tolfloat, optional

Tolerance for termination. For detailed control, use solver-specific options.

callbackcallable, optional

Called after each iteration, as callback(xk), where xk is the current parameter vector. Only available using Scipy >= 0.12.

optionsdict, optional

A dictionary of solver options. All methods accept the following generic options:

maxiterint

Maximum number of iterations to perform.

dispbool

Set to True to print convergence messages.

For method-specific options, see show_options(‘minimize’, method).

evolutionbool, optional

save history of x for each iteration. Only available using Scipy >= 0.12.

scipy.optimize.minimize
property evolution
property fopt
property message
property nfev
property nit
print_summary()
property xopt

### ParzenJointHistogram

class dipy.align.imaffine.ParzenJointHistogram(nbins)

Bases: object

Methods

 bin_index(xnorm) Bin index associated with the given normalized intensity Maps intensity x to the range covered by the moving histogram Maps intensity x to the range covered by the static histogram setup(static, moving[, smask, mmask]) Compute histogram settings to store the PDF of input images update_gradient_dense(theta, transform, ...) Computes the Gradient of the joint PDF w.r.t. update_gradient_sparse(theta, transform, ...) Computes the Gradient of the joint PDF w.r.t. update_pdfs_dense(static, moving[, smask, mmask]) Computes the Probability Density Functions of two images update_pdfs_sparse(sval, mval) Computes the Probability Density Functions from a set of samples
__init__(nbins)

Computes joint histogram and derivatives with Parzen windows

Base class to compute joint and marginal probability density functions and their derivatives with respect to a transform’s parameters. The smooth histograms are computed by using Parzen windows [Parzen62] with a cubic spline kernel, as proposed by Mattes et al. [Mattes03]. This implementation is not tied to any optimization (registration) method, the idea is that information-theoretic matching functionals (such as Mutual Information) can inherit from this class to perform the low-level computations of the joint intensity distributions and its gradient w.r.t. the transform parameters. The derived class can then compute the similarity/dissimilarity measure and gradient, and finally communicate the results to the appropriate optimizer.

Parameters:
nbinsint

the number of bins of the joint and marginal probability density functions (the actual number of bins of the joint PDF is nbins**2)

Notes

We need this class in cython to allow _joint_pdf_gradient_dense_2d and _joint_pdf_gradient_dense_3d to use a nogil Jacobian function (obtained from an instance of the Transform class), which allows us to evaluate Jacobians at all the sampling points (maybe the full grid) inside a nogil loop.

The reason we need a class is to encapsulate all the parameters related to the joint and marginal distributions.

References

[Parzen62] E. Parzen. On the estimation of a probability density

function and the mode. Annals of Mathematical Statistics, 33(3), 1065-1076, 1962.

[Mattes03] Mattes, D., Haynor, D. R., Vesselle, H., Lewellen, T. K.,

& Eubank, W. PET-CT image registration in the chest using free-form deformations. IEEE Transactions on Medical Imaging, 22(1), 120-8, 2003.

bin_index(xnorm)

Bin index associated with the given normalized intensity

The return value is an integer in [padding, nbins - 1 - padding]

Parameters:
xnormfloat

intensity value normalized to the range covered by the histogram

Returns:
binint

the bin index associated with the given normalized intensity

bin_normalize_moving(x)

Maps intensity x to the range covered by the moving histogram

If the input intensity is in [self.mmin, self.mmax] then the normalized intensity will be in [self.padding, self.nbins - self.padding]

Parameters:
xfloat

the intensity to be normalized

Returns:
xnormfloat

normalized intensity to the range covered by the moving histogram

bin_normalize_static(x)

Maps intensity x to the range covered by the static histogram

If the input intensity is in [self.smin, self.smax] then the normalized intensity will be in [self.padding, self.nbins - self.padding]

Parameters:
xfloat

the intensity to be normalized

Returns:
xnormfloat

normalized intensity to the range covered by the static histogram

Compute histogram settings to store the PDF of input images

Parameters:
staticarray

static image

movingarray

moving image

mask of static object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). If None, the behaviour is equivalent to smask=ones_like(static)

mask of moving object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). If None, the behaviour is equivalent to mmask=ones_like(static)

Computes the Gradient of the joint PDF w.r.t. transform parameters

Computes the vector of partial derivatives of the joint histogram w.r.t. each transformation parameter.

Parameters:
thetaarray, shape (n,)

parameters of the transformation to compute the gradient from

transforminstance of Transform

the transformation with respect to whose parameters the gradient must be computed

staticarray, shape (S, R, C)

static image

movingarray, shape (S, R, C)

moving image

grid2worldarray, shape (4, 4)

we assume that both images have already been sampled at a common grid. This transform must map voxel coordinates of this common grid to physical coordinates of its corresponding voxel in the moving image. For example, if the moving image was sampled on the static image’s grid (this is the typical setting) using an aligning matrix A, then

1. grid2world = A.dot(static_affine)

where static_affine is the transformation mapping static image’s grid coordinates to physical space.

mgradientarray, shape (S, R, C, 3)

the gradient of the moving image

smaskarray, shape (S, R, C), optional

mask of static object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). The default is None, indicating all voxels are considered.

mmaskarray, shape (S, R, C), optional

mask of moving object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). The default is None, indicating all voxels are considered.

Computes the Gradient of the joint PDF w.r.t. transform parameters

Computes the vector of partial derivatives of the joint histogram w.r.t. each transformation parameter.

The list of intensities sval and mval are assumed to be sampled from the static and moving images, respectively, at the same physical points. Of course, the images may not be perfectly aligned at the moment the sampling was performed. The resulting gradient corresponds to the paired intensities according to the alignment at the moment the images were sampled.

Parameters:
thetaarray, shape (n,)

parameters to compute the gradient at

transforminstance of Transform

the transformation with respect to whose parameters the gradient must be computed

svalarray, shape (m,)

sampled intensities from the static image at sampled_points

mvalarray, shape (m,)

sampled intensities from the moving image at sampled_points

sample_pointsarray, shape (m, 3)

coordinates (in physical space) of the points the images were sampled at

the gradient of the moving image at the sample points

Computes the Probability Density Functions of two images

The joint PDF is stored in self.joint. The marginal distributions corresponding to the static and moving images are computed and stored in self.smarginal and self.mmarginal, respectively.

Parameters:
staticarray, shape (S, R, C)

static image

movingarray, shape (S, R, C)

moving image

mask of static object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). If None, ones_like(static) is used as mask.

mask of moving object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). If None, ones_like(moving) is used as mask.

update_pdfs_sparse(sval, mval)

Computes the Probability Density Functions from a set of samples

The list of intensities sval and mval are assumed to be sampled from the static and moving images, respectively, at the same physical points. Of course, the images may not be perfectly aligned at the moment the sampling was performed. The resulting distributions corresponds to the paired intensities according to the alignment at the moment the images were sampled.

The joint PDF is stored in self.joint. The marginal distributions corresponding to the static and moving images are computed and stored in self.smarginal and self.mmarginal, respectively.

Parameters:
svalarray, shape (n,)

sampled intensities from the static image at sampled_points

mvalarray, shape (n,)

sampled intensities from the moving image at sampled_points

### ScaleSpace

class dipy.align.imaffine.ScaleSpace(image, num_levels, image_grid2world=None, input_spacing=None, sigma_factor=0.2, mask0=False)

Bases: object

Methods

 get_affine(level) Voxel-to-space transformation at a given level. get_affine_inv(level) Space-to-voxel transformation at a given level. get_domain_shape(level) Shape the sub-sampled image must have at a particular level. get_expand_factors(from_level, to_level) Ratio of voxel size from pyramid level from_level to to_level. get_image(level) Smoothed image at a given level. get_scaling(level) Adjustment factor for input-spacing to reflect voxel sizes at level. get_sigmas(level) Smoothing parameters used at a given level. get_spacing(level) Spacings the sub-sampled image must have at a particular level. print_level(level) Prints properties of a pyramid level.
__init__(image, num_levels, image_grid2world=None, input_spacing=None, sigma_factor=0.2, mask0=False)

ScaleSpace.

Computes the Scale Space representation of an image. The scale space is simply a list of images produced by smoothing the input image with a Gaussian kernel with increasing smoothing parameter. If the image’s voxels are isotropic, the smoothing will be the same along all directions: at level L = 0, 1, …, the sigma is given by $$s * ( 2^L - 1 )$$. If the voxel dimensions are not isotropic, then the smoothing is weaker along low resolution directions.

Parameters:
imagearray, shape (r,c) or (s, r, c) where s is the number of

slices, r is the number of rows and c is the number of columns of the input image.

num_levelsint

the desired number of levels (resolutions) of the scale space

image_grid2worldarray, shape (dim + 1, dim + 1), optional

the grid-to-space transform of the image grid. The default is the identity matrix

input_spacingarray, shape (dim,), optional

the spacing (voxel size) between voxels in physical space. The default is 1.0 along all axes

sigma_factorfloat, optional

the smoothing factor to be used in the construction of the scale space. The default is 0.2

if True, all smoothed images will be zero at all voxels that are zero in the input image. The default is False.

get_affine(level)

Voxel-to-space transformation at a given level.

Returns the voxel-to-space transformation associated with the sub-sampled image at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the sub-sampled images must have).

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get affine transform from

Returns:
the affine (voxel-to-space) transform at the requested resolution
or None if an invalid level was requested
get_affine_inv(level)

Space-to-voxel transformation at a given level.

Returns the space-to-voxel transformation associated with the sub-sampled image at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the sub-sampled images must have).

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the inverse transform from

Returns:
the inverse (space-to-voxel) transform at the requested resolution or
None if an invalid level was requested
get_domain_shape(level)

Shape the sub-sampled image must have at a particular level.

Returns the shape the sub-sampled image must have at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the sub-sampled images must have).

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the sub-sampled shape from

Returns:
the sub-sampled shape at the requested resolution or None if an
invalid level was requested
get_expand_factors(from_level, to_level)

Ratio of voxel size from pyramid level from_level to to_level.

Given two scale space resolutions a = from_level, b = to_level, returns the ratio of voxels size at level b to voxel size at level a (the factor that must be used to multiply voxels at level a to ‘expand’ them to level b).

Parameters:
from_levelint, 0 <= from_level < L, (L = number of resolutions)

the resolution to expand voxels from

to_levelint, 0 <= to_level < from_level

the resolution to expand voxels to

Returns:
factorsarray, shape (k,), k = 2, 3

the expand factors (a scalar for each voxel dimension)

get_image(level)

Smoothed image at a given level.

Returns the smoothed image at the requested level in the Scale Space.

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the smooth image from

Returns:
the smooth image at the requested resolution or None if an invalid
level was requested
get_scaling(level)

Adjustment factor for input-spacing to reflect voxel sizes at level.

Returns the scaling factor that needs to be applied to the input spacing (the voxel sizes of the image at level 0 of the scale space) to transform them to voxel sizes at the requested level.

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the scalings from

Returns:
the scaling factors from the original spacing to the spacings at the
requested level
get_sigmas(level)

Smoothing parameters used at a given level.

Returns the smoothing parameters (a scalar for each axis) used at the requested level of the scale space

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the smoothing parameters from

Returns:
the smoothing parameters at the requested level
get_spacing(level)

Spacings the sub-sampled image must have at a particular level.

Returns the spacings (voxel sizes) the sub-sampled image must have at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the sub-sampled images must have).

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the sub-sampled shape from

Returns:
the spacings (voxel sizes) at the requested resolution or None if an
invalid level was requested
print_level(level)

Prints properties of a pyramid level.

Prints the properties of a level of this scale space to standard output

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to be printed

### compute_parzen_mi

dipy.align.imaffine.compute_parzen_mi()

Computes the mutual information and its gradient (if requested)

Parameters:
jointarray, shape (nbins, nbins)

the joint intensity distribution

the gradient of the joint distribution w.r.t. the transformation parameters

smarginalarray, shape (nbins,)

the marginal intensity distribution of the static image

mmarginalarray, shape (nbins,)

the marginal intensity distribution of the moving image

the buffer in which to write the gradient of the mutual information. If None, the gradient is not computed

### deprecated_params

dipy.align.imaffine.deprecated_params(old_name, new_name=None, since='', until='', version_comparator=<function cmp_pkg_version>, arg_in_kwargs=False, warn_class=<class 'dipy.utils.deprecator.ArgsDeprecationWarning'>, error_class=<class 'dipy.utils.deprecator.ExpiredDeprecationError'>, alternative='')

Deprecate a renamed or removed function argument.

The decorator assumes that the argument with the old_name was removed from the function signature and the new_name replaced it at the same position in the signature. If the old_name argument is given when calling the decorated function the decorator will catch it and issue a deprecation warning and pass it on as new_name argument.

Parameters:
old_namestr or list/tuple thereof

The old name of the argument.

new_namestr or list/tuple thereof or None, optional

The new name of the argument. Set this to None to remove the argument old_name instead of renaming it.

sincestr or number or list/tuple thereof, optional

The release at which the old argument became deprecated.

untilstr or number or list/tuple thereof, optional

Last released version at which this function will still raise a deprecation warning. Versions higher than this will raise an error.

version_comparatorcallable

Callable accepting string as argument, and return 1 if string represents a higher version than encoded in the version_comparator, 0 if the version is equal, and -1 if the version is lower. For example, the version_comparator may compare the input version string to the current package version string.

arg_in_kwargsbool or list/tuple thereof, optional

If the argument is not a named argument (for example it was meant to be consumed by **kwargs) set this to True. Otherwise the decorator will throw an Exception if the new_name cannot be found in the signature of the decorated function. Default is False.

warn_classwarning, optional

Warning to be issued.

error_classException, optional

Error to be issued

alternativestr, optional

An alternative function or class name that the user may use in place of the deprecated object if new_name is None. The deprecation warning will tell the user about this alternative if provided.

Raises:
TypeError

If the new argument name cannot be found in the function signature and arg_in_kwargs was False or if it is used to deprecate the name of the *args-, **kwargs-like arguments. At runtime such an Error is raised if both the new_name and old_name were specified when calling the function and “relax=False”.

Notes

This function is based on the Astropy (major modification). https://github.com/astropy/astropy. See COPYING file distributed along with the astropy package for the copyright and license terms.

Examples

The deprecation warnings are not shown in the following examples. To deprecate a positional or keyword argument:: >>> from dipy.utils.deprecator import deprecated_params >>> @deprecated_params(‘sig’, ‘sigma’, ‘0.3’) … def test(sigma): … return sigma >>> test(2) 2 >>> test(sigma=2) 2 >>> test(sig=2) # doctest: +SKIP 2

It is also possible to replace multiple arguments. The old_name, new_name and since have to be tuple or list and contain the same number of entries:: >>> @deprecated_params([‘a’, ‘b’], [‘alpha’, ‘beta’], … [‘0.2’, 0.4]) … def test(alpha, beta): … return alpha, beta >>> test(a=2, b=3) # doctest: +SKIP (2, 3)

### get_direction_and_spacings

dipy.align.imaffine.get_direction_and_spacings(affine, dim)

Extracts the rotational and spacing components from a matrix

Extracts the rotational and spacing (voxel dimensions) components from a matrix. An image gradient represents the local variation of the image’s gray values per voxel. Since we are iterating on the physical space, we need to compute the gradients as variation per millimeter, so we need to divide each gradient’s component by the voxel size along the corresponding axis, that’s what the spacings are used for. Since the image’s gradients are oriented along the grid axes, we also need to re-orient the gradients to be given in physical space coordinates.

Parameters:
affinearray, shape (k, k), k = 3, 4

the matrix transforming grid coordinates to physical space.

Returns:
directionarray, shape (k-1, k-1)

the rotational component of the input matrix

spacingsarray, shape (k-1,)

the scaling component (voxel size) of the matrix

### interpolate_scalar_2d

dipy.align.imaffine.interpolate_scalar_2d(image, locations)

Bilinear interpolation of a 2D scalar image

Interpolates the 2D image at the given locations. This function is a wrapper for _interpolate_scalar_2d for testing purposes, it is equivalent to scipy.ndimage.map_coordinates with bilinear interpolation

Parameters:
fieldarray, shape (S, R)

the 2D image to be interpolated

locationsarray, shape (n, 2)

(locations[i,0], locations[i,1]), 0<=i<n must contain the row and column coordinates to interpolate the image at

Returns:
outarray, shape (n,)

out[i], 0<=i<n will be the interpolated scalar at coordinates locations[i,:], or 0 if locations[i,:] is outside the image

insidearray, (n,)

if locations[i:] is inside the image then inside[i]=1, else inside[i]=0

### interpolate_scalar_3d

dipy.align.imaffine.interpolate_scalar_3d(image, locations)

Trilinear interpolation of a 3D scalar image

Interpolates the 3D image at the given locations. This function is a wrapper for _interpolate_scalar_3d for testing purposes, it is equivalent to scipy.ndimage.map_coordinates with trilinear interpolation

Parameters:
fieldarray, shape (S, R, C)

the 3D image to be interpolated

locationsarray, shape (n, 3)

(locations[i,0], locations[i,1], locations[i,2), 0<=i<n must contain the coordinates to interpolate the image at

Returns:
outarray, shape (n,)

out[i], 0<=i<n will be the interpolated scalar at coordinates locations[i,:], or 0 if locations[i,:] is outside the image

insidearray, (n,)

if locations[i,:] is inside the image then inside[i]=1, else inside[i]=0

### sample_domain_regular

dipy.align.imaffine.sample_domain_regular()

Take floor(total_voxels/k) samples from a (2D or 3D) grid

The sampling is made by taking all pixels whose index (in lexicographical order) is a multiple of k. Each selected point is slightly perturbed by adding a realization of a normally distributed random variable and then mapped to physical space by the given grid-to-space transform.

The lexicographical order of a pixels in a grid of shape (a, b, c) is defined by assigning to each voxel position (i, j, k) the integer index

F((i, j, k)) = i * (b * c) + j * (c) + k

and sorting increasingly by this index.

Parameters:
kint

the sampling rate, as described before

shapearray, shape (dim,)

the shape of the grid to be sampled

grid2worldarray, shape (dim+1, dim+1)

the grid-to-space transform

sigmafloat

the standard deviation of the Normal random distortion to be applied to the sampled points

Returns:
samplesarray, shape (total_pixels//k, dim)

the matrix whose rows are the sampled points

Examples

>>> from dipy.align.parzenhist import sample_domain_regular
>>> import dipy.align.vector_fields as vf
>>> shape = np.array((10, 10), dtype=np.int32)
>>> sigma = 0
>>> dim = len(shape)
>>> grid2world = np.eye(dim+1)
>>> n = shape[0]*shape[1]
>>> k = 2
>>> samples = sample_domain_regular(k, shape, grid2world, sigma)
>>> (samples.shape[0], samples.shape[1]) == (n//k, dim)
True
>>> isamples = np.array(samples, dtype=np.int32)
>>> indices = (isamples[:, 0] * shape[1] + isamples[:, 1])
>>> len(set(indices)) == len(indices)
True
>>> (indices%k).sum()
0


### transform_centers_of_mass

dipy.align.imaffine.transform_centers_of_mass(static, static_grid2world, moving, moving_grid2world)

Transformation to align the center of mass of the input images.

Parameters:
staticarray, shape (S, R, C)

static image

static_grid2worldarray, shape (dim+1, dim+1)

the voxel-to-space transformation of the static image

movingarray, shape (S, R, C)

moving image

moving_grid2worldarray, shape (dim+1, dim+1)

the voxel-to-space transformation of the moving image

Returns:
affine_mapinstance of AffineMap

the affine transformation (translation only, in this case) aligning the center of mass of the moving image towards the one of the static image

### transform_geometric_centers

dipy.align.imaffine.transform_geometric_centers(static, static_grid2world, moving, moving_grid2world)

Transformation to align the geometric center of the input images.

With “geometric center” of a volume we mean the physical coordinates of its central voxel

Parameters:
staticarray, shape (S, R, C)

static image

static_grid2worldarray, shape (dim+1, dim+1)

the voxel-to-space transformation of the static image

movingarray, shape (S, R, C)

moving image

moving_grid2worldarray, shape (dim+1, dim+1)

the voxel-to-space transformation of the moving image

Returns:
affine_mapinstance of AffineMap

the affine transformation (translation only, in this case) aligning the geometric center of the moving image towards the one of the static image

### transform_origins

dipy.align.imaffine.transform_origins(static, static_grid2world, moving, moving_grid2world)

Transformation to align the origins of the input images.

With “origin” of a volume we mean the physical coordinates of voxel (0,0,0)

Parameters:
staticarray, shape (S, R, C)

static image

static_grid2worldarray, shape (dim+1, dim+1)

the voxel-to-space transformation of the static image

movingarray, shape (S, R, C)

moving image

moving_grid2worldarray, shape (dim+1, dim+1)

the voxel-to-space transformation of the moving image

Returns:
affine_mapinstance of AffineMap

the affine transformation (translation only, in this case) aligning the origin of the moving image towards the one of the static image

### warn

dipy.align.imaffine.warn(/, message, category=None, stacklevel=1, source=None)

Issue a warning, or maybe ignore it or raise an exception.

### Bunch

class dipy.align.imwarp.Bunch(**kwds)

Bases: object

__init__(**kwds)

A ‘bunch’ of values (a replacement of Enum)

This is a temporary replacement of Enum, which is not available on all versions of Python 2

### DiffeomorphicMap

class dipy.align.imwarp.DiffeomorphicMap(dim, disp_shape, disp_grid2world=None, domain_shape=None, domain_grid2world=None, codomain_shape=None, codomain_grid2world=None, prealign=None)

Bases: object

Methods

 Creates a zero displacement field Inversion error of the displacement fields expand_fields(expand_factors, new_shape) Expands the displacement fields from current shape to new_shape Deformation field to transform an image in the backward direction Deformation field to transform an image in the forward direction Constructs a simplified version of this Diffeomorhic Map Try to interpret obj as a matrix Inverse of this DiffeomorphicMap instance Shallow copy of this DiffeomorphicMap instance transform(image[, interpolation, ...]) Warps an image in the forward direction transform_inverse(image[, interpolation, ...]) Warps an image in the backward direction Composition of this DiffeomorphicMap with a given endomorphism
 transform_points transform_points_inverse
__init__(dim, disp_shape, disp_grid2world=None, domain_shape=None, domain_grid2world=None, codomain_shape=None, codomain_grid2world=None, prealign=None)

DiffeomorphicMap

Implements a diffeomorphic transformation on the physical space. The deformation fields encoding the direct and inverse transformations share the same domain discretization (both the discretization grid shape and voxel-to-space matrix). The input coordinates (physical coordinates) are first aligned using prealign, and then displaced using the corresponding vector field interpolated at the aligned coordinates.

Parameters:
dimint, 2 or 3

the transformation’s dimension

disp_shapearray, shape (dim,)

the number of slices (if 3D), rows and columns of the deformation field’s discretization

disp_grid2worldthe voxel-to-space transform between the def. fields

grid and space

domain_shapearray, shape (dim,)

the number of slices (if 3D), rows and columns of the default discretization of this map’s domain

domain_grid2worldarray, shape (dim+1, dim+1)

the default voxel-to-space transformation between this map’s discretization and physical space

codomain_shapearray, shape (dim,)

the number of slices (if 3D), rows and columns of the images that are ‘normally’ warped using this transformation in the forward direction (this will provide default transformation parameters to warp images under this transformation). By default, we assume that the inverse transformation is ‘normally’ used to warp images with the same discretization and voxel-to-space transformation as the deformation field grid.

codomain_grid2worldarray, shape (dim+1, dim+1)

the voxel-to-space transformation of images that are ‘normally’ warped using this transformation (in the forward direction).

prealignarray, shape (dim+1, dim+1)

the linear transformation to be applied to align input images to the reference space before warping under the deformation field.

allocate()

Creates a zero displacement field

Creates a zero displacement field (the identity transformation).

compute_inversion_error()

Inversion error of the displacement fields

Estimates the inversion error of the displacement fields by computing statistics of the residual vectors obtained after composing the forward and backward displacement fields.

Returns:
residualarray, shape (R, C) or (S, R, C)

the displacement field resulting from composing the forward and backward displacement fields of this transformation (the residual should be zero for a perfect diffeomorphism)

statsarray, shape (3,)

statistics from the norms of the vectors of the residual displacement field: maximum, mean and standard deviation

Notes

Since the forward and backward displacement fields have the same discretization, the final composition is given by

comp[i] = forward[ i + Dinv * backward[i]]

where Dinv is the space-to-grid transformation of the displacement fields

expand_fields(expand_factors, new_shape)

Expands the displacement fields from current shape to new_shape

Up-samples the discretization of the displacement fields to be of new_shape shape.

Parameters:
expand_factorsarray, shape (dim,)

the factors scaling current spacings (voxel sizes) to spacings in the expanded discretization.

new_shapearray, shape (dim,)

the shape of the arrays holding the up-sampled discretization

get_backward_field()

Deformation field to transform an image in the backward direction

Returns the deformation field that must be used to warp an image under this transformation in the backward direction (note the ‘is_inverse’ flag).

get_forward_field()

Deformation field to transform an image in the forward direction

Returns the deformation field that must be used to warp an image under this transformation in the forward direction (note the ‘is_inverse’ flag).

get_simplified_transform()

Constructs a simplified version of this Diffeomorhic Map

The simplified version incorporates the pre-align transform, as well as the domain and codomain affine transforms into the displacement field. The resulting transformation may be regarded as operating on the image spaces given by the domain and codomain discretization. As a result, self.prealign, self.disp_grid2world, self.domain_grid2world and self.codomain affine will be None (denoting Identity) in the resulting diffeomorphic map.

interpret_matrix(obj)

Try to interpret obj as a matrix

Some operations are performed faster if we know in advance if a matrix is the identity (so we can skip the actual matrix-vector multiplication). This function returns None if the given object is None or the ‘identity’ string. It returns the same object if it is a numpy array. It raises an exception otherwise.

Parameters:
objobject

any object

Returns:
objobject

the same object given as argument if obj is None or a numpy array. None if obj is the ‘identity’ string.

inverse()

Inverse of this DiffeomorphicMap instance

Returns a diffeomorphic map object representing the inverse of this transformation. The internal arrays are not copied but just referenced.

Returns:
invDiffeomorphicMap object

the inverse of this diffeomorphic map.

shallow_copy()

Shallow copy of this DiffeomorphicMap instance

Creates a shallow copy of this diffeomorphic map (the arrays are not copied but just referenced)

Returns:
new_mapDiffeomorphicMap object

the shallow copy of this diffeomorphic map

transform(image, interpolation='linear', image_world2grid=None, out_shape=None, out_grid2world=None)

Warps an image in the forward direction

Transforms the input image under this transformation in the forward direction. It uses the “is_inverse” flag to switch between “forward” and “backward” (if is_inverse is False, then transform(…) warps the image forwards, else it warps the image backwards).

Parameters:
imagearray, shape (s, r, c) if dim = 3 or (r, c) if dim = 2

the image to be warped under this transformation in the forward direction

interpolationstring, either ‘linear’ or ‘nearest’

the type of interpolation to be used for warping, either ‘linear’ (for k-linear interpolation) or ‘nearest’ for nearest neighbor

image_world2gridarray, shape (dim+1, dim+1)

the transformation bringing world (space) coordinates to voxel coordinates of the image given as input

out_shapearray, shape (dim,)

the number of slices, rows and columns of the desired warped image

out_grid2worldthe transformation bringing voxel coordinates of the

warped image to physical space

Returns:
warpedarray, shape = out_shape or self.codomain_shape if None

the warped image under this transformation in the forward direction

Notes

See _warp_forward and _warp_backward documentation for further information.

transform_inverse(image, interpolation='linear', image_world2grid=None, out_shape=None, out_grid2world=None)

Warps an image in the backward direction

Transforms the input image under this transformation in the backward direction. It uses the “is_inverse” flag to switch between “forward” and “backward” (if is_inverse is False, then transform_inverse(…) warps the image backwards, else it warps the image forwards)

Parameters:
imagearray, shape (s, r, c) if dim = 3 or (r, c) if dim = 2

the image to be warped under this transformation in the forward direction

interpolationstring, either ‘linear’ or ‘nearest’

the type of interpolation to be used for warping, either ‘linear’ (for k-linear interpolation) or ‘nearest’ for nearest neighbor

image_world2gridarray, shape (dim+1, dim+1)

the transformation bringing world (space) coordinates to voxel coordinates of the image given as input

out_shapearray, shape (dim,)

the number of slices, rows, and columns of the desired warped image

out_grid2worldthe transformation bringing voxel coordinates of the

warped image to physical space

Returns:
warpedarray, shape = out_shape or self.codomain_shape if None

warped image under this transformation in the backward direction

Notes

See _warp_forward and _warp_backward documentation for further information.

transform_points(points, coord2world=None, world2coord=None)
transform_points_inverse(points, coord2world=None, world2coord=None)
warp_endomorphism(phi)

Composition of this DiffeomorphicMap with a given endomorphism

Creates a new DiffeomorphicMap C with the same properties as self and composes its displacement fields with phi’s corresponding fields. The resulting diffeomorphism is of the form C(x) = phi(self(x)) with inverse C^{-1}(y) = self^{-1}(phi^{-1}(y)). We assume that phi is an endomorphism with the same discretization and domain affine as self to ensure that the composition inherits self’s properties (we also assume that the pre-aligning matrix of phi is None or identity).

Parameters:
phiDiffeomorphicMap object

the endomorphism to be warped by this diffeomorphic map

Returns:
compositionthe composition of this diffeomorphic map with the

endomorphism given as input

Notes

The problem with our current representation of a DiffeomorphicMap is that the set of Diffeomorphism that can be represented this way (a pre-aligning matrix followed by a non-linear endomorphism given as a displacement field) is not closed under the composition operation.

Supporting a general DiffeomorphicMap class, closed under composition, may be extremely costly computationally, and the kind of transformations we actually need for Avants’ mid-point algorithm (SyN) are much simpler.

### DiffeomorphicRegistration

class dipy.align.imwarp.DiffeomorphicRegistration(metric=None)

Bases: object

Methods

 Returns the resulting diffeomorphic map after optimization Starts the metric optimization set_level_iters(level_iters) Sets the number of iterations at each pyramid level
__init__(metric=None)

Diffeomorphic Registration

This abstract class defines the interface to be implemented by any optimization algorithm for diffeomorphic registration.

Parameters:
metricSimilarityMetric object

the object measuring the similarity of the two images. The registration algorithm will minimize (or maximize) the provided similarity.

abstract get_map()

Returns the resulting diffeomorphic map after optimization

abstract optimize()

Starts the metric optimization

This is the main function each specialized class derived from this must implement. Upon completion, the deformation field must be available from the forward transformation model.

set_level_iters(level_iters)

Sets the number of iterations at each pyramid level

Establishes the maximum number of iterations to be performed at each level of the Gaussian pyramid, similar to ANTS.

Parameters:
level_iterslist

the number of iterations at each level of the Gaussian pyramid. level_iters[0] corresponds to the finest level, level_iters[n-1] the coarsest, where n is the length of the list

### ScaleSpace

class dipy.align.imwarp.ScaleSpace(image, num_levels, image_grid2world=None, input_spacing=None, sigma_factor=0.2, mask0=False)

Bases: object

Methods

 get_affine(level) Voxel-to-space transformation at a given level. get_affine_inv(level) Space-to-voxel transformation at a given level. get_domain_shape(level) Shape the sub-sampled image must have at a particular level. get_expand_factors(from_level, to_level) Ratio of voxel size from pyramid level from_level to to_level. get_image(level) Smoothed image at a given level. get_scaling(level) Adjustment factor for input-spacing to reflect voxel sizes at level. get_sigmas(level) Smoothing parameters used at a given level. get_spacing(level) Spacings the sub-sampled image must have at a particular level. print_level(level) Prints properties of a pyramid level.
__init__(image, num_levels, image_grid2world=None, input_spacing=None, sigma_factor=0.2, mask0=False)

ScaleSpace.

Computes the Scale Space representation of an image. The scale space is simply a list of images produced by smoothing the input image with a Gaussian kernel with increasing smoothing parameter. If the image’s voxels are isotropic, the smoothing will be the same along all directions: at level L = 0, 1, …, the sigma is given by $$s * ( 2^L - 1 )$$. If the voxel dimensions are not isotropic, then the smoothing is weaker along low resolution directions.

Parameters:
imagearray, shape (r,c) or (s, r, c) where s is the number of

slices, r is the number of rows and c is the number of columns of the input image.

num_levelsint

the desired number of levels (resolutions) of the scale space

image_grid2worldarray, shape (dim + 1, dim + 1), optional

the grid-to-space transform of the image grid. The default is the identity matrix

input_spacingarray, shape (dim,), optional

the spacing (voxel size) between voxels in physical space. The default is 1.0 along all axes

sigma_factorfloat, optional

the smoothing factor to be used in the construction of the scale space. The default is 0.2

if True, all smoothed images will be zero at all voxels that are zero in the input image. The default is False.

get_affine(level)

Voxel-to-space transformation at a given level.

Returns the voxel-to-space transformation associated with the sub-sampled image at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the sub-sampled images must have).

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get affine transform from

Returns:
the affine (voxel-to-space) transform at the requested resolution
or None if an invalid level was requested
get_affine_inv(level)

Space-to-voxel transformation at a given level.

Returns the space-to-voxel transformation associated with the sub-sampled image at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the sub-sampled images must have).

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the inverse transform from

Returns:
the inverse (space-to-voxel) transform at the requested resolution or
None if an invalid level was requested
get_domain_shape(level)

Shape the sub-sampled image must have at a particular level.

Returns the shape the sub-sampled image must have at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the sub-sampled images must have).

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the sub-sampled shape from

Returns:
the sub-sampled shape at the requested resolution or None if an
invalid level was requested
get_expand_factors(from_level, to_level)

Ratio of voxel size from pyramid level from_level to to_level.

Given two scale space resolutions a = from_level, b = to_level, returns the ratio of voxels size at level b to voxel size at level a (the factor that must be used to multiply voxels at level a to ‘expand’ them to level b).

Parameters:
from_levelint, 0 <= from_level < L, (L = number of resolutions)

the resolution to expand voxels from

to_levelint, 0 <= to_level < from_level

the resolution to expand voxels to

Returns:
factorsarray, shape (k,), k = 2, 3

the expand factors (a scalar for each voxel dimension)

get_image(level)

Smoothed image at a given level.

Returns the smoothed image at the requested level in the Scale Space.

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the smooth image from

Returns:
the smooth image at the requested resolution or None if an invalid
level was requested
get_scaling(level)

Adjustment factor for input-spacing to reflect voxel sizes at level.

Returns the scaling factor that needs to be applied to the input spacing (the voxel sizes of the image at level 0 of the scale space) to transform them to voxel sizes at the requested level.

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the scalings from

Returns:
the scaling factors from the original spacing to the spacings at the
requested level
get_sigmas(level)

Smoothing parameters used at a given level.

Returns the smoothing parameters (a scalar for each axis) used at the requested level of the scale space

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the smoothing parameters from

Returns:
the smoothing parameters at the requested level
get_spacing(level)

Spacings the sub-sampled image must have at a particular level.

Returns the spacings (voxel sizes) the sub-sampled image must have at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the sub-sampled images must have).

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the sub-sampled shape from

Returns:
the spacings (voxel sizes) at the requested resolution or None if an
invalid level was requested
print_level(level)

Prints properties of a pyramid level.

Prints the properties of a level of this scale space to standard output

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to be printed

### Streamlines

dipy.align.imwarp.Streamlines

alias of ArraySequence

### SymmetricDiffeomorphicRegistration

class dipy.align.imwarp.SymmetricDiffeomorphicRegistration(metric, level_iters=None, step_length=0.25, ss_sigma_factor=0.2, opt_tol=1e-05, inv_iter=20, inv_tol=0.001, callback=None)

Methods

 Return the resulting diffeomorphic map. optimize(static, moving[, ...]) Starts the optimization set_level_iters(level_iters) Sets the number of iterations at each pyramid level update(current_displacement, ...) Composition of the current displacement field with the given field
__init__(metric, level_iters=None, step_length=0.25, ss_sigma_factor=0.2, opt_tol=1e-05, inv_iter=20, inv_tol=0.001, callback=None)

Symmetric Diffeomorphic Registration (SyN) Algorithm

Performs the multi-resolution optimization algorithm for non-linear registration using a given similarity metric.

Parameters:
metricSimilarityMetric object

the metric to be optimized

level_iterslist of int

the number of iterations at each level of the Gaussian Pyramid (the length of the list defines the number of pyramid levels to be used)

opt_tolfloat

the optimization will stop when the estimated derivative of the energy profile w.r.t. time falls below this threshold

inv_iterint

the number of iterations to be performed by the displacement field inversion algorithm

step_lengthfloat

the length of the maximum displacement vector of the update displacement field at each iteration

ss_sigma_factorfloat

parameter of the scale-space smoothing kernel. For example, the std. dev. of the kernel will be factor*(2^i) in the isotropic case where i = 0, 1, …, n_scales is the scale

inv_tolfloat

the displacement field inversion algorithm will stop iterating when the inversion error falls below this threshold

callbackfunction(SymmetricDiffeomorphicRegistration)

a function receiving a SymmetricDiffeomorphicRegistration object to be called after each iteration (this optimizer will call this function passing self as parameter)

get_map()

Return the resulting diffeomorphic map.

Returns the DiffeomorphicMap registering the moving image towards the static image.

optimize(static, moving, static_grid2world=None, moving_grid2world=None, prealign=None)

Starts the optimization

Parameters:
staticarray, shape (S, R, C) or (R, C)

the image to be used as reference during optimization. The displacement fields will have the same discretization as the static image.

movingarray, shape (S, R, C) or (R, C)

the image to be used as “moving” during optimization. Since the deformation fields’ discretization is the same as the static image, it is necessary to pre-align the moving image to ensure its domain lies inside the domain of the deformation fields. This is assumed to be accomplished by “pre-aligning” the moving image towards the static using an affine transformation given by the ‘prealign’ matrix

static_grid2worldarray, shape (dim+1, dim+1)

the voxel-to-space transformation associated to the static image

moving_grid2worldarray, shape (dim+1, dim+1)

the voxel-to-space transformation associated to the moving image

prealignarray, shape (dim+1, dim+1)

the affine transformation (operating on the physical space) pre-aligning the moving image towards the static

Returns:
static_to_refDiffeomorphicMap object

the diffeomorphic map that brings the moving image towards the static one in the forward direction (i.e. by calling static_to_ref.transform) and the static image towards the moving one in the backward direction (i.e. by calling static_to_ref.transform_inverse).

update(current_displacement, new_displacement, disp_world2grid, time_scaling)

Composition of the current displacement field with the given field

Interpolates new displacement at the locations defined by current_displacement. Equivalently, computes the composition C of the given displacement fields as C(x) = B(A(x)), where A is current_displacement and B is new_displacement. This function is intended to be used with deformation fields of the same sampling (e.g. to be called by a registration algorithm).

Parameters:
current_displacementarray, shape (R’, C’, 2) or (S’, R’, C’, 3)

the displacement field defining where to interpolate new_displacement

new_displacementarray, shape (R, C, 2) or (S, R, C, 3)

the displacement field to be warped by current_displacement

disp_world2gridarray, shape (dim+1, dim+1)

the space-to-grid transform associated with the displacements’ grid (we assume that both displacements are discretized over the same grid)

time_scalingfloat

scaling factor applied to d2. The effect may be interpreted as moving d1 displacements along a factor (time_scaling) of d2.

Returns:
updatedarray, shape (the same as new_displacement)

the warped displacement field

mean_normthe mean norm of all vectors in current_displacement

### floating

dipy.align.imwarp.floating

alias of float32

### get_direction_and_spacings

dipy.align.imwarp.get_direction_and_spacings(affine, dim)

Extracts the rotational and spacing components from a matrix

Extracts the rotational and spacing (voxel dimensions) components from a matrix. An image gradient represents the local variation of the image’s gray values per voxel. Since we are iterating on the physical space, we need to compute the gradients as variation per millimeter, so we need to divide each gradient’s component by the voxel size along the corresponding axis, that’s what the spacings are used for. Since the image’s gradients are oriented along the grid axes, we also need to re-orient the gradients to be given in physical space coordinates.

Parameters:
affinearray, shape (k, k), k = 3, 4

the matrix transforming grid coordinates to physical space.

Returns:
directionarray, shape (k-1, k-1)

the rotational component of the input matrix

spacingsarray, shape (k-1,)

the scaling component (voxel size) of the matrix

### mult_aff

dipy.align.imwarp.mult_aff(A, B)

Returns the matrix product A.dot(B) considering None as the identity

Parameters:
Aarray, shape (n,k)
Barray, shape (k,m)
Returns:
The matrix product A.dot(B). If any of the input matrices is None, it is
treated as the identity matrix. If both matrices are None, None is returned

### CCMetric

Methods

 Computes one step bringing the static image towards the moving. Computes one step bringing the moving image towards the static. Frees the resources allocated during initialization Numerical value assigned by this metric to the current image pair Prepares the metric to compute one displacement field iteration. set_levels_above(levels) Informs the metric how many pyramid levels are above the current one set_levels_below(levels) Informs the metric how many pyramid levels are below the current one set_moving_image(moving_image, ...) Sets the moving image being compared against the static one. set_static_image(static_image, ...) Sets the static image being compared against the moving one. use_moving_image_dynamics(...) This is called by the optimizer just after setting the moving image use_static_image_dynamics(...) This is called by the optimizer just after setting the static image.

Normalized Cross-Correlation Similarity metric.

Parameters:
dimint (either 2 or 3)

the dimension of the image domain

sigma_diffthe standard deviation of the Gaussian smoothing kernel to

be applied to the update field at each iteration

the radius of the squared (cubic) neighborhood at each voxel to be considered to compute the cross correlation

compute_backward()

Computes one step bringing the static image towards the moving.

Computes the update displacement field to be used for registration of the static image towards the moving image

compute_forward()

Computes one step bringing the moving image towards the static.

Computes the update displacement field to be used for registration of the moving image towards the static image

free_iteration()

Frees the resources allocated during initialization

get_energy()

Numerical value assigned by this metric to the current image pair

Returns the Cross Correlation (data term) energy computed at the largest iteration

initialize_iteration()

Prepares the metric to compute one displacement field iteration.

Pre-computes the cross-correlation factors for efficient computation of the gradient of the Cross Correlation w.r.t. the displacement field. It also pre-computes the image gradients in the physical space by re-orienting the gradients in the voxel space using the corresponding affine transformations.

### EMMetric

class dipy.align.metrics.EMMetric(dim, smooth=1.0, inner_iter=5, q_levels=256, double_gradient=True, step_type='gauss_newton')

Methods

 Computes one step bringing the static image towards the moving. compute_demons_step([forward_step]) Demons step for EM metric Computes one step bringing the reference image towards the static. compute_gauss_newton_step([forward_step]) Computes the Gauss-Newton energy minimization step Frees the resources allocated during initialization The numerical value assigned by this metric to the current image pair Prepares the metric to compute one displacement field iteration. set_levels_above(levels) Informs the metric how many pyramid levels are above the current one set_levels_below(levels) Informs the metric how many pyramid levels are below the current one set_moving_image(moving_image, ...) Sets the moving image being compared against the static one. set_static_image(static_image, ...) Sets the static image being compared against the moving one. This is called by the optimizer just after setting the moving image. This is called by the optimizer just after setting the static image.
__init__(dim, smooth=1.0, inner_iter=5, q_levels=256, double_gradient=True, step_type='gauss_newton')

Expectation-Maximization Metric

Similarity metric based on the Expectation-Maximization algorithm to handle multi-modal images. The transfer function is modeled as a set of hidden random variables that are estimated at each iteration of the algorithm.

Parameters:
dimint (either 2 or 3)

the dimension of the image domain

smoothfloat

smoothness parameter, the larger the value the smoother the deformation field

inner_iterint

number of iterations to be performed at each level of the multi- resolution Gauss-Seidel optimization algorithm (this is not the number of steps per Gaussian Pyramid level, that parameter must be set for the optimizer, not the metric)

q_levelsnumber of quantization levels (equal to the number of hidden

variables in the EM algorithm)

if True, the gradient of the expected static image under the moving modality will be added to the gradient of the moving image, similarly, the gradient of the expected moving image under the static modality will be added to the gradient of the static image.

step_typestring (‘gauss_newton’, ‘demons’)

the optimization schedule to be used in the multi-resolution Gauss-Seidel optimization algorithm (not used if Demons Step is selected)

compute_backward()

Computes one step bringing the static image towards the moving.

Computes the update displacement field to be used for registration of the static image towards the moving image

compute_demons_step(forward_step=True)

Demons step for EM metric

Parameters:
forward_stepboolean

if True, computes the Demons step in the forward direction (warping the moving towards the static image). If False, computes the backward step (warping the static image to the moving image)

Returns:
displacementarray, shape (R, C, 2) or (S, R, C, 3)

the Demons step

compute_forward()

Computes one step bringing the reference image towards the static.

Computes the forward update field to register the moving image towards the static image in a gradient-based optimization algorithm

compute_gauss_newton_step(forward_step=True)

Computes the Gauss-Newton energy minimization step

Computes the Newton step to minimize this energy, i.e., minimizes the linearized energy function with respect to the regularized displacement field (this step does not require post-smoothing, as opposed to the demons step, which does not include regularization). To accelerate convergence we use the multi-grid Gauss-Seidel algorithm proposed by Bruhn and Weickert et al [Bruhn05]

Parameters:
forward_stepboolean

if True, computes the Newton step in the forward direction (warping the moving towards the static image). If False, computes the backward step (warping the static image to the moving image)

Returns:
displacementarray, shape (R, C, 2) or (S, R, C, 3)

the Newton step

References

[Bruhn05] Andres Bruhn and Joachim Weickert, “Towards ultimate motion

estimation: combining highest accuracy with real-time performance”, 10th IEEE International Conference on Computer Vision, 2005. ICCV 2005.

free_iteration()

Frees the resources allocated during initialization

get_energy()

The numerical value assigned by this metric to the current image pair

Returns the EM (data term) energy computed at the largest iteration

initialize_iteration()

Prepares the metric to compute one displacement field iteration.

Pre-computes the transfer functions (hidden random variables) and variances of the estimators. Also pre-computes the gradient of both input images. Note that once the images are transformed to the opposite modality, the gradient of the transformed images can be used with the gradient of the corresponding modality in the same fashion as diff-demons does for mono-modality images. If the flag self.use_double_gradient is True these gradients are averaged.

use_moving_image_dynamics(original_moving_image, transformation)

This is called by the optimizer just after setting the moving image.

EMMetric takes advantage of the image dynamics by computing the current moving image mask from the original_moving_image mask (warped by nearest neighbor interpolation)

Parameters:
original_moving_imagearray, shape (R, C) or (S, R, C)

the original moving image from which the current moving image was generated, the current moving image is the one that was provided via ‘set_moving_image(…)’, which may not be the same as the original moving image but a warped version of it.

transformationDiffeomorphicMap object

the transformation that was applied to the original_moving_image to generate the current moving image

use_static_image_dynamics(original_static_image, transformation)

This is called by the optimizer just after setting the static image.

EMMetric takes advantage of the image dynamics by computing the current static image mask from the originalstaticImage mask (warped by nearest neighbor interpolation)

Parameters:
original_static_imagearray, shape (R, C) or (S, R, C)

the original static image from which the current static image was generated, the current static image is the one that was provided via ‘set_static_image(…)’, which may not be the same as the original static image but a warped version of it (even the static image changes during Symmetric Normalization, not only the moving one).

transformationDiffeomorphicMap object

the transformation that was applied to the original_static_image to generate the current static image

### SSDMetric

class dipy.align.metrics.SSDMetric(dim, smooth=4, inner_iter=10, step_type='demons')

Methods

 Computes one step bringing the static image towards the moving. compute_demons_step([forward_step]) Demons step for SSD metric Computes one step bringing the reference image towards the static. compute_gauss_newton_step([forward_step]) Computes the Gauss-Newton energy minimization step Nothing to free for the SSD metric The numerical value assigned by this metric to the current image pair Prepares the metric to compute one displacement field iteration. set_levels_above(levels) Informs the metric how many pyramid levels are above the current one set_levels_below(levels) Informs the metric how many pyramid levels are below the current one set_moving_image(moving_image, ...) Sets the moving image being compared against the static one. set_static_image(static_image, ...) Sets the static image being compared against the moving one. use_moving_image_dynamics(...) This is called by the optimizer just after setting the moving image use_static_image_dynamics(...) This is called by the optimizer just after setting the static image.
__init__(dim, smooth=4, inner_iter=10, step_type='demons')

Sum of Squared Differences (SSD) Metric

Similarity metric for (mono-modal) nonlinear image registration defined by the sum of squared differences (SSD)

Parameters:
dimint (either 2 or 3)

the dimension of the image domain

smoothfloat

smoothness parameter, the larger the value the smoother the deformation field

inner_iterint

number of iterations to be performed at each level of the multi- resolution Gauss-Seidel optimization algorithm (this is not the number of steps per Gaussian Pyramid level, that parameter must be set for the optimizer, not the metric)

step_typestring

the displacement field step to be computed when ‘compute_forward’ and ‘compute_backward’ are called. Either ‘demons’ or ‘gauss_newton’

compute_backward()

Computes one step bringing the static image towards the moving.

Computes the updated displacement field to be used for registration of the static image towards the moving image

compute_demons_step(forward_step=True)

Demons step for SSD metric

Computes the demons step proposed by Vercauteren et al.[Vercauteren09] for the SSD metric.

Parameters:
forward_stepboolean

if True, computes the Demons step in the forward direction (warping the moving towards the static image). If False, computes the backward step (warping the static image to the moving image)

Returns:
displacementarray, shape (R, C, 2) or (S, R, C, 3)

the Demons step

References

[Vercauteren09] Tom Vercauteren, Xavier Pennec, Aymeric Perchant,

Nicholas Ayache, “Diffeomorphic Demons: Efficient Non-parametric Image Registration”, Neuroimage 2009

compute_forward()

Computes one step bringing the reference image towards the static.

Computes the update displacement field to be used for registration of the moving image towards the static image

compute_gauss_newton_step(forward_step=True)

Computes the Gauss-Newton energy minimization step

Minimizes the linearized energy function (Newton step) defined by the sum of squared differences of corresponding pixels of the input images with respect to the displacement field.

Parameters:
forward_stepboolean

if True, computes the Newton step in the forward direction (warping the moving towards the static image). If False, computes the backward step (warping the static image to the moving image)

Returns:
displacementarray, shape = static_image.shape + (3,)

if forward_step==True, the forward SSD Gauss-Newton step, else, the backward step

free_iteration()

Nothing to free for the SSD metric

get_energy()

The numerical value assigned by this metric to the current image pair

Returns the Sum of Squared Differences (data term) energy computed at the largest iteration

initialize_iteration()

Prepares the metric to compute one displacement field iteration.

Pre-computes the gradient of the input images to be used in the computation of the forward and backward steps.

### SimilarityMetric

class dipy.align.metrics.SimilarityMetric(dim)

Bases: object

Methods

 Computes one step bringing the static image towards the moving. Computes one step bringing the reference image towards the static. Releases the resources no longer needed by the metric Numerical value assigned by this metric to the current image pair Prepares the metric to compute one displacement field iteration. set_levels_above(levels) Informs the metric how many pyramid levels are above the current one set_levels_below(levels) Informs the metric how many pyramid levels are below the current one set_moving_image(moving_image, ...) Sets the moving image being compared against the static one. set_static_image(static_image, ...) Sets the static image being compared against the moving one. This is called by the optimizer just after setting the moving image This is called by the optimizer just after setting the static image.
__init__(dim)

Similarity Metric abstract class

A similarity metric is in charge of keeping track of the numerical value of the similarity (or distance) between the two given images. It also computes the update field for the forward and inverse displacement fields to be used in a gradient-based optimization algorithm. Note that this metric does not depend on any transformation (affine or non-linear) so it assumes the static and moving images are already warped

Parameters:
dimint (either 2 or 3)

the dimension of the image domain

abstract compute_backward()

Computes one step bringing the static image towards the moving.

Computes the backward update field to register the static image towards the moving image in a gradient-based optimization algorithm

abstract compute_forward()

Computes one step bringing the reference image towards the static.

Computes the forward update field to register the moving image towards the static image in a gradient-based optimization algorithm

abstract free_iteration()

Releases the resources no longer needed by the metric

This method is called by the RegistrationOptimizer after the required iterations have been computed (forward and / or backward) so that the SimilarityMetric can safely delete any data it computed as part of the initialization

abstract get_energy()

Numerical value assigned by this metric to the current image pair

Must return the numeric value of the similarity between the given static and moving images

abstract initialize_iteration()

Prepares the metric to compute one displacement field iteration.

This method will be called before any compute_forward or compute_backward call, this allows the Metric to pre-compute any useful information for speeding up the update computations. This initialization was needed in ANTS because the updates are called once per voxel. In Python this is unpractical, though.

set_levels_above(levels)

Informs the metric how many pyramid levels are above the current one

Informs this metric the number of pyramid levels above the current one. The metric may change its behavior (e.g. number of inner iterations) accordingly

Parameters:
levelsint

the number of levels above the current Gaussian Pyramid level

set_levels_below(levels)

Informs the metric how many pyramid levels are below the current one

Informs this metric the number of pyramid levels below the current one. The metric may change its behavior (e.g. number of inner iterations) accordingly

Parameters:
levelsint

the number of levels below the current Gaussian Pyramid level

set_moving_image(moving_image, moving_affine, moving_spacing, moving_direction)

Sets the moving image being compared against the static one.

Sets the moving image. The default behavior (of this abstract class) is simply to assign the reference to an attribute, but generalizations of the metric may need to perform other operations

Parameters:
moving_imagearray, shape (R, C) or (S, R, C)

the moving image

set_static_image(static_image, static_affine, static_spacing, static_direction)

Sets the static image being compared against the moving one.

Sets the static image. The default behavior (of this abstract class) is simply to assign the reference to an attribute, but generalizations of the metric may need to perform other operations

Parameters:
static_imagearray, shape (R, C) or (S, R, C)

the static image

use_moving_image_dynamics(original_moving_image, transformation)

This is called by the optimizer just after setting the moving image

This method allows the metric to compute any useful information from knowing how the current static image was generated (as the transformation of an original static image). This method is called by the optimizer just after it sets the static image. Transformation will be an instance of DiffeomorficMap or None if the original_moving_image equals self.moving_image.

Parameters:
original_moving_imagearray, shape (R, C) or (S, R, C)

original image from which the current moving image was generated

transformationDiffeomorphicMap object

the transformation that was applied to the original image to generate the current moving image

use_static_image_dynamics(original_static_image, transformation)

This is called by the optimizer just after setting the static image.

This method allows the metric to compute any useful information from knowing how the current static image was generated (as the transformation of an original static image). This method is called by the optimizer just after it sets the static image. Transformation will be an instance of DiffeomorficMap or None if the original_static_image equals self.moving_image.

Parameters:
original_static_imagearray, shape (R, C) or (S, R, C)

original image from which the current static image was generated

transformationDiffeomorphicMap object

the transformation that was applied to original image to generate the current static image

### floating

dipy.align.metrics.floating

alias of float32

Return the gradient of an N-dimensional array.

The gradient is computed using second order accurate central differences in the interior points and either first or second order accurate one-sides (forward or backwards) differences at the boundaries. The returned gradient hence has the same shape as the input array.

Parameters:
farray_like

An N-dimensional array containing samples of a scalar function.

varargslist of scalar or array, optional

Spacing between f values. Default unitary spacing for all dimensions. Spacing can be specified using:

1. single scalar to specify a sample distance for all dimensions.

2. N scalars to specify a constant sample distance for each dimension. i.e. dx, dy, dz, …

3. N arrays to specify the coordinates of the values along each dimension of F. The length of the array must match the size of the corresponding dimension

4. Any combination of N scalars/arrays with the meaning of 2. and 3.

If axis is given, the number of varargs must equal the number of axes. Default: 1.

edge_order{1, 2}, optional

Gradient is calculated using N-th order accurate differences at the boundaries. Default: 1.

New in version 1.9.1.

axisNone or int or tuple of ints, optional

Gradient is calculated only along the given axis or axes The default (axis = None) is to calculate the gradient for all the axes of the input array. axis may be negative, in which case it counts from the last to the first axis.

New in version 1.11.0.

Returns:

A list of ndarrays (or a single ndarray if there is only one dimension) corresponding to the derivatives of f with respect to each dimension. Each derivative has the same shape as f.

Notes

Assuming that $$f\in C^{3}$$ (i.e., $$f$$ has at least 3 continuous derivatives) and let $$h_{*}$$ be a non-homogeneous stepsize, we minimize the “consistency error” $$\eta_{i}$$ between the true gradient and its estimate from a linear combination of the neighboring grid-points:

$\eta_{i} = f_{i}^{\left(1\right)} - \left[ \alpha f\left(x_{i}\right) + \beta f\left(x_{i} + h_{d}\right) + \gamma f\left(x_{i}-h_{s}\right) \right]$

By substituting $$f(x_{i} + h_{d})$$ and $$f(x_{i} - h_{s})$$ with their Taylor series expansion, this translates into solving the following the linear system:

$\begin{split}\left\{ \begin{array}{r} \alpha+\beta+\gamma=0 \\ \beta h_{d}-\gamma h_{s}=1 \\ \beta h_{d}^{2}+\gamma h_{s}^{2}=0 \end{array} \right.\end{split}$

The resulting approximation of $$f_{i}^{(1)}$$ is the following:

$\hat f_{i}^{(1)} = \frac{ h_{s}^{2}f\left(x_{i} + h_{d}\right) + \left(h_{d}^{2} - h_{s}^{2}\right)f\left(x_{i}\right) - h_{d}^{2}f\left(x_{i}-h_{s}\right)} { h_{s}h_{d}\left(h_{d} + h_{s}\right)} + \mathcal{O}\left(\frac{h_{d}h_{s}^{2} + h_{s}h_{d}^{2}}{h_{d} + h_{s}}\right)$

It is worth noting that if $$h_{s}=h_{d}$$ (i.e., data are evenly spaced) we find the standard second order approximation:

$\hat f_{i}^{(1)}= \frac{f\left(x_{i+1}\right) - f\left(x_{i-1}\right)}{2h} + \mathcal{O}\left(h^{2}\right)$

With a similar procedure the forward/backward approximations used for boundaries can be derived.

References

[1]

Quarteroni A., Sacco R., Saleri F. (2007) Numerical Mathematics (Texts in Applied Mathematics). New York: Springer.

[2]

Durran D. R. (1999) Numerical Methods for Wave Equations in Geophysical Fluid Dynamics. New York: Springer.

[3]

Fornberg B. (1988) Generation of Finite Difference Formulas on Arbitrarily Spaced Grids, Mathematics of Computation 51, no. 184 : 699-706. PDF.

Examples

>>> f = np.array([1, 2, 4, 7, 11, 16], dtype=float)
array([1. , 1.5, 2.5, 3.5, 4.5, 5. ])
array([0.5 ,  0.75,  1.25,  1.75,  2.25,  2.5 ])


Spacing can be also specified with an array that represents the coordinates of the values F along the dimensions. For instance a uniform spacing:

>>> x = np.arange(f.size)
array([1. ,  1.5,  2.5,  3.5,  4.5,  5. ])


Or a non uniform one:

>>> x = np.array([0., 1., 1.5, 3.5, 4., 6.], dtype=float)
array([1. ,  3. ,  3.5,  6.7,  6.9,  2.5])


For two dimensional arrays, the return will be two arrays ordered by axis. In this example the first array stands for the gradient in rows and the second one in columns direction:

>>> np.gradient(np.array([[1, 2, 6], [3, 4, 5]], dtype=float))
[array([[ 2.,  2., -1.],
[ 2.,  2., -1.]]), array([[1. , 2.5, 4. ],
[1. , 1. , 1. ]])]


In this example the spacing is also specified: uniform for axis=0 and non uniform for axis=1

>>> dx = 2.
>>> y = [1., 1.5, 3.5]
>>> np.gradient(np.array([[1, 2, 6], [3, 4, 5]], dtype=float), dx, y)
[array([[ 1. ,  1. , -0.5],
[ 1. ,  1. , -0.5]]), array([[2. , 2. , 2. ],
[2. , 1.7, 0.5]])]


It is possible to specify how boundaries are treated using edge_order

>>> x = np.array([0, 1, 2, 3, 4])
>>> f = x**2
array([1.,  2.,  4.,  6.,  7.])
array([0., 2., 4., 6., 8.])


The axis keyword can be used to specify a subset of axes of which the gradient is calculated

>>> np.gradient(np.array([[1, 2, 6], [3, 4, 5]], dtype=float), axis=0)
array([[ 2.,  2., -1.],
[ 2.,  2., -1.]])


### v_cycle_2d

dipy.align.metrics.v_cycle_2d(n, k, delta_field, sigma_sq_field, gradient_field, target, lambda_param, displacement, depth=0)

Multi-resolution Gauss-Seidel solver using V-type cycles

Multi-resolution Gauss-Seidel solver: solves the Gauss-Newton linear system by first filtering (GS-iterate) the current level, then solves for the residual at a coarser resolution and finally refines the solution at the current resolution. This scheme corresponds to the V-cycle proposed by Bruhn and Weickert[Bruhn05].

Parameters:
nint

number of levels of the multi-resolution algorithm (it will be called recursively until level n == 0)

kint

the number of iterations at each multi-resolution level

delta_fieldarray, shape (R, C)

the difference between the static and moving image (the ‘derivative w.r.t. time’ in the optical flow model)

sigma_sq_fieldarray, shape (R, C)

the variance of the gray level value at each voxel, according to the EM model (for SSD, it is 1 for all voxels). Inf and 0 values are processed specially to support infinite and zero variance.

the gradient of the moving image

targetarray, shape (R, C, 2)

right-hand side of the linear system to be solved in the Weickert’s multi-resolution algorithm

lambda_paramfloat

smoothness parameter, the larger its value the smoother the displacement field

displacementarray, shape (R, C, 2)

the displacement field to start the optimization from

Returns:
energythe energy of the EM (or SSD if sigmafield[…]==1) metric at this

iteration

References

[Bruhn05] Andres Bruhn and Joachim Weickert, “Towards ultimate motion

estimation: combining the highest accuracy with real-time performance”, 10th IEEE International Conference on Computer Vision, 2005. ICCV 2005.

### v_cycle_3d

dipy.align.metrics.v_cycle_3d(n, k, delta_field, sigma_sq_field, gradient_field, target, lambda_param, displacement, depth=0)

Multi-resolution Gauss-Seidel solver using V-type cycles

Multi-resolution Gauss-Seidel solver: solves the linear system by first filtering (GS-iterate) the current level, then solves for the residual at a coarser resolution and finally refines the solution at the current resolution. This scheme corresponds to the V-cycle proposed by Bruhn and Weickert[1]. [1] Andres Bruhn and Joachim Weickert, “Towards ultimate motion estimation:

combining highest accuracy with real-time performance”, 10th IEEE International Conference on Computer Vision, 2005. ICCV 2005.

Parameters:
nint

number of levels of the multi-resolution algorithm (it will be called recursively until level n == 0)

kint

the number of iterations at each multi-resolution level

delta_fieldarray, shape (S, R, C)

the difference between the static and moving image (the ‘derivative w.r.t. time’ in the optical flow model)

sigma_sq_fieldarray, shape (S, R, C)

the variance of the gray level value at each voxel, according to the EM model (for SSD, it is 1 for all voxels). Inf and 0 values are processed specially to support infinite and zero variance.

gradient_fieldarray, shape (S, R, C, 3)

the gradient of the moving image

targetarray, shape (S, R, C, 3)

right-hand side of the linear system to be solved in the Weickert’s multi-resolution algorithm

lambda_paramfloat

smoothness parameter, the larger its value the smoother the displacement field

displacementarray, shape (S, R, C, 3)

the displacement field to start the optimization from

Returns:
energythe energy of the EM (or SSD if sigmafield[…]==1) metric at this

iteration

### ParzenJointHistogram

class dipy.align.parzenhist.ParzenJointHistogram(nbins)

Bases: object

Methods

 bin_index(xnorm) Bin index associated with the given normalized intensity Maps intensity x to the range covered by the moving histogram Maps intensity x to the range covered by the static histogram setup(static, moving[, smask, mmask]) Compute histogram settings to store the PDF of input images update_gradient_dense(theta, transform, ...) Computes the Gradient of the joint PDF w.r.t. update_gradient_sparse(theta, transform, ...) Computes the Gradient of the joint PDF w.r.t. update_pdfs_dense(static, moving[, smask, mmask]) Computes the Probability Density Functions of two images update_pdfs_sparse(sval, mval) Computes the Probability Density Functions from a set of samples
__init__(nbins)

Computes joint histogram and derivatives with Parzen windows

Base class to compute joint and marginal probability density functions and their derivatives with respect to a transform’s parameters. The smooth histograms are computed by using Parzen windows [Parzen62] with a cubic spline kernel, as proposed by Mattes et al. [Mattes03]. This implementation is not tied to any optimization (registration) method, the idea is that information-theoretic matching functionals (such as Mutual Information) can inherit from this class to perform the low-level computations of the joint intensity distributions and its gradient w.r.t. the transform parameters. The derived class can then compute the similarity/dissimilarity measure and gradient, and finally communicate the results to the appropriate optimizer.

Parameters:
nbinsint

the number of bins of the joint and marginal probability density functions (the actual number of bins of the joint PDF is nbins**2)

Notes

We need this class in cython to allow _joint_pdf_gradient_dense_2d and _joint_pdf_gradient_dense_3d to use a nogil Jacobian function (obtained from an instance of the Transform class), which allows us to evaluate Jacobians at all the sampling points (maybe the full grid) inside a nogil loop.

The reason we need a class is to encapsulate all the parameters related to the joint and marginal distributions.

References

[Parzen62] E. Parzen. On the estimation of a probability density

function and the mode. Annals of Mathematical Statistics, 33(3), 1065-1076, 1962.

[Mattes03] Mattes, D., Haynor, D. R., Vesselle, H., Lewellen, T. K.,

& Eubank, W. PET-CT image registration in the chest using free-form deformations. IEEE Transactions on Medical Imaging, 22(1), 120-8, 2003.

bin_index(xnorm)

Bin index associated with the given normalized intensity

The return value is an integer in [padding, nbins - 1 - padding]

Parameters:
xnormfloat

intensity value normalized to the range covered by the histogram

Returns:
binint

the bin index associated with the given normalized intensity

bin_normalize_moving(x)

Maps intensity x to the range covered by the moving histogram

If the input intensity is in [self.mmin, self.mmax] then the normalized intensity will be in [self.padding, self.nbins - self.padding]

Parameters:
xfloat

the intensity to be normalized

Returns:
xnormfloat

normalized intensity to the range covered by the moving histogram

bin_normalize_static(x)

Maps intensity x to the range covered by the static histogram

If the input intensity is in [self.smin, self.smax] then the normalized intensity will be in [self.padding, self.nbins - self.padding]

Parameters:
xfloat

the intensity to be normalized

Returns:
xnormfloat

normalized intensity to the range covered by the static histogram

Compute histogram settings to store the PDF of input images

Parameters:
staticarray

static image

movingarray

moving image

mask of static object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). If None, the behaviour is equivalent to smask=ones_like(static)

mask of moving object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). If None, the behaviour is equivalent to mmask=ones_like(static)

Computes the Gradient of the joint PDF w.r.t. transform parameters

Computes the vector of partial derivatives of the joint histogram w.r.t. each transformation parameter.

Parameters:
thetaarray, shape (n,)

parameters of the transformation to compute the gradient from

transforminstance of Transform

the transformation with respect to whose parameters the gradient must be computed

staticarray, shape (S, R, C)

static image

movingarray, shape (S, R, C)

moving image

grid2worldarray, shape (4, 4)

we assume that both images have already been sampled at a common grid. This transform must map voxel coordinates of this common grid to physical coordinates of its corresponding voxel in the moving image. For example, if the moving image was sampled on the static image’s grid (this is the typical setting) using an aligning matrix A, then

1. grid2world = A.dot(static_affine)

where static_affine is the transformation mapping static image’s grid coordinates to physical space.

mgradientarray, shape (S, R, C, 3)

the gradient of the moving image

smaskarray, shape (S, R, C), optional

mask of static object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). The default is None, indicating all voxels are considered.

mmaskarray, shape (S, R, C), optional

mask of moving object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). The default is None, indicating all voxels are considered.

Computes the Gradient of the joint PDF w.r.t. transform parameters

Computes the vector of partial derivatives of the joint histogram w.r.t. each transformation parameter.

The list of intensities sval and mval are assumed to be sampled from the static and moving images, respectively, at the same physical points. Of course, the images may not be perfectly aligned at the moment the sampling was performed. The resulting gradient corresponds to the paired intensities according to the alignment at the moment the images were sampled.

Parameters:
thetaarray, shape (n,)

parameters to compute the gradient at

transforminstance of Transform

the transformation with respect to whose parameters the gradient must be computed

svalarray, shape (m,)

sampled intensities from the static image at sampled_points

mvalarray, shape (m,)

sampled intensities from the moving image at sampled_points

sample_pointsarray, shape (m, 3)

coordinates (in physical space) of the points the images were sampled at

the gradient of the moving image at the sample points

Computes the Probability Density Functions of two images

The joint PDF is stored in self.joint. The marginal distributions corresponding to the static and moving images are computed and stored in self.smarginal and self.mmarginal, respectively.

Parameters:
staticarray, shape (S, R, C)

static image

movingarray, shape (S, R, C)

moving image

mask of static object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). If None, ones_like(static) is used as mask.

mask of moving object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). If None, ones_like(moving) is used as mask.

update_pdfs_sparse(sval, mval)

Computes the Probability Density Functions from a set of samples

The list of intensities sval and mval are assumed to be sampled from the static and moving images, respectively, at the same physical points. Of course, the images may not be perfectly aligned at the moment the sampling was performed. The resulting distributions corresponds to the paired intensities according to the alignment at the moment the images were sampled.

The joint PDF is stored in self.joint. The marginal distributions corresponding to the static and moving images are computed and stored in self.smarginal and self.mmarginal, respectively.

Parameters:
svalarray, shape (n,)

sampled intensities from the static image at sampled_points

mvalarray, shape (n,)

sampled intensities from the moving image at sampled_points

### compute_parzen_mi

dipy.align.parzenhist.compute_parzen_mi()

Computes the mutual information and its gradient (if requested)

Parameters:
jointarray, shape (nbins, nbins)

the joint intensity distribution

the gradient of the joint distribution w.r.t. the transformation parameters

smarginalarray, shape (nbins,)

the marginal intensity distribution of the static image

mmarginalarray, shape (nbins,)

the marginal intensity distribution of the moving image

the buffer in which to write the gradient of the mutual information. If None, the gradient is not computed

### cubic_spline

dipy.align.parzenhist.cubic_spline()

Evaluates the cubic spline at a set of values

Parameters:
xarray, shape (n)

input values

### cubic_spline_derivative

dipy.align.parzenhist.cubic_spline_derivative()

Evaluates the cubic spline derivative at a set of values

Parameters:
xarray, shape (n)

input values

### sample_domain_regular

dipy.align.parzenhist.sample_domain_regular()

Take floor(total_voxels/k) samples from a (2D or 3D) grid

The sampling is made by taking all pixels whose index (in lexicographical order) is a multiple of k. Each selected point is slightly perturbed by adding a realization of a normally distributed random variable and then mapped to physical space by the given grid-to-space transform.

The lexicographical order of a pixels in a grid of shape (a, b, c) is defined by assigning to each voxel position (i, j, k) the integer index

F((i, j, k)) = i * (b * c) + j * (c) + k

and sorting increasingly by this index.

Parameters:
kint

the sampling rate, as described before

shapearray, shape (dim,)

the shape of the grid to be sampled

grid2worldarray, shape (dim+1, dim+1)

the grid-to-space transform

sigmafloat

the standard deviation of the Normal random distortion to be applied to the sampled points

Returns:
samplesarray, shape (total_pixels//k, dim)

the matrix whose rows are the sampled points

Examples

>>> from dipy.align.parzenhist import sample_domain_regular
>>> import dipy.align.vector_fields as vf
>>> shape = np.array((10, 10), dtype=np.int32)
>>> sigma = 0
>>> dim = len(shape)
>>> grid2world = np.eye(dim+1)
>>> n = shape[0]*shape[1]
>>> k = 2
>>> samples = sample_domain_regular(k, shape, grid2world, sigma)
>>> (samples.shape[0], samples.shape[1]) == (n//k, dim)
True
>>> isamples = np.array(samples, dtype=np.int32)
>>> indices = (isamples[:, 0] * shape[1] + isamples[:, 1])
>>> len(set(indices)) == len(indices)
True
>>> (indices%k).sum()
0


### Pool

Returns a process pool object

### affine_transform

dipy.align.reslice.affine_transform(input, matrix, offset=0.0, output_shape=None, output=None, order=3, mode='constant', cval=0.0, prefilter=True)

Apply an affine transformation.

Given an output image pixel index vector o, the pixel value is determined from the input image at position np.dot(matrix, o) + offset.

This does ‘pull’ (or ‘backward’) resampling, transforming the output space to the input to locate data. Affine transformations are often described in the ‘push’ (or ‘forward’) direction, transforming input to output. If you have a matrix for the ‘push’ transformation, use its inverse (numpy.linalg.inv()) in this function.

Parameters:
inputarray_like

The input array.

matrixndarray

The inverse coordinate transformation matrix, mapping output coordinates to input coordinates. If ndim is the number of dimensions of input, the given matrix must have one of the following shapes:

• (ndim, ndim): the linear transformation matrix for each output coordinate.

• (ndim,): assume that the 2-D transformation matrix is diagonal, with the diagonal specified by the given value. A more efficient algorithm is then used that exploits the separability of the problem.

• (ndim + 1, ndim + 1): assume that the transformation is specified using homogeneous coordinates [1]. In this case, any value passed to offset is ignored.

• (ndim, ndim + 1): as above, but the bottom row of a homogeneous transformation matrix is always [0, 0, ..., 1], and may be omitted.

offsetfloat or sequence, optional

The offset into the array where the transform is applied. If a float, offset is the same for each axis. If a sequence, offset should contain one value for each axis.

output_shapetuple of ints, optional

Shape tuple.

outputarray or dtype, optional

The array in which to place the output, or the dtype of the returned array. By default an array of the same dtype as input will be created.

orderint, optional

The order of the spline interpolation, default is 3. The order has to be in the range 0-5.

mode{‘reflect’, ‘grid-mirror’, ‘constant’, ‘grid-constant’, ‘nearest’, ‘mirror’, ‘grid-wrap’, ‘wrap’}, optional

The mode parameter determines how the input array is extended beyond its boundaries. Default is ‘constant’. Behavior for each valid value is as follows (see additional plots and details on boundary modes):

‘reflect’ (d c b a | a b c d | d c b a)

The input is extended by reflecting about the edge of the last pixel. This mode is also sometimes referred to as half-sample symmetric.

‘grid-mirror’

This is a synonym for ‘reflect’.

‘constant’ (k k k k | a b c d | k k k k)

The input is extended by filling all values beyond the edge with the same constant value, defined by the cval parameter. No interpolation is performed beyond the edges of the input.

‘grid-constant’ (k k k k | a b c d | k k k k)

The input is extended by filling all values beyond the edge with the same constant value, defined by the cval parameter. Interpolation occurs for samples outside the input’s extent as well.

‘nearest’ (a a a a | a b c d | d d d d)

The input is extended by replicating the last pixel.

‘mirror’ (d c b | a b c d | c b a)

The input is extended by reflecting about the center of the last pixel. This mode is also sometimes referred to as whole-sample symmetric.

‘grid-wrap’ (a b c d | a b c d | a b c d)

The input is extended by wrapping around to the opposite edge.

‘wrap’ (d b c d | a b c d | b c a b)

The input is extended by wrapping around to the opposite edge, but in a way such that the last point and initial point exactly overlap. In this case it is not well defined which sample will be chosen at the point of overlap.

cvalscalar, optional

Value to fill past edges of input if mode is ‘constant’. Default is 0.0.

prefilterbool, optional

Determines if the input array is prefiltered with spline_filter before interpolation. The default is True, which will create a temporary float64 array of filtered values if order > 1. If setting this to False, the output will be slightly blurred if order > 1, unless the input is prefiltered, i.e. it is the result of calling spline_filter on the original input.

Returns:
affine_transformndarray

The transformed input.

Notes

The given matrix and offset are used to find for each point in the output the corresponding coordinates in the input by an affine transformation. The value of the input at those coordinates is determined by spline interpolation of the requested order. Points outside the boundaries of the input are filled according to the given mode.

Changed in version 0.18.0: Previously, the exact interpretation of the affine transformation depended on whether the matrix was supplied as a 1-D or a 2-D array. If a 1-D array was supplied to the matrix parameter, the output pixel value at index o was determined from the input image at position matrix * (o + offset).

For complex-valued input, this function transforms the real and imaginary components independently.

New in version 1.6.0: Complex-valued support added.

References

### determine_num_processes

dipy.align.reslice.determine_num_processes(num_processes)

Determine the effective number of processes for parallelization.

• For num_processes = None return the maximum number of cores retrieved

by cpu_count().

• For num_processes > 0, return this value.

• For num_processes < 0, return the maximal number of cores minus

num_processes + 1. In particular num_processes = -1 will use as many cores as possible.

• For num_processes = 0 a ValueError is raised.

Parameters:
num_processesint or None

Desired number of processes to be used.

### reslice

dipy.align.reslice.reslice(data, affine, zooms, new_zooms, order=1, mode='constant', cval=0, num_processes=1)

Reslice data with new voxel resolution defined by new_zooms.

Parameters:
dataarray, shape (I,J,K) or (I,J,K,N)

3d volume or 4d volume with datasets

affinearray, shape (4,4)

mapping from voxel coordinates to world coordinates

zoomstuple, shape (3,)

voxel size for (i,j,k) dimensions

new_zoomstuple, shape (3,)

new voxel size for (i,j,k) after resampling

orderint, from 0 to 5

order of interpolation for resampling/reslicing, 0 nearest interpolation, 1 trilinear etc.. if you don’t want any smoothing 0 is the option you need.

modestring (‘constant’, ‘nearest’, ‘reflect’ or ‘wrap’)

Points outside the boundaries of the input are filled according to the given mode.

cvalfloat

Value used for points outside the boundaries of the input if mode=’constant’.

num_processesint, optional

Split the calculation to a pool of children processes. This only applies to 4D data arrays. Default is 1. If < 0 the maximal number of cores minus num_processes + 1 is used (enter -1 to use as many cores as possible). 0 raises an error.

Returns:
data2array, shape (I,J,K) or (I,J,K,N)

datasets resampled into isotropic voxel size

affine2array, shape (4,4)

new affine for the resampled image

Examples

>>> from dipy.io.image import load_nifti
>>> from dipy.align.reslice import reslice
>>> from dipy.data import get_fnames
>>> f_name = get_fnames('aniso_vox')
>>> data, affine, zooms = load_nifti(f_name, return_voxsize=True)
>>> data.shape == (58, 58, 24)
True
>>> zooms
(4.0, 4.0, 5.0)
>>> new_zooms = (3.,3.,3.)
>>> new_zooms
(3.0, 3.0, 3.0)
>>> data2, affine2 = reslice(data, affine, zooms, new_zooms)
>>> data2.shape == (77, 77, 40)
True


### IsotropicScaleSpace

class dipy.align.scalespace.IsotropicScaleSpace(image, factors, sigmas, image_grid2world=None, input_spacing=None, mask0=False)

Bases: ScaleSpace

Methods

 get_affine(level) Voxel-to-space transformation at a given level. get_affine_inv(level) Space-to-voxel transformation at a given level. get_domain_shape(level) Shape the sub-sampled image must have at a particular level. get_expand_factors(from_level, to_level) Ratio of voxel size from pyramid level from_level to to_level. get_image(level) Smoothed image at a given level. get_scaling(level) Adjustment factor for input-spacing to reflect voxel sizes at level. get_sigmas(level) Smoothing parameters used at a given level. get_spacing(level) Spacings the sub-sampled image must have at a particular level. print_level(level) Prints properties of a pyramid level.
__init__(image, factors, sigmas, image_grid2world=None, input_spacing=None, mask0=False)

IsotropicScaleSpace.

Computes the Scale Space representation of an image using isotropic smoothing kernels for all scales. The scale space is simply a list of images produced by smoothing the input image with a Gaussian kernel with different smoothing parameters.

This specialization of ScaleSpace allows the user to provide custom scale and smoothing factors for all scales.

Parameters:
imagearray, shape (r,c) or (s, r, c) where s is the number of

slices, r is the number of rows and c is the number of columns of the input image.

factorslist of floats

custom scale factors to build the scale space (one factor for each scale).

sigmaslist of floats

custom smoothing parameter to build the scale space (one parameter for each scale).

image_grid2worldarray, shape (dim + 1, dim + 1), optional

the grid-to-space transform of the image grid. The default is the identity matrix.

input_spacingarray, shape (dim,), optional

the spacing (voxel size) between voxels in physical space. The default if 1.0 along all axes.

if True, all smoothed images will be zero at all voxels that are zero in the input image. The default is False.

### ScaleSpace

class dipy.align.scalespace.ScaleSpace(image, num_levels, image_grid2world=None, input_spacing=None, sigma_factor=0.2, mask0=False)

Bases: object

Methods

 get_affine(level) Voxel-to-space transformation at a given level. get_affine_inv(level) Space-to-voxel transformation at a given level. get_domain_shape(level) Shape the sub-sampled image must have at a particular level. get_expand_factors(from_level, to_level) Ratio of voxel size from pyramid level from_level to to_level. get_image(level) Smoothed image at a given level. get_scaling(level) Adjustment factor for input-spacing to reflect voxel sizes at level. get_sigmas(level) Smoothing parameters used at a given level. get_spacing(level) Spacings the sub-sampled image must have at a particular level. print_level(level) Prints properties of a pyramid level.
__init__(image, num_levels, image_grid2world=None, input_spacing=None, sigma_factor=0.2, mask0=False)

ScaleSpace.

Computes the Scale Space representation of an image. The scale space is simply a list of images produced by smoothing the input image with a Gaussian kernel with increasing smoothing parameter. If the image’s voxels are isotropic, the smoothing will be the same along all directions: at level L = 0, 1, …, the sigma is given by $$s * ( 2^L - 1 )$$. If the voxel dimensions are not isotropic, then the smoothing is weaker along low resolution directions.

Parameters:
imagearray, shape (r,c) or (s, r, c) where s is the number of

slices, r is the number of rows and c is the number of columns of the input image.

num_levelsint

the desired number of levels (resolutions) of the scale space

image_grid2worldarray, shape (dim + 1, dim + 1), optional

the grid-to-space transform of the image grid. The default is the identity matrix

input_spacingarray, shape (dim,), optional

the spacing (voxel size) between voxels in physical space. The default is 1.0 along all axes

sigma_factorfloat, optional

the smoothing factor to be used in the construction of the scale space. The default is 0.2

if True, all smoothed images will be zero at all voxels that are zero in the input image. The default is False.

get_affine(level)

Voxel-to-space transformation at a given level.

Returns the voxel-to-space transformation associated with the sub-sampled image at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the sub-sampled images must have).

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get affine transform from

Returns:
the affine (voxel-to-space) transform at the requested resolution
or None if an invalid level was requested
get_affine_inv(level)

Space-to-voxel transformation at a given level.

Returns the space-to-voxel transformation associated with the sub-sampled image at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the sub-sampled images must have).

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the inverse transform from

Returns:
the inverse (space-to-voxel) transform at the requested resolution or
None if an invalid level was requested
get_domain_shape(level)

Shape the sub-sampled image must have at a particular level.

Returns the shape the sub-sampled image must have at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the sub-sampled images must have).

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the sub-sampled shape from

Returns:
the sub-sampled shape at the requested resolution or None if an
invalid level was requested
get_expand_factors(from_level, to_level)

Ratio of voxel size from pyramid level from_level to to_level.

Given two scale space resolutions a = from_level, b = to_level, returns the ratio of voxels size at level b to voxel size at level a (the factor that must be used to multiply voxels at level a to ‘expand’ them to level b).

Parameters:
from_levelint, 0 <= from_level < L, (L = number of resolutions)

the resolution to expand voxels from

to_levelint, 0 <= to_level < from_level

the resolution to expand voxels to

Returns:
factorsarray, shape (k,), k = 2, 3

the expand factors (a scalar for each voxel dimension)

get_image(level)

Smoothed image at a given level.

Returns the smoothed image at the requested level in the Scale Space.

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the smooth image from

Returns:
the smooth image at the requested resolution or None if an invalid
level was requested
get_scaling(level)

Adjustment factor for input-spacing to reflect voxel sizes at level.

Returns the scaling factor that needs to be applied to the input spacing (the voxel sizes of the image at level 0 of the scale space) to transform them to voxel sizes at the requested level.

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the scalings from

Returns:
the scaling factors from the original spacing to the spacings at the
requested level
get_sigmas(level)

Smoothing parameters used at a given level.

Returns the smoothing parameters (a scalar for each axis) used at the requested level of the scale space

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the smoothing parameters from

Returns:
the smoothing parameters at the requested level
get_spacing(level)

Spacings the sub-sampled image must have at a particular level.

Returns the spacings (voxel sizes) the sub-sampled image must have at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the sub-sampled images must have).

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to get the sub-sampled shape from

Returns:
the spacings (voxel sizes) at the requested resolution or None if an
invalid level was requested
print_level(level)

Prints properties of a pyramid level.

Prints the properties of a level of this scale space to standard output

Parameters:
levelint, 0 <= from_level < L, (L = number of resolutions)

the scale space level to be printed

### floating

dipy.align.scalespace.floating

alias of float32

### gaussian_filter

dipy.align.scalespace.gaussian_filter(input, sigma, order=0, output=None, mode='reflect', cval=0.0, truncate=4.0, *, radius=None)

Multidimensional Gaussian filter.

Parameters:
inputarray_like

The input array.

sigmascalar or sequence of scalars

Standard deviation for Gaussian kernel. The standard deviations of the Gaussian filter are given for each axis as a sequence, or as a single number, in which case it is equal for all axes.

orderint or sequence of ints, optional

The order of the filter along each axis is given as a sequence of integers, or as a single number. An order of 0 corresponds to convolution with a Gaussian kernel. A positive order corresponds to convolution with that derivative of a Gaussian.

outputarray or dtype, optional

The array in which to place the output, or the dtype of the returned array. By default an array of the same dtype as input will be created.

modestr or sequence, optional

The mode parameter determines how the input array is extended when the filter overlaps a border. By passing a sequence of modes with length equal to the number of dimensions of the input array, different modes can be specified along each axis. Default value is ‘reflect’. The valid values and their behavior is as follows:

‘reflect’ (d c b a | a b c d | d c b a)

The input is extended by reflecting about the edge of the last pixel. This mode is also sometimes referred to as half-sample symmetric.

‘constant’ (k k k k | a b c d | k k k k)

The input is extended by filling all values beyond the edge with the same constant value, defined by the cval parameter.

‘nearest’ (a a a a | a b c d | d d d d)

The input is extended by replicating the last pixel.

‘mirror’ (d c b | a b c d | c b a)

The input is extended by reflecting about the center of the last pixel. This mode is also sometimes referred to as whole-sample symmetric.

‘wrap’ (a b c d | a b c d | a b c d)

The input is extended by wrapping around to the opposite edge.

For consistency with the interpolation functions, the following mode names can also be used:

‘grid-constant’

This is a synonym for ‘constant’.

‘grid-mirror’

This is a synonym for ‘reflect’.

‘grid-wrap’

This is a synonym for ‘wrap’.

cvalscalar, optional

Value to fill past edges of input if mode is ‘constant’. Default is 0.0.

truncatefloat, optional

Truncate the filter at this many standard deviations. Default is 4.0.

radiusNone or int or sequence of ints, optional

Radius of the Gaussian kernel. The radius are given for each axis as a sequence, or as a single number, in which case it is equal for all axes. If specified, the size of the kernel along each axis will be 2*radius + 1, and truncate is ignored. Default is None.

Returns:
gaussian_filterndarray

Returned array of same shape as input.

Notes

The multidimensional filter is implemented as a sequence of 1-D convolution filters. The intermediate arrays are stored in the same data type as the output. Therefore, for output types with a limited precision, the results may be imprecise because intermediate results may be stored with insufficient precision.

The Gaussian kernel will have size 2*radius + 1 along each axis where radius = round(truncate * sigma).

Examples

>>> from scipy.ndimage import gaussian_filter
>>> import numpy as np
>>> a = np.arange(50, step=2).reshape((5,5))
>>> a
array([[ 0,  2,  4,  6,  8],
[10, 12, 14, 16, 18],
[20, 22, 24, 26, 28],
[30, 32, 34, 36, 38],
[40, 42, 44, 46, 48]])
>>> gaussian_filter(a, sigma=1)
array([[ 4,  6,  8,  9, 11],
[10, 12, 14, 15, 17],
[20, 22, 24, 25, 27],
[29, 31, 33, 34, 36],
[35, 37, 39, 40, 42]])

>>> from scipy import datasets
>>> import matplotlib.pyplot as plt
>>> fig = plt.figure()
>>> plt.gray()  # show the filtered result in grayscale
>>> ax1 = fig.add_subplot(121)  # left side
>>> ax2 = fig.add_subplot(122)  # right side
>>> ascent = datasets.ascent()
>>> result = gaussian_filter(ascent, sigma=5)
>>> ax1.imshow(ascent)
>>> ax2.imshow(result)
>>> plt.show()


### BundleMinDistanceAsymmetricMetric

Asymmetric Bundle-based Minimum distance.

This is a cost function that can be used by the StreamlineLinearRegistration class.

Methods

 distance(xopt) Distance calculated from this Metric. setup(static, moving) Setup static and moving sets of streamlines.

An abstract class for the metric used for streamline registration.

If the two sets of streamlines match exactly then method distance of this object should be minimum.

Parameters:

Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus |num_threads + 1| is used (enter -1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.

distance(xopt)

Distance calculated from this Metric.

Parameters:
xoptsequence

List of affine parameters as an 1D vector

### BundleMinDistanceMatrixMetric

Bundle-based Minimum Distance aka BMD

This is the cost function used by the StreamlineLinearRegistration

Notes

The difference with BundleMinDistanceMetric is that this creates the entire distance matrix and therefore requires more memory.

Methods

 setup(static, moving) distance(xopt)

An abstract class for the metric used for streamline registration.

If the two sets of streamlines match exactly then method distance of this object should be minimum.

Parameters:

Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus |num_threads + 1| is used (enter -1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.

distance(xopt)

Distance calculated from this Metric.

Parameters:
xoptsequence

List of affine parameters as an 1D vector

setup(static, moving)

Setup static and moving sets of streamlines.

Parameters:
staticstreamlines

Fixed or reference set of streamlines.

movingstreamlines

Moving streamlines.

Notes

Call this after the object is initiated and before distance.

Num_threads is not used in this class. Use BundleMinDistanceMetric for a faster, threaded and less memory hungry metric

### BundleMinDistanceMetric

Bundle-based Minimum Distance aka BMD

This is the cost function used by the StreamlineLinearRegistration.

References

[Garyfallidis14]

Garyfallidis et al., “Direct native-space fiber bundle alignment for group comparisons”, ISMRM, 2014.

Methods

 setup(static, moving) distance(xopt)

An abstract class for the metric used for streamline registration.

If the two sets of streamlines match exactly then method distance of this object should be minimum.

Parameters:

Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus |num_threads + 1| is used (enter -1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.

distance(xopt)

Distance calculated from this Metric.

Parameters:
xoptsequence

List of affine parameters as an 1D vector,

setup(static, moving)

Setup static and moving sets of streamlines.

Parameters:
staticstreamlines

Fixed or reference set of streamlines.

movingstreamlines

Moving streamlines.

Notes

Call this after the object is initiated and before distance.

### BundleSumDistanceMatrixMetric

Bundle-based Sum Distance aka BMD

This is a cost function that can be used by the StreamlineLinearRegistration class.

Notes

The difference with BundleMinDistanceMatrixMetric is that it uses uses the sum of the distance matrix and not the sum of mins.

Methods

 setup(static, moving) distance(xopt)

An abstract class for the metric used for streamline registration.

If the two sets of streamlines match exactly then method distance of this object should be minimum.

Parameters:

Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus |num_threads + 1| is used (enter -1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.

distance(xopt)

Distance calculated from this Metric

Parameters:
xoptsequence

List of affine parameters as an 1D vector

### JointBundleMinDistanceMetric

Bundle-based Minimum Distance for joint optimization.

This cost function is used by the StreamlineLinearRegistration class when running halfway streamline linear registration for unbiased groupwise bundle registration and atlasing.

It computes the BMD distance after moving both static and moving bundles to a halfway space in between both.

Notes

In this metric both static and moving bundles are treated equally (i.e., there is no static reference bundle as both are intended to move). The naming convention is kept for consistency.

Methods

 setup(static, moving) distance(xopt)

An abstract class for the metric used for streamline registration.

If the two sets of streamlines match exactly then method distance of this object should be minimum.

Parameters:

Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus |num_threads + 1| is used (enter -1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.

distance(xopt)

Distance calculated from this Metric.

Parameters:
xoptsequence

List of affine parameters as an 1D vector. These affine parameters are used to derive the corresponding halfway transformation parameters for each bundle.

setup(static, moving)

Setup static and moving sets of streamlines.

Parameters:
staticstreamlines

Set of streamlines

movingstreamlines

Set of streamlines

Notes

Call this after the object is initiated and before distance. Num_threads is not used in this class.

### JointStreamlineRegistrationMap

class dipy.align.streamlinear.JointStreamlineRegistrationMap(xopt, fopt, matopt_history, funcs, iterations)

Bases: object

Methods

 transform(static, moving) Transform both static and moving bundles to the halfway space.
__init__(xopt, fopt, matopt_history, funcs, iterations)

A map holding the optimum affine matrices for halfway streamline linear registration and some other parameters of the optimization.

xopt is optimized by StreamlineLinearRegistration using the JointBundleMinDistanceMetric. In that case the mat argument of the optimize method needs to be np.eye(4) to avoid streamline centering.

This constructor derives and stores the transformations to move both static and moving bundles to the halfway space.

Parameters:
xoptarray

1d array with the parameters of the transformation.

foptfloat

Final value of the metric.

matopt_historyarray

All transformation matrices created during the optimization.

funcsint

Number of function evaluations of the optimizer.

iterationsint

Number of iterations of the optimizer.

transform(static, moving)

Transform both static and moving bundles to the halfway space.

All this does is apply self.matrix1 and self.matrix2 to the static and moving bundles, respectively.

Parameters:
staticstreamlines
movingstreamlines
Returns:
staticstreamlines
movingstreamlines

### Optimizer

class dipy.align.streamlinear.Optimizer(fun, x0, args=(), method='L-BFGS-B', jac=None, hess=None, hessp=None, bounds=None, constraints=(), tol=None, callback=None, options=None, evolution=False)

Bases: object

Attributes:
evolution
fopt
message
nfev
nit
xopt

Methods

 print_summary
__init__(fun, x0, args=(), method='L-BFGS-B', jac=None, hess=None, hessp=None, bounds=None, constraints=(), tol=None, callback=None, options=None, evolution=False)

A class for handling minimization of scalar function of one or more variables.

Parameters:
funcallable

Objective function.

x0ndarray

Initial guess.

argstuple, optional

Extra arguments passed to the objective function and its derivatives (Jacobian, Hessian).

methodstr, optional

Type of solver. Should be one of

• ‘Powell’

• ‘CG’

• ‘BFGS’

• ‘Newton-CG’

• ‘Anneal’

• ‘L-BFGS-B’

• ‘TNC’

• ‘COBYLA’

• ‘SLSQP’

• ‘dogleg’

• ‘trust-ncg’

jacbool or callable, optional

Jacobian of objective function. Only for CG, BFGS, Newton-CG, dogleg, trust-ncg. If jac is a Boolean and is True, fun is assumed to return the value of Jacobian along with the objective function. If False, the Jacobian will be estimated numerically. jac can also be a callable returning the Jacobian of the objective. In this case, it must accept the same arguments as fun.

hess, hesspcallable, optional

Hessian of objective function or Hessian of objective function times an arbitrary vector p. Only for Newton-CG, dogleg, trust-ncg. Only one of hessp or hess needs to be given. If hess is provided, then hessp will be ignored. If neither hess nor hessp is provided, then the hessian product will be approximated using finite differences on jac. hessp must compute the Hessian times an arbitrary vector.

boundssequence, optional

Bounds for variables (only for L-BFGS-B, TNC and SLSQP). (min, max) pairs for each element in x, defining the bounds on that parameter. Use None for one of min or max when there is no bound in that direction.

constraintsdict or sequence of dict, optional

Constraints definition (only for COBYLA and SLSQP). Each constraint is defined in a dictionary with fields:

typestr

Constraint type: ‘eq’ for equality, ‘ineq’ for inequality.

funcallable

The function defining the constraint.

jaccallable, optional

The Jacobian of fun (only for SLSQP).

argssequence, optional

Extra arguments to be passed to the function and Jacobian.

Equality constraint means that the constraint function result is to be zero whereas inequality means that it is to be non-negative. Note that COBYLA only supports inequality constraints.

tolfloat, optional

Tolerance for termination. For detailed control, use solver-specific options.

callbackcallable, optional

Called after each iteration, as callback(xk), where xk is the current parameter vector. Only available using Scipy >= 0.12.

optionsdict, optional

A dictionary of solver options. All methods accept the following generic options:

maxiterint

Maximum number of iterations to perform.

dispbool

Set to True to print convergence messages.

For method-specific options, see show_options(‘minimize’, method).

evolutionbool, optional

save history of x for each iteration. Only available using Scipy >= 0.12.

scipy.optimize.minimize
property evolution
property fopt
property message
property nfev
property nit
print_summary()
property xopt

### StreamlineDistanceMetric

Bases: object

Methods

 distance(xopt) calculate distance for current set of parameters.
 setup

An abstract class for the metric used for streamline registration.

If the two sets of streamlines match exactly then method distance of this object should be minimum.

Parameters:

Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus |num_threads + 1| is used (enter -1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.

abstract distance(xopt)

calculate distance for current set of parameters.

abstract setup(static, moving)

### StreamlineLinearRegistration

class dipy.align.streamlinear.StreamlineLinearRegistration(metric=None, x0='rigid', method='L-BFGS-B', bounds=None, verbose=False, options=None, evolution=False, num_threads=None)

Bases: object

Methods

 optimize(static, moving[, mat]) Find the minimum of the provided metric.
__init__(metric=None, x0='rigid', method='L-BFGS-B', bounds=None, verbose=False, options=None, evolution=False, num_threads=None)

Linear registration of 2 sets of streamlines [Garyfallidis15].

Parameters:
metricStreamlineDistanceMetric,

If None and fast is False then the BMD distance is used. If fast is True then a faster implementation of BMD is used. Otherwise, use the given distance metric.

x0array or int or str

Initial parametrization for the optimization.

If 1D array with:

a) 6 elements then only rigid registration is performed with the 3 first elements for translation and 3 for rotation. b) 7 elements also isotropic scaling is performed (similarity). c) 12 elements then translation, rotation (in degrees), scaling and shearing is performed (affine).

Here is an example of x0 with 12 elements: x0=np.array([0, 10, 0, 40, 0, 0, 2., 1.5, 1, 0.1, -0.5, 0])

This has translation (0, 10, 0), rotation (40, 0, 0) in degrees, scaling (2., 1.5, 1) and shearing (0.1, -0.5, 0).

If int:
1. 6

x0 = np.array([0, 0, 0, 0, 0, 0])

2. 7

x0 = np.array([0, 0, 0, 0, 0, 0, 1.])

3. 12

x0 = np.array([0, 0, 0, 0, 0, 0, 1., 1., 1, 0, 0, 0])

If str:
1. “rigid”

x0 = np.array([0, 0, 0, 0, 0, 0])

2. “similarity”

x0 = np.array([0, 0, 0, 0, 0, 0, 1.])

3. “affine”

x0 = np.array([0, 0, 0, 0, 0, 0, 1., 1., 1, 0, 0, 0])

methodstr,

‘L_BFGS_B’ or ‘Powell’ optimizers can be used. Default is ‘L_BFGS_B’.

boundslist of tuples or None,

If method == ‘L_BFGS_B’ then we can use bounded optimization. For example for the six parameters of rigid rotation we can set the bounds = [(-30, 30), (-30, 30), (-30, 30),

(-45, 45), (-45, 45), (-45, 45)]

That means that we have set the bounds for the three translations and three rotation axes (in degrees).

verbosebool, optional.

If True, if True then information about the optimization is shown. Default: False.

optionsNone or dict,

Extra options to be used with the selected method.

evolutionboolean

If True save the transformation for each iteration of the optimizer. Default is False. Supported only with Scipy >= 0.11.

Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus |num_threads + 1| is used (enter -1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.

References

Garyfallidis et al. “Robust and efficient linear registration of white-matter fascicles in the space of streamlines”, NeuroImage, 117, 124–140, 2015

[Garyfallidis14]

Garyfallidis et al., “Direct native-space fiber bundle alignment for group comparisons”, ISMRM, 2014.

[Garyfallidis17]

Garyfallidis et al. Recognition of white matter bundles using local and global streamline-based registration and clustering, Neuroimage, 2017.

optimize(static, moving, mat=None)

Find the minimum of the provided metric.

Parameters:
staticstreamlines

Reference or fixed set of streamlines.

movingstreamlines

Moving set of streamlines.

matarray

Transformation (4, 4) matrix to start the registration. mat is applied to moving. Default value None which means that initial transformation will be generated by shifting the centers of moving and static sets of streamlines to the origin.

Returns:
mapStreamlineRegistrationMap

### StreamlineRegistrationMap

class dipy.align.streamlinear.StreamlineRegistrationMap(matopt, xopt, fopt, matopt_history, funcs, iterations)

Bases: object

Methods

 transform(moving) Transform moving streamlines to the static.
__init__(matopt, xopt, fopt, matopt_history, funcs, iterations)

A map holding the optimum affine matrix and some other parameters of the optimization

Parameters:
matrixarray,

4x4 affine matrix which transforms the moving to the static streamlines

xoptarray,

1d array with the parameters of the transformation after centering

foptfloat,

final value of the metric

matrix_historyarray

All transformation matrices created during the optimization

funcsint,

Number of function evaluations of the optimizer

iterationsint

Number of iterations of the optimizer

transform(moving)

Transform moving streamlines to the static.

Parameters:
movingstreamlines
Returns:
movedstreamlines

Notes

All this does is apply self.matrix to the input streamlines.

### Streamlines

dipy.align.streamlinear.Streamlines

alias of ArraySequence

### combinations

class dipy.align.streamlinear.combinations(iterable, r)

Bases: object

Return successive r-length combinations of elements in the iterable.

combinations(range(4), 3) –> (0,1,2), (0,1,3), (0,2,3), (1,2,3)

__init__(*args, **kwargs)

### bundle_min_distance

dipy.align.streamlinear.bundle_min_distance(t, static, moving)

MDF-based pairwise distance optimization function (MIN).

We minimize the distance between moving streamlines as they align with the static streamlines.

Parameters:
tndarray

t is a vector of affine transformation parameters with size at least 6. If size is 6, t is interpreted as translation + rotation. If size is 7, t is interpreted as translation + rotation + isotropic scaling. If size is 12, t is interpreted as translation + rotation + scaling + shearing.

staticlist

Static streamlines

movinglist

Moving streamlines.

Returns:
cost: float

### bundle_min_distance_asymmetric_fast

dipy.align.streamlinear.bundle_min_distance_asymmetric_fast(t, static, moving, block_size)

MDF-based pairwise distance optimization function (MIN).

We minimize the distance between moving streamlines as they align with the static streamlines.

Parameters:
tarray

1D array. t is a vector of affine transformation parameters with size at least 6. If the size is 6, t is interpreted as translation + rotation. If the size is 7, t is interpreted as translation + rotation + isotropic scaling. If size is 12, t is interpreted as translation + rotation + scaling + shearing.

staticarray

N*M x 3 array. All the points of the static streamlines. With order of streamlines intact. Where N is the number of streamlines and M is the number of points per streamline.

movingarray

K*M x 3 array. All the points of the moving streamlines. With order of streamlines intact. Where K is the number of streamlines and M is the number of points per streamline.

block_sizeint

Number of points per streamline. All streamlines in static and moving should have the same number of points M.

Returns:
cost: float

### bundle_min_distance_fast

MDF-based pairwise distance optimization function (MIN).

We minimize the distance between moving streamlines as they align with the static streamlines.

Parameters:
tarray

1D array. t is a vector of affine transformation parameters with size at least 6. If the size is 6, t is interpreted as translation + rotation. If the size is 7, t is interpreted as translation + rotation + isotropic scaling. If size is 12, t is interpreted as translation + rotation + scaling + shearing.

staticarray

N*M x 3 array. All the points of the static streamlines. With order of streamlines intact. Where N is the number of streamlines and M is the number of points per streamline.

movingarray

K*M x 3 array. All the points of the moving streamlines. With order of streamlines intact. Where K is the number of streamlines and M is the number of points per streamline.

block_sizeint

Number of points per streamline. All streamlines in static and moving should have the same number of points M.

Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus |num_threads + 1| is used (enter -1 to use as many threads as possible). 0 raises an error.

Returns:
cost: float

Notes

This is a faster implementation of bundle_min_distance, which requires that all the points of each streamline are allocated into an ndarray (of shape N*M by 3, with N the number of points per streamline and M the number of streamlines). This can be done by calling dipy.tracking.streamlines.unlist_streamlines.

### bundle_sum_distance

MDF distance optimization function (SUM).

We minimize the distance between moving streamlines as they align with the static streamlines.

Parameters:
tndarray

t is a vector of affine transformation parameters with size at least 6. If the size is 6, t is interpreted as translation + rotation. If the size is 7, t is interpreted as translation + rotation + isotropic scaling. If size is 12, t is interpreted as translation + rotation + scaling + shearing.

staticlist

Static streamlines

movinglist

Moving streamlines. These will be transformed to align with the static streamlines

Number of threads. If -1 then all available threads will be used.

Returns:
cost: float

### center_streamlines

dipy.align.streamlinear.center_streamlines(streamlines)

Move streamlines to the origin

Parameters:
streamlineslist

List of 2D ndarrays of shape[-1]==3

Returns:
new_streamlineslist

List of 2D ndarrays of shape[-1]==3

inv_shiftndarray

Translation in x,y,z to go back in the initial position

### compose_matrix

dipy.align.streamlinear.compose_matrix(scale=None, shear=None, angles=None, translate=None, perspective=None)

Return 4x4 transformation matrix from sequence of transformations.

Code modified from the work of Christoph Gohlke link provided here http://www.lfd.uci.edu/~gohlke/code/transformations.py.html

This is the inverse of the decompose_matrix function.

Parameters:
scale(3,) array_like

Scaling factors.

sheararray_like

Shear factors for x-y, x-z, y-z axes.

anglesarray_like

Euler angles about static x, y, z axes.

translatearray_like

Translation vector along x, y, z axes.

perspectivearray_like

Perspective partition of matrix.

Returns:
matrix4x4 array

Examples

>>> import math
>>> import numpy as np
>>> import dipy.core.geometry as gm
>>> scale = np.random.random(3) - 0.5
>>> shear = np.random.random(3) - 0.5
>>> angles = (np.random.random(3) - 0.5) * (2*math.pi)
>>> trans = np.random.random(3) - 0.5
>>> persp = np.random.random(4) - 0.5
>>> M0 = gm.compose_matrix(scale, shear, angles, trans, persp)


### compose_matrix44

dipy.align.streamlinear.compose_matrix44(t, dtype=<class 'numpy.float64'>)

Compose a 4x4 transformation matrix.

Parameters:
tndarray

This is a 1D vector of affine transformation parameters with size at least 3. If the size is 3, t is interpreted as translation. If the size is 6, t is interpreted as translation + rotation. If the size is 7, t is interpreted as translation + rotation + isotropic scaling. If the size is 9, t is interpreted as translation + rotation + anisotropic scaling. If size is 12, t is interpreted as translation + rotation + scaling + shearing.

Returns:
Tndarray

Homogeneous transformation matrix of size 4x4.

### compose_transformations

dipy.align.streamlinear.compose_transformations(*mats)

Compose multiple 4x4 affine transformations in one 4x4 matrix

Parameters:
mat1array, (4, 4)
mat2array, (4, 4)
matNarray, (4, 4)
Returns:
matN x … x mat2 x mat1array, (4, 4)

### decompose_matrix

dipy.align.streamlinear.decompose_matrix(matrix)

Return sequence of transformations from transformation matrix.

Code modified from the excellent work of Christoph Gohlke link provided here: http://www.lfd.uci.edu/~gohlke/code/transformations.py.html

Parameters:
matrixarray_like

Non-degenerate homogeneous transformation matrix

Returns:
scale(3,) ndarray

Three scaling factors.

shear(3,) ndarray

Shear factors for x-y, x-z, y-z axes.

angles(3,) ndarray

Euler angles about static x, y, z axes.

translate(3,) ndarray

Translation vector along x, y, z axes.

perspectivendarray

Perspective partition of matrix.

Raises:
ValueError

If matrix is of wrong type or degenerate.

Examples

>>> import numpy as np
>>> T0=np.diag([2,1,1,1])
>>> scale, shear, angles, trans, persp = decompose_matrix(T0)


### decompose_matrix44

dipy.align.streamlinear.decompose_matrix44(mat, size=12)

Given a 4x4 homogeneous matrix return the parameter vector.

Parameters:
matarray

Homogeneous 4x4 transformation matrix

sizeint

Size of the output vector. 3, for translation, 6 for rigid, 7 for similarity, 9 for scaling and 12 for affine. Default is 12.

Returns:
tndarray

One dimensional ndarray of 3, 6, 7, 9 or 12 affine parameters.

### distance_matrix_mdf

dipy.align.streamlinear.distance_matrix_mdf()

Minimum direct flipped distance matrix between two streamline sets

All streamlines need to have the same number of points

Parameters:
streamlines_asequence

of streamlines as arrays, [(N, 3) .. (N, 3)]

streamlines_bsequence

of streamlines as arrays, [(N, 3) .. (N, 3)]

Returns:
DMarray, shape (len(streamlines_a), len(streamlines_b))

distance matrix

### get_unique_pairs

dipy.align.streamlinear.get_unique_pairs(n_bundle, pairs=None)

Make unique pairs from n_bundle bundles.

The function allows to input a previous pairs asignment so that the new pairs are different.

Parameters:
n_bundleint

Number of bundles to be matched in pairs.

pairsarray, optional

array containing the indexes of previous pairs.

### groupwise_slr

dipy.align.streamlinear.groupwise_slr(bundles, x0='affine', tol=0, max_iter=20, qbx_thr=[4], nb_pts=20, select_random=10000, verbose=False, rng=None)

Function to perform unbiased groupwise bundle registration.

All bundles are moved to the same space by iteratively applying halfway streamline linear registration in pairs. With each iteration, bundles get closer to each other until the procedure converges and there is no more improvement.

Parameters:
bundleslist

List with streamlines of the bundles to be registered.

x0str, optional

rigid, similarity or affine transformation model. Default: affine.

tolfloat, optional

Tolerance value to be used to assume convergence. Default: 0.

max_iterint, optional

Maximum number of iterations. Depending on the number of bundles to be registered this may need to be larger. Default: 20.

qbx_thrvariable int, optional

Thresholds for Quickbundles used for clustering streamlines and reduce computational time. If None, no clustering is performed. Higher values cluster streamlines into a smaller number of centroids. Default: [4].

nb_ptsint, optional

Number of points for discretizing each streamline. Default: 20.

select_randomint, optional

Maximum number of streamlines for each bundle. If None, all the streamlines are used. Deafult: 10000.

verbosebool, optional

If True, logs information. Default: False.

rngRandomState

If None, creates RandomState in function. Default: None.

References

[Garyfallidis15]

Garyfallidis et al. “Robust and efficient linear

registration of white-matter fascicles in the space of streamlines”, NeuroImage, 117, 124–140, 2015 .. [Rd50d670f3909-Garyfallidis14] Garyfallidis et al., “Direct native-space fiber

bundle alignment for group comparisons”, ISMRM, 2014.

[Garyfallidis17]

Garyfallidis et al. Recognition of white matter

bundles using local and global streamline-based registration and clustering, Neuroimage, 2017.

### length

dipy.align.streamlinear.length()

Euclidean length of streamlines

Length is in mm only if streamlines are expressed in world coordinates.

Parameters:
streamlinesndarray or a list or dipy.tracking.Streamlines

If ndarray, must have shape (N,3) where N is the number of points of the streamline. If list, each item must be ndarray shape (Ni,3) where Ni is the number of points of streamline i. If dipy.tracking.Streamlines, its common_shape must be 3.

Returns:
lengthsscalar or ndarray shape (N,)

If there is only one streamline, a scalar representing the length of the streamline. If there are several streamlines, ndarray containing the length of every streamline.

Examples

>>> from dipy.tracking.streamline import length
>>> import numpy as np
>>> streamline = np.array([[1, 1, 1], [2, 3, 4], [0, 0, 0]])
>>> expected_length = np.sqrt([1+2**2+3**2, 2**2+3**2+4**2]).sum()
>>> length(streamline) == expected_length
True
>>> streamlines = [streamline, np.vstack([streamline, streamline[::-1]])]
>>> expected_lengths = [expected_length, 2*expected_length]
>>> lengths = [length(streamlines[0]), length(streamlines[1])]
>>> np.allclose(lengths, expected_lengths)
True
>>> length([])
0.0
>>> length(np.array([[1, 2, 3]]))
0.0


### progressive_slr

dipy.align.streamlinear.progressive_slr(static, moving, metric, x0, bounds, method='L-BFGS-B', verbose=False, num_threads=None)

Progressive SLR.

This is an utility function that allows for example to do affine registration using Streamline-based Linear Registration (SLR) [Garyfallidis15] by starting with translation first, then rigid, then similarity, scaling and finally affine.

Similarly, if for example, you want to perform rigid then you start with translation first. This progressive strategy can helps with finding the optimal parameters of the final transformation.

Parameters:
staticStreamlines
movingStreamlines
metricStreamlineDistanceMetric
x0string

Could be any of ‘translation’, ‘rigid’, ‘similarity’, ‘scaling’, ‘affine’

boundsarray

Boundaries of registration parameters. See variable DEFAULT_BOUNDS for example.

methodstring

L_BFGS_B’ or ‘Powell’ optimizers can be used. Default is ‘L_BFGS_B’.

verbose<