matched shore computation to brainsuite

[mattcieslak]

Hi Dipy community,

@dvya and I went through the Dipy 3dSHORE implementation and compared it to the BrainSuite implementation. We found a difference in the number of bases generated in Dipy and Brainsuite. For example with radial order set to 6 there were 50 bases in Dipy and 72 in BrainSuite.

In the Merlet 2013 paper right after eq 18 they write

For SHORE basis, the radial order is n = 6 and because l is bounded by n, this gives 72 atoms.

However, the documentation for 3dSHORE in Dipy says that a modified version is being used, so maybe the difference is due to the modification?

[skoudoro]

Thank you @mattcieslak and @mattcieslak!

I think this is @maurozucchelli who implemented this modified version of SHORE so it will be good if you could answer the @mattcieslak's question. Thank you!

[mattcieslak]

This definitely breaks some of the unit tests.

I see the differences from Merlet are described in the appendix of @maurozucchelli's 2016 paper. Is this just a case of differences in the interpretation of what radial_order should mean?

[arokem]

I don't understand the issue deeply, but if it is indeed a difference in interpretation, I think that it would be great to have both options possible. I would make the current implementation the default, and then add a kwarg to the model object that allows to use this interpretation instead.

[arokem]

Also, tiny detail: is there a reason you forgo the preallocations that were there before?

[mattcieslak]

the n_c was evaluating to the truncated basis set (50) instead of the full set (72). I can change the calculation so it matches the new loops. Funny thing is the brainsuite code used n_c too, but since it's matlab it just kept appending as more than n_c values were calculated.

Do you think it makes sense to make a separate class? Like a BrainSuiteShore or MerletShore or something? It would need different methods to calculate all the rtop/rtap/etc and adding lots of ifs in the functions makes life messier.

[arokem]

No - I think this is better handled through key-word arguments to one class that implements either.

[Garyfallidis]

@mattcieslak please also report speed difference between the two versions. Thank you! Really appreciate your help. :)

[arokem]

Hey @mattcieslak: would be good to pick this back up. How do you feel about providing both options? This would be similar to what we do with DTI, so that people can fit it with different options: least-squares, weighted least-squares, RESTORE, and so on.

[skoudoro]

Hi @mattcieslak, it would be good to have this feature for DIPY 1.0 release (mid July). Any chance to have an update with both options and a small test? Let us know if you need any help to make it happen!

[mattcieslak]

I think this is doable. I don't think it should be an option in the currently-implemented 3dSHORE class, though. Unlike DTI, where the many ways to fit always result in 6 coefficients that are used the same way for FA/MD/etc, the BrainSuiteShore functions for rtop/rtap/odf/pdf/etc would all need to be changed to check if the instance is using the Zucchelli bases or the BrainSuite bases. It will be easier to maintain and read if it's just a separate class. Would that work for you?

[arokem]

Yes. That makes sense. It seems more similar to the distinction between fwdti and dti than between different ways of fitting dti.

[pep8speaks]

Hello @mattcieslak, Thank you for updating !

• In the file dipy/reconst/brainsuite_shore.py:

Line 19:1: E302 expected 2 blank lines, found 1 Line 612:32: E201 whitespace after '(' Line 613:21: E128 continuation line under-indented for visual indent Line 613:46: E502 the backslash is redundant between brackets Line 614:21: E128 continuation line under-indented for visual indent

• In the file dipy/reconst/shm.py:

Line 295:17: E127 continuation line over-indented for visual indent

Comment last updated at 2019-07-14 22:45:37 UTC

[codecov-io]

Codecov Report

Merging #1639 into master will decrease coverage by 1.27%. The diff coverage is 17.03%.

@@            Coverage Diff             @@
##           master    #1639      +/-   ##
==========================================
- Coverage   85.31%   84.04%   -1.28%     
==========================================
  Files         118      119       +1     
  Lines       14221    14490     +269     
  Branches     2236     2273      +37     
==========================================
+ Hits        12133    12178      +45     
- Misses       1577     1800     +223     
- Partials      511      512       +1

Impacted Files	Coverage Δ
dipy/reconst/brainsuite_shore.py	17.95% <17.95%> (ø)
dipy/reconst/shm.py	80.28% <8%> (-5.7%)	:arrow_down:

[arokem]

Hey @mattcieslak : is this ready for a review or are you still working on it?

[mattcieslak]

It's almost there. I am working on the tests and getting rtop/rtap/etc. Should be ready in a couple days

[skoudoro]

Hi @mattcieslak, just for your information, we plan to release this Wednesday (31rst). Let us know if you find time to complete this PR for a review. Thanks!

[mattcieslak]

Sorry @skoudoro I'm still waiting to get the equations for rtop, rtap, etc for these bases. I'm not sure I'll have it ready to go in time :(

[skoudoro]

ok, no worries and no problem! I will tag your PR for the release 1.1 in October.

[Borda]

@skoudoro maybe you may remove tag "need review" since there already a couple of reviewers :)

[skoudoro]

The tag has already been removed but Thank you for the notice! :)

[skoudoro]

Hi @mattcieslak,

I'm still waiting to get the equations for rtop, rtap, etc for these bases

Any updates concerning this PR? Did you get the equations? How can we help?

[mattcieslak]

@skoudoro do you know any of the brainsuite developers? Their message board was a dead end

[skoudoro]

Hi @mattcieslak,

We plan to talk about your PR during the DIPY online meeting tomorrow (October 7th, 2021 | 1pm EST / 7pm CET / 10am PT).

It would be great if you could attend this meeting to provide us a short update.

More information on https://github.com/dipy/dipy/discussions/2398

Thank you!