Skip to content

Disordered Tracing Modules

Generates disordered traces (pruned skeletons) and metrics.

disorderedTrace

Calculate disordered traces for a DNA molecule and calculates statistics from those traces.

Parameters:

Name Type Description Default
image NDArray

Cropped image, typically padded beyond the bounding box.

 required
mask NDArray

Labelled mask for the grain, typically padded beyond the bounding box.

 required
filename str

Filename being processed.

 required
pixel_to_nm_scaling float

Pixel to nm scaling.

 required
min_skeleton_size int

Minimum skeleton size below which tracing statistics are not calculated.

 10
mask_smoothing_params dict

Dictionary of parameters to smooth the grain mask for better quality skeletonisation results. Contains a gaussian 'sigma' and number of dilation iterations.

 None
skeletonisation_params dict

Skeletonisation Parameters. Method of skeletonisation to use 'topostats' is the original TopoStats method. Three methods from scikit-image are available 'zhang', 'lee' and 'thin'.

 None
pruning_params dict

Dictionary of pruning parameters. Contains 'method', 'max_length', 'height_threshold', 'method_values' 'method_outlier' and 'only_height_prune_endpoints'.

 None
n_grain int

Grain number being processed (only used in logging).

 None
Source code in topostats\tracing\disordered_tracing.py 
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
class disorderedTrace:  # pylint: disable=too-many-instance-attributes
    """
    Calculate disordered traces for a DNA molecule and calculates statistics from those traces.

    Parameters
    ----------
    image : npt.NDArray
        Cropped image, typically padded beyond the bounding box.
    mask : npt.NDArray
        Labelled mask for the grain, typically padded beyond the bounding box.
    filename : str
        Filename being processed.
    pixel_to_nm_scaling : float
        Pixel to nm scaling.
    min_skeleton_size : int
        Minimum skeleton size below which tracing statistics are not calculated.
    mask_smoothing_params : dict
        Dictionary of parameters to smooth the grain mask for better quality skeletonisation results. Contains
        a gaussian 'sigma' and number of dilation iterations.
    skeletonisation_params : dict
        Skeletonisation Parameters. Method of skeletonisation to use 'topostats' is the original TopoStats
        method. Three methods from scikit-image are available 'zhang', 'lee' and 'thin'.
    pruning_params : dict
        Dictionary of pruning parameters. Contains 'method', 'max_length', 'height_threshold', 'method_values'
        'method_outlier' and 'only_height_prune_endpoints'.
    n_grain : int
        Grain number being processed (only  used in logging).
    """

    def __init__(  # pylint: disable=too-many-arguments
        self,
        image: npt.NDArray,
        mask: npt.NDArray,
        filename: str,
        pixel_to_nm_scaling: float,
        min_skeleton_size: int = 10,
        mask_smoothing_params: dict | None = None,
        skeletonisation_params: dict | None = None,
        pruning_params: dict | None = None,
        n_grain: int = None,
    ):
        """
        Calculate disordered traces for a DNA molecule and calculates statistics from those traces.

        Parameters
        ----------
        image : npt.NDArray
            Cropped image, typically padded beyond the bounding box.
        mask : npt.NDArray
            Labelled mask for the grain, typically padded beyond the bounding box.
        filename : str
            Filename being processed.
        pixel_to_nm_scaling : float
            Pixel to nm scaling.
        min_skeleton_size : int
            Minimum skeleton size below which tracing statistics are not calculated.
        mask_smoothing_params : dict
            Dictionary of parameters to smooth the grain mask for better quality skeletonisation results. Contains
            a gaussian 'sigma' and number of dilation iterations.
        skeletonisation_params : dict
            Skeletonisation Parameters. Method of skeletonisation to use 'topostats' is the original TopoStats
            method. Three methods from scikit-image are available 'zhang', 'lee' and 'thin'.
        pruning_params : dict
            Dictionary of pruning parameters. Contains 'method', 'max_length', 'height_threshold', 'method_values',
            'method_outlier' and 'only_height_prune_endpoints'.
        n_grain : int
            Grain number being processed (only  used in logging).
        """
        self.image = image
        self.mask = mask
        self.filename = filename
        self.pixel_to_nm_scaling = pixel_to_nm_scaling
        self.min_skeleton_size = min_skeleton_size
        self.mask_smoothing_params = mask_smoothing_params
        self.skeletonisation_params = (
            skeletonisation_params if skeletonisation_params is not None else {"method": "zhang"}
        )
        self.pruning_params = pruning_params if pruning_params is not None else {"method": "topostats"}
        self.n_grain = n_grain
        # Images
        self.smoothed_mask = np.zeros_like(self.image)
        self.skeleton = np.zeros_like(self.image)
        self.pruned_skeleton = np.zeros_like(self.image)
        # Trace
        self.disordered_trace = None

        # suppresses scipy splining warnings
        warnings.filterwarnings("ignore")

        LOGGER.debug(f"[{self.filename}] Performing Disordered Tracing")

    def trace_dna(self):
        """Perform the DNA skeletonisation and cleaning pipeline."""
        self.smoothed_mask = self.smooth_mask(self.mask, **self.mask_smoothing_params)
        if check_pixel_touching_edge(self.smoothed_mask):
            LOGGER.warning(
                f"[{self.filename}] : Grain {self.n_grain} skipped as padding is too small. Please consider "
                "increasing the value of grains.grain_crop_padding in your config file and try again."
            )
            self.disordered_trace = None
            return
        self.skeleton = getSkeleton(
            self.image,
            self.smoothed_mask,
            method=self.skeletonisation_params["method"],
            height_bias=self.skeletonisation_params["height_bias"],
        ).get_skeleton()
        self.pruned_skeleton = prune_skeleton(
            self.image, self.skeleton, self.pixel_to_nm_scaling, **self.pruning_params.copy()
        )
        self.pruned_skeleton = self.remove_touching_edge(self.pruned_skeleton)
        self.disordered_trace = np.argwhere(self.pruned_skeleton == 1)

        if self.disordered_trace is None:
            LOGGER.warning(f"[{self.filename}] : Grain {self.n_grain} failed to Skeletonise.")
            self.disordered_trace = None
        elif len(self.disordered_trace) < self.min_skeleton_size:
            LOGGER.warning(f"[{self.filename}] : Grain {self.n_grain} skeleton < {self.min_skeleton_size}, skipping.")
            self.disordered_trace = None

    def re_add_holes(
        self,
        orig_mask: npt.NDArray,
        smoothed_mask: npt.NDArray,
        holearea_min_max: tuple[float | int | None] = (2, None),
    ) -> npt.NDArray:
        """
        Restore holes in masks that were occluded by dilation.

        As Gaussian dilation smoothing methods can close holes in the original mask, this function obtains those holes
        (based on the general background being the first due to padding) and adds them back into the smoothed mask. When
        paired with ``smooth_mask``, this essentially just smooths the outer edge of the mask.

        Parameters
        ----------
        orig_mask : npt.NDArray
            Original mask.
        smoothed_mask : npt.NDArray
            Original mask but with inner and outer edged smoothed. The smoothing operation may have closed up important
            holes in the mask.
        holearea_min_max : tuple[float | int | None]
            Tuple of minimum and maximum hole area (in nanometers) to replace from the original mask into the smoothed
            mask.

        Returns
        -------
        npt.NDArray
            Smoothed mask with holes restored.
        """
        # handle none's
        if set(holearea_min_max) == {None}:
            return smoothed_mask
        if None in holearea_min_max:
            none_index = holearea_min_max.index(None)
            holearea_min_max[none_index] = 0 if none_index == 0 else np.inf

        # obtain px holesizes
        holesize_min_px = holearea_min_max[0] / ((self.pixel_to_nm_scaling) ** 2)
        holesize_max_px = holearea_min_max[1] / ((self.pixel_to_nm_scaling) ** 2)

        # obtain a hole mask
        holes = 1 - orig_mask
        holes = label(holes)
        hole_sizes = [holes[holes == i].size for i in range(1, holes.max() + 1)]
        holes[holes == 1] = 0  # set background to 0 assuming it is the first hole seen (from top left)

        # remove too small or too big holes from mask
        for i, hole_size in enumerate(hole_sizes):
            if hole_size < holesize_min_px or hole_size > holesize_max_px:  # small holes may be fake are left out
                holes[holes == i + 1] = 0
        holes[holes != 0] = 1  # set correct sixe holes to 1

        # replace correct sized holes
        return np.where(holes == 1, 0, smoothed_mask)

    @staticmethod
    def remove_touching_edge(skeleton: npt.NDArray) -> npt.NDArray:
        """
        Remove any skeleton points touching the border (to prevent errors later).

        Parameters
        ----------
        skeleton : npt.NDArray
            A binary array where touching clusters of 1's become 0's if touching the edge of the array.

        Returns
        -------
        npt.NDArray
            Skeleton without points touching the border.
        """
        for edge in [skeleton[0, :-1], skeleton[:-1, -1], skeleton[-1, 1:], skeleton[1:, 0]]:
            uniques = np.unique(edge)
            for i in uniques:
                skeleton[skeleton == i] = 0
        return skeleton

    def smooth_mask(
        self,
        grain: npt.NDArray,
        dilation_iterations: int = 2,
        gaussian_sigma: float | int = 2,
        holearea_min_max: tuple[int | float | None] = (0, None),
    ) -> npt.NDArray:
        """
        Smooth a grain mask based on the lower number of binary pixels added from dilation or gaussian.

        This method ensures gaussian smoothing isn't too aggressive and covers / creates gaps in the mask.

        Parameters
        ----------
        grain : npt.NDArray
            Numpy array of the grain mask.
        dilation_iterations : int
            Number of times to dilate the grain to smooth it. Default is 2.
        gaussian_sigma : float | None
            Gaussian sigma value to smooth the grains after an Otsu threshold. If None, defaults to 2.
        holearea_min_max : tuple[float | int | None]
            Tuple of minimum and maximum hole area (in nanometers) to replace from the original mask into the smoothed
            mask.

        Returns
        -------
        npt.NDArray
            Numpy array of smoothed image.
        """
        # Option to disable the smoothing (i.e. U-Net masks are already smooth)
        if dilation_iterations is None and gaussian_sigma is None:
            LOGGER.debug(f"[{self.filename}] : no grain smoothing done")
            return grain

        # Option to only do gaussian or dilation
        if dilation_iterations is not None:
            dilation = ndimage.binary_dilation(grain, iterations=dilation_iterations).astype(np.int32)
        else:
            gauss = filters.gaussian(grain, sigma=gaussian_sigma)
            gauss = np.where(gauss > filters.threshold_otsu(gauss) * 1.3, 1, 0)
            gauss = gauss.astype(np.int32)
            LOGGER.debug(f"[{self.filename}] : smoothing done by gaussian {gaussian_sigma}")
            return self.re_add_holes(grain, gauss, holearea_min_max)
        if gaussian_sigma is not None:
            gauss = filters.gaussian(grain, sigma=gaussian_sigma)
            gauss = np.where(gauss > filters.threshold_otsu(gauss) * 1.3, 1, 0)
            gauss = gauss.astype(np.int32)
        else:
            LOGGER.debug(f"[{self.filename}] : smoothing done by dilation {dilation_iterations}")
            return self.re_add_holes(grain, dilation, holearea_min_max)

        # Competition option between dilation and gaussian mask differences wrt original grains
        if abs(dilation.sum() - grain.sum()) > abs(gauss.sum() - grain.sum()):
            LOGGER.debug(f"[{self.filename}] : smoothing done by gaussian {gaussian_sigma}")
            return self.re_add_holes(grain, gauss, holearea_min_max)
        LOGGER.debug(f"[{self.filename}] : smoothing done by dilation {dilation_iterations}")
        return self.re_add_holes(grain, dilation, holearea_min_max)

    @staticmethod
    def calculate_dna_width(
        smoothed_mask: npt.NDArray, pruned_skeleton: npt.NDArray, pixel_to_nm_scaling: float = 1
    ) -> float:
        """
        Calculate the mean width in metres of the DNA using the trace and mask.

        Parameters
        ----------
        smoothed_mask : npt.NDArray
            Smoothed mask to be measured.
        pruned_skeleton : npt.NDArray
            Pruned skeleton.
        pixel_to_nm_scaling : float
            Scaling of pixels to nanometres.

        Returns
        -------
        float
            Width of grain in metres.
        """
        dist_trans = ndimage.distance_transform_edt(smoothed_mask)
        comb = np.where(pruned_skeleton == 1, dist_trans, 0)

        return comb[comb != 0].mean() * 2 * pixel_to_nm_scaling

