WIP: Improve the axonal water fraction estimation.

[jhlegarreta]

This patch set proposes to use the more accurate estimation of the axonal
water fraction (AWF) given by eq. (8) in Fieremans' paper. The current
implementation assumes that there is a direction for which the
intra-axonal space diffusivity Da,i for the i-th voxel is 0 for some
direction. This assumption may not be met in some circumstances, and hence
the expression of the current implementation (eq. (7)) is a lower bound
for the AWF.

[jhlegarreta]

@arokem @RafaelNH

A couple of issues discussed in a separate thread with Ariel and Rafael:

I haven't find the way to get the Di, and Ki i=1..3 components. Fieremans' paper says they are the diffusion/kurtosis components along the axis directions of a chosen reference frame. If I'm not mistaken, all components are indeed computed at dki_micro.py:157 and dki_micro.py:158, but I still do not see which components we are interested in.

Any help is appreciated.

As for the fit of the method within the current dki framework, not sure whether the method would replace the current implementation. I thought that we could may be have a boolean to indicate whether we want the exact or the approximate version, or else whether we should check the Da,i=0 assumption to decide for one or the other.

As for the test, I did not think in which ways we could get a ground truth. I've seen that the current test assumes a single fiber; I don't know whether a non-zero Da,i condition can be effectively simulated if the axonal radius is increased beyond 5 micrometer.

According to Rafael, using the more general AWF estimation function increases substantially the complexity of the implementation (e.g. this will decrease computing speed and might also add singularities that impact the robustness of the fit). Moreover, in practice using the more general approach is challenging because a non-zero estimate is hard to obtain (or even impossible using current MRI acquisitions). Given the later, checking the Di,a=0 condition automatically will not be feasible as well. So the only way to implement this is giving Di,a as an additional input.

To test the code, you can simulate a well-aligned fiber where Di,a is different to zero. The fitted model parameters should match the ground truth values if the ground truth non-zero Di,a value is given as input.

Still missing to include a call in the dki_reconst.py example.

[arokem]

Hey @jhlegarreta : do you have any more thoughts here? I'm going to close this for now, but if you want to keep going here, please feel free to reopen.

[jhlegarreta]

@arokem Thanks for keeping an eye on this. It's been a long time, but I would definitely like to take up this again under your and Rafael's guidance.

Just to establish some priorities or finish some contributions that would close a couple of issues, and that should be doable in a shorter time frame than the current one, I'd like to get some feedback on the following two, do the necessary changes, and hopefully get them ready to be merged:
#1709
#1711

Once that done, I think it is reasonable to come back to this challenge.

Does this look sensible to you?