rejection map

[bpinsard]

To add to a whishlist. It would be great to add to tckgen an option to output a map that counts, per voxel, the number of times a streamline have been rejected arriving there. This could help debug 5tt or choose tracking parameters and understand why we have difficulty finding some tracks in some subjects, sometimes due to faulty anatomical priors or noisy dwi data. Do you think it's feasible? Thanks.

[thijsdhollander]

Hi @bpinsard,

That sounds like a really useful suggestion. I'm not an expert on the beast that is our tracking code, but I think that should technically be possible to be plugged in somewhere (without any major redesign). The map could even be coded (4D or something) to provide information on what kind of rejection it was (e.g. due to FOD amplitude, leaving a mask and still being too short, hitting an exclusion region, violating anatomical constraints, etc...). Thanks for suggesting!

Cheers, Thijs

[Lestropie]

@bpinsard: Try un-commenting this line (delete the two leading forward slashes) and re-compiling: tckgen will then write a number of images to the working directory, indicating where streamlines have been terminated for various reasons. This is slightly different to the way streamlines are rejected: the latter usually results from certain types of the former, but they are technically different concepts (and I'm therefore very picky about the naming :-P ).

[bpinsard]

Great, that should be very useful. Thanks.

[arnaudbore]

@Lestropie Quick question following the debug mode to get these maps. Trying to track a specific bundle, I get a lot of Calibrate Fail because of max_val == 0 iFOD2. Can you explain what max_val refer to exactly and why I get these zeros ? I can send you the images and command line so you can replicate if you need to.

[Lestropie]

That's a pretty common termination mechanism; I probably shouldn't have called it a "failure"...

I'll try to give a brief explanation of what's going on, since a full explanation with demonstration should be in either an iFOD2 paper or (more likely) the eventual MRtrix3 paper.

Before extending the streamline from a particular point, the algorithm projects a pre-determined set of potential pathways (arcs), and samples the 'path probabilities' (~ FOD amplitudes) along them. Based on a prior 'calibration' step (which is done before tracking starts at all), it then estimates the maximal likely path probability within the angular cone in which the streamline is free to track. By knowing this upper bound, subsequent sampling of the FODs to determine the actual tracking direction is more efficient.

This termination occurs when, for all calibration paths, at some point along the path, the FOD amplitude is below threshold.

My personal criticism of the particulars of this approach is that even if all calibration paths have at least one point that is below the FOD amplitude threshold (and hence the streamline is terminated), there may still be a pathway within the arc of curvature that is above the FOD amplitude threshold along its entire length. Perhaps instead, it should be the rejection sampling ceiling that is tested against the FOD amplitude threshold (with appropriate transformation). I suspect a lot of this code pre-dates extensive use of ACT, where thresholding by FOD amplitude is the most important mechanism for terminating streamlines. Or maybe I just lost the argument... @jdtournier Any thoughts?

[jdtournier]

That's a pretty common termination mechanism; I probably shouldn't have called it a "failure"...

Yes, I'd agree it's a bit unfair to call it a failure - it's basically saying there's no valid FOD amplitude within the calibration set (calibration is performed at each step, it's required to set up the rejection sampler), which means it's never going to find a valid next step from here anyway (with the one caveat raised by @Lestropie, see below). Any chance we can relabel that as a termination?

My personal criticism of the particulars of this approach is that even if all calibration paths have at least one point that is below the FOD amplitude threshold (and hence the streamline is terminated), there may still be a pathway within the arc of curvature that is above the FOD amplitude threshold along its entire length.

That's probably fair, and we could fix that relatively easily, I think. But it would involve disabling the explicit "return zero if below threshold at any point" aspect during calibration only, which would then allow the calibrator to at least produce an estimate of the maximum amplitude, and the decision to terminate could be taken then. But note that this is only going to matter for very borderline cases anyway, where the amplitudes for the best path in the calibration set are close to threshold. So I don't expect it'll have a noticeable impact overall, these should be very low probability steps in any case.

Perhaps instead, it should be the rejection sampling ceiling that is tested against the FOD amplitude threshold (with appropriate transformation). I suspect a lot of this code pre-dates extensive use of ACT, where thresholding by FOD amplitude is the most important mechanism for terminating streamlines.

Yes, but we'd need to be careful to do this for the initial calibration only. The reason for that explicit "set to zero" line in the path probability calculation is to prevent paths where the tangent direction points along a below-threshold direction at any of the sample points along the arc - something that we wouldn't allow using other methods, and which brings the potential of allowing non-sensible paths (for example, switching between two dominant directions in a crossing region, if the sample points happen to have below-threshold, but above-zero FOD amplitudes - the large dominant FOD amplitudes at the other sample points can then compensate for the low ones). I don't think we'd want to disable that, it's not so much about termination as just ensuring the next step matches with expectation - and that's relevant whether using ACT or any other form of termination.

And yes, iFOD2 does predates ACT (ISMRM abstracts in 2010 & 2012 respectively)...

[Lestropie]

Any chance we can relabel that as a termination?

Done (in tag_0.3.16).

That's probably fair, and we could fix that relatively easily, I think.

Yeah it's easy, we just need to test how the behaviour changes both with and without ACT (I suspect it'll be much better for ACT, less erroneous terminations in deep WM crossings, I just never got around to trying it). There's already a branch tckgen_mods on GitLab for playing with such things, so keep an eye there if you're interested.

Yes, but we'd need to be careful to do this for the initial calibration only.

Yeah I'm on top of all of that. Easy to do: Just make path_prob() a template <bool>.

And yes, iFOD2 does predates ACT (ISMRM abstracts in 2010 & 2012 respectively)...

But was the calibrator implemented at the same time as iFOD2, or did it come later?

[Lestropie]

OK, quick tests (since this is more fun than grant writing):

Without ACT:

                                 Current        Proposed
Time (100k tracks)               1m36           3m35
Generated                        215857         148146
Calibrator terminations          98.3%          73.4%
No valid path found              0.015%         20.5%
DWI brain mask terminations      1.4%           5.7%

With ACT:

                                 Current        Proposed
Time (100k tracks)               3m58           3m45
Generated                        601460         256540
Calibrator terminations          54.4%          20.8%
No valid path found              0.008%         7.9%

So basically what I predicted:

• Without ACT, iFOD2 currently relies heavily on the calibrator to terminate streamlines. When the calibrator doesn't use a direct FOD threshold along the path ('proposed'), streamlines can propagate further (more reach the edge of the brain mask), but tracking also takes much longer because the algorithm is frequently generating 1000 paths from the endpoint before terminating.

• With ACT, there's far less tracking terminations / failures in the middle of the white matter. Still sometimes the calibrator perhaps should have caused a termination but didn't; but the computational expense of that is offset by the fact that far fewer streamlines need to be generated overall.

P.S. I'm also now convinced that the default iFOD2 power should be 0.33, not 0.25 :-/

[jdtournier]

Right, well that's fairly convincing. The effect is already larger that I'd have anticipated without ACT, but basically what it does is reduce the number of terminations due to FOD thresholding (whether at the calibrator stage or during sampling proper) from 98.3% to 93.4% - so in this case ~5% of terminations were premature.

I'm struggling to understand why the effect would be so much larger with ACT though, the corresponding change is from 54.4% down to 28.7% - pretty dramatic. And from what you're saying, these are predominantly in the middle of white matter? It'll take me a while to wrap my head around this one...

But was the calibrator implemented at the same time as iFOD2, or did it come later?

More or less: first commit dates from June 2010. But that's not really the issue: I wouldn't have expected such a large interaction with ACT, so even if they had been developed at the same time, there's nothing to suggest things wouldn't have ended up the same. The logic at the time was to use the path probability calculation for both iFOD1 & iFOD2 and run it over a sufficiently dense set of points to figure out a decent estimate of what the maximum value might be. So by that rationale, we should use exactly the same path probability calculation as is actually used in tracking. The problem is that the thresholding applied in the calculation for iFOD2 interferes with that calibration in ways I hadn't anticipated - or certainly hadn't anticipated making such a large difference...

P.S. I'm also now convinced that the default iFOD2 power should be 0.33, not 0.25 :-/

You may very well be right - I also think it's not quite right, and should be a touch higher, very much in that same ballpark as you suggest. We'll need to nut this one out properly, convince ourselves that it's correct, and make the change if required. I'll see if I can find some time to construct a reasoned argument around this today...

[Lestropie]

... so in this case ~5% of terminations were premature.

Well, not necessarily: most of that difference is going to streamlines leaving the DWI brain mask, which it could be argued are not terminating early enough. But then, ill-posed terminations is what ACT is there for.

I'm struggling to understand why the effect would be so much larger with ACT though, the corresponding change is from 54.4% down to 28.7% - pretty dramatic. And from what you're saying, these are predominantly in the middle of white matter? It'll take me a while to wrap my head around this one...

I've always had my suspicions about the calibrator, but had mostly been looking at the number of calibration directions. You can see it in the ACT paper, Figure 6 bottom row: All of the terminations in WM are due to the calibrator. In fact IIRC this is what led me to create back-tracking more so than other biologically-erroneous terminations.