__init__(image: npt.NDArray, mask: npt.NDArray, filename: str, pixel_to_nm_scaling: float, min_skeleton_size: int = 10, mask_smoothing_params: dict | None = None, skeletonisation_params: dict | None = None, pruning_params: dict | None = None, n_grain: int = None)

Calculate disordered traces for a DNA molecule and calculates statistics from those traces.

Parameters:

Name Type Description Default
image NDArray

Cropped image, typically padded beyond the bounding box.

 required
mask NDArray

Labelled mask for the grain, typically padded beyond the bounding box.

 required
filename str

Filename being processed.

 required
pixel_to_nm_scaling float

Pixel to nm scaling.

 required
min_skeleton_size int

Minimum skeleton size below which tracing statistics are not calculated.

 10
mask_smoothing_params dict

Dictionary of parameters to smooth the grain mask for better quality skeletonisation results. Contains a gaussian 'sigma' and number of dilation iterations.

 None
skeletonisation_params dict

Skeletonisation Parameters. Method of skeletonisation to use 'topostats' is the original TopoStats method. Three methods from scikit-image are available 'zhang', 'lee' and 'thin'.

 None
pruning_params dict

Dictionary of pruning parameters. Contains 'method', 'max_length', 'height_threshold', 'method_values', 'method_outlier' and 'only_height_prune_endpoints'.

 None
