Cosine-Pruned Medial Axis: A New Method for Isometric Equivariant and Noise-Free Medial Axis Extraction

We present the CPMA, a new method for medial axis pruning with noise robustness and equivariance to isometric transformations. The CPMA leverages the discrete cosine transform to create smooth versions of a shape <inline-formula> <tex-math notation="LaTeX">$\Omega $ </tex-math></inline-formula>. We use the smooth shapes to compute a score function <inline-formula> <tex-math notation="LaTeX">$\mathcal {F}_{\Omega }$ </tex-math></inline-formula> that filters out spurious branches from the medial axis of the original shape <inline-formula> <tex-math notation="LaTeX">$\Omega $ </tex-math></inline-formula>. Our method generalizes to <inline-formula> <tex-math notation="LaTeX">$n$ </tex-math></inline-formula>-dimensional shapes given the properties of the Discrete Cosine Transform. We extensively compare with state-of-the-art pruning methods to highlight the CPMA’s noise robustness and isometric equivariance. We conducted experiments using two 2D datasets — Kimia216 and Animal2000 — and one 3D dataset — the Groningen benchmark. We found that our pruning approach achieves competitive results and yields stable medial axes even in scenarios with significant contour perturbations.


Introduction
Shape analysis arises naturally in computer vision applications where geometric information plays an essential role.The shape of an object is a useful tool in fields such as: non-destructive reconstruction of archaeology and cultural heritage [49,51]; object classification and retrieval from large collections of images [23,39]; human action and pose recognition for gaming and entertainment [12,24]; environment sensing in robot navigation and planning [34,25]; and industry for automatic visual quality inspection of product defects [37,5].
We visually perceive shape as the collections of all the features that constitute an object.However, to perform computer-based shape analysis, one must rely on an accurate discrete mathematical representation of an object's shape.This representation should exhibit the same geometric and topological properties inherent to the shape itself.Accordingly, we can think of shape representation as a way to store the shape's information in a different format, which benefits speed, compactness, and efficiency.
Many authors have proposed a variety of shape representations such as voxel/pixel grids, point clouds, triangular meshes, medial axes, or signed distance functions [45,50,28,18,54].These representations differ greatly in their formulation, and aim to provide a method for extracting descriptive features from objects, while also preserving their appearance and geometric properties [9,1,48,17,29].However, these methods also have disadvantages that limit their application.For example, medial axis representations are highly sensitive to contour noise; voxel/pixel grids are inaccurate after isometric transformations; signed distance functions and triangular meshes are memory-consuming representations when high-frequency details of the shape want to be stored.
We focus this study on the medial axis, also called the topological skeleton.The medial axis represents shapes as a collection of one-dimensional curves that define the central axis (or backbone) of an object.It provides dimensionality reduction of the amount of data needed to represent an entire shape while preserving its topological structure.Moreover, the medial axis is a rotation equivariant shape representation because the medial axis of a rotated object should ideally be the rotated medial axis of the same object.The medial axis is also robust to small deformation, such as articulation, because of its graph-like structure.For instance, a human-like shape moving only its arm will not affect all of the points in the medial axis, only the connections between the arm's nodes.
Despite its advantages, the medial representation is extremely sensitive to noise on the object's contour [40,36].Even small amounts of noise can cause erroneous sections of the skeleton called spurious branches.Consequently, many medial axis extraction algorithms are equipped with a mechanism to avoid or remove these spurious branches.
There are two main strategies reported in the state-of-the-art to deal with this problem: prior smoothing of the curve representing the object's boundary, and pruning the spurious branches after the medial axis' computation.In the former, the smoothed boundary is obtained by removing small structures along the curve or surface.It is interesting to note that smoothing curves does not always result in a simplified skeleton [53,6].
Effective pruning techniques focus instead on criteria to evaluate the significance of individual medial axis branches.However, pruning often requires user-defined parameters that depend on the size and complexity of the object [43,7,44], making the pruning method domain-dependent.Moreover, some pruning strategies result in a violation of the equivariant property [36,40].As a result, medial axis pruning is still an open problem in computer vision, and this problem is in need of noise-robust methods that concurrently preserve the isometric equivariance of the medial axis.This paper presents a new method for medial axis pruning that employs mechanisms from the two aforementioned branch-removal strategies.Our method works by computing a controlled set of smoothed versions of the original shape via the discrete cosine transform.We combine these smoothed shapes' medial axis to create a score function that rates points and branches by their degree of importance.We use our score function to prune spurious branches while preserving the medial axis' ability to reconstruct the original object.Our method is robust to boundary noise and exhibits isometric equivariance.
We benchmark our approach on three datasets of 2D and 3D segmented objects.
We use the Kimia216 [42] and the Animal dataset [8] to evaluate 2D medial axis extraction.These two datasets provide a method to assess 2D medial axis extraction in the presence of intra-class variation.We also use the University of Groningen Benchmark [46,47,13] to evaluate our approach on 3D objects.Our results show that our approach achieves competitive results on isometric equivariance and noise robustness compared to the state-of-the-art.
The main contributions of this paper are summarized as follows: • We define a novel method to compute medial axes that are robust to several degrees of boundary noise without losing the capacity to reconstruct the original object.
• Our computation pipeline guarantees that the isometric equivariance of the medial axis is preserved.
• The definition of our score function allows for a medial axis pruning that is efficiently computed in parallel.

