A new fractal index to classify forest fragmentation and disorder

Context Forest loss and fragmentation pose extreme threats to biodiversity. Their efficient characterization from remotely sensed data therefore has strong practical implications. Data are often separately analyzed for spatial fragmentation and disorder


Introduction
Forested landscapes have complex spatial patterns, and quantification of those patterns has long been of interest for scientists in diverse fields in ecology and beyond (Turner et al. 1989;Wu and David 2002;Newman et al. 2019).Forest loss poses an extreme threat to biodiversity (Haddad et al. 2015;Wilson et al. 2016;Liu et al. 2019).The resulting forest fragmentation produces changes in spatial patterns leading to isolation of patches, habitat degradation, and other ecological processes that further biodiversity loss (Fahrig 2003).Forests support human society through a wide array of functions and processes, each of which may be affected negatively by deforestation, habitat loss, and ecological disruption following fragmentation of habitat.This includes the maintenance of biodiversity (Thompson et al. 2009), stabilizing soils and topography (Mandal and Mondal 2019;Song et al. 2019), regulating stream flows (Cecílio et al. 2019), land surface temperatures (Arroyo-Rodríguez et al. 2017), and acting as carbon sinks mitigating climate change (Bonan 2008;Thompson et al. 2009).Because of forests' many functions, regional patterns and spatial dynamics of forest disturbance and land-use change have global conservation relevance.
Forest fragmentation and land-use change are often assessed by remote sensing (Wickham and Riitters 2019;Fischer et al. 2021;Batar et al. 2021).Although satellite images can be used to easily track the fraction of land usage conversion over time, there is a growing need to quantify temporal and spatial patterns of forest loss.These quantifications are necessary to assess habitat loss and degradation, as well as to predict outcomes such as loss of habitat complexity, biodiversity, or carbon stocks.The fraction of forest cover is one of the variables that can link spatial Vol.: (0123456789) patterns of land use to their ecological importance, whereas fractal fragmentation analysis and Information Dimension are mathematical tools that have demonstrated utility in analysis of satellite images at multiple spatial scales.
Fractal analysis of fragmentation and information dimension measures of disorder Fractal analysis is a set of methods that quantifies the complexity of mathematical or physical objects.Connecting ecosystem and biodiversity complexity measures continues to be challenging (Jin et al. 1995;Loke and Chisholm 2022), as is selecting the scale of analysis, when different scales of forest fragmentation are relevant to different organisms and their interactions (Kareiva 1987;Wiens et al. 1993;Boyce et al. 2017).Fractal metrics have the advantages of ease of calculation, as well as scale-invariance (Ma et al. 2009).
Some of the earliest studies that quantified landscape pattern have used forms of fractal analysis, including fractal neutral landscape models (Krummel et al. 1987;Gardner et al. 1987;Milne 1988;Gardner and O'Neill 1991;With andKing 1997, 1999), the midpoint displacement method (Barnsley et al. 1988), lacunarity (a fractal method that relies on box counting at different scales) (Plotnick et al. 1993(Plotnick et al. , 1996) ) and hierarchically structured random maps (Plotnick et al. 1993).Fractal dimension (Mandelbrot 1982;Gao et al. 2019) of the neutral landscape is the foundation of a great variety of models.These include models for source-sink relationships (Milne 1992;With 1997), invasive spread of exotic species (Lavorel and Chesson 1995;With 2002With , 2004)), edge effects (With and King 1997), and landscape fragmentation and its effects on population and community dynamics (With 2002;With and King 2004).Fractal dimension has also been useful in identifying critical thresholds for correlation of ecological processes with landscape models (Homan et al. 2004;Groffman et al. 2006), characterizing disorder of landscapes (Wimberly 2006) and analysis of dispersion patterns (Walters 2007), including spatial patterns of disturbance (Coops et al. 2018).As useful as fractal analysis has proven, ecological patterns are not truly fractal in the mathematical definition, but are strictly multiscale (Halley et al. 2004).Furthermore, fractal equations tend to be relatively simple and contain very little information about the complex systems they model (Halley et al. 2004;Newman et al. 2019).While fractals were among the earliest algorithms applied in landscape ecology, indices and algorithms developed in recent years that have broader applicability to real ecosystems have largely replaced them (Gustafson 2019;Gao et al. 2019;Turner and Gardner 2015).
Information entropy is a non-fractal approach to quantify the amount of spatial structure or disorder in a matrix of objects.Information entropy has also been used widely in landscape ecology (Nowosad and Stepinski 2019), for example in the study of landslides (Yufeng and Fengxiang 2009;Park 2015), in generating models of species distribution (Phillips et al. 2004;Elith and Leathwick 2009), as an indicator of fragmented landscapes (Joshi et al. 2006), for identification of types of mature woody forest ecosystems (Weber 2011), deriving macroecological patterns (Shipley et al. 2006;Harte and Newman 2014), and for quantification of the effects of disturbance on biodiversity (Supp and Ernest 2014;Newman et al. 2020).The most used form is Shannon information entropy, which measures the average information content in cases when the value of the random variable is not known (Shannon 1948).In this study, we use a generalization of Shannon information called the Rényi information dimension (or "Information Dimension") (Rényi 1970).The Information Dimension is a general way to quantify the diversity, uncertainty, or randomness of a system.The Information Dimension, denoted D 1 , measures the fractal dimension of a probability distribution, and relates the rate of growth of Shannon information to how finely the system under study is discretized.This method links scale-dependence and fractal structure to information entropy metrics.

Previous work and the introduction of the FFDI
In this study, we build on the previously developed fractal fragmentation index (FFI) (Andronache et al. 2016), by introducing the fractal fragmentation and disorder index (FFDI).The FFI is an efficient, fractalbased index that determines the degree of fragmentation or compaction of objects based on their shape (Andronache et al. 2016), which is not a feature of classical fractal analysis.The FFI is simply calculated as the dimension of the area minus the dimension of the perimeter of a shape, representing a measure of occupancy of space ("area" of a fractal can also be called "mass" or "core").We improved the FFI performance by incorporating the Information Dimension, which further quantifies the spatial distribution of these objects.We observed that FFI does not differentiate patterns where the component objects are extremely small such that no borders can be extracted, or where objects completely fill the space (in which case, FFI = 0, because D area = D perimeter ).The FFDI was also designed to overcome these limitations by additionally quantifying informational entropy in the image using the Information Dimension (D 1 ; see Table 1).
The FFI has previously been applied in the analysis of deforestation effects at the county level, territorial administrative units, and mountain groups of Romania (Drăghici et al. 2017;Pintilii et al. 2017;Andronache et al. 2017Andronache et al. , 2019;;Diaconu et al. 2019), indicating in all cases that fractal fragmentation increases as a result of deforestation.However, deforestation and total habitat loss is not the same as habitat fragmentation per se, so the degree of fragmentation of a landscape only measures one aspect of remaining habitat on landscapes after deforestation (Fahrig 2017).The spatial clustering and organization of the remaining fragments, captured by the Information Dimension, correlates with the distance among patches, which in turn captures ecological processes such as dispersal among patches, or patch isolation (Plotnick et al. 1993).The FFI was also previously tested on synthetically generated hierarchically structured random map (HRM) images, which are characterized by increasing fragmentation and disorder but decreasing average size over multiple scales.
In this paper, we address the limitations of the FFI by developing the fractal fragmentation and disorder index (FFDI), which quantifies the organization and arrangement of objects, thereby providing an estimate of spatial complexity.It can be used for the relative categorization of images with the same size, that is, extent and resolution.We validated the FFDI using several sets of binarized fractal and nonfractal test images, as well as synthetically generated hierarchically structured random maps (HRMs), which are binary maps with known properties, including a specified amount of cover and a hierarchical, nested property of patches (Lavorel et al. 1993;Plotnick et al. 1993).We show the added value of FFDI by demonstrating the complementarity in information content between the measures of fragmentation and disorder (definitions of abbreviations are provided in Table 1).We applied the index to images of forest cover from the Carpathian Mountains of Romania, from the Global Land Analysis and Discovery database (Hansen et al. 2013) (https:// glad.earth engine.app/ view/ global-forest-change; accessed 2020-2022).Finally, we discuss the utility of FFDI for

Methods
Defining the fractal fragmentation and disorder index (FFDI) Because spatial disorder and fragmentation of landscape patterns are often correlated but not perfectly so, a mixed index can quantify these two characteristics in a single value.Spatial disorder refers to the way in which the objects analyzed are organized in space: whether they are arranged with a high level of pattern and order, or alternatively show disorder.Fragmentation, on the other hand, refers to how much the objects analyzed break up a space.If all the space is occupied entirely by a single object, then the image has no fragmentation at all.However, the objects that break up, or "perforate" the image, indicate fragmentation.Fragmentation tends to correlate with coverage (Fig. 1), such that, for example, images with 10% cover are more fragmented, and images with 90% cover are more compact.In the case of a single object, if it is Euclidean, then it is compact, and the more irregular it is, the more fragmented it will be.
The fractal fragmentation and disorder index, or FFDI, is derived from both the fractal fragmentation index (FFI) and the information dimension (D 1 ) (Fig. 2).Like the FFI, it differentiates spatial organization patterns for processed images that represent two categories (target and non-target).Initial work and limitations of the FFI led us to the testing of a model which includes both FFI and D 1 to explain more than one kind of empirical pattern.We hypothesized that a functional form of FFI multiplied by D 1 would be better able to differentiate compaction/fragmentation and disorder patterns than either metric on its own.In ecology, these patterns of compaction, fragmentation, and disorder could correspond to, for example, natural patchiness of forest habitat, natural patterns of vegetation in response to aridity, or deforested landscapes due to road building and resource extraction.A more informative metric that contains measures of fragmentation and disorder may be better able to quantify these different ecological processes, and be applied in areas where on-the-ground information is not readily obtained.
The Information Dimension is a measure of the amount of "disorder" presented in the image (Bianciardi et al. 2014;de Souza Lins Borba et al. 2016;Borowska et al. 2017).Here, D 1 is derived from the multiscale Rényi entropies.The calculation of the Information Dimension (D 1 ) is efficiently estimated with a box counting method (Loke and Chisholm 2022), and is initialized with the information entropy (Fig. 2).
The information entropy I( ) for a set of N( ) boxes of linear size ε is defined as: The variable m i is defined as where M i is the number of points in the ith box, and M is the total number of points in the fractal object.Following from these definitions, D 1 is given by the equation (from Baker and Gollub 1996): (1) (2) If we consider a set of points distributed evenly in the two-dimensional plane, we will have approximately N number of points: The area is normalized to unity.Equation (1) therefore becomes: For a set of points situated along a smooth line, we get And finally, the information dimension can be defined in relation to information entropy as: (3) (5) Information entropy I(ε) is calculated using Eq. ( 7).If the set is fractal, then the plot I(ε) versus log(ε) will follow a straight line with a negative slope equal to − D 1 .The information dimension differs from the box-counting dimension in that it weights the boxes according to the number of points in the boxes, according to Eqs. ( 1) and ( 8).
The information dimension effectively assigns greater weights if a greater number of points are present (Kunicki et al. 2009).
The FFI is an index that describes the degree of fragmentation of objects, but can also be used as a compaction index (Andronache et al. 2016).A smaller FFI indicates more fragmentation and a larger FFI indicates more compaction.The FFI derives from the box-counting fractal analysis for the area of the analyzed objects and perimeter of the same objects, and is calculated: where FFI is the fractal fragmentation index; D A is the fractal dimension of the summed areas and D P is the fractal dimension of the summed perimeters; ε represents the size of the box; N( ) represents the number of contiguous and non-overlapping boxes necessary to cover the area of the object; and N � ( ) represents the number of contiguous and non-overlapping boxes needed to cover only the perimeter of the object (Andronache et al. 2016).
The FFI has the following properties: (1) FFI = 0 is the case when the analyzed objects are very small (in the 1-4-pixel range), such that their contour cannot be extracted, such that D A = D P = 0; (2) FFI = 1 is calculated when analyzing a 100% compact Euclidean object, without any discontinuity (D P = 1 and D A = 2); (3) when the objects are compact, the FFI value approaches 1, and when they are more fragmented, the FFI value approaches 0. Therefore, the new index FFDI is defined by: where m i= M i M as defined above, and ε is the size of the box.The slope of the log-log relationship provides the fractal dimension.In this definition, 1-FFI is used because FFI = 0 indicates fragmentation and FFI = 1 indicates lack of fragmentation (Andronache et al. 2016).
While the 1-FFI values range between 0 and 1 and D 1 between 1 and 2, the FFDI range is between 0 and 2, thus allowing a clearer differentiation than the FFI.The maximum value of FFDI approaches 2 when objects are strongly disordered and fragmented, whereas the lowest value approaches 0 when objects are weakly disordered and compact.
To calculate D 1 for the mass of objects, we used built-in tools for the extraction of the primary objects, i.e., the boundary image and mass of the object in the image (Fig. 1).We then employed a built-in algorithm for measuring the global Information Dimension from the Rényi family of fractal dimensions (10) , Fig. 3 Examples of hierarchically structured random map (HRM) images generated: a high s 1 set (p 1 = 0.5 and p 3 = p 2 = 1), b middle s 2 set (p 2 = 0.5 and p 1 = p 3 = 1), c low s 3 set (p 3 = 0.5 and p 1 = p 2 = 1), and several combined sets: d mixed s 4 set (p 1 = p 2 = p 3 = 0.5), e high + middle s 1 + s 2 set (p 1 = p 2 = 0.5 and p 3 = 1), f high + low s 1 + s 3 set (p 1 = p 3 = 0.5 and p 2 = 1) and g middle + low s 1 + s 3 set (p 2 = p 3 = 0.5 and p 1 = 1).For the s 1 images with p 2 = p 3 = 1 and p 1 = 0.5, 50% of the space is randomly occupied by HRM objects sized 100 × 100 pixels.For the s 2 images with p 1 = p 3 = 1 and p 2 = 0.5, 50% of the space is randomly occupied by HRM objects sized 10 × 10 pixels.For the s 3 images with p 1 = p 2 = 1 and p 3 = 0.5, 50% of the space is randomly occupied by HRM objects sized 1 × 1 pixel.In the case of s 4 images, p 1 = p 2 = p 3 = 0.5, which represent images where 50% of the space contains HRM objects sized 100 × 100 pixels (s 1 ), 10 × 10 pixels (s 2 ), and 1 × 1 pixel (s 3 ) available in the open-source plugin ComsystanJ 2.0.1 (Ahammer 2011) for the Java-based open-source software Fiji.The ComsystanJ box-counting dimension method was used to calculate both the FFDI and FFI.
Generation of hierarchically structured random map images for testing FFDI HRMs are maps of pixels that are randomly set to 1 or 0, with probability distributions of 1's and 0's that vary by region.HRMs can simulate landscapes that have different levels of clustering and fragmentation, which is accomplished by curdling and random trema.Two surfaces can have the same total amount of filled space, but differ in their level of curdling, that is, the clustering or nucleation of the filled space.The nonfilled space is generated by random trema, or randomly placed holes of a given shape in the mathematical surface (Mandelbrot 1982;Gardner et al. 1987;Milne 1992).For this study, HRM images were generated by the Fiji plugin ComsystanJ using a recursive algorithm derived from the curdling and random trema generation method employed by Plotnick (Mandelbrot 1982;Plotnick et al. 1993).Note that for studies that compare images, all images should have the same resolution (pixel size) and extent, in order to avoid introducing bias from automatic resizing and other sources of error (Loke and Chisholm 2022).This applies both to images generated algorithmically and to images processed from empirical data.
The HRM images were generated for our tests with numbers of rows and columns (M) fixed at 1000 (1,000,000 pixels; or 1000 × 1000-pixel resolution).Variables p 1 , p 2 and p 3 represent the percentage of space occupied by the HRM objects at different grain sizes (e.g.0.5 = 50%).These were generated by creating arrays with L rows and L columns and randomly setting the L 2 elements to 1 with a probability p′, with p 1 , p 2 and p 3 , representing the percentage of space occupied by the HRM objects at different scales.For p 1 , the size of the generated HRM objects is at 100 × 100 pixels grain; for p 2 it is 10 × 10 pixels; and for p 3 it is 1 × 1 pixel.In each case, p′, the total filled area of the image, is 0.5 (i.e.50%) for each type of HRM object generated.Due to their random distribution, the images are evenly occupied by HRM objects (Fig. 3).
We generated seven sets of 2D HRMs synthetic images with different values of p 1 , p 2 and p 3 , labeled s 1, s 2, s 3, s 4, s 5 = s 1+2, s 6 = s 1+3 and s 7 = s 2+3 .Following the p 1 -p 3 designations, the s 1 -s 3 sets are each characterized by a single "grain size": s 1 being the "high" set of 100 × 100 pixels grain size, s 2 being the "medium" set of 10 × 10 pixels, and s 3 being the "low" set of 1 × 1 pixel grain size.The remaining s 4 -s 7 sets were "mixed" sets with combinations of the high, medium, and low sets.We then calculated the FFI, D 1 , and FFDI for these images to compare the performance of the different metrics.
Although the exact levels of fragmentation and disorder vary, even among images with the same levels of p 1 , p 2 and p 3 , we can say general things about the image classes before analyzing them, based on the definitions of fragmentation and disorder.Images in the s 1 class are characterized by large, compact HRM objects with less disordered spatial distribution; s 2 and s 1 + s 2 images by medium-sized HRM objects with moderate fragmentation and spatial disorder; while images in the s 3 set are very fragmented and their distribution highly disordered.The s 4 set images include s 1 , s 2 , and s 3 HRM objects and are obtained by filling the space with HRM objects from the s 1 and s 2 image sets and fractal spaces from the s 2 with s 3 HRM objects.Therefore, s 4 is fragmented and disordered.The most fragmented and disordered images are the s 1 + s 3 and s 2 + s 3 sets of images.
For FFDI validation, we generated 20 horizontal sinusoidal images with an increasing number of iterations, from 1-20; 20 vertical sinusoidal images generated with 1-20 iterations; 20 radial sinusoidal images generated with 1-20 iterations; 10 IFS Sierpiński Gasket images generated with 0-9 iterations; and 10 IFS Menger Carpet images generated with 0-9 iterations (80 images total).All testing and validation images generated for this study are provided in the Supplementary Information 1-4.

Case study: dynamics of forest loss in the Romanian Carpathians
To test the proposed FFDI in a real-world setting, we analyzed the dynamics of forest loss and the impact f Parâng Mountains Group; g Retezat-Godeanu Mountains Group; h Banat Mountains; i Poiana Ruscă Mountains and j Apuseni Mountains.Compass points are approximate over the regional extent given the map projection of forest dynamics in the Romanian Carpathian Mountains.The study area constitutes a major landscape unit with mountainous features, which forms a circular arc from the border with Ukraine in the north of the country to the border with Serbia in the south (Fig. 4).
Geological and geomorphological features divide the Carpathians into three major sectors (the Eastern Carpathians, the Southern Carpathians and the Western Carpathians).Each sector is divided into 10 more-or-less coherent mountain groups, depending on the degree of forest fragmentation, controlled by the geology and the cross-sectional valley corridors (Fig. 4).The Eastern Carpathians include the Northern, Central and Curvatures Groups.The Bucegi, Făgăraș, Parâng and Retezat-Godeanu Mountains belong to the Southern Carpathians, whereas the Apuseni, Poiana Ruscă and Banat Mountains are a part of the Western Carpathians.Their average elevation is ~ 2000 m, and the absolute maximum elevation is the Moldoveanu Peak in the Făgăraș Mountains (2544 m) (Balteanu et al. 1998).
In comparison to the situation at the global level (Song et al. 2018), the deforestation in Romania shows continuous growth (Andronache et al. 2019).The 2017 Report on the State of the Forests of Romania states that in that year, forests of Romania occupied an area of 6,565,000 ha, or 27.5% of the country's land area.Although deforestation has increased in total during recent decades, areas of forest loss have been partly compensated for by forest regeneration, reforestation, refurbishment of forested pastures and the inclusion of degraded lands in afforestation efforts, aligned with the provisions of Law no.46/2008 (Ministry of Water and Forests 2017).

Forest image preprocessing
The Global Forest Change database used in the current research is based on Landsat satellite time series data, provided by the Geographical Sciences Department at the University of Maryland, College Park.The data set provided is based on a change detection analysis, which highlights year-to-year forest dynamics since the year 2000 CE.We therefore chose the year 2000 to be used as the base year for which the rest of the analysis is reported.Forest cover images for Romania were downloaded with 30 m spatial resolution in GeoTIFF format from the Global Forest Change database.Total tree cover for the year 2000 was used as a baseline comparison for forest loss for each year from 2001 to 2021.Subsets of forest loss area images for each of the 21 years were generated for each of the mountain groups, for subsequent binarization and analysis.
We calculated differences among the images of tree cover for the year of interest and the base year, 2000, to determine area of forest loss for individual years and the cumulative forest loss (using the ImageJ 1.52 Process and Image Calculator functions).Subsequent processing of satellite imagery for extraction of annual images and preparation for the fractal analysis was performed with ArcGIS software (ESRI 2020) (Fig. 5).The grayscale GeoTIFF images were converted to 8-bit TIFF binary images, compatible with binary fractal analysis.Binarization was done using the Threshold operator in the Image-Adjust menu in ImageJ, with the thresholds being 1 and 255.Two classes of pixels resulted: black pixels, with value 0 representing the background and white pixels with value 255, representing the foreground and the objects to be analyzed respectively.In our case, these correspond to patches of forest and deforested areas.

Statistics for images and metrics
In order to apply the FFDI, we first took into account the sample size and the general rule that box counting using at least ten box sizes should be obtained over the scaling range (Kenkel 2013), where the number of boxes was the number of distinct scales of pixels differing by a power of 2 (because "fractal analysis" for empirical images is technically multiscale analysis).Thus, for our fractal analysis, the image resolution allowed us to set the number of boxes to 11, with the regression ranging from 1 (min) to 11 (max).
Second, we checked that the images were sufficiently fractal-like that we could calculate the metrics D A , D P , and D 1 .According to Benguigui et al. (2000), the quality of the fractal dimension estimation is controlled by using the Pearson correlation coefficient r, and estimations can be accepted when r exceeds the value of 0.99.If the fit between the empirical curve and the estimated curve is low, it can be concluded that either the pattern being studied is not fractal, or that it is multi-fractal (Tannier and Pumain 2005).For the images of the Romanian Carpathians as well as all the test images, we performed preliminary analyses (in ComsystanJ) and found that the coefficient of determination, R 2 , was greater than 0.99, and in the case of a linear regression, R 2 = r 2 .This constrains the Pearson correlation coefficient to be greater than 0.99, so we were able to perform the fractal analyses on the images.
We note that high R 2 values do not always and in general correspond to high linearity.We therefore tested the linear regressions for several box size ranges, from 2 0 -2 6 up to 2 0 -2 10 (measured in pixels), and compared the gained results to each other.This was done for all images s 1 , s 2 , s 3 , s 4 , s 1 + s 2 , s 1 + s 3 , s 2 + s 3 , Gaussian, Random, FFT, and MPD (with FD ranging from 2.0 to 3.0).The correlations between these results were very high (e.g.R 2 = 0.996 between 2 0 -2 9 and 2 0 -2 11 , and R 2 = 0.958 between 2 0 -2 6 and 2 0 -2 11 ), with only the absolute values deceasing slightly for smaller box size ranges.This decrease was negligible and did not affect the outcome of this study.Differences between HRM, FFT and MPD images were maintained.Final results were obtained with a range of 2 0 -2 11 , although a smaller range would have been acceptable as well.
Statistical analyses were carried out with SPSS 27 (IBM Corp 2017).Data distributions were tested for normality using Shapiro-Wilk tests (p < 0.05).Not all distributions were normally distributed.Thus, nonparametric Kruskal-Wallis tests including Bonferroni corrected post hoc tests were applied.

Evaluation of the FFDI using HRM images
We calculated the FFDI, 1-FFI, and D 1 for all the HRM image groups to compare the outcomes of each metric, and their ability to differentiate among image classes.Objects that have high fragmentation are usually spatially disordered, and the most compact objects are usually uniformly spatially distributed.Moreover, when fragmentation is very high, FFI cannot identify patterns (Andronache et al. 2016(Andronache et al. , 2017)).As stated previously, a very high level of fragmentation (FFI near 0) is achieved when images have objects of size 1, or they have very few pixels such that no border can be extracted.In the first stage of fractal analysis, the three algorithms used (1-FFI, D 1 and FFDI) were tested in the three-level HRM images (Fig. 3).
The distributions of 1-FFI, D 1 and FFDI values were compared for each of the 9 × 200 generated images in the high, medium, low, mixed, highmedium, high-low and medium-low HRM sets (Fig. 6).In each case, every metric yields a clear differentiation among the distinct HRM image sets with a sample size of n = 1800, with FFDI offering the best differentiation when considering all comparisons, as assessed by Bonferroni-corrected Kruskal-Wallis significance tests.The number of statistically different pairs varied for each of the three metrics.From 21 possible pairs, the FFDI yielded 20 highly significant pairs and one significant pair (Fig. 6a, d).D 1 yielded 19 highly significant pairs and one significant pair (Fig. 6b, e).1-FFI yielded 14 highly significant pairs and again one significant pair (Fig. 6c, f).
A linear regression with all 12,600 HRM images shows a very low correlation between disorder and fragmentation (R 2 = 0.075) in the synthetic images.The FFDI shows a very low correlation with D 1 of R 2 = 0.017 and a strong correlation with 1-FFI of R 2 = 0.837 in the synthetic images.Due to higher overall discrimination ability, FFDI is expected to work better than the comparison metrics at lower sample sizes.
Evaluation of the FFDI using non-fractal test images Additional test images were generated by Midpoint displacement (MPD), fast Fourier transform (FFT), Random, and Gaussian algorithms, with different fractal dimensions (FD).As expected, the midpoint displacement and fast Fourier transform algorithms yielded nearly identical images, and thus showed high correlation.On the other hand, the random and Gaussian images had substantially different patterns, without significant correlation.Statistical analyses of discrimination among test images yield very high significances between the distinct fractal dimensions of MPD and FFT generated and binarized images.From 276 possible pairs, the FFDI yields 212 highly significant pairs and 17 significant pairs (Fig. 7).D 1 yields 174 highly significant pairs and 20 significant pairs (Fig. 8).1-FFI yields 213 highly significant pairs and again 12 significant pairs (Fig. 9).Thus, there are again more significant pairs for FFDI compared to D 1 and 1-FFI.
Several other sets of test images were analyzed as well: horizontal sinusoidal, vertical sinusoidal and radial sinusoidal images with 1-20 iterations; and Menger Carpet and Sierpiński Gasket images with 0-9 iterations.For these images, the value of the FFDI generally increased with the number of iterations, indicating that both the resolution and depth of scale affected the calculated FFDI values (Supplementary Information 5).
Using FFDI to analyze tree cover and forest loss in the Romanian Carpathians As a practical application, we analyzed and differentiated forest characteristics concerning disorder, and fragmentation of forest cover area (Fig. 10a), total area of loss in the year analyzed (Fig. 10b), and for the cumulative loss area in the Romanian Carpathian Mountains from 2000 to 2021 (Fig. 10c).The FFDI differentiated among three patterns of mountain groups, based on their fragmentation and disorder: the highest values of FFDI (Central, Northern, and Parâng mountain groups) coincided with high fragmentation in the images, which also corresponded to areas of intense deforestation (Fig. 10a).The natural fragmentation and compaction of this landscape lead to the construction of access roads into the forest for resources, which intensifies deforestation pressures.These characteristics and values of fractal measures were very similar to the s 4 set of HRM generated images (Fig. 3).Intermediate values of FFDI (Apuseni, Bucegi, Făgăraș, Poiana Ruscă, and Godeanu mountain groups) were specific for mountain groups with moderate fragmentation and disorder.The lower values of FFDI (Banat, Retezat-Godeanu, and Southern Group mountains) were specific to mountain groups characterized by a higher degree of mountain terrain continuity and lower fragmentation of terrain (relief), implicitly favoring a stronger compaction of forests.This pattern resembled the s 2 "middle" set of HRM images.These patterns agree with previously described patterns of deforestation and fragmentation (Andronache et al. 2019).Figure 10c shows a clear increase in FFDI values associated with the areas of cumulative loss over time, which is reflected in the different range of values for the ten mountain groups.In these cases, FFDI increased with increasing loss of forests, because the cumulative loss areas became fragmented and disordered as they expand.We note that for the Făgăraș and Parâng Mountains Groups the cumulative forest loss is increasingly disordered, even though the forest loss areas are not as naturally fragmented as the other three mountain groups from the highly fragmented pattern.Low FFDI values were found for the third pattern of cumulative loss areas due to deforestation patches being very small but numerous, with less disordered spatial distribution (Ciobotaru et al. 2019).
Interestingly, due to very low FFDI values and similar fragmentation of forest loss areas in all 10 mountain groups (Fig. 5), Information Dimension was found to be the measure that could differentiate among these three structural types of mountain groups, similar to the s 4 and s 3 sets of HRM images.With further investigation, the influence of FFI on the FFDI was very small when similar images of mountain groups or structural types are analyzed, similar to s 3 or s 4 HRM images when FFI < 0.006, but becomes important when FFI exceeds 0.1, as is the case for s 1 and s 2 HRM images.

Discussion
Quantifying the spatial structure of landscapes is the key to understanding the ecological effects of landscape patterns and habitat fragmentation (Costanza et al. 2019).Metrics based on fractal analysis are commonly used because of ease of their calculation and multiscale properties.However, another ecologically relevant pattern of the arrangement of habitat fragments, or "landscape disorder," had not been previously measured in fractal-based landscape metrics.Consequently, we formulated a new metric called the Fractal Fragmentation and Disorder Index, or FFDI, that can identify and also classify patterns and subpatterns of image spatial fragmentation and disorder.
In this study, we showed that FFDI outperforms various other metrics in quantifying both fragmentation and disorder of images and solves a previous limitation of the FFI specific to images with highly fragmented, spatially distributed objects.We demonstrated the ability of FFDI to better differentiate and classify HRM patterns compared to fractal fragmentation or Information Dimension, and that this index could be used in the future as a versatile and accurate tool for pattern differentiation and sub-patterns in any field that uses binary imaging.The FFDI performed well in the classification of spatial characteristics of artificially synthesized HRM objects as well as images of real forest patches.We recognize a limitation in the method applied here that arises as a result of binarization of images to be analyzed with the FFDI.Depending on the chosen threshold, binarization may introduce analytical bias (Eguiraun et al. 2014); however, this effect can be investigated further and might be rectified with the incorporation of a scaling relationship, or use of the FFDI as a relative index for comparison of different system states, rather than as an absolute measure of fragmentation and disorder.Even with this limitation, FFDI offers a clearer differentiation in categorizing synthetic images and in real applications.
Compared to D 1 and FFI, the FFDI achieved a clearer differentiation of forest tree cover, cumulative forest loss and annual forest loss patterns for the 10 mountain groups in the Romanian Carpathians.It is therefore a promising approach for differentiating patterns of deforestation and regrowth and will be useful in environmental monitoring and protection.
In conclusion, although FFDI was developed with the aim to improve the understanding of forest fragmentation, it might be applied generally to landscape fragmentation resulting from urban development, wildfires, and other natural and anthropogenic phenomena.FFDI can also be used for comparative analyses across different regions.We demonstrated the use of the FFDI in estimation fractal dimension as a tool for relative classification of images of the same size and resolution, based on their fragmentation and disorder.
We expect a broad applicability of the FFDI in various disciplines of medicine or biology.For example, in histological slides, the distribution of cells or cell nuclei shows a fragmented pattern which is often irregular or random due to the underlying nonlinear dynamics of gene regulation and cell differentiation (Fabrizii et al. 2014).Improved quantification of such patterns could improve diagnosis, prognosis and prediction in cancer patients.(Macaluso et al. 2015).

Fig. 1
Fig. 1 Schematic of the correlation between disorder and fragmentation.Increasing fragmentation is shown from left to right, and increasing disorder is shown from top to bottom.Black indicates objects, and white indicates background

Fig. 2
Fig. 2 Workflow diagram for producing both the fractal fragmentation index (FFI) and the fractal fragmentation and disorder index (FFDI)

Fig. 4
Fig. 4 The geographic location of the study area, indicating a Northern group of Oriental Carpathians; b Central group of Oriental Carpathians; c Southern group of Oriental Carpathians; d Bucegi Mountains Group; e Făgăraș Mountains Group;

Fig. 5
Fig. 5 Examples of tree cover binary images (from 2000 CE): a Northern group of Oriental Carpathians; b Central group of Oriental Carpathians; c Southern group of Oriental Carpathians; d Bucegi Mountains Group; e Făgăraș Mountains Group; f Parâng Mountains Group; g Retezat-Godeanu Moun-

Fig. 6
Fig. 6 Statistical analysis of discrimination among HRM images for seven test image sets s 1 -s 7 , across three metrics: a the Fractal Fragmentation and Disorder Index (FFDI); b Rényi information dimension (D 1 ) (note: x-axis does not cross at 0); and c 1-the Fractal Fragmentation Index (1-FFI).Values on

Fig. 7
Fig. 7 Statistical analysis of discrimination among binarized test images generated by Midpoint displacement (MPD), fast Fourier transform (FFT), Random, and Gaussian algorithms, with a range of fractal dimensions (FD), for the fractal fragmentation and disorder index (FFDI).a For FFDI, values on

Fig. 8 Fig. 9
Fig. 8 Statistical analysis of discrimination among binarized test images generated by Midpoint displacement (MPD), fast Fourier transform (FFT), Random, and Gaussian algorithms, with different fractal dimensions (FD), for the Rényi information dimension (D 1 ) (note: x-axis does not cross at 0). a For D 1 , values on y-axis are dimensionless.b Comparisons are demonstrated with the p < 0.01 Kruskal-Wallis test significance matrix, with Bonferroni corrected pairwise post hoc tests (*p < 0.05, **p < 0.001) ◂

Fig. 10
Fig. 10 The compaction and disorder patterns of forest cover and loss in the Romanian Carpathians as represented by the fractal fragmentation and disorder index (FFDI): a forest cover areas, b forest loss areas, and c cumulative loss areas.The highly fragmented and disordered patterns are depicted with dashed lines (Central, Northern, and Parâng mountain groups); the moderately fragmented and disordered patterns are rep-