n_grain int

Grain number being processed (only used in logging).

 None
Source code in topostats\tracing\disordered_tracing.py 
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
def __init__(  # pylint: disable=too-many-arguments
    self,
    image: npt.NDArray,
    mask: npt.NDArray,
    filename: str,
    pixel_to_nm_scaling: float,
    min_skeleton_size: int = 10,
    mask_smoothing_params: dict | None = None,
    skeletonisation_params: dict | None = None,
    pruning_params: dict | None = None,
    n_grain: int = None,
):
    """
    Calculate disordered traces for a DNA molecule and calculates statistics from those traces.

    Parameters
    ----------
    image : npt.NDArray
        Cropped image, typically padded beyond the bounding box.
    mask : npt.NDArray
        Labelled mask for the grain, typically padded beyond the bounding box.
    filename : str
        Filename being processed.
    pixel_to_nm_scaling : float
        Pixel to nm scaling.
    min_skeleton_size : int
        Minimum skeleton size below which tracing statistics are not calculated.
    mask_smoothing_params : dict
        Dictionary of parameters to smooth the grain mask for better quality skeletonisation results. Contains
        a gaussian 'sigma' and number of dilation iterations.
    skeletonisation_params : dict
        Skeletonisation Parameters. Method of skeletonisation to use 'topostats' is the original TopoStats
        method. Three methods from scikit-image are available 'zhang', 'lee' and 'thin'.
    pruning_params : dict
        Dictionary of pruning parameters. Contains 'method', 'max_length', 'height_threshold', 'method_values',
        'method_outlier' and 'only_height_prune_endpoints'.
    n_grain : int
        Grain number being processed (only  used in logging).
    """
    self.image = image
    self.mask = mask
    self.filename = filename
    self.pixel_to_nm_scaling = pixel_to_nm_scaling
    self.min_skeleton_size = min_skeleton_size
    self.mask_smoothing_params = mask_smoothing_params
    self.skeletonisation_params = (
        skeletonisation_params if skeletonisation_params is not None else {"method": "zhang"}
    )
    self.pruning_params = pruning_params if pruning_params is not None else {"method": "topostats"}
    self.n_grain = n_grain
    # Images
    self.smoothed_mask = np.zeros_like(self.image)
    self.skeleton = np.zeros_like(self.image)
    self.pruned_skeleton = np.zeros_like(self.image)
    # Trace
    self.disordered_trace = None

    # suppresses scipy splining warnings
    warnings.filterwarnings("ignore")

    LOGGER.debug(f"[{self.filename}] Performing Disordered Tracing")

calculate_dna_width(smoothed_mask: npt.NDArray, pruned_skeleton: npt.NDArray, pixel_to_nm_scaling: float = 1) -> float staticmethod

Calculate the mean width in metres of the DNA using the trace and mask.

Parameters:

Name Type Description Default
smoothed_mask NDArray

Smoothed mask to be measured.

 required
pruned_skeleton NDArray

Pruned skeleton.

 required
pixel_to_nm_scaling float

Scaling of pixels to nanometres.

 1

Returns:

Type Description
float

Width of grain in metres.
Source code in topostats\tracing\disordered_tracing.py 
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
@staticmethod
def calculate_dna_width(
    smoothed_mask: npt.NDArray, pruned_skeleton: npt.NDArray, pixel_to_nm_scaling: float = 1
) -> float:
    """
    Calculate the mean width in metres of the DNA using the trace and mask.

    Parameters
    ----------
    smoothed_mask : npt.NDArray
        Smoothed mask to be measured.
    pruned_skeleton : npt.NDArray
        Pruned skeleton.
    pixel_to_nm_scaling : float
        Scaling of pixels to nanometres.

    Returns
    -------
    float
        Width of grain in metres.
    """
    dist_trans = ndimage.distance_transform_edt(smoothed_mask)
    comb = np.where(pruned_skeleton == 1, dist_trans, 0)

    return comb[comb != 0].mean() * 2 * pixel_to_nm_scaling

re_add_holes(orig_mask: npt.NDArray, smoothed_mask: npt.NDArray, holearea_min_max: tuple[float | int | None] = (2, None)) -> npt.NDArray

Restore holes in masks that were occluded by dilation.

As Gaussian dilation smoothing methods can close holes in the original mask, this function obtains those holes (based on the general background being the first due to padding) and adds them back into the smoothed mask. When paired with smooth_mask, this essentially just smooths the outer edge of the mask.

Parameters:

Name Type Description Default
orig_mask NDArray

Original mask.

 required
smoothed_mask NDArray

Original mask but with inner and outer edged smoothed. The smoothing operation may have closed up important holes in the mask.

 required
holearea_min_max tuple[float | int | None]

Tuple of minimum and maximum hole area (in nanometers) to replace from the original mask into the smoothed mask.

 (2, None)

Returns:

Type Description
NDArray

Smoothed mask with holes restored.
Source code in topostats\tracing\disordered_tracing.py 
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
def re_add_holes(
    self,
    orig_mask: npt.NDArray,
    smoothed_mask: npt.NDArray,
    holearea_min_max: tuple[float | int | None] = (2, None),
) -> npt.NDArray:
    """
    Restore holes in masks that were occluded by dilation.

    As Gaussian dilation smoothing methods can close holes in the original mask, this function obtains those holes
    (based on the general background being the first due to padding) and adds them back into the smoothed mask. When
    paired with ``smooth_mask``, this essentially just smooths the outer edge of the mask.

    Parameters
    ----------
    orig_mask : npt.NDArray
        Original mask.
    smoothed_mask : npt.NDArray
        Original mask but with inner and outer edged smoothed. The smoothing operation may have closed up important
        holes in the mask.
    holearea_min_max : tuple[float | int | None]
        Tuple of minimum and maximum hole area (in nanometers) to replace from the original mask into the smoothed
        mask.

    Returns
    -------
    npt.NDArray
        Smoothed mask with holes restored.
    """
    # handle none's
    if set(holearea_min_max) == {None}:
        return smoothed_mask
    if None in holearea_min_max:
        none_index = holearea_min_max.index(None)
        holearea_min_max[none_index] = 0 if none_index == 0 else np.inf

    # obtain px holesizes
    holesize_min_px = holearea_min_max[0] / ((self.pixel_to_nm_scaling) ** 2)
    holesize_max_px = holearea_min_max[1] / ((self.pixel_to_nm_scaling) ** 2)

    # obtain a hole mask
    holes = 1 - orig_mask
    holes = label(holes)
    hole_sizes = [holes[holes == i].size for i in range(1, holes.max() + 1)]
    holes[holes == 1] = 0  # set background to 0 assuming it is the first hole seen (from top left)

    # remove too small or too big holes from mask
    for i, hole_size in enumerate(hole_sizes):
        if hole_size < holesize_min_px or hole_size > holesize_max_px:  # small holes may be fake are left out
            holes[holes == i + 1] = 0
    holes[holes != 0] = 1  # set correct sixe holes to 1

    # replace correct sized holes
    return np.where(holes == 1, 0, smoothed_mask)