Related work
Many algorithms and strategies exist to extract the medial axis and simplify it when affected by contour noise.This section briefly reviews the most representative algorithms for medial axis computation and discusses their key advantages and disadvantages.

The Medial Axis
Blum [11] first introduced the medial axis as an analogy of a fire propagating with uniform velocity on a grass field.The field is assumed to have the form of a given shape.If one "sets fire" at all boundary points, the medial axis is the set of quench points.There are other equivalent definitions of the medial axis.In this work we use a geometric definition as follows: Definition 1. Medial Axis.Let Ω be a connected bounded domain in R n , and x, x two points such that x, x ∈ Ω.The medial axis of Ω is defined as all the points x where x is the center of a maximal ball B r of radius r that is inscribed inside Ω. Formally, The medial axis, together with the associated radius of the maximally inscribed ball, is called the medial axis transform (MAT(Ω)).The medial axis transform is a complete shape descriptor, meaning that it can be used to reconstruct the shape of the original domain.In some work, MA and MAT are also referred to as shape skeletonization.Figure 1 shows an example of a 2D shape and its medial axis as the center of maximal discs.In R 3 definition 1 may result in a 2-dimensional medial axis sometimes called the medial surface.We will restrict our examples and analysis to only one-dimensional medial axes.

Medial Axis Computation
There are three primary mechanisms to compute the MA: 1) layer by layer morphological erosion also called thinning methods, 2) extraction of the medial axis from the edges of the Voronoi diagram generated by the boundary points, and 3) detection of ridges in distance map of the boundary points.In digital spaces, only an approximation to the "true medial axis" can be extracted.
When using thinning methods [16,27,26,32], points which belong to Ω are deleted from the outer boundary first.Later, the deletion proceeds iteratively inside until it results in a single-pixel wide medial axis.The medial axis by thinning can be approximated in terms of erosion and opening morphological operations [55].Thinning methods are easy to implement in a discrete setting, but they are not robust to isometric transformations.
The most well-known algorithm for thinning skeletonization is perhaps the Zhang Suen [55] algorithm.However, other approaches have been developed using similar principles [52,26,32], mainly focused on parallel computation.
Another way to estimate the medial axis works by computing the Voronoi Diagram (VD) of a polygonal approximation of the object's contour.The contour is expressed as line segments in 2D or a polygonal mesh in 3D.The seed points for the VD are the vertices of the polygonal representation.The medial axis is then computed as the union of all of the edges e i j of the VD, such that the points i, and j are not neighbors in the polygonal approximation [33].
Voronoi skeletonization methods preserve the topology of Ω.However, a suitable polygonal approximation of an object is crucial to guarantee the medial axis' accu- racy.Noise in the boundary forms convex areas in the contour, which induce spurious branches on the medial axis.In general, the better the polygonal approximation, the closer the Voronoi skeleton will be to the real medial axis.Nevertheless, this is an expensive process, particularly for large and complex objects [36].
The most common methods to extract the medial axis are those based on the distance transform.Within these methods, the medial axis is computed as the ridges of the distance transform inside the object [3,40,2,22,14,35,41].This interpretation of the medial axis follows definition 1, because the centers of the maximal balls are located on points x along the ridgeline of the distance transform, and the radius of the balls correspond to the distance value at x.
When computed in a discrete framework, distance-based approaches suffer from the same isometric robustness limitations as thinning and Voronoi methods [36].Moreover, noise in the contour produced by a low discretization resolution directly affects the medial axis' computation by introducing artificial ridges in the distance transform and, consequently, spurious branches.

