ConTrack comprises three essential parts. First, scoring procedures encode knowledge about pathway likelihood; the scoring algorithm is a key component in finding anatomically correct pathways. Second, the vast number of possible pathways makes it impossible to score every pathway, so ConTrack includes a sampling procedure to select potential pathways for scoring. Third, the algorithm specifies a simple inference method that defines how to use sampled and scored pathways to select the most likely pathways that connect two regions. illustrates the ConTrack algorithm by showing example results from each of the algorithm’s three stages: ROI specification, pathway sampling, and pathway scoring and selection.

Scoring procedures

Principles Pathway scoring encapsulates our knowledge about the likely properties of white matter pathways. The scoring used in ConTrack also incorporates two basic principles that have not been adopted by other algorithms: symmetry and independence ().

Symmetry means that the pathway score measured from *R*_{1} to *R*_{2} should equal the score of the same pathway measured from *R*_{2} to *R*_{1}. Independence means that the score of a pathway between *R*_{1} and *R*_{2} depends only on the data along that pathway. Therefore, the score along the green pathway between *R*_{1} and *R*_{2} should be the same in the presence or absence of another pathway. For example, if the orange tensors in become isotropic due to a disease, thereby lowering the score of a pathway traversing them, the score along the *R*_{1}–*R*_{2} pathway should remain unchanged.

Notice also that the pathway from *R*_{1} to *R*_{2} is unlikely to be identified by a deterministic algorithm with typical stopping rules based on a fractional anisotropy (FA) threshold; the FA of the green tensor in the fourth column is very low and would stop the growth. But if we know *R*_{1} is connected to *R*_{2}, the green pathway is the most likely route.

Algorithm Consider a single sampled pathway with endpoints in the two regions of interest (ROIs) of . We call the sample points along the pathway “nodes.” illustrates a short pathway composed of two linear segments connecting three nodes (*s*_{i}). The pathway score is determined by the diffusion data at each node and the angles between the linear segments.

We define the score, *Q*(*s*), as,

where

*s* is a pathway and

*D* are the data along the pathway. The scoring function is defined in this way to separate the data-dependent portion of the score,

*p*(

*D* |

*s*), from the data independent portion,

*p*(

*s*). Both the

*p*(

*D* |

*s*) and

*p*(

*s*) components of the pathway score are a product of the local scores along all the nodes of the pathway, as described in the following paragraphs.

Scoring the pathway fit to the data

The function *p*(*D* | *s*) measures how well the existence of pathway *s* is supported by the tensors *D* derived from diffusion measurements. The pathway score aggregates local terms multiplicatively, i.e.,

where

*v*_{1},

*v*_{2}, and

*v*_{3} are the three ordered eigenvectors of the

*i*th diffusion tensor,

*D*_{i};

*C*(

*σ*_{3},

*σ*_{2}) is a normalizing constant so that the function integrates to one over the sphere; and

*σ*_{i} is the decay rate in the direction of the eigenvector,

*v*_{i}. The two decay rate parameters are determined by summing two dispersion parameters that are based on (a) uncertainty in the data (

*σ*_{m}) and (b) the ellipsoid shape (

*σ*_{i}*). That is,

*σ*_{i} =

*σ*_{m} +

*σ*_{i}*.

We calculate one dispersion parameter (

*σ*_{m}) from repeated measures of the principal direction of diffusion or PDD (

*v*_{1}). The repeatability of the data sets a lower bound on fiber direction uncertainty. We estimate

*σ*_{m} using a bootstrap procedure (

Efron, 1979): a tensor is fit to 1000 permutations of the raw DWI. The bootstrap procedure estimates the mean PDD and the dispersion about that mean. Other bootstrap procedures, such as the wild bootstrap (

Liu, 1988), can be used if the DWI data lacks an adequate number of repetitions (

Whitcher, Tuch, & L., 2005). We summarize the bootstrapped dispersion,

*σ*_{m}, using the Watson distribution on the 3D sphere (

Mardia & Jupp, 2000;

Schwartzman, Dougherty, & Taylor, 2005). The Watson distribution is identical to a Bingham distribution with radial symmetry, i.e., when

*σ*_{2} =

*σ*_{3}.

The bootstrap sometimes produces an unreasonably low dispersion parameter. Hence, we set a minimum dispersion to 4°—the mean dispersion of all brain voxels with a high linearity index (

*C*_{L} >0.4) in all of our data sets. The linearity index is a measure of anisotropy (

Peled, Gudbjartsson, Westin, Kikinis, & Jolesz, 1998) and is the positive difference between the largest two eigenvalues of the diffusion tensor divided by the sum of its eigenvalues. In our algorithm, increasing the linearity criterion does not change the minimum dispersion value; decreasing the linearity criterion increases the minimum dispersion to a value larger than 4°. The precise value of this parameter is not critical to the performance of ConTrack.

We derive a second dispersion term, *σ*_{i}*, from the diffusion shape. This dispersion term is necessary because the PDD of an oblate or nearly isotropic diffusion ellipsoid might be highly reliable, yet such voxels clearly have high fiber direction uncertainty. Reliable measurements of oblate and spherical ellipsoids occur when voxels include many crossing fibers, such as where the tapetum interdigitates with the inferior longitudinal fasciculus.

