Landscape pattern changes over 25 years across a hotspot zone in southern Brazil

Forest fragmentation caused by human activities has many implications for natural landscapes, such as habitat reduction and the loss of biodiversity. This study investigated the temporal fragmentation process of forest remnants in a strongly agro-industrialised region in southern Brazil over 25 years. The studied watershed area hosts two important typologies of the Atlantic Forest biome as well as Quaternary remnants of the Brazilian Savanna biome, which are considered hotspots of biodiversity and reflect the intense process of forest fragmentation caused by Brazilian urban and agro-industrial development. Thus, studies encompassing multitemporal scales are paramount to understanding changes in forest patterns and are fundamental for trend predictions of landscape dynamics. To perform the calculation of the mean normalised difference vegetation index, Landsat satellite images from 1991 to 2016 were processed using SPRING® software. Subsequently, FRAGSTATS® software was used to calculate landscape metrics. A reduction in the number of forest fragments since 1991 was observed, with a maximum amount of 5 243 fragments in 1993 that declined to 4 015 fragments in 2016. Although the number of fragments in the watershed decreased, the mean area increased by 72.9% and the mean of the shape index increased from 1.3 in 1991 to 1.5 in 2016. In addition, there was a 64.7% increase in the edge density and a reduction of 35.6 m in the isolation between the nearest neighbours. The degree of isolation of the fragments underwent a process of expansion and reduction when compared to 1991, presenting results that support the hypothesis that the Atlantic Forest is in a process of stabilisation and forest restoration.