Medial Axis Simplification
The medial axis' sensitivity to boundary noise limits its applications to real-life problems [10].Even negligible boundary noise can cause spurious branches, as shown in Figure 2.
One strategy for removal of spurious branches consists of computing MA(Ω) as the medial axis of a smoother version of the shape, Ω [38,31].The main disadvantage of this approach is that, in most cases, the resulting medial axis is not a good approximation.Additionally, the smoothing of Ω can potentially change its topology.Miklos et al. [30] introduced a slightly different approach they call Scaled Axis Transform (SAT).
The SAT involves scaling the distance transform and computing the medial axis of the original un-scaled shape as the medial axis of the scaled one.However, in [35], the authors show that the SAT is not necessarily a subset of the medial axis of the original shape.Instead, they propose a solution that guarantees a better approximation by including additional constraints on the scaled distance transform.
Another method to overcome the noise sensitivity limitation of the medial axis is spurious branch pruning.Pruning methods are a family of regularization processes incorporated in most medial axis extraction algorithms [43].Effective pruning techniques focus on different criteria to evaluate the significance of medial axis branches.
Thus, the decision is whether to remove the branch (and its points) or not.We can say that pruning methods are adequate if the resulting MA is noticeably simplified while preserving the topology and its ability to reconstruct the original object.Most pruning methods rely on ad hoc heuristic rules, which are invented and often reinvented in a variety of equivalent application-driven formulations [43].Some authors apply these rules while computing all medial axis points.Others do so by removing branches that are considered useless after the computation [36,19,22].
One of the most popular pruning methods is Couprie et al. [14].They consider the angle θ formed by a point x ∈ Ω, and its two closest boundary points denoted by the set Π Ω (x).This solution removes points from MA for which θ is lower than a constant threshold.This criterion allows different scales within a shape but generally leads to an unconnected medial axis.Another pruning method found in the state-ofthe-art is the work of Hesselink et al. [22].They introduce the Gamma Integer Medial Axis (GIMA), where a point belongs to the medial axis if the distance between its two closest boundary points is at least equal to γ.
Many pruning methods rely on the distance transform D Ω (x).The distance transform acts as a generator function for the medial axis, such that points x ∈ MA if and only if they satisfy some constraint involving their distance to the boundary.However, other authors have proposed alternative generator functions in their pruning strategies [21,4].For example, [21] and [4] introduce what they call Poisson skeletons by approximating D Ω (x) as the solution of the Poisson equation with constant source function.
Poisson skeletons rely on a solid mathematical formulation.Among other concepts, they use the local minimums and maximums of the curvature of δΩ.However, when such methodology is applied in a discrete environment, many spurious branches appear due to the need to define the length of a kernel size to estimate these local extreme points.