A user-specified function relates ordered tensor eigenvalues (*λ*_{i}) to the shape dispersion term.

The value *δ* is a dispersion factor calculated from *C*_{L} the tensor linearity index, and a user parameter *η*, discussed below.

The constant, 100°, is maximum total dispersion. A perfect spherical tensor, where *σ*_{m} = 4° and *λ*_{1} = *λ*_{2} = *λ*_{3}, would have a uniform distribution about the sphere (*σ*_{2} = *σ*_{3} = 54°).

The constant 0.015 is chosen so that *δ* varies between the maximum added dispersion and zero within ~0.15 units of *C*_{L}. Adjusting these constants did not alter the results significantly.

Four examples of local distribution functions are shown in . When the linearity index is small, the shape of the tensor adds significant uncertainty to the local distribution function. This uncertainty may be concentrated along one axis for more oblate-shaped ellipsoids or distributed uniformly across the sphere for more spherical ellipsoids. When the ellipsoid has a very high linearity index (i.e., a prolate ellipsoid), the total fiber direction uncertainty equals the bootstrap dispersion estimate.

Scoring the pathway shape independent of the data

The function

*p*(

*s*) (

Equation 1) measures how well the pathway represents prior knowledge about neuronal connections. This term currently encodes assumptions such as the smoothness of white matter pathways, length preference, as well as the possible locations for endpoints of pathways. We expect that over the years our knowledge will increase based on studies of the neuroanatomy and the function can evolve to incorporate this knowledge.

Using the angles of the pathway as it enters and exits each node (*θ*_{i}, see ), the location of pathway nodes, and the number of nodes, we assess the score of a pathway (independent of the data) as

where we set 0 ≤ |

*θ*_{i} | ≤

*π*/2, i.e., only directions on the hemisphere are possible,

*σ*_{c} is the angular dispersion parameter in degrees, and

*C*(

*σ*_{c}) is a normalizing constant. The Watson distribution of axes on a sphere is easily transformed into a distribution of directions on a hemisphere and this was a motivation for selecting this distribution. The possible directions are distributed on the hemisphere because it is assumed that the pathways have a node spacing small enough such that an intra-pathway angle >90° would be too sharp to reconstruct the desired anatomy. When

*σ*_{c} is less than 54° (uniform dispersion), a pathway is penalized for having many large angular deviations.

The two functions, *p*_{length} and *p*_{end} are defined on the pathway nodes.

where

*λ is* a manually chosen length scoring parameter and the white matter mask is defined as those voxels that are within the brain mask and have an FA >0.15 and either mean diffusivity (MD) less than 1.1

*μ*m

^{2}/ms or FA >0.4. This mask is dilated by one voxel to allow for partial voluming at the edges of the white matter. Using FA and MD to create the white matter mask is a common method (e.g.,

Basser et al., 2000;

Behrens, Woolrich, et al., 2003;

Mori, Crain, Chacko, & van Zijl, 1999;

Parker et al., 2003). As with other algorithms, the mask prevents pathways from entering CSF or gray matter. The accuracy of the white matter mask was ensured by visual inspection.

Pathway sampling

The pathway sampling algorithm identifies a large set of potential pathways that have endpoints in the two ROIs. The pathway sampling procedure differs from procedures used in other probabilistic tractography algorithms (

Behrens, Woolrich, et al., 2003;

Friman & Westin, 2005;

Parker & Alexander, 2005). Rather than sampling strictly proportionally to estimated anatomical validity (likelihood), the ConTrack sampler includes pathways with a broad range of scores.

The algorithm begins by placing a seed point randomly within one of the voxels in an ROI. A pathway is grown from this seed position by iteratively choosing step directions and moving a unit step distance along that direction. A step direction is chosen in one of two ways; (a) if this is the first step or if the fiber direction uncertainty is small (*σ*_{3} < 14°), then the next step direction is chosen according to *p*(*D*_{i} | *t*_{i}); (b) otherwise the step direction is drawn from the *p*_{curve}(*θ*_{i}) distribution. The step distance is currently 1 mm. The pathway is terminated if the current node is outside of the white matter mask or within an ROI voxel. Only pathways with endpoints in both ROIs are retained. The process repeats until a desired number (~10^{5}) of pathway samples is collected.

While the generation of a single pathway must start in one ROI, pathways are seeded equally from each of the two ROIs. Hence, the population of pathways is symmetric and not biased to either ROI.

Inference methodology

ConTrack scores the likelihood of each sampled pathway. We infer that the highest scoring pathways are the most likely to exist; the selection criterion for anatomical validity is based on a combination of visual inspection of the pathways and comparison of the anatomical validity scores for pathways identified by STT. We derive various quantitative measures from the highest scoring pathways, including the pathway destination (projection zone), mean fractional anisotropy, pathway length, and so forth.

This inference—that the highest scoring pathways are the most likely to exist—is a simple and perhaps obvious design choice for our application; yet, the reader should note that the inference methodology used in ConTrack differs significantly from that used by probabilistic tractography algorithms, which measure the connection probability between two regions. In the Discussion section, we describe why the combination of separating sampling from scoring and the difference in inference strategy enables ConTrack to find pathways that other algorithms miss.