remove_touching_edge(skeleton: npt.NDArray) -> npt.NDArray staticmethod

Remove any skeleton points touching the border (to prevent errors later).

Parameters:

Name Type Description Default
skeleton NDArray

A binary array where touching clusters of 1's become 0's if touching the edge of the array.

 required

Returns:

Type Description
NDArray

Skeleton without points touching the border.
Source code in topostats\tracing\disordered_tracing.py 
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
@staticmethod
def remove_touching_edge(skeleton: npt.NDArray) -> npt.NDArray:
    """
    Remove any skeleton points touching the border (to prevent errors later).

    Parameters
    ----------
    skeleton : npt.NDArray
        A binary array where touching clusters of 1's become 0's if touching the edge of the array.

    Returns
    -------
    npt.NDArray
        Skeleton without points touching the border.
    """
    for edge in [skeleton[0, :-1], skeleton[:-1, -1], skeleton[-1, 1:], skeleton[1:, 0]]:
        uniques = np.unique(edge)
        for i in uniques:
            skeleton[skeleton == i] = 0
    return skeleton

smooth_mask(grain: npt.NDArray, dilation_iterations: int = 2, gaussian_sigma: float | int = 2, holearea_min_max: tuple[int | float | None] = (0, None)) -> npt.NDArray

Smooth a grain mask based on the lower number of binary pixels added from dilation or gaussian.

This method ensures gaussian smoothing isn't too aggressive and covers / creates gaps in the mask.

Parameters:

Name Type Description Default
grain NDArray

Numpy array of the grain mask.

 required
dilation_iterations int

Number of times to dilate the grain to smooth it. Default is 2.

 2
gaussian_sigma float | None

Gaussian sigma value to smooth the grains after an Otsu threshold. If None, defaults to 2.

 2
holearea_min_max tuple[float | int | None]

Tuple of minimum and maximum hole area (in nanometers) to replace from the original mask into the smoothed mask.

 (0, None)

Returns:

Type Description
NDArray

Numpy array of smoothed image.
Source code in topostats\tracing\disordered_tracing.py 
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
def smooth_mask(
    self,
    grain: npt.NDArray,
    dilation_iterations: int = 2,
    gaussian_sigma: float | int = 2,
    holearea_min_max: tuple[int | float | None] = (0, None),
) -> npt.NDArray:
    """
    Smooth a grain mask based on the lower number of binary pixels added from dilation or gaussian.

    This method ensures gaussian smoothing isn't too aggressive and covers / creates gaps in the mask.

    Parameters
    ----------
    grain : npt.NDArray
        Numpy array of the grain mask.
    dilation_iterations : int
        Number of times to dilate the grain to smooth it. Default is 2.
    gaussian_sigma : float | None
        Gaussian sigma value to smooth the grains after an Otsu threshold. If None, defaults to 2.
    holearea_min_max : tuple[float | int | None]
        Tuple of minimum and maximum hole area (in nanometers) to replace from the original mask into the smoothed
        mask.

    Returns
    -------
    npt.NDArray
        Numpy array of smoothed image.
    """
    # Option to disable the smoothing (i.e. U-Net masks are already smooth)
    if dilation_iterations is None and gaussian_sigma is None:
        LOGGER.debug(f"[{self.filename}] : no grain smoothing done")
        return grain

    # Option to only do gaussian or dilation
    if dilation_iterations is not None:
        dilation = ndimage.binary_dilation(grain, iterations=dilation_iterations).astype(np.int32)
    else:
        gauss = filters.gaussian(grain, sigma=gaussian_sigma)
        gauss = np.where(gauss > filters.threshold_otsu(gauss) * 1.3, 1, 0)
        gauss = gauss.astype(np.int32)
        LOGGER.debug(f"[{self.filename}] : smoothing done by gaussian {gaussian_sigma}")
        return self.re_add_holes(grain, gauss, holearea_min_max)
    if gaussian_sigma is not None:
        gauss = filters.gaussian(grain, sigma=gaussian_sigma)
        gauss = np.where(gauss > filters.threshold_otsu(gauss) * 1.3, 1, 0)
        gauss = gauss.astype(np.int32)
    else:
        LOGGER.debug(f"[{self.filename}] : smoothing done by dilation {dilation_iterations}")
        return self.re_add_holes(grain, dilation, holearea_min_max)

    # Competition option between dilation and gaussian mask differences wrt original grains
    if abs(dilation.sum() - grain.sum()) > abs(gauss.sum() - grain.sum()):
        LOGGER.debug(f"[{self.filename}] : smoothing done by gaussian {gaussian_sigma}")
        return self.re_add_holes(grain, gauss, holearea_min_max)
    LOGGER.debug(f"[{self.filename}] : smoothing done by dilation {dilation_iterations}")
    return self.re_add_holes(grain, dilation, holearea_min_max)

trace_dna()

Perform the DNA skeletonisation and cleaning pipeline.

Source code in topostats\tracing\disordered_tracing.py 
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
def trace_dna(self):
    """Perform the DNA skeletonisation and cleaning pipeline."""
    self.smoothed_mask = self.smooth_mask(self.mask, **self.mask_smoothing_params)
    if check_pixel_touching_edge(self.smoothed_mask):
        LOGGER.warning(
            f"[{self.filename}] : Grain {self.n_grain} skipped as padding is too small. Please consider "
            "increasing the value of grains.grain_crop_padding in your config file and try again."
        )
        self.disordered_trace = None
        return
    self.skeleton = getSkeleton(
        self.image,
        self.smoothed_mask,
        method=self.skeletonisation_params["method"],
        height_bias=self.skeletonisation_params["height_bias"],
    ).get_skeleton()
    self.pruned_skeleton = prune_skeleton(
        self.image, self.skeleton, self.pixel_to_nm_scaling, **self.pruning_params.copy()
    )
    self.pruned_skeleton = self.remove_touching_edge(self.pruned_skeleton)
    self.disordered_trace = np.argwhere(self.pruned_skeleton == 1)

    if self.disordered_trace is None:
        LOGGER.warning(f"[{self.filename}] : Grain {self.n_grain} failed to Skeletonise.")
        self.disordered_trace = None
    elif len(self.disordered_trace) < self.min_skeleton_size:
        LOGGER.warning(f"[{self.filename}] : Grain {self.n_grain} skeleton < {self.min_skeleton_size}, skipping.")
        self.disordered_trace = None