Method
We propose a new pruning approach to remove spurious branches from the medial axis of a n-dimensional closed shape Ω.We call our method the Cosine-Pruned Medial Axis (CPMA).The CPMA works by filtering out points from the Euclidean Medial Axis MAT(Ω) with a score function F Ω (x) : N n → [0, 1].We define the function F Ω by aggregating the medial axis of controlled smoothed versions of Ω.Our formulation of F Ω must yield high values at points x that belong to the real medial axis and low values at points that belong to spurious branches.Additionally, we require F Ω to be equivariant to isometric transformations.

The Cosine-Pruned Medial Axis
Let us represent Ω as a square binary image I : N 2 → {0, 1} with a resolution of M × M pixels.We start the computation of the CPMA by estimating a set of smoothed versions of the I via the Discrete Cosine Transform (DCT) and its inverse (IDCT): where (u, v) are coordinates in the frequency domain, and (x, y) are the spatial coordinates of the Euclidean space where Ω is defined.The values of C u and C v are determined by: The DCT is closely related to the discrete Fourier transform of real valued-functions.
However, it has better energy compaction properties with only a few of the transform coefficients representing the majority of the energy.Multidimensional variants of the various DCT types follow directly from the one-dimensional definition.They are simply a separable product along each dimension.
Let us now denote by I (i) with i = 1, 2, ..., M the reconstructions of I using only the first i frequencies as per equation 2. We seek to obtain a F Ω acting as a sort of probability indicating how likely it is for a point x to be in the medial axis of Ω. Points on relevant branches will appear regularly in the smoothed shapes' medial axis, resulting in high score function values.In contrast, spurious branches will only appear occasionally, resulting in low values.
Definition 2. Score function.Let I : N n → {0, 1} be a square binary image such that Let I (i) also be the i-frequency reconstruction of I via the IDCT.
We define F Ω (x) : N n → [0, 1] as the per pixel average over a set of estimations of the MAT on the smoothed shapes I (i) .
The score function is defined for all x in the domain of I. Notice that we represent MAT as another binary image where MAT(x) = 1 only when x belongs to the medial axis.Finally, with F Ω , we have all the elements to present our definition of the CPMA.
Definition 3. Cosine-Pruned Medial Axis Given a binary image I : representing a shape Ω, the CPMA(Ω) consist of all the pairs (x, r) ∈ MAT(Ω) such that F Ω (x) is greater than a threshold τ: We empirically set the value of the threshold to τ = 0.47.However, we conducted an additional experiment to show that this value is consistent across different shapes.
Although the CPMA results in a noise-free MA, there is no restriction in its formulation to force the CPMA to create a connected medial axis.We solve this by finding all the disconnected pieces of the CPMA.Later, we connect them using a minimum distance criterion g(x i , x j ), where x i and x j are endpoints of two distinct pieces.However, neither the Euclidean distance nor the geodesic distance are suitable criteria because they lead to connections between nodes that do not follow the medial axis (See figure 4).We instead connect x i and x j with a minimum energy path using an energy function E Ω .We must guarantee that E Ω (x) has high values when x is close to δΩ and low values when x is close to the medial axis.This way, we enforce the paths to be close to the MAT.We call the result of this strategy the Connected CPMA (C-CPMA).In section 3.3, we provide details of E Ω computation.

Isometric Equivariance of the CPMA
The distance transform-based medial axis depends only on the shape Ω, not on the position or size in the embedding Euclidean space.Therefore the medial axis should be equivariant under isometric transformations.Let us now define y = R(x).Applying R to every element of MAT(Ω), we have: Moreover, the CPMA depends primarily on F Ω , which also holds the isometric equivariant property.
Corollary 1.Let F Ω be the score function of Ω as per definition 3, and let R be an isometric transformation.We say that Proof.Using the results from proposition 1 and recalling that R is a linear transformation, we have that: However, in a discrete domain, this equivariance is only an approximation because points on both Ω and MAT(Ω) are constrained to be on a fixed, regular grid.In a continuous domain, it is easy to demonstrate that the cosine transform has exact isometric equivariance.

