dipy logo

Site Navigation

NIPY Community

Previous topic

dipy.reconst.dsi

dipy.reconst.dti

Classes and functions for fitting tensors

class dipy.reconst.dti.TensorModel(gtab, fit_method='WLS', *args, **kwargs)

Diffusion Tensor

Methods

fit(data[, mask]) Fit method of the DTI model class
fit(data, mask=None)

Fit method of the DTI model class

Parameters:

data : array

The measured signal from one voxel.

mask : array

A boolean array used to mark the coordinates in the data that should be analyzed that has the shape data.shape[-1]

dipy.reconst.dti.apparent_diffusion_coef(q_form, sphere)

Calculate the apparent diffusion coefficient (ADC) in each direction of a sphere.

Parameters:

q_form : ndarray

The quadratic form of a tensor, or an array with quadratic forms of tensors. Should be of shape (..., 3, 3)

sphere : a Sphere class instance

The ADC will be calculated for each of the vertices in the sphere

Notes

The calculation of ADC, relies on the following relationship:

ADC = \vec{b} Q \vec{b}^T

Where Q is the quadratic form of the tensor.

dipy.reconst.dti.axial_diffusivity(evals, axis=-1)

Axial Diffusivity (AD) of a diffusion tensor. Also called parallel diffusivity.

Parameters:

evals : array-like

Eigenvalues of a diffusion tensor, must be sorted in descending order along axis.

axis : int

Axis of evals which contains 3 eigenvalues.

Returns:

ad : array

Calculated AD.

Notes

AD is calculated with the following equation:

AD = \lambda_1

dipy.reconst.dti.color_fa(fa, evecs)

Color fractional anisotropy of diffusion tensor

Parameters:

fa : array-like

Array of the fractional anisotropy (can be 1D, 2D or 3D)

evecs : array-like

eigen vectors from the tensor model

Returns:

rgb : Array with 3 channels for each color as the last dimension.

Colormap of the FA with red for the x value, y for the green value and z for the blue value.

dipy.reconst.dti.decompose_tensor(tensor, min_diffusivity=0)

Returns eigenvalues and eigenvectors given a diffusion tensor

Computes tensor eigen decomposition to calculate eigenvalues and eigenvectors (Basser et al., 1994a).

Parameters:

tensor : array (3, 3)

Hermitian matrix representing a diffusion tensor.

min_diffusivity : float

Because negative eigenvalues are not physical and small eigenvalues, much smaller than the diffusion weighting, cause quite a lot of noise in metrics such as fa, diffusivity values smaller than min_diffusivity are replaced with min_diffusivity.

Returns:

eigvals : array (3,)

Eigenvalues from eigen decomposition of the tensor. Negative eigenvalues are replaced by zero. Sorted from largest to smallest.

eigvecs : array (3, 3)

Associated eigenvectors from eigen decomposition of the tensor. Eigenvectors are columnar (e.g. eigvecs[:,j] is associated with eigvals[j])

dipy.reconst.dti.design_matrix(gtab, dtype=None)

Constructs design matrix for DTI weighted least squares or least squares fitting. (Basser et al., 1994a)

Parameters:

gtab : A GradientTable class instance

dtype : string

Parameter to control the dtype of returned designed matrix

Returns:

design_matrix : array (g,7)

Design matrix or B matrix assuming Gaussian distributed tensor model design_matrix[j, :] = (Bxx, Byy, Bzz, Bxy, Bxz, Byz, dummy)

dipy.reconst.dti.determinant(q_form)

The determinant of a tensor, given in quadratic form

Parameters:

q_form : ndarray

The quadratic form of a tensor, or an array with quadratic forms of tensors. Should be of shape (x, y, z, 3, 3) or (n, 3, 3) or (3, 3).

Returns:

det : array

The determinant of the tensor in each spatial coordinate

dipy.reconst.dti.deviatoric(q_form)

Calculate the deviatoric (anisotropic) part of the tensor [R7].

Parameters:

q_form : ndarray

The quadratic form of a tensor, or an array with quadratic forms of tensors. Should be of shape (x,y,z,3,3) or (n, 3, 3) or (3,3).

Returns:

A_squiggle : ndarray

The deviatoric part of the tensor in each spatial coordinate.

Notes

The deviatoric part of the tensor is defined as (equations 3-5 in [R7]):

\widetilde{A} = A - \bar{A}

Where A is the tensor quadratic form and \bar{A} is the anisotropic part of the tensor.