check_pixel_touching_edge(mask: npt.NDArray) -> bool

Check if any pixels in a mask touch the edge of the image.

Parameters:

Name Type Description Default
mask NDArray

Numpy array of mask to be checked.

 required

Returns:

Type Description
bool

True or False if a pixel is found on the edge of the image.
Source code in topostats\tracing\disordered_tracing.py 
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
def check_pixel_touching_edge(mask: npt.NDArray) -> bool:
    """
    Check if any pixels in a mask touch the edge of the image.

    Parameters
    ----------
    mask : npt.NDArray
        Numpy array of mask to be checked.

    Returns
    -------
    bool
        True or False if a pixel is found on the edge of the image.
    """
    return mask[:, 0].any() or mask[:, -1].any() or mask[0, :].any() or mask[-1, :].any()

crop_array(array: npt.NDArray, bounding_box: tuple, pad_width: int = 0) -> npt.NDArray

Crop an array.

Ideally we pad the array that is being cropped so that we have heights outside of the grains bounding box. However, in some cases, if a grain is near the edge of the image scan this results in requesting indexes outside of the existing image. In which case we get as much of the image padded as possible.

Parameters:

Name Type Description Default
array NDArray

2D Numpy array to be cropped.

 required
bounding_box Tuple

Tuple of coordinates to crop, should be of form (min_row, min_col, max_row, max_col).

 required
pad_width int

Padding to apply to bounding box.

 0

Returns:

Type Description
NDArray()

Cropped array.
Source code in topostats\tracing\disordered_tracing.py 
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
def crop_array(array: npt.NDArray, bounding_box: tuple, pad_width: int = 0) -> npt.NDArray:
    """
    Crop an array.

    Ideally we pad the array that is being cropped so that we have heights outside of the grains bounding box. However,
    in some cases, if a grain is near the edge of the image scan this results in requesting indexes outside of the
    existing image. In which case we get as much of the image padded as possible.

    Parameters
    ----------
    array : npt.NDArray
        2D Numpy array to be cropped.
    bounding_box : Tuple
        Tuple of coordinates to crop, should be of form (min_row, min_col, max_row, max_col).
    pad_width : int
        Padding to apply to bounding box.

    Returns
    -------
    npt.NDArray()
        Cropped array.
    """
    bounding_box = list(bounding_box)
    bounding_box = pad_bounding_box(array.shape, bounding_box, pad_width)
    return array[
        bounding_box[0] : bounding_box[2],
        bounding_box[1] : bounding_box[3],
    ]

disordered_trace_grain(cropped_image: npt.NDArray, cropped_mask: npt.NDArray, pixel_to_nm_scaling: float, mask_smoothing_params: dict, skeletonisation_params: dict, pruning_params: dict, filename: str = None, min_skeleton_size: int = 10, n_grain: int = None) -> dict

Trace an individual grain.

Tracing involves multiple steps...

  1. Skeletonisation
  2. Pruning of side branches (artefacts from skeletonisation).
  3. Ordering of the skeleton.

Parameters:

Name Type Description Default
cropped_image NDArray

Cropped array from the original image defined as the bounding box from the labelled mask.

 required
cropped_mask NDArray

Cropped array from the labelled image defined as the bounding box from the labelled mask. This should have been converted to a binary mask.

 required
pixel_to_nm_scaling float

Pixel to nm scaling.

 required
mask_smoothing_params dict

Dictionary of parameters to smooth the grain mask for better quality skeletonisation results. Contains a gaussian 'sigma' and number of dilation iterations.

 required
skeletonisation_params dict

Dictionary of skeletonisation parameters, options are 'zhang' (scikit-image) / 'lee' (scikit-image) / 'thin' (scikitimage) or 'topostats' (original TopoStats method).

 required
pruning_params dict

Dictionary of pruning parameters.

 required
filename str

File being processed.

 None
min_skeleton_size int

Minimum size of grain in pixels after skeletonisation.

 10
n_grain int

Grain number being processed.

 None

Returns:

Type Description
dict

Dictionary of the contour length, whether the image is circular or linear, the end-to-end distance and an array of coordinates.
Source code in topostats\tracing\disordered_tracing.py 
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
def disordered_trace_grain(  # pylint: disable=too-many-arguments
    cropped_image: npt.NDArray,
    cropped_mask: npt.NDArray,
    pixel_to_nm_scaling: float,
    mask_smoothing_params: dict,
    skeletonisation_params: dict,
    pruning_params: dict,
    filename: str = None,
    min_skeleton_size: int = 10,
    n_grain: int = None,
) -> dict:
    """
    Trace an individual grain.

    Tracing involves multiple steps...

    1. Skeletonisation
    2. Pruning of side branches (artefacts from skeletonisation).
    3. Ordering of the skeleton.

    Parameters
    ----------
    cropped_image : npt.NDArray
        Cropped array from the original image defined as the bounding box from the labelled mask.
    cropped_mask : npt.NDArray
        Cropped array from the labelled image defined as the bounding box from the labelled mask. This should have been
        converted to a binary mask.
    pixel_to_nm_scaling : float
        Pixel to nm scaling.
    mask_smoothing_params : dict
        Dictionary of parameters to smooth the grain mask for better quality skeletonisation results. Contains
        a gaussian 'sigma' and number of dilation iterations.
    skeletonisation_params : dict
        Dictionary of skeletonisation parameters, options are 'zhang' (scikit-image) / 'lee' (scikit-image) / 'thin'
        (scikitimage) or 'topostats' (original TopoStats method).
    pruning_params : dict
        Dictionary of pruning parameters.
    filename : str
        File being processed.
    min_skeleton_size : int
        Minimum size of grain in pixels after skeletonisation.
    n_grain : int
        Grain number being processed.

    Returns
    -------
    dict
        Dictionary of the contour length, whether the image is circular or linear, the end-to-end distance and an array
        of coordinates.
    """
    disorderedtrace = disorderedTrace(
        image=cropped_image,
        mask=cropped_mask,
        filename=filename,
        pixel_to_nm_scaling=pixel_to_nm_scaling,
        min_skeleton_size=min_skeleton_size,
        mask_smoothing_params=mask_smoothing_params,
        skeletonisation_params=skeletonisation_params,
        pruning_params=pruning_params,
        n_grain=n_grain,
    )

    disorderedtrace.trace_dna()

    if disorderedtrace.disordered_trace is None:
        return None

    return {
        "image": cropped_image,
        "grain": cropped_mask,
        "smoothed_mask": disorderedtrace.smoothed_mask,
        "skeleton": disorderedtrace.skeleton,
        "pruned_skeleton": disorderedtrace.pruned_skeleton,
        "branch_types": get_skan_image(cropped_image, disorderedtrace.pruned_skeleton, "branch_type"),
        "branch_indexes": get_skan_image(cropped_image, disorderedtrace.pruned_skeleton, "node_id_src"),
    }