Implementation Details
To compute the CPMA enforcing the connectivity, we create a lattice graph G = (V, E).A point p in the domain of I is a node of G, if and only if p ∈ Ω.The node p shares an edge with all its neighbors in the lattice only if the neighbors are also inside Ω.We used an 8-connectivity neighborhood in 2D and a 26-connectivity neighborhood in 3D.
In order to determine the minimum energy path between pairs of pixels/voxels, we compute the minimum distance path inside G using Dijkstra's algorithm.We assign weights to every edge with values extracted from E Ω .Given (x, y) ∈ E, we compute the energy of every edge as follows: This method guarantees connectivity, but it is inefficient because of the minimum energy path's iterative computation.We sacrifice performance in favor of connectivity.
We include the pseudo-code to compute the CPMA and the C-CPMA in algorithms 1 and 2.
Algorithm 1: Cosine-Pruned Medial Axis (CPMA) Input: I: N-dimensional binary array representing the object M: number of frequencies of I to be used in the computation Output: CPMA: Cosine-Pruned Medial Axis τ ← 0.5 // Reconstructs I using only the first i frequencies // The final F Ω is the average of all reconstructions The CPMA only relies on one parameter, τ.The value of τ is the threshold that determines whether a point of F Ω is a medial axis point.We empirically set the value of τ = 0.47.However, in section 5, we present the result of an additional experiment to show how sensitive the CPMA is to different threshold values.
Another essential consideration when computing the CPMA is the maximum frequency used to reconstruct the original shape through the IDCT.Using less than M frequencies enables a faster computation of the CPMA without losing accuracy.We found that using frequencies greater than M 2 does not yield significant improvement for

Experiments
In this section, we describe experiments used to evaluate our approach compared to state-of-the-art medial axis pruning methods.

Datasets
We chose three extensively used datasets of 2D and 3D segmented objects to evaluate our methodology on medial axis extraction robustness.These datasets are part of the accepted benchmarks in literature, enabling us to compare our results.
Kimia216 dataset [42].It consists of 18 classes of different shapes with 12 samples in each class.The dataset's images are a collection of slightly different views of a set of shapes with varying topology.Figure 5a shows two samples from each class.Contour noise and random rotations are also present in some images in the dataset.Kimia216 has been largely used to test a wide range of medial axis extraction algorithms.Because of the large variety of shapes, we assume that this benchmark ensured a fair comparison with the state-of-the-art.the contour/surface.We are also interested in assessing how stable the CPMA responds in the presence of rotations of Ω to test its isometric equivariance.
We employ the Hausdorff distance (d H ), and Dubuisson-Jain dissimilarity (d D ) as metrics between shapes.The Dubuisson-Jain similarity is a normalization of the Hausdorff distance [15], which aims to overcome d H sensitivity to outliers.The Dubuisson-Jain similarity between point sets X and Y is defined as: with We must first choose a strategy to induce noise to the boundary to evaluate the noise sensitivity.We use a stochastic approach denoted by E, where a random number of points p ∈ δΩ are deformed by a vector v in the direction orthogonal to the boundary, with a deformation magnitude that is normally distributed, |v| ∼ N(p, 1).This noise model is recursively applied n times to every shape in our datasets.We denote as MAT(E(Ω, k)) the medial axis of a shape Ω after applying the noise model k times.In our experiments, we used k = 1, ..., 20 noise levels.
For every object in each dataset, we compared the medial axis MAT(Ω) with the noisy versions MAT(E(Ω, k)) to determine how sensitive a method is to boundary noise.
The medial axis is ideally an isometric equivariant shape representation so that MA(R(Ω)) = R(MA(Ω)).Due to sampling factors, this relationship is only an approximation.However, we can measure the equivariance by comparing MAT(R(Ω)) with

R(MAT(Ω))
. The more similar they are, the more equivariant the method.
Because the translation equivariance is trivial, we evaluate isometric equivariance only with rotations of Ω.We do not use reflections because the properties of rotation matrices we use in this study can be extended to reflection matrices.In our experiments, 2D rotations are counterclockwise in the range [0, 90] degrees around the origin.We use 30 rotations at regular intervals.In 3D, we use a combination of azimuthal (θ ∈ [0, 90]) and elevation (φ ∈ [0, 90]) rotations around the origin.The angles θ and φ take values each 18 degrees.