[R7](1, 2) Daniel B. Ennis and G. Kindlmann, “Orthogonal Tensor Invariants and the Analysis of Diffusion Tensor Magnetic Resonance Images”, Magnetic Resonance in Medicine, vol. 55, no. 1, pp. 136-146, 2006.
dipy.reconst.dti.eig_from_lo_tri(data)

Calculates parameters for creating a Tensor instance

Calculates tensor parameters from the six unique tensor elements. This function can be passed to the Tensor class as a fit_method for creating a Tensor instance from tensors stored in a nifti file.

Parameters:

data : array_like (..., 6)

diffusion tensors elements stored in lower triangular order

Returns:

dti_params :

Eigen values and vectors, used by the Tensor class to create an instance

dipy.reconst.dti.fractional_anisotropy(evals, axis=-1)

Fractional anisotropy (FA) of a diffusion tensor.

Parameters:

evals : array-like

Eigenvalues of a diffusion tensor.

axis : int

Axis of evals which contains 3 eigenvalues.

Returns:

fa : array

Calculated FA. Range is 0 <= FA <= 1.

Notes

FA is calculated using the following equation:

FA = \sqrt{\frac{1}{2}\frac{(\lambda_1-\lambda_2)^2+(\lambda_1- \lambda_3)^2+(\lambda_2-\lambda_3)^2}{\lambda_1^2+ \lambda_2^2+\lambda_3^2}}

dipy.reconst.dti.from_lower_triangular(D)

Returns a tensor given the six unique tensor elements

Given the six unique tensor elments (in the order: Dxx, Dxy, Dyy, Dxz, Dyz, Dzz) returns a 3 by 3 tensor. All elements after the sixth are ignored.

Parameters:

D : array_like, (..., >6)

Unique elements of the tensors

Returns:

tensor : ndarray (..., 3, 3)

3 by 3 tensors

dipy.reconst.dti.isotropic(q_form)

Calculate the isotropic part of the tensor [R8].

Parameters:

q_form : ndarray

The quadratic form of a tensor, or an array with quadratic forms of tensors. Should be of shape (x,y,z,3,3) or (n, 3, 3) or (3,3).

Returns:

A_hat: ndarray :

The isotropic part of the tensor in each spatial coordinate

Notes

The isotropic part of a tensor is defined as (equations 3-5 of [R8]):

\bar{A} = \frac{1}{2} tr(A) I

[R8](1, 2) Daniel B. Ennis and G. Kindlmann, “Orthogonal Tensor Invariants and the Analysis of Diffusion Tensor Magnetic Resonance Images”, Magnetic Resonance in Medicine, vol. 55, no. 1, pp. 136-146, 2006.
dipy.reconst.dti.linearity(evals, axis=-1)

The linearity of the tensor [1]

Parameters:

evals : array-like

Eigenvalues of a diffusion tensor.

axis : int

Axis of evals which contains 3 eigenvalues.

Returns:

linearity : array

Calculated linearity of the diffusion tensor.

Notes

[1] Westin C.-F., Peled S., Gubjartsson H., Kikinis R., Jolesz F.,
“Geometrical diffusion measures for MRI from tensor basis analysis” in Proc. 5th Annual ISMRM, 1997.
dipy.reconst.dti.lower_triangular(tensor, b0=None)

Returns the six lower triangular values of the tensor and a dummy variable if b0 is not None

Parameters:

tensor : array_like (..., 3, 3)

a collection of 3, 3 diffusion tensors

b0 : float

if b0 is not none log(b0) is returned as the dummy variable

Returns:

D : ndarray

If b0 is none, then the shape will be (..., 6) otherwise (..., 7)

dipy.reconst.dti.mean_diffusivity(evals, axis=-1)

Mean Diffusivity (MD) of a diffusion tensor.

Parameters:

evals : array-like

Eigenvalues of a diffusion tensor.

axis : int

Axis of evals which contains 3 eigenvalues.

Returns:

md : array

Calculated MD.

Notes

MD is calculated with the following equation:

MD = \frac{\lambda_1 + \lambda_2 + \lambda_3}{3}

dipy.reconst.dti.mode(q_form)

Mode (MO) of a diffusion tensor [R9].

Parameters:

q_form : ndarray

The quadratic form of a tensor, or an array with quadratic forms of tensors. Should be of shape (x, y, z, 3, 3) or (n, 3, 3) or (3, 3).

Returns:

mode : array

Calculated tensor mode in each spatial coordinate.

Notes

Mode ranges between -1 (linear anisotropy) and +1 (planar anisotropy) with 0 representing orthotropy. Mode is calculated with the following equation (equation 9 in [R9]):