find_connections(row: pd.Series, skan_df: pd.DataFrame) -> str

Compile the neighbouring branch indexes of the row.

Parameters:

Name Type Description Default
row Series

A row from the Skan summarize dataframe.

 required
skan_df DataFrame

The statistics DataFrame produced by Skan's summarize function.

 required

Returns:

Type Description
str

A string representation of a list of matching row indices where the node src and dst columns match that of the rows. String is needed for csv compatibility since csvs can't hold lists.
Source code in topostats\tracing\disordered_tracing.py 
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
def find_connections(row: pd.Series, skan_df: pd.DataFrame) -> str:
    """
    Compile the neighbouring branch indexes of the row.

    Parameters
    ----------
    row : pd.Series
        A row from the Skan summarize dataframe.
    skan_df : pd.DataFrame
        The statistics DataFrame produced by Skan's `summarize` function.

    Returns
    -------
    str
        A string representation of a list of matching row indices where the node src and dst
        columns match that of the rows.
        String is needed for csv compatibility since csvs can't hold lists.
    """
    connections = skan_df[
        (skan_df["node_id_src"] == row["node_id_src"])
        | (skan_df["node_id_dst"] == row["node_id_dst"])
        | (skan_df["node_id_src"] == row["node_id_dst"])
        | (skan_df["node_id_dst"] == row["node_id_src"])
    ].index.tolist()

    # Remove the index of the current row itself from the list of connections
    connections.remove(row.name)
    return str(connections)

get_skan_image(original_image: npt.NDArray, pruned_skeleton: npt.NDArray, skan_column: str) -> npt.NDArray

Label each branch with it's Skan branch type label.

Branch types (+1 compared to Skan docs) are defined as: 1 = Endpoint-to-endpoint (isolated branch) 2 = Junction-to-endpoint 3 = Junction-to-junction 4 = Isolated cycle

Parameters:

Name Type Description Default
original_image NDArray

Height image from which the pruned skeleton is derived from.

 required
pruned_skeleton NDArray

Single pixel thick skeleton mask.

 required
skan_column str

A column from Skan's summarize function to colour the branch segments with.

 required

Returns:

Type Description
NDArray

2D array where the background is 0, and skeleton branches label as their Skan branch type.
Source code in topostats\tracing\disordered_tracing.py 
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
def get_skan_image(original_image: npt.NDArray, pruned_skeleton: npt.NDArray, skan_column: str) -> npt.NDArray:
    """
    Label each branch with it's Skan branch type label.

    Branch types (+1 compared to Skan docs) are defined as:
    1 = Endpoint-to-endpoint (isolated branch)
    2 = Junction-to-endpoint
    3 = Junction-to-junction
    4 = Isolated cycle

    Parameters
    ----------
    original_image : npt.NDArray
        Height image from which the pruned skeleton is derived from.
    pruned_skeleton : npt.NDArray
        Single pixel thick skeleton mask.
    skan_column : str
        A column from Skan's summarize function to colour the branch segments with.

    Returns
    -------
    npt.NDArray
        2D array where the background is 0, and skeleton branches label as their Skan branch type.
    """
    branch_field_image = np.zeros_like(original_image)
    skeleton_image = np.where(pruned_skeleton == 1, original_image, 0)
    try:
        skan_skeleton = skan.Skeleton(skeleton_image, spacing=1e-9, value_is_height=True)
        res = skan.summarize(skan_skeleton, separator="_")
        for i, branch_field in enumerate(res[skan_column]):
            path_coords = skan_skeleton.path_coordinates(i)
            if skan_column == "node_id_src":
                branch_field = i
            branch_field_image[path_coords[:, 0], path_coords[:, 1]] = branch_field + 1
    except ValueError:  # when no skeleton to skan
        LOGGER.warning("Skeleton has been pruned out of existence.")

    return branch_field_image

grain_anchor(array_shape: tuple, bounding_box: list, pad_width: int) -> list

Extract anchor (min_row, min_col) from labelled regions and align individual traces over the original image.

Parameters:

Name Type Description Default
array_shape tuple

Shape of original array.

 required
bounding_box list

A list of region properties returned by 'skimage.measure.regionprops()'.

 required
pad_width int

Padding for image.

 required

Returns:

Type Description
list(Tuple)

A list of tuples of the min_row, min_col of each bounding box.
Source code in topostats\tracing\disordered_tracing.py 
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
def grain_anchor(array_shape: tuple, bounding_box: list, pad_width: int) -> list:
    """
    Extract anchor (min_row, min_col) from labelled regions and align individual traces over the original image.

    Parameters
    ----------
    array_shape : tuple
        Shape of original array.
    bounding_box : list
        A list of region properties returned by 'skimage.measure.regionprops()'.
    pad_width : int
        Padding for image.

    Returns
    -------
    list(Tuple)
        A list of tuples of the min_row, min_col of each bounding box.
    """
    bounding_coordinates = pad_bounding_box(array_shape, bounding_box, pad_width)
    return (bounding_coordinates[0], bounding_coordinates[1])

pad_bounding_box(array_shape: tuple, bounding_box: list, pad_width: int) -> list

Pad coordinates, if they extend beyond image boundaries stop at boundary.

Parameters:

Name Type Description Default
array_shape tuple

Shape of original image (row, columns).

 required
bounding_box list

List of coordinates 'min_row', 'min_col', 'max_row', 'max_col'.

 required
pad_width int

Cells to pad arrays by.

 required

Returns:

Type Description
list

