Skip to content

Create ae script#1302

Open
EmmaRenauld wants to merge 7 commits intoscilus:masterfrom
EmmaRenauld:create_ae_script
Open

Create ae script#1302
EmmaRenauld wants to merge 7 commits intoscilus:masterfrom
EmmaRenauld:create_ae_script

Conversation

@EmmaRenauld
Copy link
Contributor

@EmmaRenauld EmmaRenauld commented Jan 28, 2026

Quick description

New script to compute the local angular error at each point of the streamline, compared to underlying peaks. For each point, the following segment direction is compared with the peaks in this voxel (nearest neighbor). If peaks come from DTI; simply compared to the vector e1. If peaks come from fODF, compared to the most aligned peak. The output is an angle in degree, saved in DPP. Value at the last coordinate of streamlines is always 0. Value if the voxel has no peak is also 0.

Optionally, the AE can also be saved as color, or projected to a map (mean of segments in that voxel).

Then, I also wanted a script to extract the worst N streamlines based on that value, so here it is at the same time.

I will add official tests soon, but this can already be tested.

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)

Using on ISMRM here

  • The ground truth streamlines
  • The DWI that I preprocessed using Tractoinferno: fodf peaks.
scil_tractogram_compute_ae INPUTS_AE/GT_ismrm.tck INPUTS_AE/ISMRM_fodf_peaks.nii.gz \
       TESTS_AE/test_ismrm_ae_fodf_GT.trk --no_bbox --save_mean_map TESTS_AE/test_ismrm_ae_fodf_GT_map.nii.gz \
      --save_as_color --processes 6 -v
scil_tractogram_extract_streamlines TESTS_AE/test_ismrm_ae_fodf_GT.trk TESTS_AE/test_ismrm_ae_fodf_GT_WorstNb500.trk \
      --from_dpp AE --top --nb 500 -v 

Results:

The tractogram:
image

With "Toggle 2D fiber mappers"
image

The projected map:
image

Extracted worst 500 streamlines:
image

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 ----> todo!
  • New and existing unit tests pass locally with my changes

@EmmaRenauld EmmaRenauld requested a review from mdesco January 28, 2026 21:18
@codecov
Copy link

codecov bot commented Jan 29, 2026

Codecov Report

❌ Patch coverage is 75.59809% with 51 lines in your changes missing coverage. Please review.
✅ Project coverage is 72.55%. Comparing base (0ba0902) to head (5e1f493).

Additional details and impacted files
@@            Coverage Diff             @@
##           master    #1302      +/-   ##
==========================================
+ Coverage   72.51%   72.55%   +0.04%     
==========================================
  Files         295      297       +2     
  Lines       25442    25651     +209     
  Branches     3565     3586      +21     
==========================================
+ Hits        18449    18611     +162     
- Misses       5488     5523      +35     
- Partials     1505     1517      +12     
Flag Coverage Δ
smoketests 69.77% <75.59%> (+0.06%) ⬆️
unittests 13.13% <0.00%> (-0.11%) ⬇️

Flags with carried forward coverage won't be shown. Click here to find out more.

Components Coverage Δ
Scripts 75.45% <82.11%> (+0.07%) ⬆️
Library 69.06% <58.62%> (-0.02%) ⬇️
🚀 New features to boost your workflow:
  • ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

Copy link
Contributor

@mdesco mdesco left a comment

Choose a reason for hiding this comment

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

This looks great. I tested and reproduced the outputs. It would be nice to have a formal technical review from @karanphil or @CHrlS98 that have played with local fixels in the past PR.

Compute the angular error (AE) for each segment of the streamlines.

For each segment of each streamline, the direction is compared with the
underlying peak (for single peak files like DTI) or with the closest peak
Copy link
Contributor

Choose a reason for hiding this comment

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

si possible, ajouter le terme eigen_vector pour le peak ici.

@arnaudbore arnaudbore requested a review from CHrlS98 February 6, 2026 01:56
help="Saves the streamlines in the top / lowest percentile.\n"
"Default if set: The top / bottom 5%%")
gg.add_argument(
'--mean_std', type=int, const=3, nargs='?', dest='std',
Copy link
Contributor

Choose a reason for hiding this comment

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

Not obvious that this argument takes the number of standard deviations as parameter.

- WPC: wrong path connections, streamlines connecting correct ROIs but not
respecting the other criteria for that bundle. Such streamlines always
exist but they are only saved separately if specified in the options.
exist but they are only saved separately if specified{} in the options.
Copy link
Contributor

Choose a reason for hiding this comment

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

This may have been added by mistake?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

oups thanks

dirs = dirs / np.linalg.norm(dirs, axis=-1, keepdims=True)

# Concatenating streamlines for faster processing + nearest neighbor
coords = np.floor(np.vstack([s[1:] for s in sft.streamlines])).astype(int)
Copy link
Contributor

Choose a reason for hiding this comment

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

this won't be the intersection between the streamline and all voxels (if the step is greater than the voxel size or if there was streamlines compression). Is this what you want? Maybe, since you modify the dpp property.

But, if you wanted the intersection with all voxels, you could use subdivide_streamlines_at_voxel_faces in tractanalysis.voxel_boundary_intersection.


# Using the abs value because vectors are undirected.
cos_theta = np.abs(np.dot(current_peaks, dirs[i]))
cos_theta = np.clip(cos_theta, -1.0, 1.0) # numerical safety
Copy link
Contributor

Choose a reason for hiding this comment

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

Since it's an absolute value you can clip between 0-1 no?

Copy link
Contributor

Choose a reason for hiding this comment

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

could we name the script something like scil_tractogram_filter_by_dps? Feels like a filtering script.

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Ok, filtering is good. But it would be scil_tractogram_filter_by_dps_or_dpp. Too long, or still ok?

Copy link
Contributor

Choose a reason for hiding this comment

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

Not that long but it is a bit weird having "or" in a script name, filter_by_streamlines_properties? maybe

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.

3 participants