Comparative Studies
We chose seven of the most relevant methods in the scientific literature to compare with CPMA extraction results.Each method was selected based on a careful review of the state-of-the-art.These methods illustrate the variety of approaches that authors employ to prune the medial axis.Notice that the first method we used in our study is the un-pruned MAT.
Table 1 summarizes all of the methods in our comparative study.In many cases, the performance of a pruning method depended on its parametrization.We selected parameters for all of the methods following the best performance parametrization reported in the state-of-the-art.

3D
Thinning [55] Zhang-Suen Algorithm N/A TEASAR [41] Tree-structure skeleton extraction This table shows name and parameter description for each method.The point x ∈ MAT is an element of the MAT that can be potentially pruned.Π Ω (x) refers to the set of closest boundary points of x.All the elements of Π Ω (x) have the same distance from x.

Results
In this section, we present and discuss the results of our approach on medial axis pruning.We focus on two properties: 1) robustness to noise of the contour and 2) isometric equivariance.

Stability Under Boundary Noise
We compared the stability of the CPMA under boundary noise against other approaches in Table 1.We conducted our experiments on Kimia216 and the Animal2000 Dataset for 2D images.Additionally, we used a set of three-dimensional triangular meshes from the Groningen Benchmark for 3D experimentation.
For our noise sensitivity experiments, we applied 20 times the noise model E(Ω, k) to every object of each dataset.We then computed their MAT using every method in Table 1 with different parameters.Finally, each MAT(E(Ω, k)) was compared with MAT(Ω) using both the Hausdorff distance and Dubuisson-Jain dissimilarity.We report the per method average of each metric over all the elements of each dataset.
First, we tested our medial axis pruning method on the Kimia216 dataset and present the results in Table 2.The table shows that the CPMA and the C-CPMA are competitive against state-of-the-art methods for medial axis extraction.Our results show similar performance to two state-of-the-art methods: the GIMA and SFEMA.
The CPMA and C-CPMA also performed better than Poisson Skeletons, SAT, topological thinning, and the MAT itself.For visual comparison, Figure 6  It is clear that the CPMA and the C-CPMA are among the three best results when we use the Dubuisson-Jain dissimilarity metric.However, we observe a decrease in performance when we the Hausdorff distance metric.We believe this occurs because of Hausdorff's distance sensitivity to outliers.
The Animal2000 dataset contains nearly ten times more shapes than Kimia216.
This leads to more variation between shapes, and therefore a more challenging setting.Table 3 shows similar results compared to Kimia216, confirming that the noise invariant properties of the CPMA still hold in a more robust dataset.The GIMA and the SFEMA are still the best methods measured with the Dubuisson-Jain dissimilarity, closely followed by both the CPMA and the C-CPMA.Results of using the Dubuisson-Jain dissimilarity as a metric show that the CPMA is close to methods such as BEMA The table shows the average Hausdorff distance and Dubuisson-Jain dissimilarity for different noise levels (5)(6)(7)(8)(9)(10)(11)(12)(13)(14)(15)(16)(17)(18)(19)(20) over each element of the dataset.
and SFEMA.However, the results are not as good as when using the Hausdorff distance metric.Figure 6 (middle row) depicts the best performance for every method in comparison to ours.Our experiment's results suggest that the CPMA noise invariant properties generalize across different datasets.
For our 3D experiments, we selected 14 objects from the Groningen Benchmark, reflecting the most common shapes used in the literature.Each object was voxelized to a binary voxel grid with resolution 150 × 150 × 150.This resolution offered sufficient details as well as a sufficiently low computational cost.In contrast to the 2D case, we applied E(Ω, k) only 10 times to the 3D object.We did this for two reasons: 1) to reduce computational complexity and 2) noise tends to be more extreme in 3D at the chosen resolution.The results on the Groningen dataset are shown in Table 4 and among the other methods when compared with the dissimilarity measure.These results show that our methodology has noise-invariance properties, and it is stable in the presence of small surface deformation.However, the results show unusual patterns when compared with the Hausdorff distance.In fact, for some methods, the metric decreases  when the noise level increases.We attribute this behavior to the outlier sensibility of the Hausdorff distance.The table shows the average Hausdorff distance and Dubuisson-Jain dissimilarity for different noise levels (5)(6)(7)(8)(9)(10)(11)(12)(13)(14)(15)(16)(17)(18)(19)(20) over each element of the dataset.
We complete the noise stability analysis showing examples of the MAT computed with our methodology, and compared to the other methods in this study, Figure 7.

