<no title>

class designer.fitting.dwipy.DWI(imPath, nthreads=- 1)¶

Bases: object

The DWI object handles tensor estimation and parameter extraction of dwiffusion weighted images

Variables

hdr (class) – Nibabel class object of input nifti file
img (ndarray) – 3D or 4D input image array
grad (ndarray) – [N x 4] gradient and bvalue scheme, where the first three columns are the X, Y, and Z gradient vectors respectively, and the fourth column is B-values
mask (ndarray(dtype=bool)) – 3D array of brain mask
MaskStatus (bool) – True if brain_mask.nii is present, False otherwise
workers (int) – Number of CPU workers to use in processing

akccorrect(akc_out, window=3, connectivity='face')¶

Applies AKC outlier map to DT to replace outliers with a moving median. Run this only after tensor fitting and akc outlier detection.

Parameters

akc_out (ndarray(dtype=bool)) – 3D map containing outliers from DWI.akcoutliers
window (int, optional) – Width of square matrix filter (Default: 5)
connectivity (str, {‘face’, ‘all’}, optional) – Specifies what kind of connected-component connectivity to use for median determination

Examples

dwi.akccorrect(akc_out), where dwi is the DWI class object

akcoutliers(iter=10)¶

Uses 100,000 direction in chunks of 10 to iteratively find outliers. Returns a mask of locations where said violations occur. Multiprocessing is disabled because this is a memory-intensive task. To be run only after tensor fitting.

Parameters: iter (int, optional) – number of iterations to perform out of 10. Reduce this number if your computer does not have sufficient RAM. (Default: 10)
Returns: akc_out – 3D map containing outliers where AKC falls fails the inequality test -2 < AKC < 10
Return type: ndarray(dtype=bool)

Examples

akc_out = dwi.akoutliers(), where dwi is the DWI class object

createConstraints(constraints=[0, 1, 0])¶

Generates constraint array for constrained minimization quadratic programming

Parameters: constraints (array_like(dtype=int)) – [1 X 3] logical vector indicating which constraints out of three to enable (Default: [0, 1, 0]) C1 is Dapp > 0 C1 is Kapp > 0 C3 is Kapp < 3/(b*Dapp)
Returns: C – Array containing constraints to consider during minimization, C is shaped [number of constraints enforced * number of directions, 22]
Return type: ndarray(dtype=float)

Examples

C = dwi.createConstraints([0, 1, 0])

createTensorOrder(order=None)¶

Creates tensor order array and indices

Parameters

order (2 or 4 (int or None)) – Tensor order number, 2 for diffusion and 4 for kurtosis. Default: None; auto-detect

Returns

cnt (ndarray(dtype=int)) – Count of number of times a tensor appears
ind (ndarray(dtype=int)) – Indices of count

Examples

(cnt, ind) = dwi.createTensorOrder(order)

Notes

The tensors for this pipeline are based on NYU’s designer layout as depicted in the table below. This will soon be depreciated and updated with MRTRIX3’s layout.

~~~~~~D~~~~~~
|    D11
|    D12
|    D13
|    D22
|    D23
|    D33
~~~~~~K~~~~~~
|   W1111
|   W1112
|   W1113
|   W1122
|   W1123
|   W1133
|   W1222
|   W1223
|   W1233
|   W1333
|   W2222
|   W2223
|   W2233
|   W2333
|   W3333

diffusionCoeff(dt, dir)¶

Computes apparent diffusion coefficient (ADC)

Parameters

dt (ndarray(dtype=float)) – [21 x nvoxel] array containing diffusion tensor
dir (ndarray(dtype=float)) – [n x 3] array containing gradient directions

Returns

adc – Array containing apparent diffusion coefficient

Return type

ndarray(dtype=float)

Examples

adc = dwi.diffusionCoeff(dt, dir)

dkiTensorParams(v1, dt)¶

Uses average directional statistics to approximate axial kurtosis(AK) and radial kurtosis (RK)