List of padded coordinates.
Source code in topostats\tracing\disordered_tracing.py 
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
def pad_bounding_box(array_shape: tuple, bounding_box: list, pad_width: int) -> list:
    """
    Pad coordinates, if they extend beyond image boundaries stop at boundary.

    Parameters
    ----------
    array_shape : tuple
        Shape of original image (row, columns).
    bounding_box : list
        List of coordinates 'min_row', 'min_col', 'max_row', 'max_col'.
    pad_width : int
        Cells to pad arrays by.

    Returns
    -------
    list
       List of padded coordinates.
    """
    # Top Row : Make this the first column if too close
    bounding_box[0] = 0 if bounding_box[0] - pad_width < 0 else bounding_box[0] - pad_width
    # Left Column : Make this the first column if too close
    bounding_box[1] = 0 if bounding_box[1] - pad_width < 0 else bounding_box[1] - pad_width
    # Bottom Row : Make this the last row if too close
    bounding_box[2] = array_shape[0] if bounding_box[2] + pad_width > array_shape[0] else bounding_box[2] + pad_width
    # Right Column : Make this the last column if too close
    bounding_box[3] = array_shape[1] if bounding_box[3] + pad_width > array_shape[1] else bounding_box[3] + pad_width
    return bounding_box

segment_heights(row: pd.Series, skan_skeleton: skan.Skeleton, image: npt.NDArray) -> npt.NDArray

Obtain an ordered list of heights from the skan defined skeleton segment.

Parameters:

Name Type Description Default
row Series

A row from the Skan summarize dataframe.

 required
skan_skeleton Skeleton

The graphical representation of the skeleton produced by Skan.

 required
image NDArray

The image the skeleton was produced from.

 required

Returns:

Type Description
NDArray

Heights along the segment, naturally ordered by Skan.
Source code in topostats\tracing\disordered_tracing.py 
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
def segment_heights(row: pd.Series, skan_skeleton: skan.Skeleton, image: npt.NDArray) -> npt.NDArray:
    """
    Obtain an ordered list of heights from the skan defined skeleton segment.

    Parameters
    ----------
    row : pd.Series
        A row from the Skan summarize dataframe.
    skan_skeleton : skan.Skeleton
        The graphical representation of the skeleton produced by Skan.
    image : npt.NDArray
        The image the skeleton was produced from.

    Returns
    -------
    npt.NDArray
        Heights along the segment, naturally ordered by Skan.
    """
    coords = skan_skeleton.path_coordinates(row.name)
    return image[coords[:, 0], coords[:, 1]]

segment_middles(row: pd.Series, skan_skeleton: skan.csr.Skeleton, image: npt.NDArray) -> float

Obtain the pixel value in the middle of the ordered segment.

Parameters:

Name Type Description Default
row Series

A row from the Skan summarize dataframe.

 required
skan_skeleton Skeleton

The graphical representation of the skeleton produced by Skan.

 required
image NDArray

The image the skeleton was produced from.

 required

Returns:

Type Description
float

The single or mean pixel value corresponding to the middle coordinate(s) of the segment.
Source code in topostats\tracing\disordered_tracing.py 
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
def segment_middles(row: pd.Series, skan_skeleton: skan.csr.Skeleton, image: npt.NDArray) -> float:
    """
    Obtain the pixel value in the middle of the ordered segment.

    Parameters
    ----------
    row : pd.Series
        A row from the Skan summarize dataframe.
    skan_skeleton : skan.csr.Skeleton
        The graphical representation of the skeleton produced by Skan.
    image : npt.NDArray
        The image the skeleton was produced from.

    Returns
    -------
    float
        The single or mean pixel value corresponding to the middle coordinate(s) of the segment.
    """
    heights = segment_heights(row, skan_skeleton, image)
    middle_idx, middle_remainder = (len(heights) + 1) // 2 - 1, (len(heights) + 1) % 2
    return heights[[middle_idx, middle_idx + middle_remainder]].mean()

trace_image_disordered(topostats_object: TopoStats, class_index: int, min_skeleton_size: int, mask_smoothing_params: dict, skeletonisation_params: dict, pruning_params: dict) -> None

Processor function for tracing image.

Parameters:

Name Type Description Default
topostats_object TopoStats

TopoStats object post ''Grains'' so that there are ''GrainCrops'' nested within.

 required
class_index int

Index of the class to trace.

 required
min_skeleton_size int

Minimum size of grain in pixels after skeletonisation.

 required
mask_smoothing_params dict

Dictionary of parameters to smooth the grain mask for better quality skeletonisation results. Contains a gaussian 'sigma' and number of dilation iterations.

 required
skeletonisation_params dict

Dictionary of options for skeletonisation, options are 'zhang' (scikit-image) / 'lee' (scikit-image) / 'thin' (scikitimage) or 'topostats' (original TopoStats method).

 required
pruning_params dict

Dictionary of options for pruning.

 required
