Note
Go to the end to download the full example code.
Precise voxel features
Background. A voxel feature is only useful for habitats if it gives nearly the same value when the same tumour is imaged again and when its computation settings change slightly. Precise screening measures this and keeps only the features that pass.
Purpose. You get ICC forest plots for three experiments on one demo subject, the list of kept and dropped features, and a side-by-side check of habitats clustered with all features vs only the precise ones (mean Dice and disagreement under the same simulated retest).
Key terms.
retest perturbation – a simulated repeat scan: Gaussian noise, a half-voxel shift and a 0.5 degree in-plane rotation applied to the image (Prior et al. Appendix S2); the original ROI is kept.
repeatability / reproducibility – agreement between original and retest image / between two settings (kernel radius 1 vs 3, bin width 12 vs 25).
ICC (intraclass correlation coefficient) – agreement of a feature between two measurements of the same voxels; near 1 means repeatable. HABIT uses ICC(3A,1) (absolute agreement) for repeatability and ICC(3C,1) (consistency) for reproducibility, after min-max scaling each feature map.
LCL – lower limit of the ICC’s 95% confidence interval; a feature is precise here when LCL >= 0.5 in all three experiments.
whitelist – the list of precise features; a preprocessing step that drops every other column before clustering.
Dice – overlap between two label maps (0 = none, 1 = identical), scored after matching ids (see Why habitat ids must be matched).
Decide which voxel features may define habitats, then cluster only those robust features. This is the Prior et al. precision screen (Radiol Artif Intell 2024;6(2):e230118; DOI).
Evaluating stability under perturbation
The scientific claim of precise features is stability under the same
simulated retest the method was designed for (Appendix S2: Gaussian
noise, sub-voxel translation, small in-plane rotation). This page clusters
habitats twice on one subject – once with all texture features, once
with the precise whitelist – and scores original vs perturbed maps
with habitat_stability() (mean Dice) plus the
voxel-wise disagreement panel of
plot_habitat_label_compare(). Read the printed numbers:
precise is only “more stable” on this demo when those scores improve.
An optional MONAI bspline_deform ROI edge perturbation is shown later
to inspect contour deformations; that is separate from the Appendix S2
retest used for the habitat stability comparison.
Load one demo subject. extract_voxel_texture crops to the ROI box
internally (crop_to_roi=True).
sphinx_gallery_thumbnail_number = 5
from pathlib import Path
from typing import Dict, List, Tuple
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from matplotlib.lines import Line2D
from matplotlib.patches import Patch
from habit.contracts import Cohort, Subject, cohort_from_directory
from habit.datasets import fetch_demo
from habit.kernels.habitat_label_match import adjusted_rand_index
from habit.kernels.image_perturbation import binary_mask_dice
from habit.precision import (
ImagePerturbationRegistry,
aggregate_panels,
align_habitat_map,
habitat_stability,
identify_precise_features,
perturb_image,
precision_panel,
)
from habit.recipes import Study
from habit.spec import HabitatSpec, Spec
from habit.voxel_features import extract_voxel_texture
from habit.viz import plot_habitat_label_compare, plot_intensity_slice, plot_precision_icc
from habit.viz import use_style
from habit.viz.labels import sanitize_label
DATA = fetch_demo()
MODALITIES = ("LAP",)
ROI = "LAP"
cohort = cohort_from_directory(DATA, modalities=MODALITIES, roi=ROI)[:1]
subject = cohort[0]
image = subject.image(MODALITIES[0])
mask = subject.mask(ROI)
Path("out").mkdir(exist_ok=True)
print(f"Grid shape: {image.data.shape}")
HABIT demo data (cached)
DATA (preprocessed root): C:\Users\dongm\.habit_data\demo-data-v1\preprocessed
On-disk inventory of this folder:
subjects (5): subj001, subj002, subj003, subj004, subj005
image series: LAP, PVP, delay_3min, pre_contrast
mask keys: LAP, PVP, delay_3min, pre_contrast
example image: images/subj001/delay_3min/WATER__BH_Ax_LAVA_Flex_3min_Series0012.nrrd
example mask: masks/subj001/delay_3min/WATER__BH_Ax_LAVA_Flex_10min_Series0017_mask.nrrd
Your own data must use the same folder tree (change IDs / series names):
DATA/
images/<subject_id>/<modality>/<one image file>
masks/<subject_id>/<roi>/<one mask file>
Then load it with the same call the demos use:
cohort = cohort_from_directory(DATA, modalities=("LAP",), roi="LAP")
Swap DATA / modalities / roi to match your tree. Mask key is often the
same as one image series (here LAP).
Grid shape: (200, 360, 360)
Appendix S2 retest chain on one shared RNG. Sequentially applies Gaussian noise -> translation -> rotation.
retest_rng = np.random.default_rng(7)
noisy = perturb_image(image, method="gaussian_noise", rng=retest_rng)
shifted = perturb_image(
noisy, method="translation", shift_fraction=0.5, rng=retest_rng
)
perturbed = perturb_image(
shifted, method="rotation", angle_degrees=0.5, rng=retest_rng
)
print("Appendix S2: gaussian_noise -> translation -> rotation")
fig_s2 = plot_intensity_slice(
perturbed,
before=image,
roi_mask=mask,
roi_contour=True,
title="Appendix S2 chain (original vs perturbed)",
before_label="Original",
image_label="+ noise / shift / rotation",
)
fig_s2.savefig("out/precise_features_perturb_methods.png", dpi=150, bbox_inches="tight")
plt.show()
Appendix S2: gaussian_noise -> translation -> rotation
F:\work\habit_project\habit\viz\intensity.py:616: UserWarning: Display geometry conflict: image/anatomy direction does not match mask/label direction. Using the mask/label direction so coronal/sagittal superior-up follows the labelled anatomy. Pass direction= to override.
resolved_direction, resolved_spacing = resolve_display_geometry(
Extract texture features at base R3/B12 and the two reproducibility contrasts.
FEATURE_CLASSES: Dict[str, Tuple[str, ...]] = {
"firstorder": ("Entropy", "Mean", "Variance", "Skewness", "Kurtosis"),
"glcm": (
"Contrast",
"Correlation",
"JointEntropy",
"Idm",
"DifferenceEntropy",
),
}
# Base setting R3/B12 is compared with R1 (kernel radius), B25 (bin width)
# and the perturbed image; every other setting is held fixed.
feat_r1 = extract_voxel_texture(
image, mask, kernel_radius=1, bin_width=12, feature_classes=FEATURE_CLASSES
)
feat_r3 = extract_voxel_texture(
image, mask, kernel_radius=3, bin_width=12, feature_classes=FEATURE_CLASSES
)
feat_b25 = extract_voxel_texture(
image, mask, kernel_radius=3, bin_width=25, feature_classes=FEATURE_CLASSES
)
feat_pert = extract_voxel_texture(
perturbed, mask, kernel_radius=3, bin_width=12, feature_classes=FEATURE_CLASSES
)
print(f"Texture features ({len(feat_r3.feature_names)}): {list(feat_r3.feature_names)}")
feat_r3.feature_frame().head()
voxel_radiomics: firstorder batch 1/35 0.203s peak=13MiB
voxel_radiomics: firstorder batch 2/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 3/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 4/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 5/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 6/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 7/35 0.005s peak=13MiB
voxel_radiomics: firstorder batch 8/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 9/35 0.005s peak=13MiB
voxel_radiomics: firstorder batch 10/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 11/35 0.003s peak=13MiB
voxel_radiomics: firstorder batch 12/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 13/35 0.005s peak=13MiB
voxel_radiomics: firstorder batch 14/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 15/35 0.005s peak=13MiB
voxel_radiomics: firstorder batch 16/35 0.005s peak=13MiB
voxel_radiomics: firstorder batch 17/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 18/35 0.003s peak=13MiB
voxel_radiomics: firstorder batch 19/35 0.005s peak=13MiB
voxel_radiomics: firstorder batch 20/35 0.005s peak=13MiB
voxel_radiomics: firstorder batch 21/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 22/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 23/35 0.005s peak=13MiB
voxel_radiomics: firstorder batch 24/35 0.003s peak=13MiB
voxel_radiomics: firstorder batch 25/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 26/35 0.003s peak=13MiB
voxel_radiomics: firstorder batch 27/35 0.006s peak=13MiB
voxel_radiomics: firstorder batch 28/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 29/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 30/35 0.005s peak=13MiB
voxel_radiomics: firstorder batch 31/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 32/35 0.005s peak=13MiB
voxel_radiomics: firstorder batch 33/35 0.005s peak=13MiB
voxel_radiomics: firstorder batch 34/35 0.004s peak=13MiB
voxel_radiomics: firstorder batch 35/35 0.004s peak=13MiB
voxel_radiomics: glcm batch 1/35 0.066s peak=12MiB
voxel_radiomics: glcm batch 2/35 0.023s peak=13MiB
voxel_radiomics: glcm batch 3/35 0.023s peak=13MiB
voxel_radiomics: glcm batch 4/35 0.022s peak=14MiB
voxel_radiomics: glcm batch 5/35 0.025s peak=15MiB
voxel_radiomics: glcm batch 6/35 0.027s peak=16MiB
voxel_radiomics: glcm batch 7/35 0.022s peak=17MiB
voxel_radiomics: glcm batch 8/35 0.023s peak=18MiB
voxel_radiomics: glcm batch 9/35 0.020s peak=19MiB
voxel_radiomics: glcm batch 10/35 0.020s peak=20MiB
voxel_radiomics: glcm batch 11/35 0.019s peak=21MiB
voxel_radiomics: glcm batch 12/35 0.021s peak=22MiB
voxel_radiomics: glcm batch 13/35 0.019s peak=22MiB
voxel_radiomics: glcm batch 14/35 0.020s peak=23MiB
voxel_radiomics: glcm batch 15/35 0.019s peak=24MiB
voxel_radiomics: glcm batch 16/35 0.026s peak=25MiB
voxel_radiomics: glcm batch 17/35 0.025s peak=26MiB
voxel_radiomics: glcm batch 18/35 0.024s peak=27MiB
voxel_radiomics: glcm batch 19/35 0.025s peak=28MiB
voxel_radiomics: glcm batch 20/35 0.019s peak=29MiB
voxel_radiomics: glcm batch 21/35 0.023s peak=30MiB
voxel_radiomics: glcm batch 22/35 0.022s peak=31MiB
voxel_radiomics: glcm batch 23/35 0.020s peak=31MiB
voxel_radiomics: glcm batch 24/35 0.019s peak=32MiB
voxel_radiomics: glcm batch 25/35 0.021s peak=33MiB
voxel_radiomics: glcm batch 26/35 0.021s peak=34MiB
voxel_radiomics: glcm batch 27/35 0.017s peak=35MiB
voxel_radiomics: glcm batch 28/35 0.019s peak=36MiB
voxel_radiomics: glcm batch 29/35 0.019s peak=37MiB
voxel_radiomics: glcm batch 30/35 0.018s peak=38MiB
voxel_radiomics: glcm batch 31/35 0.019s peak=39MiB
voxel_radiomics: glcm batch 32/35 0.018s peak=40MiB
voxel_radiomics: glcm batch 33/35 0.018s peak=40MiB
voxel_radiomics: glcm batch 34/35 0.020s peak=41MiB
voxel_radiomics: glcm batch 35/35 0.017s peak=42MiB
voxel_radiomics: firstorder batch 1/35 0.017s peak=48MiB
voxel_radiomics: firstorder batch 2/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 3/35 0.009s peak=48MiB
voxel_radiomics: firstorder batch 4/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 5/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 6/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 7/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 8/35 0.009s peak=48MiB
voxel_radiomics: firstorder batch 9/35 0.006s peak=48MiB
voxel_radiomics: firstorder batch 10/35 0.008s peak=48MiB
voxel_radiomics: firstorder batch 11/35 0.006s peak=48MiB
voxel_radiomics: firstorder batch 12/35 0.006s peak=48MiB
voxel_radiomics: firstorder batch 13/35 0.005s peak=48MiB
voxel_radiomics: firstorder batch 14/35 0.005s peak=48MiB
voxel_radiomics: firstorder batch 15/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 16/35 0.006s peak=48MiB
voxel_radiomics: firstorder batch 17/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 18/35 0.006s peak=48MiB
voxel_radiomics: firstorder batch 19/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 20/35 0.006s peak=48MiB
voxel_radiomics: firstorder batch 21/35 0.006s peak=48MiB
voxel_radiomics: firstorder batch 22/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 23/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 24/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 25/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 26/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 27/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 28/35 0.006s peak=48MiB
voxel_radiomics: firstorder batch 29/35 0.006s peak=48MiB
voxel_radiomics: firstorder batch 30/35 0.006s peak=48MiB
voxel_radiomics: firstorder batch 31/35 0.008s peak=48MiB
voxel_radiomics: firstorder batch 32/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 33/35 0.007s peak=48MiB
voxel_radiomics: firstorder batch 34/35 0.006s peak=48MiB
voxel_radiomics: firstorder batch 35/35 0.008s peak=47MiB
voxel_radiomics: glcm batch 1/35 0.064s peak=43MiB
voxel_radiomics: glcm batch 2/35 0.061s peak=45MiB
voxel_radiomics: glcm batch 3/35 0.061s peak=46MiB
voxel_radiomics: glcm batch 4/35 0.061s peak=47MiB
voxel_radiomics: glcm batch 5/35 0.063s peak=48MiB
voxel_radiomics: glcm batch 6/35 0.063s peak=49MiB
voxel_radiomics: glcm batch 7/35 0.066s peak=50MiB
voxel_radiomics: glcm batch 8/35 0.061s peak=51MiB
voxel_radiomics: glcm batch 9/35 0.059s peak=53MiB
voxel_radiomics: glcm batch 10/35 0.063s peak=54MiB
voxel_radiomics: glcm batch 11/35 0.062s peak=55MiB
voxel_radiomics: glcm batch 12/35 0.062s peak=56MiB
voxel_radiomics: glcm batch 13/35 0.062s peak=57MiB
voxel_radiomics: glcm batch 14/35 0.064s peak=58MiB
voxel_radiomics: glcm batch 15/35 0.064s peak=59MiB
voxel_radiomics: glcm batch 16/35 0.065s peak=61MiB
voxel_radiomics: glcm batch 17/35 0.066s peak=62MiB
voxel_radiomics: glcm batch 18/35 0.064s peak=1MiB
voxel_radiomics: glcm batch 19/35 0.063s peak=2MiB
voxel_radiomics: glcm batch 20/35 0.064s peak=3MiB
voxel_radiomics: glcm batch 21/35 0.064s peak=5MiB
voxel_radiomics: glcm batch 22/35 0.060s peak=6MiB
voxel_radiomics: glcm batch 23/35 0.059s peak=7MiB
voxel_radiomics: glcm batch 24/35 0.062s peak=8MiB
voxel_radiomics: glcm batch 25/35 0.060s peak=9MiB
voxel_radiomics: glcm batch 26/35 0.060s peak=10MiB
voxel_radiomics: glcm batch 27/35 0.061s peak=11MiB
voxel_radiomics: glcm batch 28/35 0.062s peak=13MiB
voxel_radiomics: glcm batch 29/35 0.060s peak=14MiB
voxel_radiomics: glcm batch 30/35 0.056s peak=15MiB
voxel_radiomics: glcm batch 31/35 0.055s peak=16MiB
voxel_radiomics: glcm batch 32/35 0.060s peak=17MiB
voxel_radiomics: glcm batch 33/35 0.055s peak=18MiB
voxel_radiomics: glcm batch 34/35 0.046s peak=20MiB
voxel_radiomics: glcm batch 35/35 0.032s peak=21MiB
voxel_radiomics: firstorder batch 1/35 0.010s peak=26MiB
voxel_radiomics: firstorder batch 2/35 0.005s peak=26MiB
voxel_radiomics: firstorder batch 3/35 0.005s peak=26MiB
voxel_radiomics: firstorder batch 4/35 0.005s peak=26MiB
voxel_radiomics: firstorder batch 5/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 6/35 0.005s peak=26MiB
voxel_radiomics: firstorder batch 7/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 8/35 0.005s peak=26MiB
voxel_radiomics: firstorder batch 9/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 10/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 11/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 12/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 13/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 14/35 0.005s peak=26MiB
voxel_radiomics: firstorder batch 15/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 16/35 0.005s peak=26MiB
voxel_radiomics: firstorder batch 17/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 18/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 19/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 20/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 21/35 0.007s peak=26MiB
voxel_radiomics: firstorder batch 22/35 0.007s peak=26MiB
voxel_radiomics: firstorder batch 23/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 24/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 25/35 0.005s peak=26MiB
voxel_radiomics: firstorder batch 26/35 0.005s peak=26MiB
voxel_radiomics: firstorder batch 27/35 0.005s peak=26MiB
voxel_radiomics: firstorder batch 28/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 29/35 0.006s peak=26MiB
voxel_radiomics: firstorder batch 30/35 0.005s peak=26MiB
voxel_radiomics: firstorder batch 31/35 0.005s peak=26MiB
voxel_radiomics: firstorder batch 32/35 0.005s peak=26MiB
voxel_radiomics: firstorder batch 33/35 0.005s peak=26MiB
voxel_radiomics: firstorder batch 34/35 0.005s peak=26MiB
voxel_radiomics: firstorder batch 35/35 0.004s peak=26MiB
voxel_radiomics: glcm batch 1/35 0.045s peak=22MiB
voxel_radiomics: glcm batch 2/35 0.042s peak=23MiB
voxel_radiomics: glcm batch 3/35 0.042s peak=24MiB
voxel_radiomics: glcm batch 4/35 0.041s peak=25MiB
voxel_radiomics: glcm batch 5/35 0.040s peak=26MiB
voxel_radiomics: glcm batch 6/35 0.040s peak=28MiB
voxel_radiomics: glcm batch 7/35 0.043s peak=29MiB
voxel_radiomics: glcm batch 8/35 0.041s peak=30MiB
voxel_radiomics: glcm batch 9/35 0.041s peak=31MiB
voxel_radiomics: glcm batch 10/35 0.042s peak=32MiB
voxel_radiomics: glcm batch 11/35 0.042s peak=33MiB
voxel_radiomics: glcm batch 12/35 0.040s peak=34MiB
voxel_radiomics: glcm batch 13/35 0.041s peak=36MiB
voxel_radiomics: glcm batch 14/35 0.043s peak=37MiB
voxel_radiomics: glcm batch 15/35 0.039s peak=38MiB
voxel_radiomics: glcm batch 16/35 0.042s peak=39MiB
voxel_radiomics: glcm batch 17/35 0.041s peak=40MiB
voxel_radiomics: glcm batch 18/35 0.045s peak=41MiB
voxel_radiomics: glcm batch 19/35 0.041s peak=43MiB
voxel_radiomics: glcm batch 20/35 0.045s peak=44MiB
voxel_radiomics: glcm batch 21/35 0.042s peak=45MiB
voxel_radiomics: glcm batch 22/35 0.040s peak=46MiB
voxel_radiomics: glcm batch 23/35 0.037s peak=47MiB
voxel_radiomics: glcm batch 24/35 0.039s peak=48MiB
voxel_radiomics: glcm batch 25/35 0.037s peak=49MiB
voxel_radiomics: glcm batch 26/35 0.036s peak=51MiB
voxel_radiomics: glcm batch 27/35 0.038s peak=52MiB
voxel_radiomics: glcm batch 28/35 0.037s peak=53MiB
voxel_radiomics: glcm batch 29/35 0.036s peak=54MiB
voxel_radiomics: glcm batch 30/35 0.034s peak=55MiB
voxel_radiomics: glcm batch 31/35 0.034s peak=56MiB
voxel_radiomics: glcm batch 32/35 0.034s peak=57MiB
voxel_radiomics: glcm batch 33/35 0.033s peak=59MiB
voxel_radiomics: glcm batch 34/35 0.031s peak=60MiB
voxel_radiomics: glcm batch 35/35 0.022s peak=61MiB
voxel_radiomics: firstorder batch 1/35 0.014s peak=67MiB
voxel_radiomics: firstorder batch 2/35 0.005s peak=67MiB
voxel_radiomics: firstorder batch 3/35 0.008s peak=67MiB
voxel_radiomics: firstorder batch 4/35 0.005s peak=67MiB
voxel_radiomics: firstorder batch 5/35 0.007s peak=67MiB
voxel_radiomics: firstorder batch 6/35 0.006s peak=67MiB
voxel_radiomics: firstorder batch 7/35 0.006s peak=67MiB
voxel_radiomics: firstorder batch 8/35 0.007s peak=67MiB
voxel_radiomics: firstorder batch 9/35 0.007s peak=67MiB
voxel_radiomics: firstorder batch 10/35 0.007s peak=67MiB
voxel_radiomics: firstorder batch 11/35 0.006s peak=67MiB
voxel_radiomics: firstorder batch 12/35 0.007s peak=67MiB
voxel_radiomics: firstorder batch 13/35 0.007s peak=67MiB
voxel_radiomics: firstorder batch 14/35 0.007s peak=67MiB
voxel_radiomics: firstorder batch 15/35 0.006s peak=67MiB
voxel_radiomics: firstorder batch 16/35 0.006s peak=67MiB
voxel_radiomics: firstorder batch 17/35 0.007s peak=67MiB
voxel_radiomics: firstorder batch 18/35 0.006s peak=67MiB
voxel_radiomics: firstorder batch 19/35 0.006s peak=67MiB
voxel_radiomics: firstorder batch 20/35 0.006s peak=67MiB
voxel_radiomics: firstorder batch 21/35 0.005s peak=67MiB
voxel_radiomics: firstorder batch 22/35 0.007s peak=67MiB
voxel_radiomics: firstorder batch 23/35 0.005s peak=67MiB
voxel_radiomics: firstorder batch 24/35 0.007s peak=67MiB
voxel_radiomics: firstorder batch 25/35 0.006s peak=67MiB
voxel_radiomics: firstorder batch 26/35 0.007s peak=67MiB
voxel_radiomics: firstorder batch 27/35 0.006s peak=67MiB
voxel_radiomics: firstorder batch 28/35 0.006s peak=67MiB
voxel_radiomics: firstorder batch 29/35 0.006s peak=67MiB
voxel_radiomics: firstorder batch 30/35 0.005s peak=67MiB
voxel_radiomics: firstorder batch 31/35 0.006s peak=67MiB
voxel_radiomics: firstorder batch 32/35 0.005s peak=67MiB
voxel_radiomics: firstorder batch 33/35 0.005s peak=67MiB
voxel_radiomics: firstorder batch 34/35 0.006s peak=67MiB
voxel_radiomics: firstorder batch 35/35 0.005s peak=66MiB
voxel_radiomics: glcm batch 1/35 0.056s peak=62MiB
voxel_radiomics: glcm batch 2/35 0.056s peak=63MiB
voxel_radiomics: glcm batch 3/35 0.058s peak=64MiB
voxel_radiomics: glcm batch 4/35 0.058s peak=66MiB
voxel_radiomics: glcm batch 5/35 0.063s peak=67MiB
voxel_radiomics: glcm batch 6/35 0.062s peak=68MiB
voxel_radiomics: glcm batch 7/35 0.060s peak=69MiB
voxel_radiomics: glcm batch 8/35 0.059s peak=70MiB
voxel_radiomics: glcm batch 9/35 0.060s peak=71MiB
voxel_radiomics: glcm batch 10/35 0.062s peak=72MiB
voxel_radiomics: glcm batch 11/35 0.062s peak=74MiB
voxel_radiomics: glcm batch 12/35 0.065s peak=1MiB
voxel_radiomics: glcm batch 13/35 0.061s peak=2MiB
voxel_radiomics: glcm batch 14/35 0.061s peak=3MiB
voxel_radiomics: glcm batch 15/35 0.063s peak=5MiB
voxel_radiomics: glcm batch 16/35 0.060s peak=6MiB
voxel_radiomics: glcm batch 17/35 0.060s peak=7MiB
voxel_radiomics: glcm batch 18/35 0.063s peak=8MiB
voxel_radiomics: glcm batch 19/35 0.062s peak=9MiB
voxel_radiomics: glcm batch 20/35 0.064s peak=10MiB
voxel_radiomics: glcm batch 21/35 0.060s peak=11MiB
voxel_radiomics: glcm batch 22/35 0.061s peak=13MiB
voxel_radiomics: glcm batch 23/35 0.059s peak=14MiB
voxel_radiomics: glcm batch 24/35 0.059s peak=15MiB
voxel_radiomics: glcm batch 25/35 0.056s peak=16MiB
voxel_radiomics: glcm batch 26/35 0.060s peak=17MiB
voxel_radiomics: glcm batch 27/35 0.059s peak=18MiB
voxel_radiomics: glcm batch 28/35 0.057s peak=20MiB
voxel_radiomics: glcm batch 29/35 0.057s peak=21MiB
voxel_radiomics: glcm batch 30/35 0.055s peak=22MiB
voxel_radiomics: glcm batch 31/35 0.053s peak=23MiB
voxel_radiomics: glcm batch 32/35 0.057s peak=24MiB
voxel_radiomics: glcm batch 33/35 0.054s peak=25MiB
voxel_radiomics: glcm batch 34/35 0.050s peak=26MiB
voxel_radiomics: glcm batch 35/35 0.030s peak=28MiB
Texture features (10): ['original_firstorder_Entropy-LAP', 'original_firstorder_Mean-LAP', 'original_firstorder_Variance-LAP', 'original_firstorder_Skewness-LAP', 'original_firstorder_Kurtosis-LAP', 'original_glcm_Contrast-LAP', 'original_glcm_Correlation-LAP', 'original_glcm_JointEntropy-LAP', 'original_glcm_Idm-LAP', 'original_glcm_DifferenceEntropy-LAP']
Precise screening = Lower Confidence Limit (LCL) >= 0.5 across all 3 ICC experiments.
Each precision_panel is one subject’s per-feature ICC table;
aggregate_panels takes the per-feature median over subjects (here a
single subject, so the median is that subject’s value).
precise = identify_precise_features(
{
"repeatability": aggregate_panels(
[precision_panel({"original": feat_r3, "perturbed": feat_pert}, agreement="absolute")]
),
"reproducibility_kernel_radius": aggregate_panels(
[precision_panel({"R1": feat_r1, "R3": feat_r3}, agreement="consistency")]
),
"reproducibility_bin_width": aggregate_panels(
[precision_panel({"B12": feat_r3, "B25": feat_b25}, agreement="consistency")]
),
},
lcl_threshold=0.5,
)
evidence = precise.to_frame().round(3)
kept: List[str] = list(precise.feature_names)
dropped = [n for n in feat_r3.feature_names if n not in set(kept)]
print(f"Kept features ({len(kept)}): {kept}")
print(f"Dropped features ({len(dropped)}): {dropped}")
evidence
Kept features (7): ['original_firstorder_Entropy-LAP', 'original_firstorder_Mean-LAP', 'original_firstorder_Variance-LAP', 'original_glcm_Contrast-LAP', 'original_glcm_JointEntropy-LAP', 'original_glcm_Idm-LAP', 'original_glcm_DifferenceEntropy-LAP']
Dropped features (3): ['original_firstorder_Skewness-LAP', 'original_firstorder_Kurtosis-LAP', 'original_glcm_Correlation-LAP']
Plot one ICC forest per experiment to inspect lower confidence limits.
for experiment, fname, title in (
("repeatability", "precise_features_icc_lcl.png", "Repeatability ICC"),
(
"reproducibility_kernel_radius",
"precise_features_icc_kernel.png",
"Kernel-radius reproducibility ICC",
),
(
"reproducibility_bin_width",
"precise_features_icc_bin.png",
"Bin-width reproducibility ICC",
),
):
panel = evidence.loc[evidence["experiment"] == experiment].drop(
columns=["precise"], errors="ignore"
)
fig_icc = plot_precision_icc(
panel.dropna(subset=["value", "lcl", "ucl"]),
lcl_threshold=precise.lcl_threshold,
title=title,
orientation="row",
)
fig_icc.savefig(f"out/{fname}", dpi=150, bbox_inches="tight")
plt.show()
Same subject, same Appendix S2 perturbation, same k / seed: only the feature set changes. Without precise, habitats can shift under the retest; with the precise whitelist they should agree more – check the printed mean Dice and labelled-voxel disagreement before claiming that.
texture_params = {
"imageType": {"Original": {}},
"featureClass": {k: list(v) for k, v in FEATURE_CLASSES.items()},
"setting": {"binWidth": 12.0, "normalize": False},
}
extractor_spec = Spec(
"voxel_radiomics",
{"modalities": list(MODALITIES), "kernel_radius": 3, "params": texture_params},
)
# Same fitter on both arms (n_habitats / n_init) so only the feature set differs.
fitter_spec = Spec(
"kmeans",
{"n_habitats": 3, "n_init": 3},
)
minmax_spec = Spec("minmax", {"across_features": False})
subject_pert = Subject(
subject_id=subject.subject_id,
images={MODALITIES[0]: perturbed},
masks=subject.masks,
)
demo = Cohort(subjects=(subject,))
demo_pert = Cohort(subjects=(subject_pert,))
def _labelled_disagreement(reference_map, aligned_map) -> float:
"""Fraction of labelled voxels whose ids differ after overlap matching."""
ref = np.asarray(reference_map.label_array)
mov = np.asarray(aligned_map.label_array)
labelled = (ref > 0) | (mov > 0)
if not np.any(labelled):
return float("nan")
return float(np.mean(ref[labelled] != mov[labelled]))
# --- Without precise: all texture features ---
# No whitelist preprocessor -- every extracted texture column enters clustering.
spec_all = HabitatSpec(
name="all_texture_one_step",
voxel_feature_extractor=extractor_spec,
voxel_feature_preprocessors=(minmax_spec,),
habitat_model_fitter=fitter_spec,
habitat_assigner=Spec("nearest_centroid"),
random_seed=11,
pooling="none",
)
result_all_orig = Study(spec_all).fit_predict(demo)
result_all_pert = Study(spec_all).fit_predict(demo_pert)
map_all_orig = result_all_orig.habitat_maps[0]
map_all_pert = result_all_pert.habitat_maps[0]
# habitat_stability pairs ids by voxel overlap (Prior Hungarian step) then Dice.
stab_all = habitat_stability(map_all_orig, [map_all_pert])
mean_dice_all = float(stab_all["dice"].mean())
ari_all = float(
adjusted_rand_index(
np.asarray(map_all_orig.label_array),
np.asarray(map_all_pert.label_array),
)
)
# force=True: independent one_step runs share model_id by subject+spec, not image content.
aligned_all = align_habitat_map(map_all_orig, map_all_pert, force=True)
disagree_all = _labelled_disagreement(map_all_orig, aligned_all)
print(
"Without precise (all texture features): "
f"mean Dice={mean_dice_all:.3f}, "
f"labelled-voxel disagreement={disagree_all:.3f}, "
f"ARI={ari_all:.3f}"
)
# Original vs perturbed habitats + disagreement panel (same anatomy image).
fig_cmp_all = plot_habitat_label_compare(
image,
map_all_orig,
aligned_all,
titles=(
"Without precise: original",
f"Without precise: perturbed (mean Dice={mean_dice_all:.3f})",
),
align_labels=False,
show_disagreement=True,
crop_to="labels",
)
fig_cmp_all.savefig("out/precise_features_all_orig_vs_pert.png", dpi=150, bbox_inches="tight")
plt.show()
F:\work\habit_project\examples\guide\08_precision\plot_01_precise_features.py:247: HabitDeprecationWarning: HabitatSpec named-field constructor is deprecated since version 2.0.0 and will be removed in version 4.0.0. Use HabitatSpec(..., stages=(Stage(...), ...)) instead. Named-field YAML / from_dict payloads still load and keep their historical fingerprints.
spec_all = HabitatSpec(
Cohort.map[_DefineAndLabelWithinSubject]: 0%| | 0/1 [00:00<?, ?it/s]voxel_radiomics: firstorder batch 1/35 0.029s peak=34MiB
voxel_radiomics: firstorder batch 2/35 0.013s peak=34MiB
voxel_radiomics: firstorder batch 3/35 0.033s peak=34MiB
voxel_radiomics: firstorder batch 4/35 0.022s peak=34MiB
voxel_radiomics: firstorder batch 5/35 0.013s peak=34MiB
voxel_radiomics: firstorder batch 6/35 0.015s peak=34MiB
voxel_radiomics: firstorder batch 7/35 0.016s peak=34MiB
voxel_radiomics: firstorder batch 8/35 0.016s peak=34MiB
voxel_radiomics: firstorder batch 9/35 0.014s peak=34MiB
voxel_radiomics: firstorder batch 10/35 0.012s peak=34MiB
voxel_radiomics: firstorder batch 11/35 0.014s peak=34MiB
voxel_radiomics: firstorder batch 12/35 0.011s peak=34MiB
voxel_radiomics: firstorder batch 13/35 0.015s peak=34MiB
voxel_radiomics: firstorder batch 14/35 0.013s peak=34MiB
voxel_radiomics: firstorder batch 15/35 0.021s peak=34MiB
voxel_radiomics: firstorder batch 16/35 0.016s peak=34MiB
voxel_radiomics: firstorder batch 17/35 0.016s peak=34MiB
voxel_radiomics: firstorder batch 18/35 0.011s peak=34MiB
voxel_radiomics: firstorder batch 19/35 0.007s peak=34MiB
voxel_radiomics: firstorder batch 20/35 0.007s peak=34MiB
voxel_radiomics: firstorder batch 21/35 0.008s peak=34MiB
voxel_radiomics: firstorder batch 22/35 0.010s peak=34MiB
voxel_radiomics: firstorder batch 23/35 0.007s peak=34MiB
voxel_radiomics: firstorder batch 24/35 0.009s peak=34MiB
voxel_radiomics: firstorder batch 25/35 0.007s peak=34MiB
voxel_radiomics: firstorder batch 26/35 0.008s peak=34MiB
voxel_radiomics: firstorder batch 27/35 0.006s peak=34MiB
voxel_radiomics: firstorder batch 28/35 0.006s peak=34MiB
voxel_radiomics: firstorder batch 29/35 0.008s peak=34MiB
voxel_radiomics: firstorder batch 30/35 0.007s peak=34MiB
voxel_radiomics: firstorder batch 31/35 0.006s peak=34MiB
voxel_radiomics: firstorder batch 32/35 0.008s peak=34MiB
voxel_radiomics: firstorder batch 33/35 0.007s peak=34MiB
voxel_radiomics: firstorder batch 34/35 0.008s peak=34MiB
voxel_radiomics: firstorder batch 35/35 0.006s peak=33MiB
voxel_radiomics: glcm batch 1/35 0.065s peak=29MiB
voxel_radiomics: glcm batch 2/35 0.057s peak=30MiB
voxel_radiomics: glcm batch 3/35 0.059s peak=31MiB
voxel_radiomics: glcm batch 4/35 0.059s peak=32MiB
voxel_radiomics: glcm batch 5/35 0.062s peak=33MiB
voxel_radiomics: glcm batch 6/35 0.061s peak=34MiB
voxel_radiomics: glcm batch 7/35 0.062s peak=36MiB
voxel_radiomics: glcm batch 8/35 0.060s peak=37MiB
voxel_radiomics: glcm batch 9/35 0.062s peak=38MiB
voxel_radiomics: glcm batch 10/35 0.062s peak=39MiB
voxel_radiomics: glcm batch 11/35 0.062s peak=40MiB
voxel_radiomics: glcm batch 12/35 0.063s peak=41MiB
voxel_radiomics: glcm batch 13/35 0.060s peak=43MiB
voxel_radiomics: glcm batch 14/35 0.066s peak=44MiB
voxel_radiomics: glcm batch 15/35 0.070s peak=45MiB
voxel_radiomics: glcm batch 16/35 0.061s peak=46MiB
voxel_radiomics: glcm batch 17/35 0.062s peak=47MiB
voxel_radiomics: glcm batch 18/35 0.064s peak=48MiB
voxel_radiomics: glcm batch 19/35 0.061s peak=49MiB
voxel_radiomics: glcm batch 20/35 0.060s peak=51MiB
voxel_radiomics: glcm batch 21/35 0.059s peak=52MiB
voxel_radiomics: glcm batch 22/35 0.060s peak=53MiB
voxel_radiomics: glcm batch 23/35 0.058s peak=54MiB
voxel_radiomics: glcm batch 24/35 0.059s peak=55MiB
voxel_radiomics: glcm batch 25/35 0.061s peak=56MiB
voxel_radiomics: glcm batch 26/35 0.058s peak=57MiB
voxel_radiomics: glcm batch 27/35 0.060s peak=59MiB
voxel_radiomics: glcm batch 28/35 0.059s peak=60MiB
voxel_radiomics: glcm batch 29/35 0.059s peak=61MiB
voxel_radiomics: glcm batch 30/35 0.055s peak=62MiB
voxel_radiomics: glcm batch 31/35 0.053s peak=63MiB
voxel_radiomics: glcm batch 32/35 0.055s peak=64MiB
voxel_radiomics: glcm batch 33/35 0.054s peak=66MiB
voxel_radiomics: glcm batch 34/35 0.047s peak=67MiB
voxel_radiomics: glcm batch 35/35 0.034s peak=68MiB
Cohort.map[_DefineAndLabelWithinSubject]: 100%|██████████| 1/1 [00:04<00:00, 4.49s/it]
Cohort.map[_DefineAndLabelWithinSubject]: 100%|██████████| 1/1 [00:04<00:00, 4.49s/it]
Cohort.map[_DefineAndLabelWithinSubject]: 0%| | 0/1 [00:00<?, ?it/s]voxel_radiomics: firstorder batch 1/35 0.029s peak=74MiB
voxel_radiomics: firstorder batch 2/35 0.014s peak=74MiB
voxel_radiomics: firstorder batch 3/35 0.015s peak=74MiB
voxel_radiomics: firstorder batch 4/35 0.014s peak=74MiB
voxel_radiomics: firstorder batch 5/35 0.015s peak=74MiB
voxel_radiomics: firstorder batch 6/35 0.011s peak=74MiB
voxel_radiomics: firstorder batch 7/35 0.014s peak=74MiB
voxel_radiomics: firstorder batch 8/35 0.015s peak=74MiB
voxel_radiomics: firstorder batch 9/35 0.015s peak=74MiB
voxel_radiomics: firstorder batch 10/35 0.012s peak=74MiB
voxel_radiomics: firstorder batch 11/35 0.013s peak=74MiB
voxel_radiomics: firstorder batch 12/35 0.015s peak=74MiB
voxel_radiomics: firstorder batch 13/35 0.013s peak=74MiB
voxel_radiomics: firstorder batch 14/35 0.014s peak=74MiB
voxel_radiomics: firstorder batch 15/35 0.014s peak=74MiB
voxel_radiomics: firstorder batch 16/35 0.011s peak=74MiB
voxel_radiomics: firstorder batch 17/35 0.014s peak=74MiB
voxel_radiomics: firstorder batch 18/35 0.015s peak=74MiB
voxel_radiomics: firstorder batch 19/35 0.008s peak=74MiB
voxel_radiomics: firstorder batch 20/35 0.007s peak=74MiB
voxel_radiomics: firstorder batch 21/35 0.008s peak=74MiB
voxel_radiomics: firstorder batch 22/35 0.006s peak=74MiB
voxel_radiomics: firstorder batch 23/35 0.005s peak=74MiB
voxel_radiomics: firstorder batch 24/35 0.006s peak=74MiB
voxel_radiomics: firstorder batch 25/35 0.007s peak=74MiB
voxel_radiomics: firstorder batch 26/35 0.008s peak=74MiB
voxel_radiomics: firstorder batch 27/35 0.005s peak=74MiB
voxel_radiomics: firstorder batch 28/35 0.006s peak=74MiB
voxel_radiomics: firstorder batch 29/35 0.006s peak=74MiB
voxel_radiomics: firstorder batch 30/35 0.009s peak=74MiB
voxel_radiomics: firstorder batch 31/35 0.005s peak=74MiB
voxel_radiomics: firstorder batch 32/35 0.007s peak=74MiB
voxel_radiomics: firstorder batch 33/35 0.006s peak=74MiB
voxel_radiomics: firstorder batch 34/35 0.007s peak=74MiB
voxel_radiomics: firstorder batch 35/35 0.005s peak=73MiB
voxel_radiomics: glcm batch 1/35 0.070s peak=69MiB
voxel_radiomics: glcm batch 2/35 0.057s peak=70MiB
voxel_radiomics: glcm batch 3/35 0.059s peak=71MiB
voxel_radiomics: glcm batch 4/35 0.058s peak=72MiB
voxel_radiomics: glcm batch 5/35 0.059s peak=74MiB
voxel_radiomics: glcm batch 6/35 0.060s peak=1MiB
voxel_radiomics: glcm batch 7/35 0.062s peak=2MiB
voxel_radiomics: glcm batch 8/35 0.060s peak=3MiB
voxel_radiomics: glcm batch 9/35 0.060s peak=5MiB
voxel_radiomics: glcm batch 10/35 0.059s peak=6MiB
voxel_radiomics: glcm batch 11/35 0.060s peak=7MiB
voxel_radiomics: glcm batch 12/35 0.058s peak=8MiB
voxel_radiomics: glcm batch 13/35 0.059s peak=9MiB
voxel_radiomics: glcm batch 14/35 0.063s peak=10MiB
voxel_radiomics: glcm batch 15/35 0.060s peak=11MiB
voxel_radiomics: glcm batch 16/35 0.063s peak=13MiB
voxel_radiomics: glcm batch 17/35 0.060s peak=14MiB
voxel_radiomics: glcm batch 18/35 0.063s peak=15MiB
voxel_radiomics: glcm batch 19/35 0.060s peak=16MiB
voxel_radiomics: glcm batch 20/35 0.063s peak=17MiB
voxel_radiomics: glcm batch 21/35 0.059s peak=18MiB
voxel_radiomics: glcm batch 22/35 0.059s peak=20MiB
voxel_radiomics: glcm batch 23/35 0.060s peak=21MiB
voxel_radiomics: glcm batch 24/35 0.061s peak=22MiB
voxel_radiomics: glcm batch 25/35 0.061s peak=23MiB
voxel_radiomics: glcm batch 26/35 0.056s peak=24MiB
voxel_radiomics: glcm batch 27/35 0.062s peak=25MiB
voxel_radiomics: glcm batch 28/35 0.060s peak=26MiB
voxel_radiomics: glcm batch 29/35 0.056s peak=28MiB
voxel_radiomics: glcm batch 30/35 0.055s peak=29MiB
voxel_radiomics: glcm batch 31/35 0.053s peak=30MiB
voxel_radiomics: glcm batch 32/35 0.056s peak=31MiB
voxel_radiomics: glcm batch 33/35 0.053s peak=32MiB
voxel_radiomics: glcm batch 34/35 0.047s peak=33MiB
voxel_radiomics: glcm batch 35/35 0.032s peak=34MiB
Cohort.map[_DefineAndLabelWithinSubject]: 100%|██████████| 1/1 [00:04<00:00, 4.34s/it]
Cohort.map[_DefineAndLabelWithinSubject]: 100%|██████████| 1/1 [00:04<00:00, 4.34s/it]
Without precise (all texture features): mean Dice=0.787, labelled-voxel disagreement=0.210, ARI=0.442
F:\work\habit_project\habit\viz\habitat_core.py:948: UserWarning: Display geometry conflict: image/anatomy direction does not match mask/label direction. Using the mask/label direction so coronal/sagittal superior-up follows the labelled anatomy. Pass direction= to override.
resolved_direction, resolved_spacing = resolve_display_geometry(
— With precise: whitelist only (precise on) —
if not kept:
print("No feature passed every experiment; skip precise habitats")
mean_dice_p = float("nan")
disagree_p = float("nan")
ari_p = float("nan")
else:
# precise.on: insert the ICC whitelist before minmax; fitter/seed unchanged.
whitelist = precise.preprocessor()
spec_precise = HabitatSpec(
name="precise_one_step",
voxel_feature_extractor=extractor_spec,
voxel_feature_preprocessors=(whitelist.spec, minmax_spec),
habitat_model_fitter=fitter_spec,
habitat_assigner=Spec("nearest_centroid"),
random_seed=11,
pooling="none",
)
result_precise_orig = Study(spec_precise).fit_predict(demo)
result_precise_pert = Study(spec_precise).fit_predict(demo_pert)
map_p_orig = result_precise_orig.habitat_maps[0]
map_p_pert = result_precise_pert.habitat_maps[0]
stab_p = habitat_stability(map_p_orig, [map_p_pert])
mean_dice_p = float(stab_p["dice"].mean())
ari_p = float(
adjusted_rand_index(
np.asarray(map_p_orig.label_array),
np.asarray(map_p_pert.label_array),
)
)
aligned_p = align_habitat_map(map_p_orig, map_p_pert, force=True)
disagree_p = _labelled_disagreement(map_p_orig, aligned_p)
print(
"With precise (whitelist only): "
f"mean Dice={mean_dice_p:.3f}, "
f"labelled-voxel disagreement={disagree_p:.3f}, "
f"ARI={ari_p:.3f}"
)
fig_cmp_p = plot_habitat_label_compare(
image,
map_p_orig,
aligned_p,
titles=(
"With precise: original",
f"With precise: perturbed (mean Dice={mean_dice_p:.3f})",
),
align_labels=False,
show_disagreement=True,
crop_to="labels",
)
fig_cmp_p.savefig(
"out/precise_features_precise_orig_vs_pert.png", dpi=150, bbox_inches="tight"
)
plt.show()
stability = pd.DataFrame(
[
{
"feature_set": "Without precise (all texture)",
"mean_dice": mean_dice_all,
"disagreement": disagree_all,
"ari": ari_all,
},
{
"feature_set": "With precise (whitelist)",
"mean_dice": mean_dice_p,
"disagreement": disagree_p,
"ari": ari_p,
},
]
)
print("Stability under Appendix S2 perturbation (original vs perturbed):")
print(stability.to_string(index=False))
if mean_dice_p > mean_dice_all and disagree_p < disagree_all:
print(
"On this demo subject, precise lowers disagreement "
f"({disagree_all:.3f} -> {disagree_p:.3f}) and raises mean Dice "
f"({mean_dice_all:.3f} -> {mean_dice_p:.3f})."
)
else:
print(
"On this demo subject, precise did not improve both scores; "
"report the measured mean Dice and disagreement as printed above."
)
with use_style("radiology"):
fig_stab, ax_s = plt.subplots(figsize=(5.5, 3.8), constrained_layout=True)
x_indices = np.arange(2)
bar_width = 0.35
dices = [mean_dice_all, mean_dice_p]
disagrees = [disagree_all, disagree_p]
ax_s.bar(
x_indices - bar_width / 2,
dices,
bar_width,
label="Mean Dice",
color="#0072B2",
)
ax_s.bar(
x_indices + bar_width / 2,
disagrees,
bar_width,
label="Labelled-voxel disagreement",
color="#E69F00",
)
ax_s.set_xticks(x_indices)
ax_s.set_xticklabels(["Without precise", "With precise"])
ax_s.set_ylim(0.0, 1.05)
ax_s.set_ylabel("Score")
ax_s.set_title(sanitize_label("Habitat stability under Appendix S2 perturbation"))
ax_s.legend(loc="best", frameon=True)
fig_stab.savefig("out/precise_features_stability_bar.png", dpi=150, bbox_inches="tight")
plt.show()
stability
F:\work\habit_project\examples\guide\08_precision\plot_01_precise_features.py:305: HabitDeprecationWarning: HabitatSpec named-field constructor is deprecated since version 2.0.0 and will be removed in version 4.0.0. Use HabitatSpec(..., stages=(Stage(...), ...)) instead. Named-field YAML / from_dict payloads still load and keep their historical fingerprints.
spec_precise = HabitatSpec(
Cohort.map[_DefineAndLabelWithinSubject]: 0%| | 0/1 [00:00<?, ?it/s]voxel_radiomics: firstorder batch 1/35 0.033s peak=41MiB
voxel_radiomics: firstorder batch 2/35 0.016s peak=41MiB
voxel_radiomics: firstorder batch 3/35 0.018s peak=41MiB
voxel_radiomics: firstorder batch 4/35 0.024s peak=41MiB
voxel_radiomics: firstorder batch 5/35 0.015s peak=41MiB
voxel_radiomics: firstorder batch 6/35 0.018s peak=41MiB
voxel_radiomics: firstorder batch 7/35 0.018s peak=41MiB
voxel_radiomics: firstorder batch 8/35 0.019s peak=41MiB
voxel_radiomics: firstorder batch 9/35 0.009s peak=41MiB
voxel_radiomics: firstorder batch 10/35 0.012s peak=41MiB
voxel_radiomics: firstorder batch 11/35 0.009s peak=41MiB
voxel_radiomics: firstorder batch 12/35 0.012s peak=41MiB
voxel_radiomics: firstorder batch 13/35 0.014s peak=41MiB
voxel_radiomics: firstorder batch 14/35 0.013s peak=41MiB
voxel_radiomics: firstorder batch 15/35 0.012s peak=41MiB
voxel_radiomics: firstorder batch 16/35 0.008s peak=41MiB
voxel_radiomics: firstorder batch 17/35 0.013s peak=41MiB
voxel_radiomics: firstorder batch 18/35 0.010s peak=41MiB
voxel_radiomics: firstorder batch 19/35 0.009s peak=41MiB
voxel_radiomics: firstorder batch 20/35 0.009s peak=41MiB
voxel_radiomics: firstorder batch 21/35 0.012s peak=41MiB
voxel_radiomics: firstorder batch 22/35 0.013s peak=41MiB
voxel_radiomics: firstorder batch 23/35 0.011s peak=41MiB
voxel_radiomics: firstorder batch 24/35 0.009s peak=41MiB
voxel_radiomics: firstorder batch 25/35 0.010s peak=41MiB
voxel_radiomics: firstorder batch 26/35 0.006s peak=41MiB
voxel_radiomics: firstorder batch 27/35 0.007s peak=41MiB
voxel_radiomics: firstorder batch 28/35 0.006s peak=41MiB
voxel_radiomics: firstorder batch 29/35 0.007s peak=41MiB
voxel_radiomics: firstorder batch 30/35 0.007s peak=41MiB
voxel_radiomics: firstorder batch 31/35 0.007s peak=41MiB
voxel_radiomics: firstorder batch 32/35 0.008s peak=41MiB
voxel_radiomics: firstorder batch 33/35 0.006s peak=41MiB
voxel_radiomics: firstorder batch 34/35 0.008s peak=41MiB
voxel_radiomics: firstorder batch 35/35 0.005s peak=40MiB
voxel_radiomics: glcm batch 1/35 0.068s peak=36MiB
voxel_radiomics: glcm batch 2/35 0.064s peak=37MiB
voxel_radiomics: glcm batch 3/35 0.059s peak=38MiB
voxel_radiomics: glcm batch 4/35 0.059s peak=39MiB
voxel_radiomics: glcm batch 5/35 0.061s peak=40MiB
voxel_radiomics: glcm batch 6/35 0.059s peak=41MiB
voxel_radiomics: glcm batch 7/35 0.061s peak=43MiB
voxel_radiomics: glcm batch 8/35 0.062s peak=44MiB
voxel_radiomics: glcm batch 9/35 0.062s peak=45MiB
voxel_radiomics: glcm batch 10/35 0.064s peak=46MiB
voxel_radiomics: glcm batch 11/35 0.062s peak=47MiB
voxel_radiomics: glcm batch 12/35 0.062s peak=48MiB
voxel_radiomics: glcm batch 13/35 0.061s peak=49MiB
voxel_radiomics: glcm batch 14/35 0.062s peak=51MiB
voxel_radiomics: glcm batch 15/35 0.065s peak=52MiB
voxel_radiomics: glcm batch 16/35 0.063s peak=53MiB
voxel_radiomics: glcm batch 17/35 0.063s peak=54MiB
voxel_radiomics: glcm batch 18/35 0.061s peak=55MiB
voxel_radiomics: glcm batch 19/35 0.060s peak=56MiB
voxel_radiomics: glcm batch 20/35 0.060s peak=57MiB
voxel_radiomics: glcm batch 21/35 0.060s peak=59MiB
voxel_radiomics: glcm batch 22/35 0.061s peak=60MiB
voxel_radiomics: glcm batch 23/35 0.059s peak=61MiB
voxel_radiomics: glcm batch 24/35 0.058s peak=62MiB
voxel_radiomics: glcm batch 25/35 0.059s peak=63MiB
voxel_radiomics: glcm batch 26/35 0.060s peak=64MiB
voxel_radiomics: glcm batch 27/35 0.061s peak=66MiB
voxel_radiomics: glcm batch 28/35 0.060s peak=67MiB
voxel_radiomics: glcm batch 29/35 0.058s peak=68MiB
voxel_radiomics: glcm batch 30/35 0.056s peak=69MiB
voxel_radiomics: glcm batch 31/35 0.052s peak=70MiB
voxel_radiomics: glcm batch 32/35 0.055s peak=71MiB
voxel_radiomics: glcm batch 33/35 0.051s peak=72MiB
voxel_radiomics: glcm batch 34/35 0.046s peak=74MiB
voxel_radiomics: glcm batch 35/35 0.032s peak=1MiB
Cohort.map[_DefineAndLabelWithinSubject]: 100%|██████████| 1/1 [00:04<00:00, 4.46s/it]
Cohort.map[_DefineAndLabelWithinSubject]: 100%|██████████| 1/1 [00:04<00:00, 4.46s/it]
Cohort.map[_DefineAndLabelWithinSubject]: 0%| | 0/1 [00:00<?, ?it/s]voxel_radiomics: firstorder batch 1/35 0.028s peak=7MiB
voxel_radiomics: firstorder batch 2/35 0.011s peak=7MiB
voxel_radiomics: firstorder batch 3/35 0.014s peak=7MiB
voxel_radiomics: firstorder batch 4/35 0.020s peak=7MiB
voxel_radiomics: firstorder batch 5/35 0.019s peak=7MiB
voxel_radiomics: firstorder batch 6/35 0.015s peak=7MiB
voxel_radiomics: firstorder batch 7/35 0.016s peak=7MiB
voxel_radiomics: firstorder batch 8/35 0.013s peak=7MiB
voxel_radiomics: firstorder batch 9/35 0.014s peak=7MiB
voxel_radiomics: firstorder batch 10/35 0.012s peak=7MiB
voxel_radiomics: firstorder batch 11/35 0.014s peak=7MiB
voxel_radiomics: firstorder batch 12/35 0.015s peak=7MiB
voxel_radiomics: firstorder batch 13/35 0.009s peak=7MiB
voxel_radiomics: firstorder batch 14/35 0.007s peak=7MiB
voxel_radiomics: firstorder batch 15/35 0.007s peak=7MiB
voxel_radiomics: firstorder batch 16/35 0.006s peak=7MiB
voxel_radiomics: firstorder batch 17/35 0.007s peak=7MiB
voxel_radiomics: firstorder batch 18/35 0.007s peak=7MiB
voxel_radiomics: firstorder batch 19/35 0.006s peak=7MiB
voxel_radiomics: firstorder batch 20/35 0.005s peak=7MiB
voxel_radiomics: firstorder batch 21/35 0.005s peak=7MiB
voxel_radiomics: firstorder batch 22/35 0.007s peak=7MiB
voxel_radiomics: firstorder batch 23/35 0.007s peak=7MiB
voxel_radiomics: firstorder batch 24/35 0.009s peak=7MiB
voxel_radiomics: firstorder batch 25/35 0.006s peak=7MiB
voxel_radiomics: firstorder batch 26/35 0.006s peak=7MiB
voxel_radiomics: firstorder batch 27/35 0.006s peak=7MiB
voxel_radiomics: firstorder batch 28/35 0.008s peak=7MiB
voxel_radiomics: firstorder batch 29/35 0.005s peak=7MiB
voxel_radiomics: firstorder batch 30/35 0.006s peak=7MiB
voxel_radiomics: firstorder batch 31/35 0.008s peak=7MiB
voxel_radiomics: firstorder batch 32/35 0.005s peak=7MiB
voxel_radiomics: firstorder batch 33/35 0.008s peak=7MiB
voxel_radiomics: firstorder batch 34/35 0.006s peak=7MiB
voxel_radiomics: firstorder batch 35/35 0.006s peak=6MiB
voxel_radiomics: glcm batch 1/35 0.061s peak=2MiB
voxel_radiomics: glcm batch 2/35 0.055s peak=3MiB
voxel_radiomics: glcm batch 3/35 0.056s peak=5MiB
voxel_radiomics: glcm batch 4/35 0.057s peak=6MiB
voxel_radiomics: glcm batch 5/35 0.059s peak=7MiB
voxel_radiomics: glcm batch 6/35 0.058s peak=8MiB
voxel_radiomics: glcm batch 7/35 0.060s peak=9MiB
voxel_radiomics: glcm batch 8/35 0.062s peak=10MiB
voxel_radiomics: glcm batch 9/35 0.059s peak=11MiB
voxel_radiomics: glcm batch 10/35 0.062s peak=13MiB
voxel_radiomics: glcm batch 11/35 0.059s peak=14MiB
voxel_radiomics: glcm batch 12/35 0.062s peak=15MiB
voxel_radiomics: glcm batch 13/35 0.060s peak=16MiB
voxel_radiomics: glcm batch 14/35 0.061s peak=17MiB
voxel_radiomics: glcm batch 15/35 0.061s peak=18MiB
voxel_radiomics: glcm batch 16/35 0.064s peak=20MiB
voxel_radiomics: glcm batch 17/35 0.061s peak=21MiB
voxel_radiomics: glcm batch 18/35 0.064s peak=22MiB
voxel_radiomics: glcm batch 19/35 0.059s peak=23MiB
voxel_radiomics: glcm batch 20/35 0.062s peak=24MiB
voxel_radiomics: glcm batch 21/35 0.060s peak=25MiB
voxel_radiomics: glcm batch 22/35 0.058s peak=26MiB
voxel_radiomics: glcm batch 23/35 0.065s peak=28MiB
voxel_radiomics: glcm batch 24/35 0.060s peak=29MiB
voxel_radiomics: glcm batch 25/35 0.058s peak=30MiB
voxel_radiomics: glcm batch 26/35 0.058s peak=31MiB
voxel_radiomics: glcm batch 27/35 0.059s peak=32MiB
voxel_radiomics: glcm batch 28/35 0.061s peak=33MiB
voxel_radiomics: glcm batch 29/35 0.056s peak=34MiB
voxel_radiomics: glcm batch 30/35 0.053s peak=36MiB
voxel_radiomics: glcm batch 31/35 0.053s peak=37MiB
voxel_radiomics: glcm batch 32/35 0.059s peak=38MiB
voxel_radiomics: glcm batch 33/35 0.052s peak=39MiB
voxel_radiomics: glcm batch 34/35 0.044s peak=40MiB
voxel_radiomics: glcm batch 35/35 0.031s peak=41MiB
Cohort.map[_DefineAndLabelWithinSubject]: 100%|██████████| 1/1 [00:04<00:00, 4.23s/it]
Cohort.map[_DefineAndLabelWithinSubject]: 100%|██████████| 1/1 [00:04<00:00, 4.23s/it]
With precise (whitelist only): mean Dice=0.792, labelled-voxel disagreement=0.206, ARI=0.449
F:\work\habit_project\habit\viz\habitat_core.py:948: UserWarning: Display geometry conflict: image/anatomy direction does not match mask/label direction. Using the mask/label direction so coronal/sagittal superior-up follows the labelled anatomy. Pass direction= to override.
resolved_direction, resolved_spacing = resolve_display_geometry(
Stability under Appendix S2 perturbation (original vs perturbed):
feature_set mean_dice disagreement ari
Without precise (all texture) 0.786672 0.210007 0.442492
With precise (whitelist) 0.791963 0.205799 0.448918
On this demo subject, precise lowers disagreement (0.210 -> 0.206) and raises mean Dice (0.787 -> 0.792).
MONAI elastic / B-spline deformation of image and ROI mask. A realistic displacement field (magnitude_range=(35.0, 50.0) voxels) models anatomical and contour variation across repeat acquisitions or observer differences.
deform = ImagePerturbationRegistry.create(
"bspline_deform",
sigma_range=(2.0, 4.0),
magnitude_range=(35.0, 50.0),
)
warped = deform(subject, rng=np.random.default_rng(0))
image_w = warped.image(MODALITIES[0])
mask_w = warped.mask(ROI)
ref_bin = np.asarray(mask.data) > 0
mov_bin = np.asarray(mask_w.data) > 0
n_inter = int(np.count_nonzero(ref_bin & mov_bin))
n_union = int(np.count_nonzero(ref_bin | mov_bin))
overlap = pd.DataFrame(
[
{
"metric": "dice",
"value": binary_mask_dice(ref_bin, mov_bin),
},
{
"metric": "jaccard",
"value": float(n_inter / n_union) if n_union else float("nan"),
},
{
"metric": "intersection_voxels",
"value": float(n_inter),
},
{
"metric": "union_voxels",
"value": float(n_union),
},
]
)
print("MONAI bspline_deform ROI overlap metrics:")
print(overlap.round(4).to_string(index=False))
overlap
MONAI bspline_deform ROI overlap metrics:
metric value
dice 0.9454
jaccard 0.8965
intersection_voxels 32845.0000
union_voxels 36637.0000
Anatomy slice before and after the elastic deformation.
fig_warp = plot_intensity_slice(
image_w,
before=image,
roi_mask=mask,
roi_contour=True,
title="MONAI Rand3DElastic (image + ROI share one field)",
before_label="Original",
image_label="bspline_deform",
)
fig_warp.savefig("out/precise_features_bspline_anatomy.png", dpi=150, bbox_inches="tight")
plt.show()
F:\work\habit_project\habit\viz\intensity.py:616: UserWarning: Display geometry conflict: image/anatomy direction does not match mask/label direction. Using the mask/label direction so coronal/sagittal superior-up follows the labelled anatomy. Pass direction= to override.
resolved_direction, resolved_spacing = resolve_display_geometry(
Zoomed edge perturbation figure: Original vs Deformed ROI with XOR contour difference.
counts = np.sum(ref_bin, axis=(1, 2))
z = int(np.argmax(counts)) if int(np.max(counts)) > 0 else int(ref_bin.shape[0] // 2)
grey = np.take(np.asarray(image.data), z, axis=0)
m0 = np.take(ref_bin, z, axis=0)
m1 = np.take(mov_bin, z, axis=0)
# Crop closely around the ROI on the slice so contour differences are clearly visible
union_slice = m0 | m1
rows = np.any(union_slice, axis=1)
cols = np.any(union_slice, axis=0)
ymin, ymax = np.where(rows)[0][[0, -1]]
xmin, xmax = np.where(cols)[0][[0, -1]]
pad = 20
ymin = max(0, ymin - pad)
ymax = min(grey.shape[0], ymax + pad)
xmin = max(0, xmin - pad)
xmax = min(grey.shape[1], xmax + pad)
grey_c = grey[ymin:ymax, xmin:xmax]
m0_c = m0[ymin:ymax, xmin:xmax]
m1_c = m1[ymin:ymax, xmin:xmax]
xor_c = (m0_c != m1_c)
finite = grey_c[np.isfinite(grey_c)]
vmin, vmax = np.percentile(finite, (2.0, 98.0))
with use_style("radiology"):
fig_c, ax = plt.subplots(figsize=(6, 5.5), constrained_layout=True)
ax.imshow(grey_c, cmap="gray", origin="upper", vmin=vmin, vmax=vmax)
ax.contourf(xor_c.astype(float), levels=[0.5, 1.5], colors=["#E69F00"], alpha=0.45, origin="upper")
ax.contour(m0_c.astype(float), levels=[0.5], colors=["#00E5FF"], linewidths=2.0, origin="upper")
ax.contour(m1_c.astype(float), levels=[0.5], colors=["#D55E00"], linewidths=2.0, linestyles="--", origin="upper")
ax.set_title(sanitize_label("MONAI Elastic Edge Perturbation (ROI Zoom)"))
ax.axis("off")
ax.legend(
handles=[
Line2D([0], [0], color="#00E5FF", lw=2.0, label="Original ROI"),
Line2D([0], [0], color="#D55E00", lw=2.0, ls="--", label="Deformed ROI"),
Patch(facecolor="#E69F00", edgecolor="none", alpha=0.45, label="Contour shift (XOR)"),
],
loc="lower right",
frameon=True,
)
fig_c.savefig("out/precise_features_bspline_contours.png", dpi=150, bbox_inches="tight")
plt.show()
Total running time of the script: (1 minutes 17.560 seconds)