Parameters

v1 (ndarray(dtype=float)) – Array of first eigenvectors from DWI.dtiTensorParams()
dt (ndarray(dtype=float)) – Array of diffusion tensor

Returns

rk (ndarray(dtype=float)) – Radial Kurtosis
ak (ndarray(dtype=float)) – Axial Kurtosis
kfa (ndarray(dtype=float)) – Kurtosis Fractional Anisotropy
mkt (ndarray(dtype=float)) – Mean Kurtosis Tensor

Examples

(rk, ak) = dwi.dkiTensorParams(v1, dt)

dtiTensorParams(nn)¶

Computes sorted DTI tensor eigenvalues and eigenvectors

Parameters

DT (ndarray(dtype=float)) – Diffusion tensor array

Returns

values (ndarray(dtype=float)) – Array of sorted eigenvalues
vectors (ndarray(dtype=float)) – Array pf sorted eigenvectors

Examples

(values, vectors) = dwi.dtiTensorParams(DT)

extractDKI()¶

Extract all DKI parameters from DT tensor. Warning, this can only be run after tensor fitting dwi.fit()

Returns

mk (ndarray(dtype=float)) – Mean Diffusivity
rk (ndarray(dtype=float)) – Radial Diffusivity
ak (ndarray(dtype=float)) – Axial Diffusivity
kfa (ndarray(dtype=float)) – Kurtosis Fractional Anisotropy
mkt (ndarray(dtype=float)) – Mean Kurtosis Tensor
trace (ndarray(dtype=float)) – Sum of first eigenvalues

Examples

(mk, rk, ak, fe, trace) = dwi.extractDTI(), where dwi is the DWI class object

extractDTI()¶

Extract all DTI parameters from DT tensor. Warning, this can only be run after tensor fitting dwi.fit()

Returns

md (ndarray(dtype=float)) – Mean Diffusivity
rd (ndarray(dtype=float)) – Radial Diffusivity
ad (ndarray(dtype=float)) – Axial Diffusivity
fa (ndarray(dtype=float)) – Fractional Anisotropy
fe (ndarray(dtype=float)) – First Eigenvectors
trace (ndarray(dtype=float)) – Sum of first eigenvalues

Examples

(md, rd, ad, fa) = dwi.extractDTI(), where dwi is the DWI class object

extractWMTI()¶

Returns white matter tract integrity (WMTI) parameters. Warning: this can only be run after fitting and DWI.extractDTI().

Returns

awf (ndarray(dtype=float)) – Axonal Water Fraction
eas_ad (ndarray(dtype=float)) – Extra-axonal space Axial Diffusivity
eas_rd (ndarray(dtype=float)) – Extra-axonal Space Radial Diffusivity
eas_md (ndarray(dtype=float)) – Extra-axonal Space Mean Diffusivity
eas_tort (ndarray(dtype=float)) – Extra-axonal Space Tortuosity
ias_ad (ndarray(dtype=float)) – Intra-axonal Space Axial Diffusivity
ias_rd (ndarray(dtype=float)) – Intra-axonal Space Radial Diffusivity
ias_da (ndarray(dtype=float)) – Intra-axonal Space Intrinsic Diffusivity
ias_tort (ndarray(dtype=float)) – Intra-axonal Space Tortuosity

fibonacciSphere(samples=1, randomize=True)¶

Returns evenly spaced points on a sphere

Parameters

samples (int) – Number of points to compute from sphere, must be a positive and real integer (Default: 1)
randomize (bool) – True if sampling is randomized; False otherwise (Default: True)

Returns

points – [3 x samples] array containing evenly spaced points from a sphere

Return type

ndarray(dtype=float)

Examples

dirs = dwi.fibonacciSphere(256, True)

findViols(c=[0, 1, 0])¶

Returns a 3D violation map of voxels that violate constraints.

Parameters