Mode = 3*\sqrt{6}*det(\widetilde{A}/norm(\widetilde{A}))

Where \widetilde{A} is the deviatoric part of the tensor quadratic form.

References

[R9](1, 2, 3) Daniel B. Ennis and G. Kindlmann, “Orthogonal Tensor Invariants and the Analysis of Diffusion Tensor Magnetic Resonance Images”, Magnetic Resonance in Medicine, vol. 55, no. 1, pp. 136-146, 2006.
dipy.reconst.dti.nlls_fit_tensor(design_matrix, data, min_signal=1, weighting=None, sigma=None, jac=True)

Fit the tensor params using non-linear least-squares.

Parameters:

design_matrix : array (g, 7)

Design matrix holding the covariants used to solve for the regression coefficients.

data : array ([X, Y, Z, ...], g)

Data or response variables holding the data. Note that the last dimension should contain the data. It makes no copies of data.

min_signal : float, optional

All values below min_signal are repalced with min_signal. This is done in order to avaid taking log(0) durring the tensor fitting. Default = 1

weighting: str :

the weighting scheme to use in considering the squared-error. Default behavior is to use uniform weighting. Other options: ‘sigma’ ‘gmm’

sigma: float :

If the ‘sigma’ weighting scheme is used, a value of sigma needs to be provided here. According to [Chang2005], a good value to use is 1.5267 * std(background_noise), where background_noise is estimated from some part of the image known to contain no signal (only noise).

jac : bool

Use the Jacobian? Default: True

Returns:

nlls_params: the eigen-values and eigen-vectors of the tensor in each voxel. :

dipy.reconst.dti.norm(q_form)

Calculate the Frobenius norm of a tensor quadratic form

Parameters:

q_form: ndarray :

The quadratic form of a tensor, or an array with quadratic forms of tensors. Should be of shape (x,y,z,3,3) or (n, 3, 3) or (3,3).

Returns:

norm : ndarray

The Frobenius norm of the 3,3 tensor q_form in each spatial coordinate.

See also

np.linalg.norm

Notes

The Frobenius norm is defined as:

Math:||A||_F = [sum_{i,j} abs(a_{i,j})^2]^{1/2}
dipy.reconst.dti.ols_fit_tensor(design_matrix, data, min_signal=1)

Computes ordinary least squares (OLS) fit to calculate self-diffusion tensor using a linear regression model [1].

Parameters:

design_matrix : array (g, 7)

Design matrix holding the covariants used to solve for the regression coefficients.

data : array ([X, Y, Z, ...], g)

Data or response variables holding the data. Note that the last dimension should contain the data. It makes no copies of data.

min_signal : default = 1

All values below min_signal are repalced with min_signal. This is done in order to avaid taking log(0) durring the tensor fitting.

Returns:

eigvals : array (..., 3)

Eigenvalues from eigen decomposition of the tensor.

eigvecs : array (..., 3, 3)

Associated eigenvectors from eigen decomposition of the tensor. Eigenvectors are columnar (e.g. eigvecs[:,j] is associated with eigvals[j])

See also

WLS_fit_tensor, decompose_tensor, design_matrix

Notes

y = \mathrm{data} \\ X = \mathrm{design matrix} \\ \hat{\beta}_\mathrm{OLS} = (X^T X)^{-1} X^T y

References

[1](1, 2, 3, 4) Chung, SW., Lu, Y., Henry, R.G., 2006. Comparison of bootstrap approaches for estimation of uncertainties of DTI parameters. NeuroImage 33, 531-541.
dipy.reconst.dti.planarity(evals, axis=-1)

The planarity of the tensor [1]

Parameters:

evals : array-like

Eigenvalues of a diffusion tensor.

axis : int

Axis of evals which contains 3 eigenvalues.

Returns:

linearity : array

Calculated linearity of the diffusion tensor.

Notes

[1] Westin C.-F., Peled S., Gubjartsson H., Kikinis R., Jolesz F.,
“Geometrical diffusion measures for MRI from tensor basis analysis” in Proc. 5th Annual ISMRM, 1997.
dipy.reconst.dti.quantize_evecs(evecs, odf_vertices=None)

Find the closest orientation of an evenly distributed sphere

Parameters:

evecs : ndarray

odf_vertices : None or ndarray

If None, then set vertices from symmetric362 sphere. Otherwise use passed ndarray as vertices

Returns:

IN : ndarray

