dipy logo

Site Navigation

NIPY Community

Previous topic

dipy.reconst.dsi

dipy.reconst.dti

class dipy.reconst.dti.Tensor(data, b_values, b_vectors, mask=True, thresh=None, fit_method='WLS', verbose=False, *args, **kargs)

For backwards compatibility, we continue to support this form of the Tensor fitting.

Methods

ad() Radial diffusitivity (RD) calculated from cached eigenvalues.
fa([fill_value, nonans]) Fractional anisotropy (FA) calculated from cached eigenvalues.
ind() Quantizes eigenvectors with maximum eigenvalues on an evenly distributed sphere so that the can be used for tractography.
lower_triangular([b0])
md() Mean diffusitivity (MD) calculated from cached eigenvalues.
odf(sphere)
rd() Radial diffusitivity (RD) calculated from cached eigenvalues.
trace() Trace of the tensor calculated from cached eigenvalues.
D

Calculates the 3x3 diffusion tensor for each voxel

fa(fill_value=0, nonans=True)

Fractional anisotropy (FA) calculated from cached eigenvalues.

Parameters :

fill_value : float

value of fa where self.mask == True.

nonans : Bool

When True, fa is 0 when all eigenvalues are 0, otherwise fa is nan

Returns :

fa : array (V, 1)

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

Notes

FA is calculated with 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} }

ind()

Quantizes eigenvectors with maximum eigenvalues on an evenly distributed sphere so that the can be used for tractography.

Returns :

IN : array, shape(x,y,z) integer indices for the points of the

evenly distributed sphere representing tensor eigenvectors of :

maximum eigenvalue :

md()

Mean diffusitivity (MD) calculated from cached eigenvalues.

Returns :

md : array (V, 1)

Calculated MD.

Notes

MD is calculated with the following equation:

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

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.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.

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.

Notes

it is computed from the clipped FA between 0 and 1 using the following formula

rgb = abs(max(eigen_vector)) \times fa

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, bval, dtype=None)

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

Parameters :

gtab : array with shape (3,g)

Diffusion gradient table found in DICOM header as a numpy array.

bval : array with shape (g,)

Diffusion weighting factor b for each vector in gtab.

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.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.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. Also, called Apparent diffusion coefficient (ADC)

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.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. Use design_matrix to build a valid design matrix from bvalues and a gradient table.

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

This function is offered mainly as a quick comparison to WLS.

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

References

[1]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.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.

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.tensor_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.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:

MD = \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 [R6].

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

[R6](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.