Introduction
Southern Forests is co-published by NISC (Pty) Ltd and Informa UK Limited (trading as Taylor & Francis Group) The expansion of agriculture and urban areas in Brazil has led to the fragmentation of most Brazilian habitats, which affects landscape and biodiversity dynamics (Laurance et al. 2000a;Ribeiro et al. 2009). Fragmentation is a phenomenon in which a continuous portion of native vegetation is replaced by smaller portions that are usually separated by different types of land use (Fahrig 2003;Collinge 2009;Collins et al. 2017). Because it alters the spatial configuration of the landscape, the main effect of fragmentation is continuous habitat reduction and the appearance of a mosaic composed of smaller environments. Generally, it causes the landscape to lose ecological functions. As a result, the size of a forest patch is reduced and the number of transitional environments between the core of the fragment and the surrounding matrix increases, influencing the interactions of the local diversity (Fahrig 2003;Metzger 2009).
In fact, the fragmentation of primary forests causes a loss in quantity and quality of habitats as well as smaller fragments and an increase in the degree of forest fragment isolation (Fahrig 2003). In addition, fragmentation affects the local microclimate due to increased exposure to solar incidence at the edge of the fragment, which affects moisture, evapotranspiration, and the density of helophytic tree and epiphytic species (Teixeira et al. 2009).
Natural resource depletion and the invasion of non-native species have been some of the results of forest devastation, and could be interpreted as the result of the expansion of extensive agriculture, agro-industry and urbanisation activities in the last three decades (Machado et al. 2004;Ribeiro et al. 2009). Therefore, many researchers have focused on forest fragmentation and habitat loss as consequences of the decrease of primary forest areas in the Brazilian biomes (Laurance et al. 2000a(Laurance et al. , 2007Metzger 2009;Rocha 2011).
The Atlantic Forest is considered the most fragmented Brazilian biome and is classified as a global biodiversity hotspot (Myers 1988(Myers , 1990Myers et al. 2000;Machado et al. 2004). Only 12% of this biome remains and is distributed as forest remnants, and 83.4% of these remnants consist of less than 50 ha (Ribeiro et al. 2009). Although the Brazilian deforestation scenario began with the colonisation process, recent research on Atlantic Forest fragments point to a positive and growing trend of the forested area in the last three decades and suggest a stage of transition from loss to gain of forested areas, resulting in a positive aspect for the recovery and conservation of biodiversity (Rudel et al. 2002;Silva et al. 2016;Costa et al. 2017).
Data from remote sensing technology is one way to quantify changes in landscape structure by making comparisons and identifying differences when applied to mappings over time. The Landsat series of satellite imagery has provided a large database since the 1970s and supports current research as well as retrospective analyses. An additional means of quantifying changes is to use landscape ecology metrics because they express numerical information on spatial pattern changes and landscape process functions (McGarigal and Marks 1995;. Thus, studies focused on multitemporal scales are paramount to understanding changes in forest patterns and predict trends of landscape dynamics caused by anthropic activities (Costa et al. 2017).
In this context, the present study sought to map the remnants of Atlantic Forest over 25 years  and measure the temporal changes in the landscape of the Paraná River watershed (Mourão River tributary), the second major basin in South America.

Materials and methods
We studied the Mourão River watershed, which is located in the state of Paraná, Brazil. This area was chosen due to its ecological and economic relevance because it is an ecotone zone between the Mixed Ombrophilous Forest and the Seasonal Semideciduous Forest. The studied area contains remnants of the Brazilian savanna biome dating from the Quaternary Paleoclimatic period (Maack 1948(Maack , 1968Roderjan et al. 2002;Parolin et al. 2010) (Figure 1). It also hosts one of the largest agrobusiness poles of Latin America in Brazil, corresponding to 14% of the gross nominal value of production in the state of Paraná in 2016 (IPARDES 2017).
Medium resolution satellite imagery from the Landsat 5 TM, 7 ETM+ and 8 OLI sensor databases available on the US Geological Survey's Earth Explorer platform (https://earthexplorer.usgs.gov) was used. Only images with less than 10% obstructing cloud percentages were selected.  Two scenes of the same row/column numbers (223/76 and 223/77) overlap the studied area. Because of this, a mosaic composition needed to be processed for each analysed date. The data values consisted of the effective reflectance of red, green, blue and near-infrared bands as well as the normalised difference vegetation index (NDVI).
All of the satellite data used had already been converted to apparent reflectance and effective reflectance using Top of Atmosphere (ToA) methods and an internal algorithm, respectively, as well as being radiometrically calibrated and orthorectified by the producer through terrestrial control points and digital elevation models, ensuring the correction of displacements caused by the terrain.
The studied region is characterised by the major presence of agricultural crops with short phenological cycles, such as soybeans and corn (Lumbierres et al. 2017). These types of crops show abrupt variations in their vegetation indexes along their life stages, and for growth vegetative stages, the crops may present NDVI index values similar to forest remnant areas. On the other hand, forest remnants may also present variations in their indexes due to change of vegetative vigor caused by seasonal climatic conditions. Thus, it makes it difficult to classify using the segmentation process because of the variation in the phenological cycle in agricultural crops as well as forest remnants (Valero et al. 2016;Csillik and Belgiu 2017).
The NDVI arithmetic mean allows for the stratification of the values of each class from different elements in the landscape, facilitating the detection of forest remnants in agricultural matrixes. Therefore, the NDVI annual arithmetic mean was chosen as a parameter to associate different index intervals in the different classes such as forest, water and other elements (temporary crop, silviculture, urban and built area, and others) (Valero et al. 2016;Csillik and Belgiu 2017).
To mitigate misinterpretation of the arithmetic mean caused by the seasonality effects of forest and temporary crop cultivation reflectance, we selected years with the largest numbers of photo spatial records during the dry and rainy seasons from 1984 to 2016. Thus, we selected images from 1991, 1993, 2004, 2009, 2013 and 2016. Afterwards, the mean calculation of the NDVI indices was performed using an algebraic geoprocessing spatial language (LEGAL) algorithm. As an outcome, it obtained a raster for each sampled year based on the NDVI mean information. SPRING ® 5.5.1 software was used for these procedures.
After calculating the mean NDVI, the index intervals that illustrated the classes of native forest cover, water and the other elements were identified (Table 1). The fragments were identified by photo interpretation using false colour composition images for vegetation analysis. This composition was used to contrast elements that compose the landscape because it highlights the vegetation in shades of green, the soil in shades of purple and water in shades of red. The false colour composition was obtained by using the bands 5, 4 and 3 for images of Landsat TM 5 and ETM+ 7 or bands 6, 5 and 4 for Landsat OLI 8 images, which correspond to the infrared, near infrared and red spectral bands, respectively.
Given that some species with long growth periods, such as Pinus spp. (Glufke et al. 1997), were planted near the primary forest area, it was necessary to edit the images to fix any silvicultural patches that could be misinterpreted as forest fragments.
Based on the landscape metrics available in the FRAGSTATS ® 4.2.1 software (McGarigal and Marks 1995), the metrics addressing area, core area, edge, shape and connectivity issues were calculated using the forest fragment maps obtained in the previous procedures.
A cluster analysis was performed on the mean data to identify which sampled years had the most similar data. For the calculation, all metrics that were used were grouped. Thus, the Euclidean distance was taken as the dissimilarity factor among the n calculated individuals; the smaller the distance between two locations, the more similar and grouped the outcomes will be. The hierarchical clustering dendrogram was assembled using STATISTICA ® 7.0 software (Statsoft Inc., Tulsa, OK, USA).
Although a Student's t-test is often used to compare the statistics of two supposedly different populations, it is not trivial to measure the magnitude of the spatial dependency of random samples because 'spatial data are most frequently correlated with nearby samples' (Myers 1997). Therefore, to verify whether or not the changes in the landscape metrics were statistically significant, we conducted a Wilcoxon paired rank test.
We created a regular grid of 10 km² per cell overlapping the whole study area. It should be noted that the same grid was used over all of the analysed images, maintaining the same geographical position of the cells over the time series. Hence, a cell could be paired and analysed in two different time periods. After the metrics were calculated, 140 cells were randomly sampled and were used to perform the statistical analyses using RStudio 1.1.456 software (RStudio, Boston, MA, USA) (Supplementary Table S1).

Results
The Mourão River basin is approximately 165 000 ha and is composed of a large number of forest fragments of varying sizes and shapes ( Figure 2).
The Wilcoxon paired rank test indicated that there is enough evidence to conclude that the changes in the landscape metrics were statistically significant (Table 2).
In general, the number of forest fragments in the Mourão River watershed decreased from 1991 to 2016. In 1993, the maximum number of fragments was recorded (5 243) 1991,1993,1999,2004,2009,2013   to the decrease in the number of forest fragments through time, the total forest cover increased by 9 370 ha (Figure 3). In 1991, 16 133 ha of forest remnants covered the studied watershed, representing only 9.8% of the total watershed; however, in 2013, the forested area attained its peak with 27 925 ha, and declined to 25 503 ha in 2016. We noticed that the area of the fragments increased even though the number of forest fragments decreased by 8%. The forest cover increment is directly reflected in the mean area of the forest patches, which increased from 3.7 ha to 6.4 ha (Figure 4). Thus, the Mourão River watershed had 4 368 forest fragments with a mean area of 3.7 ha in 1991, which dropped to 4 015 forest fragments with a mean area of 6.4 ha in 2016.
The majority of sampled forest fragments were less than 50 ha, which accounted for 98% of the total samples from 1991-2009, and slightly more than 97% from 2013-2016 ( Figure 5). In contrast, the percentage of fragments greater than 1 000 ha is imperceptible (0.02%) up to 2009, and becomes absent after 2009 (Figure 5e-g).
Regarding the shape complexity index, we observed that 75% of the values in 1991 were below 1.4 (with a maximum value of 6.8), whereas in 2016 75% of the indices were below 1.6 with a maximum value of 12 ( Figure 6). Overall, the mean values for the shape index increased from 1.3 in the first analysed year to 1.5 in the last (Figure 7a).
This metric indicates the irregularity degree of the shape of a fragment when compared with a circle. As the index increases, the shape of the fragment tends to become more elongated and, because of that, it may have negative effects on biodiversity by restricting the core area and becoming vulnerable to the influence of the edge effect.
We observed an increase from 3 447.3 km (in 1991) to 5 679.3 km (in 2016) in the total edge length of the fragments in this study. Likewise, the edge density followed an increasing trend throughout the years (influenced by the edge length increase) from 20.89 m ha −1 in 1991 to 34.41 m ha −1 in 2016 (Figure 7b).
Given that the edge is the transition zone between the matrix and the core area of a fragment, its effect on the dynamics of the fragment and biodiversity depends on different forest typologies and biomes, which makes it an important feature to be considered in studies on fragmented landscapes. 2 000 4 000 6 000 8 000 2 000 4 000 6 000 8 000 2 000 4 000 6 000 8 000 2 000 4 000 6 000 8 000 2 000 4 000 6 000 8 000 2 000 4 000 6 000 8 000 2 000 4 000 6 000 8 000 Another relevant index for landscape analyses is the mean distance between the nearest neighbouring fragments because it may quantify the isolation degree of the fragments in the landscape. For this study, we noted a decrease from 139.2 m in 1991 to 103.6 m in 2016 ( Figure 8).
From the landscape ecology data of this study we detected dissimilarity between the first and the last analysed years (Figure 9).
Regarding the similarity between the metrics for yearly data, Figure 9 indicates that the outcomes may be divided into two major groups: Group I, which includes the data set of 1991, 1993, 1999 and 2004, and Group II, which includes the remaining analysed data. From Group II, the data set of 2016 was more similar to 2009, followed by 2013. The data set of 1993 was more to similar to 1999 and 2004 than the others. In addition, it should be noted that 1991 had the least similar data set in comparison to 2016, which confirms the hypothesis of changing patterns of forest fragments over time.

Discussion
The landscape of the studied area has undergone constant changes over the years due to factors such as the natural behaviour of the ecosystem in response to anthropic activities for land use and actions for the conservation and recovery of biodiversity.
The forest fragments declined in number and increased in mean size during the study period, which may explain the increase in forest cover in the watershed. However, whilst the increase in area and mean size represent positive aspects and ecological gains for the landscape, most of the fragments (approximately 98% from 1991-2009 and 97% after 2013) were less than 50 ha and therefore may be considered small. According to Ribeiro et al. (2009), the occurrence and predominance of forest fragments smaller than 50 ha are common (80%) throughout the Atlantic  1991,1993,1999,2004,2009,2013  Other studies on fragments of Atlantic Forest have found the same trend of an increase in forested area over the last three decades and support the hypothesis that the biome may be in a transition stage, resulting in a positive aspect for the recovery and conservation of biodiversity (Rudel et al. 2002;Rezende et al. 2015;Silva et al. 2016;Costa et al. 2017).
The area of a forest fragment is an important aspect that must be taken into account in approaches aiming to analyse the landscape quality. It is important because the area influences directly the capacity of an environment to host greater species richness by controlling the natural resources availability and fauna and flora density (Metzger 2009;Munguía-Rosas and Montiel 2014;Gross 2017). For instance, Almeida-Gomes et al. (2016) found a greater diversity of amphibian species in larger fragments of the Atlantic Forest than in smaller and more isolated fragments. Other studies have reported a positive relationship between species density and area for frogs in the Amazon forest (Lima et al. 2015), small mammals (Vieira et al. 2009) and for birds in fragments of Atlantic Forest (Banks-Leite et al. 2012). These results are associated with fact that larger fragments implicate a higher availability of food and drinking water resources, which controls inter-and intraspecific interactions, such as competition and predation (Townsend et al. 2010).
Besides the application of fiscal incentives for the creation of conservation units in Brazil, another aspect that may have influenced the increase of forested areas since 1991 is the increase of environmental inspections due to the application of Brazilian laws for conservation areas such as riparian forests and natural reserves. Tax incentives, such as the state value-added tax on goods and services (ICMS), later known as ecological ICMS in the state of Paraná, allow for financial transfer to municipalities with conservation units and supply sources in their territories, which has encouraged the creation of several conservation units in recent decades (Loureiro 2002). Between 1992 and 2008, the state of Paraná increased from 112 to 192 municipalities that received the ecological ICMS as a result of new conservation units, which contributed to the increase of forested areas in several municipalities of the state of Paraná (Rossi 2010). The Mourão River watershed added 14 new conservation units between 1991 and 2016, including state and municipal parks, ecological stations and private natural heritage reserves (RPPN).
Another parameter we observed was the shape complexity of the forest fragments. Shape complexity was measured using an index that compares the object of study to a standard shape; in this case, it was compared with a circle (McGarigal and Marks 1995). The minimum value of the shape index is assigned to a circle and express fragments with little shape complexity, and it increases as the shape becomes more complex. This index has great relevance in the analysis of landscape ecology once it can be associated with edge effects and size of the preserved core area of the fragment (Calegari et al. 2010).
Over the years, there was an increase from 1.4 to 1.6 (or 75%) in the shape index of the fragments, suggesting that these fragments underwent alterations in the complexity of their shape even though their mean size had increased. This adjustment leads to changes in the fauna and flora dynamics because more complex patches offer better conditions for the appearance of microhabitats, leading to greater species richness (Fletcher Junior et al. 2007) or higher emigration rates (Collinge and Palmer 2002). However, due to the small size of the fragments, these microhabitats are not capable of sustaining viable populations of a given species, and thus it becomes a highly diversified fragment with high extinction rates (Fahrig 2002;Borges et al. 2004).
The shape features of the forest patches are closely related to the edge effect because the more irregular the shape of a forest patch, the larger the contact area between the core and the surrounding matrix (Pereira et al. 2006). Therefore, changes in the shape complexity added to the increase of the mean area may have led to the increase of 64.7% in edge density compared with that of 1991 ( Figure 8).
The edge may exert a negative effect on the landscape structure by changing the biotic and abiotic aspects and, consequently, triggering a series of responses of the elements that compose the forest patches and its microclimate (Harper et al. 2005;Laurance et al. 2017). The edge region receives more light than the core of the fragment, which increases the average temperature and depletes the humidity level (Didham and Lawton 1999;Hardwick et al. 2015). The increase in the solar radiation affects the floristic composition (Laurance 1991;Ewers and Didham 2006;Hardt et al. 2013) because it can reduce the abundance of the vegetal species and interfere in their growth and pollen production, affect animal reproduction (Hartfield and Prueger 2015) and reduce the predominance of small tree species (Machado et al. 2008). In addition, it is possible to provide better conditions that increase β diversity (Didham et al. 1998;Magura 2017) by interfering in the interaction of animal species (Schneider-Maunoury et al. 2016;Rocha et al. 2017).
Moreover, sun exposure leads to liana densification at the edge, causing negative impacts due to shading heliophytic species near the edge (Laurance et al. 2000b;Laurance 2001;Gross 2017). In addition, they are negatively affected by intense winds, which may kill several plant species due to mechanical damages caused to the canopy (Gross 2017).
The connectivity between the fragments directly influences the dynamics of the species in fragmented environments because the flow and genetic variability of species may be altered and interrupted between isolated patches over long periods of time (Calegari et al. 2010;Boscolo and Metzger 2011). Therefore, the degree of isolation of the watershed decreased during the sampled years and the mean distance of the nearest neighbour was reduced from 139.2 m in 1991 to 103.6 m in 2016. This smaller distance represents a positive aspect for conservation of the local biodiversity and ensures greater species circulation among the fragments.
When the mean distance to the nearest neighbour that is reported in this work is compared with the distance data reported by Ferreira et al. (2018) in the Alonzo River watershed, which has a similar forest typology and is 120 km away from the studied area of this research, we noticed that the Mourão River watershed has a degree of isolation that is considerably smaller than the aforementioned author's research. According to Ferreira et al. (2018), the fragments of the Alonzo River watershed are 344 m from each other, which means they are approximately three times further than the fragments in the Mourão River basin.
The connectivity of the small fragments that compose a landscape plays an important role in minimising the degree of isolation between native forest areas (Ribeiro et al. 2009;Ferreira et al. 2018) because they work as ecological stepping stones for fauna, mainly birds that use them for resting, feeding and reproduction (Boscolo et al. 2008;Uezu et al. 2008;Barbosa et al. 2017).
In general, the effects of fragmentation tend to have long-term effects on biodiversity, especially in terrestrial ecosystems. This process of biodiversity loss is slow and occurs gradually, and is more noticeable in landscapes with successive episodes of fragmentation that are caused by both natural and man-made pressures (Ewers and Didham 2006;Munguía-Rosas and Montiel 2014). Nonetheless, it is worth noting that the factors that have led to forest cover increase can contribute to minimising such negative effects on biodiversity.
In addition, it should be noted that this research has some limitations due to the data set used, such as the 30 m spatial resolution of the Landsat images. In addition to spatial resolution, the radiometrical resolution varies because different Landsat sensors were used. These limitations are responsible for generating estimation and classification errors in the images, which influence the final result and analysis interpretation. Thus, new studies using finer spatial and radiometric resolutions should be carried out to achieve more accurate results.

Conclusion
The use of remote sensing data, such as the NDVI vegetation index from Landsat satellite images and geoprocessing tools to calculate landscape ecology metrics enabled us to map and analyse the forest landscape structure patterns of the Mourão River watershed from 1991-2016 quickly, practically and efficiently.
Likewise, the use of landscape ecology metrics provided robust results on changes in the structure of forest fragments over the years, as well as on the current situation of the watershed fragmentation process, making a temporal comparative study of forest cover possible.
The watershed still presents a high fragmentation stage with more than 97% of the current forest cover composed of fragments smaller than 50 ha and complex shapes. On the other hand, the outcomes also indicate a process of expansion in area and decrease of the degree of isolation of the forest fragments when compared with 1991, supporting the hypothesis that although the Atlantic Forest is highly fragmented and reduced to small remnants, there is a tendency for stabilisation and reforestation.