Sensitivity to Rotations
We measured the dissimilarity between MAT(R(Ω)) and R(MAT(Ω)) for different instances, and different definitions of the medial axis transform.The lower this dissimilarity, the more stable the method is under rotation.The rotation sensitivity analysis on the Kimia216 dataset is summarized in Table 5 and Figure 8 (top row).The results show that the curves of the CPMA and the C-CPMA fall near the average of the rest of the methods achieving state-of-the-art performance.The results also surpassed several methods, such as Poisson skeleton, SAT, and thinning methods.Notice that when using the dissimilarity metric the CPMA, the GIMA, the SFEMA, and the BEMA form a subgroup that performs significantly better compared to the others.Moreover, the performance of these methods oscillates around a value of dissimilarity of around 1 pixel on average.The intuition for this result is that regardless of the rotation, skeletons computed with these methods vary only at one pixel.Consequently, we can claim that they exhibit rotation equivariance.We applied the same analysis to the Animal2000 dataset achieving similar results.
In this case, the CPMA and the C-CPMA ranked third and fourth, respectively, among all methods when we used the dissimilarity metric.The results for all methods and parameters are presented in Table 6.As before, we also present a summary with the best parametrization for each method in Figure 8 (middle row) to facilitate the interpretation.Notice that due to the larger number of objects in the Animal2000 dataset, the curves for every method appear to be smoother, highlighting stability across different rotation angles and shapes.
Finally, we conducted the rotation sensitivity analysis on the 3D dataset.the results are summarized in Figure 8 (bottom row).The image shows the four 3D medial axis extraction methods we compared in our study for combinations of azimuthal and elevation angles.This figure shows how both the Hausdorff distance and the Dubuisson-Jain dissimilarity become higher when the rotation becomes more extreme, except in the  5 and 6.
case of C-CPMA.We believe this behavior is due to the connectivity enforcement mitigating the gaps in the medial axis, and reducing the metrics.

Hyper-parameter Selection
We tested our medial axis pruning method on the Kimia216 dataset.Many medial axis pruning methods depend on hyper-parameters to accurately estimate the medial axis [22,14,20,35].These parameters usually have a physical meaning in the context of the object whose medial axis they seek to estimate.Often, the parameters are distances or angles formed between points inside the object.In other work, some authors create score functions like ours, intending to use its values as a filter parameter to remove points on spurious branches of the MAT.In most cases, however, such parameters are subject to factors like resolution or scale.Thus, we conducted another experiment to test the sensitivity of the CPMA to the pruning parameter τ at different scale factors of the input object.and branch pruning.Moreover, around this value, the standard deviation reaches its minimum value suggesting optimal performance regardless of the object.Because the value of τ is stable for different scale factors, we conclude that scale does not affect the selection of the threshold.

Conclusions and Future Work
We presented the CPMA, a new method for medial axis pruning with noise robustness and equivariance to isometric transformations.Our method leverages the discrete cosine transform to compute a score function that rates the importance of individual points and branches within the medial axis of a shape.
Our pruning approach delivers competitive results compared to the state-of-the-art.
The results of our experiments show that our method is robust to boundary noise.Additionally, it is also equivariant to isometric transformations, and it is capable of producing a stable and connected medial axis even in scenarios with significant perturbations

Figure 1 :
Figure 1: Medial Axis Transform Computation.(Left) a shape and its boundary.(Right) Medial Axis elements consisting of the centers and radius of balls inscribed in the shape [34].

Figure 2 :
Figure 2: (Left) Spurious branch in medial axis.(Right) A new branch appears in the presence of a small perturbation in the contour [43].

Figure 3 :
Figure 3: Score Function illustrative example.The rows show F Ω of an image of size on a 180 × 180.We computed the reconstructions up to only M 1 = 10 and M 2 = 62 of the first frequencies.

Figure 4 :
Figure4: Path connectivity between CPMA segments.When using the Euclidean distance (left), two nodes can connect through a path that is not contained within Ω.The geodesic distance (center) guarantees that the path is in Ω, but does not follow the center-line.The minimum energy distance (center) E Ω is a better alternative to enforce the path to follow the medial axis.

Proposition 1 .
Let MAT(Ω) be the medial axis transform of a connected bounded domain Ω embedded in R n , and let R(x) = Mx + b be an isometric transformation in R n .The square matrix M is a composition of any number or rotations and reflection matrices, and b is an n-dimensional vector.We say that MAT(Ω) is equivariant to any R such that MAT(R(Ω)) = R(MAT(Ω)).Proof.Recall that R is an isometric transformation and thus it is invertible and preserves the Euclidean distance.The previous statement implies that R is an isomorphic map between an open ball B r (x) and B r (R(x)).Consequently, if B r (R(x)) is a maximal ball in Ω then from definition 1 we have that B r (R(x)) B r (R(x )), ∀r > r.
(top row) shows both Hausdorff distance and Dubuisson-Jain dissimilarity against noise level.The figure only displays the best parametrization of every method to improve visualization.

Fig- ure 6 (
bottom row).Notice that both the CPMA and C-CPMA achieved the best results

Figure 6 :
Figure 6: Noise sensitivity results on Kimia216 dataset (top), Animal2000 dataset (middle), and Groningen Benchmark (bottom).The figure shows the Hausdorff distance (left) and the Dubuisson-Jain dissimilarity (right) for all of the methods in Tables 2, 3, and 4.Only the best parametrization of each method is depicted for better interpretation.

Figure 7 :
Figure 7: The images show the un-pruned MAT and the results of four different pruning methods.Rows one and two are objects from Kimia216 dataset.Rows three and four are objects from Animal2000.Rows five and six are objects from the Groningen Benchmark.Notice how the CPMA and the C-CPMA yield medial axes with less spurious branches while preserving the topology.

Figure 8 :
Figure 8: Rotation equivariance results on Kimia216 dataset (top), Animal2000 (middle), and Groningen Benchmark (bottom).The top row shows the Hausdorff distance and Dubuisson-Jain dissimilarity for all the methods in Tables5 and 6.

Figure 9
Figure 9 shows the average of the Jaccard index as the reconstruction metric vs the values of τ.We compared an object Ω against its reconstruction Ω over all images in Kimia216 dataset.The figure shows how high values of τ deteriorate the reconstruction, while lower values do not prune enough spurious branches.From the figure, we can infer that values around τ = 0.47 offer a good trade-off between reconstruction

Figure 9 :
Figure 9: Sensitivity analysis of threshold τ to different scales of the input images.The graph shows the average Jaccard index of the reconstructed shape w.r.t the original object for the CPMA computed with different values of τ.Higher values of the threshold lead to less spurious branches.We also show the standard deviation error bands.

Table 1 :
Pruning methods employed for the comparative study in 2D

Table 4 :
Noise sensitivity results on Groningen Benchmark.
Rotation equivariance results on Kimia216.The table shows the average Hausdorff distance and Dubuisson-Jain dissimilarity for different rotations of each element in the dataset.

Table 6 :
Rotation equivariance results on Animal2000.Rotation equivariance results on Animal2000.The table shows the average Hausdorff distance and Dubuisson-Jain dissimilarity for different rotations of each element in the dataset.