Source code in topostats\tracing\disordered_tracing.py 
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
def trace_image_disordered(  # pylint: disable=too-many-arguments,too-many-locals, too-many-statements
    topostats_object: TopoStats,
    class_index: int,
    min_skeleton_size: int,
    mask_smoothing_params: dict,
    skeletonisation_params: dict,
    pruning_params: dict,
) -> None:
    """
    Processor function for tracing image.

    Parameters
    ----------
    topostats_object : TopoStats
        TopoStats object post ''Grains'' so that there are ''GrainCrops'' nested within.
    class_index : int
        Index of the class to trace.
    min_skeleton_size : int
        Minimum size of grain in pixels after skeletonisation.
    mask_smoothing_params : dict
        Dictionary of parameters to smooth the grain mask for better quality skeletonisation results. Contains
        a gaussian 'sigma' and number of dilation iterations.
    skeletonisation_params : dict
        Dictionary of options for skeletonisation, options are 'zhang' (scikit-image) / 'lee' (scikit-image) / 'thin'
        (scikitimage) or 'topostats' (original TopoStats method).
    pruning_params : dict
        Dictionary of options for pruning.
    """
    img_base = np.zeros_like(topostats_object.image)
    disordered_trace_crop_data = {}
    grainstats_additions = {}
    disordered_tracing_stats = pd.DataFrame()
    grain_crops = topostats_object.grain_crops
    filename = topostats_object.filename
    # These are images for diagnostics, edited during tracing to show
    # various steps
    all_images = {
        "smoothed_mask": img_base.copy(),
        "skeleton": img_base.copy(),
        "pruned_skeleton": img_base.copy(),
        "branch_indexes": img_base.copy(),
        "branch_types": img_base.copy(),
    }
    # for cropped_image_index, cropped_image in cropped_images.items():
    number_of_grains = len(grain_crops)
    for grain_number, grain_crop in grain_crops.items():
        try:
            grain_crop_tensor = grain_crop.mask
            grain_crop_class_mask = grain_crop_tensor[:, :, class_index]

            disordered_trace_images: dict | None = disordered_trace_grain(
                cropped_image=grain_crop.image,
                cropped_mask=grain_crop_class_mask,
                pixel_to_nm_scaling=grain_crop.pixel_to_nm_scaling,
                filename=filename,
                mask_smoothing_params=mask_smoothing_params,
                skeletonisation_params=skeletonisation_params,
                pruning_params=pruning_params,
                min_skeleton_size=min_skeleton_size,
                n_grain=grain_number,
            )
            LOGGER.debug(f"[{grain_crop.filename}] : Disordered Traced grain {grain_number + 1} of {number_of_grains}")
            if disordered_trace_images is not None:
                # obtain segment stats
                try:
                    skan_skeleton = skan.Skeleton(
                        np.where(disordered_trace_images["pruned_skeleton"] == 1, grain_crop.image, 0),
                        spacing=grain_crop.pixel_to_nm_scaling,
                    )
                    skan_df = skan.summarize(skel=skan_skeleton, separator="_")
                    total_branch_length_nm = skan_df["branch_distance"].sum() * 1e-9
                except ValueError:
                    LOGGER.warning(
                        f"[{grain_crop.filename}] : Skeleton for grain {grain_number} has been pruned out of existence."
                    )
                    total_branch_length_nm = 0
                    skan_df = pd.DataFrame()
                disordered_tracing_stats = pd.concat((disordered_tracing_stats, skan_df))
                # obtain stats
                conv_pruned_skeleton = convolve_skeleton(disordered_trace_images["pruned_skeleton"])
                grainstats_additions[grain_number] = {
                    "image": grain_crop.filename,
                    "grain_endpoints": np.int64((conv_pruned_skeleton == 2).sum()),
                    "grain_junctions": np.int64((conv_pruned_skeleton == 3).sum()),
                    "total_branch_lengths": total_branch_length_nm,
                    "grain_width_mean": disorderedTrace.calculate_dna_width(
                        disordered_trace_images["smoothed_mask"],
                        disordered_trace_images["pruned_skeleton"],
                        grain_crop.pixel_to_nm_scaling,
                    )
                    * 1e-9,
                }
                # Instantiate a DisorderedTrace object, setting its attributes and add to the GrainCrop.disordered_trace
                # attribute in the right direction
                topostats_object.grain_crops[grain_number].disordered_trace = DisorderedTrace(
                    images=disordered_trace_images,
                    grain_endpoints=np.int64((conv_pruned_skeleton == 2).sum()),
                    grain_junctions=np.int64((conv_pruned_skeleton == 3).sum()),
                    total_branch_length=total_branch_length_nm,
                    grain_width_mean=disorderedTrace.calculate_dna_width(
                        disordered_trace_images["smoothed_mask"],
                        disordered_trace_images["pruned_skeleton"],
                        grain_crop.pixel_to_nm_scaling,
                    )
                    * 1e-9,
                )
                # Extract disordered tracing statistics using skan
                # ns-rse 2025-12-22 : feels crude to pull things out of skan_df to build a dictionary which we
                # ultimately convert back to pd.DataFrame later, should really augment skan_df with additional metrics
                topostats_object.grain_crops[grain_number].disordered_trace.stats = {}
                for index, _ in enumerate(skan_df.iterrows()):
                    topostats_object.grain_crops[grain_number].disordered_trace.stats[index] = {}
                    topostats_object.grain_crops[grain_number].disordered_trace.stats[index]["image"] = filename
                    topostats_object.grain_crops[grain_number].disordered_trace.stats[index][
                        "branch_distance"
                    ] = skan_skeleton.path_lengths()[index]
                    topostats_object.grain_crops[grain_number].disordered_trace.stats[index]["branch_type"] = np.int64(
                        skan_df["branch_type"]
                    )[index]
                    topostats_object.grain_crops[grain_number].disordered_trace.stats[index]["connected_segments"] = (
                        skan_df.apply(find_connections, axis=1, skan_df=skan_df)[index]
                    )
                    topostats_object.grain_crops[grain_number].disordered_trace.stats[index][
                        "mean_pixel_value"
                    ] = skan_skeleton.path_means()[index]
                    topostats_object.grain_crops[grain_number].disordered_trace.stats[index][
                        "stdev_pixel_value"
                    ] = skan_skeleton.path_stdev()[index]
                    # pylint: disable=cell-var-from-loop
                    topostats_object.grain_crops[grain_number].disordered_trace.stats[index]["min_value"] = (
                        skan_df.apply(
                            lambda x: segment_heights(x, skan_skeleton=skan_skeleton, image=grain_crop.image).min(),
                            axis=1,
                        )[index]
                    )
                    topostats_object.grain_crops[grain_number].disordered_trace.stats[index]["median_value"] = (
                        skan_df.apply(
                            lambda x: np.median(
                                segment_heights(x, skan_skeleton=skan_skeleton, image=grain_crop.image)
                            ),
                            axis=1,
                        )[index]
                    )
                    # pylint: enable=cell-var-from-loop
                    topostats_object.grain_crops[grain_number].disordered_trace.stats[index]["middle_value"] = (
                        skan_df.apply(segment_middles, skan_skeleton=skan_skeleton, image=grain_crop.image, axis=1)[
                            index
                        ]
                    )
                    topostats_object.grain_crops[grain_number].disordered_trace.stats[index]["basename"] = str(
                        topostats_object.img_path
                    )
                # remap the cropped images back onto the original, there are many image crops that we want to
                #  remap back onto the original image so we iterate over them, as passed by the function
                for image_name, full_diagnostic_image in all_images.items():
                    crop = disordered_trace_images[image_name]
                    bbox = grain_crop.bbox
                    full_diagnostic_image[bbox[0] : bbox[2], bbox[1] : bbox[3]] += crop
                disordered_trace_crop_data[f"grain_{grain_number}"] = disordered_trace_images
                disordered_trace_crop_data[f"grain_{grain_number}"]["bbox"] = grain_crop.bbox
                disordered_trace_crop_data[f"grain_{grain_number}"]["pad_width"] = grain_crop.padding
            else:
                # Not sure why this is required as the GrainCrop definition sets these to None by default but I
                # encountered errors complaining about the attribute note being set and this seems to solve that.
                topostats_object.grain_crops[grain_number].disordered_trace = None
            # Ensure topostats_object.full_image_plots is a dictionary so we can append full image plots, this will
            # _always_ overwrite what is there so if disordered tracing is being re-run the old images are lost from the
            # .topostats object.
            topostats_object.full_image_plots = all_images

        # when skeleton is too small, pruned to 0's, skan -> ValueError -> skipped
        except Exception as e:  # pylint: disable=broad-exception-caught
            LOGGER.error(  # pylint: disable=logging-not-lazy
                f"[{grain_crop.filename}] : Disordered tracing of grain {grain_number} failed. "
                f"Consider raising an issue on GitHub. Error: ",
                exc_info=e,
            )

        # convert stats dict to dataframe
        grainstats_additions_df = pd.DataFrame.from_dict(grainstats_additions, orient="index")
        # Set the name of the index column to be the grain number
        grainstats_additions_df.index.name = "grain_number"