img (ndarray(dtype=float)) – 3D metric array such as mk or fa
c (array_like(dtype=int)) – [3 x 1] vector that toggles which constraints to check c[0]: Check D < 0 constraint c[1]: Check K < 0 constraint (Default) c[2]: Check K > 3/(b*D) constraint

Returns

map – 3D array containing locations of voxels that incur directional violations. Voxels with values contain violaions and voxel values represent proportion of directional violations.

Return type

ndarray(dtype=bool)

Examples

map = findViols(img, [0 1 0]

findVoxelViol(adcVox, akcVox, maxB, c)¶

Returns the proportions of violations occurring at a voxel.

Parameters

img (ndarray(dtype=float)) – 3D metric array such as mk or fa
c (array_like(dtype=int)) – [3 x 1] vector that toggles which constraints to check c[0]: Check D < 0 constraint c[1]: Check K < 0 constraint (Default) c[2]: Check K > 3/(b*D) constraint

Returns

n – percentaghhe ranging from 0 to 1 that indicates proportion of violations occuring at voxel.

Return type

ndarray(dtype=float)

Examples

map = findViols(voxel, [0 1 0]

fit(constraints=None, reject=None)¶

Returns fitted diffusion or kurtosis tensor

Parameters

constraints (array_like(dtype=int)) – [1 x 3] vector that specifies which constraints to use (Default: None)
reject (ndarray(dtype=bool)) – 4D array containing information on voxels to exclude from DT estimation (Default: None)

Examples

dwi.fit() dwi.fit(constraints=[0,1,0], reject=irlls_output)

getBvals()¶

Returns a vector of b-values, requires no input arguments

Returns: Vector array of b-values
Return type: ndarray(dtype=float)

Examples

bvals = dwi.getBvals(), where dwi is the DWI class object

getBvecs()¶

Returns an array of gradient vectors, requires no input parameters

Returns: [N x 3] array of gradient vectors
Return type: ndarray(dtype=float)

Examples

bvecs = dwi.getBvecs(), where dwi is the DWI class object

getndirs()¶

Returns the number of gradient directions acquired from the scanner

Returns: number of gradient directions
Return type: ndarray

Examples

n = dwi.getndirs(), where dwi is the DWI class object

goodDirections(outliers)¶

Creates a 3D maps of good directions from IRLLS outlier map For any given voxel, a violation is computed using logical or operant for all b-values. Whether an outlier occurs at b1000 or b1000 and b2000, that voxel is still a violation unless none of the b-values have outliers.

Parameters: outliers (ndarray(dtype=bool)) – 4D maps of outliers from IRLLS
Returns: map – 3D map of number of good directions
Return type: ndarray(dtype=int)

Examples

map = dwi.goodDirections(outliers)

irlls(excludeb0=True, maxiter=25, convcrit=0.001, mode='DKI', leverage=0.85, bounds=3)¶

This functions performs outlier detection and robust parameter estimation for diffusion MRI using the iterative reweigthed linear least squares (IRLLS) approach.

Parameters

exludeb0 (bool, optional) – Exlude the b0 images when removing outliers (Default: True)
maxiter (int, optional Integer; default: 25) – Maximum number of iterations in the iterative reweighting loop (Default: 25)
convcrit (float, optional) – Fraction of L2-norm of estimated diffusion parameter vector that the L2-norm of different vector should get in order to reach convergence in the iterative reweighted loop (Default: 1e-3)
mode (str, {‘DKI’, ‘DTI’}, optional) – Specifies whether to use DTI or DKI model (Default: ‘DKI’)
leverage (float, optional) – Measurement ranging from 0 to 1 where a leverage above this threshold will not be excluded in estimation of DT after outlier (Default: 0.85)
Bounds (int, optional) – Set the threshold of the number of standard deviation that are needed to exclude a measurement (Default: 3)

Returns

outliers (ndarray(dtype=bool)) – 4D image same size as input DWI marking voxels that are outliers
dt (ndarray(dtype=float)) – IRLLS method of DT estimation

Examples

outliers = dwi.irlls()

irllsviolmask(reject)¶

Computes 3D violation mask of outliers detected from IRLLS method

Parameters: reject (ndarray(dtype=bool)) – 4D input outlier map from IRLLS
Returns: propviol – 3D mask where voxel value is the percentage of directional violations
Return type: ndarray(dtype=float)

Examples

mask = dwi.irllsviolmask(outliers)

isdki()¶

Returns logical value to answer the mystical question whether the input image is DKI

Returns: ans – True if DKI; false otherwise
Return type: bool

Examples

ans = dwi.isdki(), where dwi is the DWI class object

kurtosisCoeff(dt, dir)¶

Computes apparent kurtosis coefficient (AKC)

Parameters

dt (ndarray(dtype=float)) – [21 x nvoxel] array containing diffusion tensor
dir (ndarray(dtype=float)) – [n x 3] array containing gradient directions

Returns

adc – Array containing apparent kurtosis coefficient

Return type

ndarray(dtype=float)

Examples

adc = dwi.kurtosisCoeff(dt, dir)

maxBval()¶

Returns the maximum b-value in a dataset to determine between DTI and DKI, requires no input parameters

Returns: maximum B-value in DWI
Return type: float

Examples

a = dwi.maxBval(), where dwi is the DWI class object

multiplyMask(img)¶

Multiplies a 3D image by the brain mask

Parameters: img (ndarray(dtype=float)) – 3D image to be multiplied
Returns: multiplied image
Return type: ndarray(dtype=float)

parfindViols(c=[0, 0, 0])¶

radialSampling(dir, n)¶

Get the radial component of a metric from a set of directions

Parameters

dir (ndarray(dtype=float)) – [n x 3] input array of directions
n (int) – number of rows, n

Returns

dirs

Return type

Matrix containing radial components

Examples

grad = dwi.radiansampling(dir, number_of_dirs)

tensorReorder(dwiType)¶

Reorders tensors in DT to those of MRTRIX in accordance to the table below

Parameters

dwiType (str, {‘dti’, ‘dki’}) – Indicates whether image is DTI or DKI

Returns

DT (ndarray(dtype=float)) – 4D image containing DT tensor
KT (ndarray(dtype=float)) – 4D image containing KT tensor

Examples

dt = dwi.tensorReorder()

Notes

MRTRIX3 and Designer tensors are described below.

MRTRIX3 Tensors                     DESIGNER Tensors
~~~~~~~~~~~~~~~                     ~~~~~~~~~~~~~~~~

 D0      1   1                       1   1
 D1      2   2                       1   2
 D2      3   3                       1   3
 D3      1   2                       2   2
 D4      1   3                       2   3
 D5      2   3                       3   3

 K0      1   1   1   1               1   1   1   1
 K1      2   2   2   2               1   1   1   2
 K2      3   3   3   3               1   1   1   3
 K3      1   1   1   2               1   1   2   2
K4      1   1   1   3               1   1   2   3
K5      1   2   2   2               1   1   3   3
K6      1   3   3   3               1   2   2   2
K7      2   2   2   3               1   2   2   3
K8      2   3   3   3               1   2   3   3
K9      1   1   2   2               1   3   3   3
K10     1   1   3   3               2   2   2   2
K11     2   2   3   3               2   2   2   3
K12     1   1   2   3               2   2   3   3
K13     1   2   2   3               2   3   3   3
K14     1   2   3   3               3   3   3   3

Value Assignment
~~~~~~~~~~~~~~~~

MRTRIX3         DESIGNER
~~~~~~~         ~~~~~~~~
             0
             3
             5
             1
             2
             4

             6
             16
             20
             7
            8
            12
            15
            17
            19
            9
            11
            18
            10
            13
            14

tensorType()¶

Returns whether input image is DTI or DKI compatible, requires no input parameters

Returns: ‘dti’ or ‘dki’
Return type: str

Examples

a = dwi.tensorType(), where dwi is the DWI class object

wlls(shat, dwi, b, cons=None, warmup=None)¶

Estimates diffusion and kurtosis tenor at voxel with unconstrained Moore-Penrose pseudoinverse or constrained quadratic convex optimization. This is a helper function for dwi.fit() so a multiprocessing parallel loop can be iterated over voxels

Parameters

shat (ndarray(dtype=float)) – [ndir x 1] array of S_hat, approximated signal intensity at voxel
dwi (ndarray(dtype=float)) – [ndir x 1] array of diffusion weighted intensity values at voxel, for all b-values
b (ndarray(dtype=float)) – [ndir x 1] array of b-values vector
cons (ndarray(dtype=float)) – [(n * dir) x 22) array containing inequality constraints for fitting (Default: None)
warmup (ndarray(dtype=float)) – Estimated dt vector (22, 0) at each voxel for warm starting constrianed tensor fitting (Default: None)

Returns

dt – Diffusion tensor

Return type

ndarray(dtype=float)

Examples

dt = dwi.wlls(shat, dwi, b, constraints)

Notes

For Unconstrained Fitting: In the absence of constraints, an exact formulation in the form Cx = b is produced. This is further simplified to x_hat = C^+ * b. One can use the Moore-Penrose method to compute the pseudoinverse to approximate diffusion tensors.

For Constrained Fitting: .. code-block:: none

The equation |Cx -b|^2 expands to 0.5*x.T(C.T*A)*x -(C.T*b).T

~~~~~ ~~~~~
P q

where A is denoted by multiplier matrix (w * b) Multiplying by a positive constant (0.5) does not change the value of optimum x*. Similarly, the constant offset b.T*b does not affect x*, therefore we can leave these out.

Minimize: || C*x -b ||_2^2: subject to A*x <= b No lower or upper bounds

designer.fitting.dwipy.clipImage(img, range)¶

Clips input matrix within desired range. Min and max values are inclusive of range Classification: Function

Parameters

img (ndarray(dtype=float)) – Input 3D image array
range (array_like) – [1 x 2] vector specifying range to clip

Returns

clippedImage

Return type

clipped image; same size as img

Examples

clippedImage = clipImage(image, [0 3]) Clips input matrix in the range 0 to 3

designer.fitting.dwipy.highprecisionexp(array, maxp=1e+32)¶

Prevents overflow warning with numpy.exp by assigning overflows to a maxumum precision value Classification: Function

Parameters

array (ndarray) – Array or scalar of number to run np.exp on
maxp (float, optional) – Maximum preicison to assign if overflow (Default: 1e32)

Returns

Return type

exponent or max-precision

Examples

a = highprecisionexp(array)

designer.fitting.dwipy.vectorize(img, mask)¶

Returns vectorized image based on brain mask, requires no input parameters If the input is 1D or 2D, unpatch it to 3D or 4D using a mask If the input is 3D or 4D, vectorize it using a mask Classification: Function

Parameters

img (ndarray) – 1D, 2D, 3D or 4D image array to vectorize
mask (ndarray) – 3D image array for masking

Returns

vec – of DWI volumes

Return type

N X number_of_voxels vector or array, where N is the number

Examples

vec = vectorize(img) if there’s no mask vec = vectorize(img, mask) if there’s a mask

designer.fitting.dwipy.writeNii(map, hdr, outDir, range=None)¶

Write clipped NifTi images

Parameters

map (ndarray(dtype=float)) – 3D array to write
header (class) – Nibabel class header file containing NifTi properties
outDir (str) – Output file name
range (array_like) – [1 x 2] vector specifying range to clip, inclusive of value in range, e.g. range = [0, 1] for FA map

Returns

Return type

None; writes out file

Examples

writeNii(matrix, header, output_directory, [0, 2])

See Also clipImage(img, range) : this function is wrapped around