[Bug]: disordered tracing : mask smoothing sometimes lets mask touch edge of mask if padding is small

[SylviaWhittle]

Checklist

• Find the offending file in the output. If processing halts, re-run analysis with topostats --core 1 process.
• Describe the bug.
• Include the configuration file.
• Copy of the log-file from running with topostats --log-level debug <command>.
• The exact command that failed. This is what you typed at the command line, including any options.
• TopoStats version, this is reported by topostats --version
• Operating System and Python Version

Describe the bug

When padding is too small, the mask smoothing at the start of disordered tracing trace_dna() causes the mask to touch the border of the image.

Original mask:

Smoothed mask:

This issue is solvable by simply increasing the grain_crop_padding to, say, 10.

I don't think we should spend time trying to prevent this occurring, rather I think we should make it clear in the logs when this error occurs that it can be solved by increasing the padding.

Here is the output currently when the error is thrown:

[Thu, 13 Feb 2025 09:58:54] [ERROR   ] [topostats] [20221019_Ni_Topo_20ng_Shelterin_40nM_5min.0_00060] : Disordered tracing of grain 1 failed. Consider raising an issue on GitHub. Error: 
Traceback (most recent call last):
  File "/Users/laura/TopoStats/topostats/tracing/disordered_tracing.py", line 367, in trace_image_disordered
    disordered_trace_images: dict | None = disordered_trace_grain(
                                           ^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/laura/TopoStats/topostats/tracing/disordered_tracing.py", line 654, in disordered_trace_grain
    disorderedtrace.trace_dna()
  File "/Users/laura/TopoStats/topostats/tracing/disordered_tracing.py", line 127, in trace_dna
    ).get_skeleton()
      ^^^^^^^^^^^^^^
  File "/Users/laura/TopoStats/topostats/tracing/skeletonize.py", line 70, in get_skeleton
    return self._get_skeletonize()
           ^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/laura/TopoStats/topostats/tracing/skeletonize.py", line 90, in _get_skeletonize
    return self._skeletonize_topostats(image=self.image, mask=self.mask, height_bias=self.height_bias).astype(
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/laura/TopoStats/topostats/tracing/skeletonize.py", line 186, in _skeletonize_topostats
    return topostatsSkeletonize(image, mask, height_bias).do_skeletonising()
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/laura/TopoStats/topostats/tracing/skeletonize.py", line 246, in do_skeletonising
    self._do_skeletonising_iteration()
  File "/Users/laura/TopoStats/topostats/tracing/skeletonize.py", line 269, in _do_skeletonising_iteration
    if self._delete_pixel_subit1(point):
       ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/laura/TopoStats/topostats/tracing/skeletonize.py", line 326, in _delete_pixel_subit1
    self.p7, self.p8, self.p9, self.p6, self.p2, self.p5, self.p4, self.p3 = self.get_local_pixels_binary(
                                                                             ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/laura/TopoStats/topostats/tracing/skeletonize.py", line 546, in get_local_pixels_binary
    return np.delete(local_pixels, 4)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/laura/miniconda3/envs/topoly/lib/python3.11/site-packages/numpy/lib/function_base.py", line 5330, in delete
    raise IndexError(
IndexError: index 4 is out of bounds for axis 0 with size 0

Proposed solution:

• Add a check just after the smoothing step in trace_dna() and fail early:

Copy of the log-file from running with topostats --log-level debug <command>

[Thu, 13 Feb 2025 09:58:54] [ERROR   ] [topostats] [20221019_Ni_Topo_20ng_Shelterin_40nM_5min.0_00060] : Disordered tracing of grain 1 failed. Consider raising an issue on GitHub. Error: 
Traceback (most recent call last):
  File "/Users/laura/TopoStats/topostats/tracing/disordered_tracing.py", line 367, in trace_image_disordered
    disordered_trace_images: dict | None = disordered_trace_grain(
                                           ^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/laura/TopoStats/topostats/tracing/disordered_tracing.py", line 654, in disordered_trace_grain
    disorderedtrace.trace_dna()
  File "/Users/laura/TopoStats/topostats/tracing/disordered_tracing.py", line 127, in trace_dna
    ).get_skeleton()
      ^^^^^^^^^^^^^^
  File "/Users/laura/TopoStats/topostats/tracing/skeletonize.py", line 70, in get_skeleton
    return self._get_skeletonize()
           ^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/laura/TopoStats/topostats/tracing/skeletonize.py", line 90, in _get_skeletonize
    return self._skeletonize_topostats(image=self.image, mask=self.mask, height_bias=self.height_bias).astype(
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/laura/TopoStats/topostats/tracing/skeletonize.py", line 186, in _skeletonize_topostats
    return topostatsSkeletonize(image, mask, height_bias).do_skeletonising()
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/laura/TopoStats/topostats/tracing/skeletonize.py", line 246, in do_skeletonising
    self._do_skeletonising_iteration()
  File "/Users/laura/TopoStats/topostats/tracing/skeletonize.py", line 269, in _do_skeletonising_iteration
    if self._delete_pixel_subit1(point):
       ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/laura/TopoStats/topostats/tracing/skeletonize.py", line 326, in _delete_pixel_subit1
    self.p7, self.p8, self.p9, self.p6, self.p2, self.p5, self.p4, self.p3 = self.get_local_pixels_binary(
                                                                             ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/laura/TopoStats/topostats/tracing/skeletonize.py", line 546, in get_local_pixels_binary
    return np.delete(local_pixels, 4)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/laura/miniconda3/envs/topoly/lib/python3.11/site-packages/numpy/lib/function_base.py", line 5330, in delete
    raise IndexError(
IndexError: index 4 is out of bounds for axis 0 with size 0

Include the configuration file

# Config file generated 2025-02-11 19:02:51
# # For more information on configuration and how to use it:
# https://afm-spm.github.io/TopoStats/main/configuration.html
base_dir: ./data # Directory in which to search for data files
output_dir: ./output # Directory to output results to
log_level: info # Verbosity of output. Options: warning, error, info, debug
cores: 1 # Number of CPU cores to utilise for processing multiple files simultaneously.
file_ext: .spm # File extension of the data files.
loading:
  channel: Height # Channel to pull data from in the data files.
filter:
  run: true # Options : true, false
  row_alignment_quantile: 0.5 # lower values may improve flattening of larger features
  threshold_method: std_dev # Options : otsu, std_dev, absolute
  otsu_threshold_multiplier: 1.0
  threshold_std_dev:
    below: 10.0 # Threshold for data below the image background
    above: 1.0 # Threshold for data above the image background
  threshold_absolute:
    below: -1.0 # Threshold for data below the image background
    above: 1.0 # Threshold for data above the image background
  gaussian_size: 1.0121397464510862 # Gaussian blur intensity in px
  gaussian_mode: nearest # Mode for Gaussian blurring. Options : nearest, reflect, constant, mirror, wrap
  # Scar remvoal parameters. Be careful with editing these as making the algorithm too sensitive may
  # result in ruining legitimate data.
  remove_scars:
    run: false
    removal_iterations: 2 # Number of times to run scar removal.
    threshold_low: 0.250 # lower values make scar removal more sensitive
    threshold_high: 0.666 # lower values make scar removal more sensitive
    max_scar_width: 4 # Maximum thickness of scars in pixels.
    min_scar_length: 16 # Minimum length of scars in pixels.
grains:
  run: true # Options : true, false
  # Thresholding by height
  grain_crop_padding: 1 # Padding to apply to grains. Needs to be at least 1, more padding may help with unets.
  threshold_method: absolute # Options : std_dev, otsu, absolute, unet
  otsu_threshold_multiplier: 1.0
  threshold_std_dev:
    below: 10.0 # Threshold for grains below the image background
    above: 1.0 # Threshold for grains above the image background
  threshold_absolute:
    below: -1.0 # Threshold for grains below the image background
    above: 1.0 # Threshold for grains above the image background
  direction: above # Options: above, below, both (defines whether to look for grains above or below thresholds or both)
  # Thresholding by area
  smallest_grain_size_nm2: 50 # Size in nm^2 of tiny grains/blobs (noise) to remove, must be > 0.0
  absolute_area_threshold:
    above: [300, 30000] # above surface [Low, High] in nm^2 (also takes null)
    below: [null, null] # below surface [Low, High] in nm^2 (also takes null)
  remove_edge_intersecting_grains: true # Whether or not to remove grains that touch the image border
  unet_config:
    model_path: /Users/sylvi/topo_data/topo-test-unets/debugging-laura/new_shelterin_model.h5 # Path to a trained U-Net model
    upper_norm_bound: 8.0 # Upper bound for normalisation of input data. This should be slightly higher than the maximum desired / expected height of grains.
    lower_norm_bound: -2.0 # Lower bound for normalisation of input data. This should be slightly lower than the minimum desired / expected height of the background.
  vetting:
    class_conversion_size_thresholds: null # Class conversion size thresholds, list of tuples of 3 integers and 2 integers, ie list[tuple[tuple[int, int, int], tuple[int, int]]] eg [[[1, 2, 3], [5, 10]]] for each region of class 1 to convert to 2 if smaller than 5 nm^2 and to class 3 if larger than 10 nm^2.
    class_region_number_thresholds: null # Class region number thresholds, list of lists, ie [[class, low, high],] eg [[1, 2, 4], [2, 1, 1]] for class 1 to have 2-4 regions and class 2 to have 1 region. Can use None to not set an upper/lower bound.
    class_size_thresholds: null # Class size thresholds (nm^2), list of tuples of 3 integers, ie [[class, low, high],] eg [[1, 100, 1000], [2, 1000, None]] for class 1 to have 100-1000 nm^2 and class 2 to have 1000-any nm^2. Can use None to not set an upper/lower bound.
    nearby_conversion_classes_to_convert: null # Class conversion for nearby regions, list of tuples of two-integer tuples, eg [[[1, 2], [3, 4]]] to convert class 1 to 2 and 3 to 4 for small touching regions
    class_touching_threshold: 5 # Number of dilation steps to use for detecting touching regions
    keep_largest_labelled_regions_classes: [1] # Classes to keep the only largest regions for, list of integers eg [1, 2] to keep only the largest regions of class 1 and 2
    class_connection_point_thresholds: null # Class connection point thresholds, [[[class_1, class_2], [min, max]]] eg [[[1, 2], [1, 1]]] for class 1 to have 1 connection point with class 2
grainstats:
  run: true # Options : true, false
  edge_detection_method: binary_erosion # Options: canny, binary erosion. Do not change this unless you are sure of what this will do.
  cropped_size: -1 # Length (in nm) of square cropped images (can take -1 for grain-sized box)
  extract_height_profile: true # Extract height profiles along maximum feret of molecules
disordered_tracing:
  run: true # Options : true, false
  class_index: 1 # The class index to trace. This is the class index of the grains.
  min_skeleton_size: 10 # Minimum number of pixels in a skeleton for it to be retained.
  mask_smoothing_params:
    gaussian_sigma: 2 # Gaussian smoothing parameter 'sigma' in pixels.
    dilation_iterations: 2 # Number of dilation iterations to use for grain smoothing.
    holearea_min_max: [0, null] # Range (min, max) of a hole area in nm to refill in the smoothed masks.
  skeletonisation_params:
    method: topostats # Options : zhang | lee | thin | topostats
    height_bias: 0.6 # Percentage of lowest pixels to remove each skeletonisation iteration. 1 equates to zhang.
  pruning_params:
    method: topostats # Method to clean branches of the skeleton. Options : topostats
    max_length: 10.0 # Maximum length in nm to remove a branch containing an endpoint.
    height_threshold: # The height to remove branches below.
    method_values: mid # The method to obtain a branch's height for pruning. Options : min | median | mid.
    method_outlier: mean_abs # The method to prune branches based on height. Options : abs | mean_abs | iqr.
nodestats:
  run: true # Options : true, false
  node_joining_length: 7.0 # The distance in nanometres over which to join nearby crossing points.
  node_extend_dist: 14.0 # The distance in nanometres over which to join nearby odd-branched nodes.
  branch_pairing_length: 20.0 # The length in nanometres from the crossing point to pair and trace, obtaining FWHM's.
  pair_odd_branches: false # Whether to try and pair odd-branched nodes. Options: true and false.
ordered_tracing:
  run: true
  ordering_method: nodestats # The method of ordering the disordered traces.
splining:
  run: true # Options : true, false
  method: "rolling_window" # Options : "spline", "rolling_window"
  rolling_window_size: 20.0e-9 # size in nm of the rolling window.
  spline_step_size: 7.0e-9 # The sampling rate of the spline in metres.
  spline_linear_smoothing: 5.0 # The amount of smoothing to apply to linear features.
  spline_circular_smoothing: 5.0 # The amount of smoothing to apply to circular features.
  spline_degree: 3 # The polynomial degree of the spline.
curvature:
  run: true # Options : true, false
  colourmap_normalisation_bounds: [-0.5, 0.5] # Radians per nm to normalise the colourmap to.
plotting:
  run: true # Options : true, false
  style: topostats.mplstyle # Options : topostats.mplstyle or path to a matplotlibrc params file
  savefig_format: null # Options : null, png, svg or pdf. tif is also available although no metadata will be saved. (defaults to png) See https://matplotlib.org/stable/api/_as_gen/matplotlib.pyplot.savefig.html
  savefig_dpi: 200 # Options : null (defaults to the value in topostats/plotting_dictionary.yaml), see https://afm-spm.github.io/TopoStats/main/configuration.html#further-customisation and https://matplotlib.org/stable/api/_as_gen/matplotlib.pyplot.savefig.html
  pixel_interpolation: null # Options : https://matplotlib.org/stable/gallery/images_contours_and_fields/interpolation_methods.html
  image_set: all # Options : all, core
  zrange: [-2, 8] # low and high height range for core images (can take [null, null]). low <= high
  colorbar: true # Options : true, false
  axes: true # Options : true, false (due to off being a bool when parsed)
  num_ticks: [null, null] # Number of ticks to have along the x and y axes. Options : null (auto) or integer > 1
  cmap: null # Colormap/colourmap to use (default is 'nanoscope' which is used if null, other options are 'afmhot', 'viridis' etc.)
  mask_cmap: blue_purple_green # Options : blu, jet_r and any in matplotlib
  histogram_log_axis: false # Options : true, false
summary_stats:
  run: true # Whether to make summary plots for output data
  config: null

To Reproduce

Files here: https://drive.google.com/drive/folders/1vKTyqujkT_Iv49X-VNDpgyX5LnlDxPL8

Run topostats with above config on those files

topostats -c my-config.yaml process

TopoStats Version

Git main branch

Python Version

3.11

Operating System

MacOS M1/M2 (post-2021)

Python Packages

absl-py==2.1.0
accessible-pygments==0.0.5
AFMReader==0.0.3
aiohappyeyeballs==2.4.0
aiohttp==3.10.5
aiohttp-retry==2.8.3
aiosignal==1.3.1
alabaster==1.0.0
amqp==5.2.0
annotated-types==0.7.0
antlr4-python3-runtime==4.9.3
appdirs==1.4.4
appnope==0.1.4
art==6.3
astroid==3.2.4
asttokens==2.4.1
astunparse==1.6.3
asyncssh==2.16.0
atpublic==5.0
attrs==24.2.0
babel==2.16.0
beautifulsoup4==4.12.3
billiard==4.2.0
biopython==1.84
black==24.8.0
cachetools==5.5.0
celery==5.4.0
certifi==2024.12.14
cffi==1.17.0
cfgv==3.4.0
charset-normalizer==3.4.1
cheap_repr==0.5.2
click==8.1.7
click-didyoumean==0.3.1
click-plugins==1.1.1
click-repl==0.3.0
colorama==0.4.6
comm==0.2.2
configobj==5.0.8
contourpy==1.3.1
coverage==7.6.1
cryptography==43.0.0
cycler==0.12.1
debugpy==1.8.5
decorator==5.1.1
deepdiff==8.0.1
dictdiffer==0.9.0
dill==0.3.8
diskcache==5.6.3
distlib==0.3.9
distro==1.9.0
docstring_parser==0.16
docutils==0.21.2
dpath==2.2.0
dulwich==0.22.1
dvc==3.54.0
dvc-data==3.16.4
dvc-gdrive==3.0.1
dvc-http==2.32.0
dvc-objects==5.1.0
dvc-render==1.0.2
dvc-studio-client==0.21.0
dvc-task==0.4.0
dvclive==3.48.0
entrypoints==0.4
et_xmlfile==2.0.0
execnet==2.1.1
executing==2.0.1
filelock==3.16.1
filetype==1.2.0
flatbuffers==24.3.25
flatten-dict==0.4.2
flufl.lock==8.1.0
fonttools==4.55.5
frozenlist==1.4.1
fsspec==2024.6.1
funcy==2.0
gast==0.6.0
gitdb==4.0.11
GitPython==3.1.43
google-api-core==2.19.1
google-api-python-client==2.141.0
google-auth==2.34.0
google-auth-httplib2==0.2.0
google-pasta==0.2.0
googleapis-common-protos==1.63.2
grandalf==0.8
grpcio==1.65.5
gto==1.7.1
h5glance==0.9.0
h5py==3.11.0
htmlgen==2.0.0
httplib2==0.22.0
hydra-core==1.3.2
icdiff==2.0.7
identify==2.6.2
idna==3.10
igor2==0.5.8
imageio==2.37.0
imagesize==1.4.1
iniconfig==2.0.0
ipykernel==6.29.5
ipython==8.26.0
isort==5.13.2
iterative-telemetry==0.0.8
jedi==0.19.1
Jinja2==3.1.4
joblib==1.4.2
jupyter_client==8.6.2
jupyter_core==5.7.2
keras==3.5.0
kiwisolver==1.4.8
kombu==5.4.0
lazy_loader==0.4
libclang==18.1.1
llvmlite==0.44.0
loguru==0.7.2
magicgui==0.10.0
Markdown==3.7
markdown-it-py==3.0.0
MarkupSafe==2.1.5
matplotlib==3.9.4
matplotlib-inline==0.1.7
mccabe==0.7.0
mdit-py-plugins==0.4.2
mdurl==0.1.2
ml-dtypes==0.3.2
multidict==6.0.5
mypy-extensions==1.0.0
myst-parser==4.0.0
namex==0.0.8
nest-asyncio==1.6.0
networkx==3.4.2
nodeenv==1.9.1
numba==0.61.0
numpy==1.26.4
numpydoc==1.8.0
numpyencoder==0.3.0
oauth2client==4.1.3
omegaconf==2.3.0
opencv-python==4.10.0.84
openpyxl==3.1.5
opt-einsum==3.3.0
optree==0.12.1
orderly-set==5.2.2
orjson==3.10.7
packaging==24.2
pandas==2.2.3
pandas-stubs==2.2.2.240807
parso==0.8.4
pathspec==0.12.1
pexpect==4.9.0
pillow==11.1.0
platformdirs==4.3.6
pluggy==1.5.0
pockets==0.9.1
pprintpp==0.4.0
pre_commit==4.0.1
prompt_toolkit==3.0.47
proto-plus==1.24.0
protobuf==4.25.4
psutil==5.9.8
psygnal==0.11.1
ptyprocess==0.7.0
pure_eval==0.2.3
pyasn1==0.6.0
pyasn1_modules==0.4.0
pyconify==0.2
pycparser==2.22
pydantic==2.8.2
pydantic_core==2.20.1
pydata-sphinx-theme==0.15.4
pydot==3.0.1
PyDrive2==1.20.0
pyFFTW==0.14.0
pyfiglet==1.0.2
pygit2==1.15.1
Pygments==2.19.1
pygtrie==2.5.0
pylint==3.2.6
pynvml==11.5.3
pyOpenSSL==24.2.1
pyparsing==3.2.1
pyspm==0.6.2
pytest==7.4.4
pytest-cov==5.0.0
pytest-durations==1.2.0
pytest-github-actions-annotate-failures==0.2.0
pytest-icdiff==0.9
pytest-lazy-fixture==0.6.3
pytest-mpl==0.17.0
pytest-regtest==2.3.1
pytest-testmon==2.1.1
pytest-xdist==3.6.1
python-dateutil==2.9.0.post0
pytz==2024.2
pyupgrade==3.17.0
PyYAML==6.0.2
pyzmq==26.1.1
QtPy==2.4.2
requests==2.32.3
rich==13.7.1
rsa==4.9
ruamel.yaml==0.18.6
ruamel.yaml.clib==0.2.8
schema==0.7.7
scikit-image==0.25.0
scikit-learn==1.5.1
scipy==1.15.1
scmrepo==3.3.7
seaborn==0.13.2
semver==3.0.2
shellingham==1.5.4
shortuuid==1.0.13
shtab==1.7.1
six==1.17.0
skan==0.12.2
smmap==5.0.1
snoop==0.4.3
snowballstemmer==2.2.0
soupsieve==2.6
Sphinx==8.1.3
sphinx-autoapi==3.1.0b0
sphinx-autodoc-typehints==2.5.0
sphinx-markdown-tables==0.0.17
sphinx-multiversion==0.2.4
sphinx-rtd-theme==3.0.1
sphinxcontrib-applehelp==2.0.0
sphinxcontrib-devhelp==2.0.0
sphinxcontrib-htmlhelp==2.1.0
sphinxcontrib-jquery==4.1
sphinxcontrib-jsmath==1.0.1
sphinxcontrib-mermaid==1.0.0
sphinxcontrib-napoleon==0.7
sphinxcontrib-qthelp==2.0.0
sphinxcontrib-serializinghtml==2.0.0
sqltrie==0.11.1
stack-data==0.6.3
superqt==0.7.1
tabulate==0.9.0
tensorboard==2.16.2
tensorboard-data-server==0.7.2
tensorflow==2.16.2
tensorflow-io-gcs-filesystem==0.37.1
tensorflow-macos==2.16.2
tensorflow-metal==1.1.0
termcolor==2.4.0
threadpoolctl==3.5.0
tifffile==2025.1.10
tokenize-rt==6.0.0
tomlkit==0.13.2
toolz==1.0.0
topoly==1.0.4
-e git+ssh://[email protected]/AFM-SPM/TopoStats.git@77e68376d5d75ea8d190ed7101e83d2a4f75cfae#egg=topostats
tornado==6.4.1
tqdm==4.67.1
traitlets==5.14.3
typer==0.12.4
types-protobuf==5.27.0.20240626
types-pytz==2024.1.0.20240417
types-requests==2.32.0.20240712
typing_extensions==4.12.2
tzdata==2025.1
uritemplate==4.1.1
urllib3==2.3.0
vine==5.1.0
virtualenv==20.27.1
voluptuous==0.15.2
wcwidth==0.2.13
Werkzeug==3.0.3
wrapt==1.16.0
yarl==1.9.4
zc.lockfile==3.0.post1