Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add script to validate gradient tables #1110

Open
wants to merge 14 commits into
base: master
Choose a base branch
from

Conversation

karanphil
Copy link
Contributor

Quick description

I wrote a script that validates the distribution of directions of b-vectors. To do so, I used the already existing functions in scilpy.gradients.gen_gradient_sampling to compute the electrostatic-like repulsion energy of inputed b-vectors and compare its energy to the energy of "optimal" b-vectors. These optimal b-vectors are computed on-the-fly with generate_gradient_sampling, the Caruyer method we use in scil_gradients_generate_sampling.py, with the same architecture of number of bvecs per shell as the input.

The script starts by looking for b0s, to remove them from the analysis. Then, it looks for duplicate b-vectors. Finally, both energies are computed and compared as the ratio between the inputed b-vectors' energy and the optimal b-vectors' energy (input_energy/optimal_energy). Above a given maximum ratio value, the script raises a warning.

The user might want to use the -v verbose option to see the computed energies. The --visualize option displays both the inputed and optimal b-vectors on a single shell.

For the default value of max_ratio (1.25), I looked at various databases like HCP, Penthera, ADNI3 and converged at this number. Multi-shell gradient tables often brake this, because many gtabs have similar bvecs across shells and it raises the energy a lot. I could add an option "split_shells".

Would be curious to know what @mdesco thinks about this!

...

Type of change

Check the relevant options.

  • Bug fix (non-breaking change which fixes an issue)
  • New feature (non-breaking change which adds functionality)
  • Breaking change (fix or feature that would cause existing functionality to not work as expected)
  • This change requires a documentation update

Provide data, screenshots, command line to test (if relevant)

...

Checklist

  • My code follows the style guidelines of this project (run autopep8)
  • I added relevant citations to scripts, modules and functions docstrings and descriptions
  • I have performed a self-review of my code
  • I have commented my code, particularly in hard-to-understand areas
  • I have made corresponding changes to the documentation
  • My changes generate no new warnings
  • I moved all functions from the script file (except the argparser and main) to scilpy modules
  • I have added tests that prove my fix is effective or that my feature works
  • New and existing unit tests pass locally with my changes

Copy link

codecov bot commented Dec 16, 2024

Codecov Report

Attention: Patch coverage is 68.62745% with 32 lines in your changes missing coverage. Please review.

Project coverage is 69.55%. Comparing base (0ac4461) to head (c3a17b5).

Additional details and impacted files
@@            Coverage Diff             @@
##           master    #1110      +/-   ##
==========================================
- Coverage   69.56%   69.55%   -0.02%     
==========================================
  Files         448      449       +1     
  Lines       24078    24176      +98     
  Branches     3290     3305      +15     
==========================================
+ Hits        16751    16816      +65     
- Misses       5934     5962      +28     
- Partials     1393     1398       +5     
Components Coverage Δ
Scripts 70.26% <71.66%> (+<0.01%) ⬆️
Library 68.62% <64.28%> (-0.04%) ⬇️

EmmaRenauld
EmmaRenauld previously approved these changes Jan 9, 2025
Copy link
Contributor

@EmmaRenauld EmmaRenauld left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Looks really clean and well-written! (Although I don't understand very well the energy section, might require an additional reviewer).

Let's not forget to update after PR #1113 is merged.

nb_shells = len(ubvals)
nb_dir_per_shell = [np.sum(shell_idx == idx) for idx in range(nb_shells)]

ubvecs = np.unique(bvecs, axis=0)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should you also check if some bvecs are the same as another inversed one (-bvecs)?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants