Mandelbrot popularised fractal geometry in the 1980s26. Since then, various fractal analyses have been described in the literature. Amongst them, the box counting method has been implemented to determine the fractal dimension of shorelines, vegetation, river networks and fracture patterns36 and to calculate the area of irregular cartographic features37,38. Given its applicability with equal effectiveness to point sets, linear features, areas, and volumes, the box counting method is widely used for determining fractal dimensions36. We are aware of only one publication where fractal dimension has been used for palatal rugae39, however only for the purpose of image registration and not for the estimation of the complexity of rugae in a systematic and well-defined manner.
Studies that have assessed rugae have used methods which are descriptive of the size (described with ordinal variables), shape, branching or direction of individual rugae, and not the entire set. In contrast, the presented method is capable of measuring the complexity of either a single ruga, a selection of rugae or the entire set of a palate’s rugae, depending on which structures are preserved at the steps of Initial and Final cropping.
As far as box-counting fractal analysis is concerned, all recommendations by Kenkel40 were taken into account for the calculation of the rugae’s fractal dimensions. First, the outlines of the rugae were converted to point coordinates, instead of lines or pixels, and the number of points was well above the minimal recommended limit of 2,50040. In order to achieve a reliable estimation of the box count dimension, the minimum box size was set at 0.4 mm, ensuring that the box count was lower than M/10, where M: the number of points. Likewise, regarding the maximum box size, it was defined in such a way that all box sizes r, for which N(r) = (1/r)2, were excluded. Furthermore, 15 different box sizes were used, in accordance with Kenkel’s proposal that at least 10 should be employed, and the log–log plots were visually checked for departures from linearity40. Finally, in an effort to minimise the quantisation error25, more than 10,000 different search positions were employed (specifically 10,125 search positions; 15 translations per axis and 45 rotations).
The proposed method is a 2D-3D hybrid. It starts with cropping of the 3D palatal surface followed by the steps of construction of a BP surface, distance mapping and creation of rugae’s contour lines. The “Flattening” stage is the turning point from 3 to 2D. It could be omitted if the method was exclusively 3D and subsequently, the final step would be a 3D fractal dimension analysis, namely 3D box-counting (also known as cube-counting) fractal analysis41. However, if fractal dimensions were calculated with cube-counting, the shape of the palatal vault would influence the final outcome. As our objective was to evaluate the complexity of palatal rugae, independently of other parameters, flattening was considered necessary.
Flattening a 3D surface to a plane necessitates, in general, some distortion. Flattening algorithms can be conformal, minimising angular distortion, equiareal, minimising distortion of areas, or can minimise some combination of angle and area distortion42. Conformal flattening is advantageous for the study of rugae. Nevertheless, existing methods either offer little direct control over the shape of the flattened surface or demand significant nonlinear optimization34. In the present method, Boundary First Flattening (BFF) was used34. This is a linear method for conformal parameterisation which is faster than traditional linear methods. It “unfolds” and flattens the rugae area of the palate automatically and offers accurate preservation of sharp corners34.
Inter- and intra-rater repeatabilities were high due to minimal human interaction and operator calibration. Only two of the steps require subjective action by the operators. The step of “Initial cropping” determines the area that will be flattened and therefore affects the mesh distortion that will occur. The “Final cropping” step determines which of the structures should be considered as rugae. This stage is inevitable in all existing methods to distinguish which structures are to be measured or evaluated according to their shape, direction or branching. Contrary to previous qualitative methods, where the raters need to be calibrated in regard to the qualitative indices, in the present method, calibration is necessary only for the distinction of the contours that correspond to rugae.
As expected, the limits of agreement for intra-rater reliability were narrower than those for inter-rater agreement. Also, both inter- and intra-observer errors increased with fractal dimension (Figs. 3 and 4). A logical explanation is that as the complexity/fractal dimension increases, it becomes harder for the operator to differentiate rugae from other adjacent anatomical structures, such as incisive papilla and palatal gingivae, however, without having a great impact on the final result.
Advantages of the proposed method include a comprehensive evaluation of the complexity of all rugae (regardless of their size), with fractal analysis and a complete set of information about their contours and heights, with distance mapping. Minimum user intervention is needed, rendering the results repeatable and objective. In addition, fractal dimensions are quantitative variables and this makes them comparable.
Limitations of the method include sensitivity to parameters, such as ball radius and selection area. The size of the ball radius could over- or underestimate the size of rugae, and therefore their fractal dimensions, as well. Moreover, distortion is unavoidable after mesh flattening. The selected area could affect the extent of distortion, and rugae’s fractal dimension. Therefore, the operators need to be calibrated as per area selection.
Another drawback is the use of fractal dimension as a single measure of complexity. Despite its advantages, this method does not provide any information about the qualitative characteristics of rugae, such as shape, position, direction or branching or their quantitative variables, such as width, length or number.
Applications of this method could be helpful in forensics and dental research. With regard to the former, the procedure of ante- or post-mortem identification can be automated by means of distance mapping. Distance mapping is an intermediate stage of our method, as fractal dimensions are calculated in the end. However, in the case of identification, this step may be sufficient. Our method, which is a 2D-3D hybrid, could outweigh other exclusively 2D imaging methods19,20,21,22,23, as it provides information about the height of all rugae and details difficult to be detected otherwise. Undoubtedly, fractal analysis could be also used for the purpose of identification, by comparing rugae’s fractal dimensions in different sets of records. Both distance mapping and fractal analysis can be also applied partially, if the entire palate is not available.
As far as dental research is concerned, this method could be applied in various topics of research. First of all, it could be used to investigate a potential association of the complexity of palatal rugae and other anatomical structures, such as the size and number of teeth or the shape of the palate43, which may have common genetic background with rugae44. It could also detect potential correlations between specific genes and the complexity of the rugae, or help identify genetic polymorphisms or developmental stress—that could lead to asymmetries45—in the formation of rugae and provide insights into the hypothesis that palatal rugae are formed by a reaction–diffusion mechanism4. Furthermore, it has been proposed that rugae may change following various interventions (maxillary expansion, tooth extractions, etc.)46,47,48,49,50. Fractal dimension could be employed to exhibit the robustness of this hypothesis. Besides, since rugae participate in mastication and deglutition7, this method could be used to disclose potential relationships. Moreover, distance-mapped rugae could be used for improving superimpositions of stable regions of palate for measuring tooth movements (e.g. before and after orthodontic treatment).
Finally, this study was performed on digitised scanned plaster models, as the objective was to develop the methodology and to examine its reliability. Should it be implemented for one of the aforementioned purposes, it is recommended that it be used on direct scans of the rugae in order to minimise potential errors and problems related to dimensional changes of the impression materials and gypsum51,52.