dipy.reconst.dti.radial_diffusivity(evals, axis=-1)

Radial Diffusivity (RD) of a diffusion tensor. Also called perpendicular diffusivity.

Parameters:

evals : array-like

Eigenvalues of a diffusion tensor, must be sorted in descending order along axis.

axis : int

Axis of evals which contains 3 eigenvalues.

Returns:

rd : array

Calculated RD.

Notes

RD is calculated with the following equation:

RD = \frac{\lambda_2 + \lambda_3}{2}

dipy.reconst.dti.restore_fit_tensor(design_matrix, data, min_signal=1.0, sigma=None, jac=True)

Use the RESTORE algorithm [Chang2005] to calculate a robust tensor fit

Parameters:

design_matrix : array of shape (g, 7)

Design matrix holding the covariants used to solve for the regression coefficients.

data : array of shape ([X, Y, Z, n_directions], g)

Data or response variables holding the data. Note that the last dimension should contain the data. It makes no copies of data.

min_signal : float, optional

All values below min_signal are repalced with min_signal. This is done in order to avaid taking log(0) durring the tensor fitting. Default = 1

sigma : float

An estimate of the variance. [Chang2005] recommend to use 1.5267 * std(background_noise), where background_noise is estimated from some part of the image known to contain no signal (only noise).

jac : bool, optional

Whether to use the Jacobian of the tensor to speed the non-linear optimization procedure used to fit the tensor paramters (see also nlls_fit_tensor()). Default: True

Returns:

restore_params : an estimate of the tensor parameters in each voxel.

dipy.reconst.dti.sphericity(evals, axis=-1)

The sphericity of the tensor [1]

Parameters:

evals : array-like

Eigenvalues of a diffusion tensor.

axis : int

Axis of evals which contains 3 eigenvalues.

Returns:

sphericity : array

Calculated sphericity of the diffusion tensor.

Notes

[1] Westin C.-F., Peled S., Gubjartsson H., Kikinis R., Jolesz F.,
“Geometrical diffusion measures for MRI from tensor basis analysis” in Proc. 5th Annual ISMRM, 1997.
dipy.reconst.dti.trace(evals, axis=-1)

Trace of a diffusion tensor.

Parameters:

evals : array-like

Eigenvalues of a diffusion tensor.

axis : int

Axis of evals which contains 3 eigenvalues.

Returns:

trace : array

Calculated trace of the diffusion tensor.

Notes

Trace is calculated with the following equation:

Trace = \lambda_1 + \lambda_2 + \lambda_3

dipy.reconst.dti.wls_fit_tensor(design_matrix, data, min_signal=1)

Computes weighted least squares (WLS) fit to calculate self-diffusion tensor using a linear regression model [R10].

Parameters:

design_matrix : array (g, 7)

Design matrix holding the covariants used to solve for the regression coefficients.

data : array ([X, Y, Z, ...], g)

Data or response variables holding the data. Note that the last dimension should contain the data. It makes no copies of data.

min_signal : default = 1

All values below min_signal are repalced with min_signal. This is done in order to avaid taking log(0) durring the tensor fitting.

Returns:

eigvals : array (..., 3)

Eigenvalues from eigen decomposition of the tensor.

eigvecs : array (..., 3, 3)

Associated eigenvectors from eigen decomposition of the tensor. Eigenvectors are columnar (e.g. eigvecs[:,j] is associated with eigvals[j])

See also

decompose_tensor

Notes

In Chung, et al. 2006, the regression of the WLS fit needed an unbiased preliminary estimate of the weights and therefore the ordinary least squares (OLS) estimates were used. A “two pass” method was implemented:

  1. calculate OLS estimates of the data
  2. apply the OLS estimates as weights to the WLS fit of the data

This ensured heteroscadasticity could be properly modeled for various types of bootstrap resampling (namely residual bootstrap).

y = \mathrm{data} \\ X = \mathrm{design matrix} \\ \hat{\beta}_\mathrm{WLS} = \mathrm{desired regression coefficients (e.g. tensor)}\\ \\ \hat{\beta}_\mathrm{WLS} = (X^T W X)^{-1} X^T W y \\ \\ W = \mathrm{diag}((X \hat{\beta}_\mathrm{OLS})^2), \mathrm{where} \hat{\beta}_\mathrm{OLS} = (X^T X)^{-1} X^T y

References

[R10](1, 2) Chung, SW., Lu, Y., Henry, R.G., 2006. Comparison of bootstrap approaches for estimation of uncertainties of DTI parameters. NeuroImage 33, 531-541.