Spatial analysis of degradation in the Reserve for San Rafael National Park, period 2005-2019

The goal of this study was to analyze the forest degradation in the Reserve for San Rafael National Park, Paraguay, during the period 2005-2019. This Reserve is one of the most important forest remnants of the Upper Paraná Atlantic Forest Ecoregion. A multitemporal analysis of degradation was carried out due to the occurrence of three disturbances: forest fires, a twister and illicit crops, using the Continuous Degradation Detection (CODED) algorithm, for which 3 factors were considered: variations due to pixel in the NDFI index values before, during and after every disturbance registered. In this context, the phenomenon with the greatest impact in terms of magnitude of degradation were the forest fires of 2005, being that year at the same time, the one that reported the highest degradation values. Secondly, there are the illicit crops established until the first semester of 2019, and lastly, the twister that occurred in 2017. Our findings demonstrate that CODED algorithm can detect multi-temporal degradation events in a Subtropical Broadleaf Forest, and the postdisturbance regeneration process after every disturbance tends to occur immediately. The response in terms of degradation-regeneration is highly variable, depending of the nature and severity of each disturbance and the vegetation recovery dynamics.


Introduction
The NDFI (Normalized Difference Fraction Index) is a spectral index of degradation, and it has the ability to synthesize information from satellite imagery, allowing the mapping of logging and damage caused to the forest canopy. Spectral mixing models are like a linear combination of spectral signatures called "endmembers". With the choice of the "endmembers", the coefficients of the linear combination are interpretable as fractions of the land cover [1,2,3].
The CODED (Continuous Degradation Detection) algorithm is an open source algorithm executable on the Google Earth Engine platform. It was designed to use Landsat images of 30 meters of resolution. Detection of change (deforestation or degradation), is made using regression coefficients to predict multitemporal observations of the NDFI index. In this way, the algorithm is executed online, which means that the change is controlled sequentially in time using data available in the platform itself [4].
On a global scale, the NDFI spectral index has been used to detect changes in the structure of tropical forests and detect degradation patterns, either with Landsat or MODIS images [5,6,7]. In Latin America, it has been used in investigations to monitor forest degradation and deforestation, mainly in the Brazilian Amazon, and to detect selective extraction activities in areas with forest management plans. Landsat images are the most used due to their availability, allowing multitemporal studies to be carried out [8,9].
In Paraguay, there has been an increase in forest degradation during the last decades, affecting the quality of the forests and producing structural changes in them [10]. This is mainly due to the overexploitation of wood to obtain firewood and charcoal, which are the main sources of energy on which more than 51% of Paraguayan households depends [11].
The San Rafael National Park Reserve in Paraguay is considered one of the main conservation priorities, due it constitutes a representative sample of the biodiversity of the Upper Parana Atlantic Forest (BAAPA). It is floristic richness is among the highest in the country [12], and it is estimated that the different plant formations all present some degree of degradation due to human activities. In the absence of definitions, defined assessment methods and limited data, forest degradation is a complex problem to address, so it is important to generate solid knowledge about this phenomenon.
In that context, just a few countries have the capacity to report the forest areas affected by this phenomenon and the degree of degradation of their resources. The information available is still limited, which is why implementing methodologies to monitoring the historical and current levels of degradation of forests and their remnants, through the use of remote sensors, is a priority.
The applications of the CODED algorithm allow to analyze retrospective and current geospatial data about the state of the forests that are within the study area. The identification and prioritization of factors and phenomena that contribute to forest degradation, and that are considered "critical", make it possible to evaluate the effectiveness of forest management and conservation.
The general goal of the study was to analyze the forest degradation in the Reserve for San Rafael National Park, Paraguay, during period 2005-2019. Specifically, the goals were to: 1) describe the multitemporal behavior of NDFI index values with respect to 3 registered forest disturbances; and 2) analyze forest disturbances in terms of magnitude of degradation.

Study area
The research was carried out in San Rafael National Park Reserve (RPNSR), which is located in the Eastern Region of the Paraguayan Republic at 380 km from the City of Asunción (Capital). To the north of Itapúa Department covers the districts of Alto Vera, Itapúa Poty and San Pedro del Paraná, and to the south of the Department of Caazapá covers the district of Tavai ( Figure 1). The main access road corresponds to number 1 "Mariscal Francisco Solano López".
The RPNSR is one of the largest Protected Wild Areas of the BAAPA Ecoregion, and is part of one of the most important freshwater reservoir basins, called the Guaraní aquifer, with several streams that irrigate the area. Its central location coordinates are 55° 42'58.43 "west longitude and 26° 26'57.17" south latitude and occupies an approximate area of 66.575 hectares [13,14]. It is made up of two main natural communities: grasslands and forest. The first belongs to the ecoregion called Mesopotamian Grasslands, and the forest is part of the socalled Upper Paraná Atlantic Forest (BAAPA) [15].
The forest is known as "Subtropical Broadleaf Forest". It is the predominant type of vegetation, occupying approximately 80% of the RPNSR. It includes a variety of habitats: primary forests, with 20 meters high trees, closed canopy and open undergrowth; degraded forests, with trees between 15-19 meters high, open canopy and closed undergrowth, gallery forests and natural forest islands. It is composed of three types of plant communities, which, according to the Rapid Ecological Assessment [16], are: Subhumid Forest, Forest Seasonally Saturated, and Forest Seasonally Flooded.

Data
The Landsat L1T image collection available on the Google Earth Engine platform was used. These have a correction based on control points, digital terrain models and a topographic correction for the displacement of the terrain by the relief [17,18].
Subsequently, the images were filtered by date ".filterDate", region ".filterBounds" and percentage of clouds less than 20% ".filterMetadata". The ".FeatureCollection" function was used to delimit the area of interest, defining as a mask the official RPNSR shapefile available in the Geoportal of the Ministry of Environment and Sustainable Development (2019), taking into account that this is the institution in charge of the application of on Protected Wild Areas [19].

Selection and analysis of disturbances
The three disturbances considered were: forest fires that occurred in 2005, a twister that occurred in 2017, and illicit crops established until the first semester of 2019. Every disturbance was spatially delimited with the ArcGIS 10.3 ® software. In the next stage was selected the "Time Series Explorer" tool, implemented on the Google Earth Engine platform. It allows selecting a geographic point for the temporal analysis of the NDFI Index and its components: Green Vegetation (GV), Non-photosynthetic Vegetation (NPV), Soil and Shadow [9].
A pixel covering a total area of 900 square meters was selected for each disturbance that occurred in the RPNSR and a multitemporal analysis was carried out, for which 3 factors were considered: the variations in the index values before, during and after registered disturbance.
For this study, it is considered that the vegetation has suffered damage to the canopy and is degraded within the range of NDFI values between 0 and 0.75. This range of NDFI values was empirically defined using canopy damage data collected by conducting forest inventories and their correlation with values reported by the NDFI index [3]. A time period was selected for the analysis of each disturbance (Table 1).

Magnitude analysis
The interpretation tool and algorithm code used were adapted from github.com/bullocke/coded. For the analysis of the magnitude of affectation resulting from the occurrence of each disturbance, training samples were taken from the coverage observed in each affected area, and the parameters for the execution of the CODED algorithm were defined in the user interface, which below are cited with the respective assigned values ( Table 2). Unlike the multitemporal monitoring of the NDFI index (by point), to determine the magnitude of each disturbance, the analysis was carried out by regions. Were considered the forest fires in the northeast area of the Reserve, the tornado in its entirety in the southern area, and illicit crops in the north zone.

Multitemporal behavior of the NDFI index respect to forest disturbances
The time series of the NDFI index values respect to the forest fires that occurred in 2005 is shown in Figure 2. This behavior, supported by satellite observations of the area, indicates a process of regrowth of forest vegetation, mainly during the last 2 years. This point out that forests pass through multiple stable states depending on the region's natural dynamics and the severity of the change with which they have evolved [20]. As forest fires are dynamic and frequent disturbances within the RPNSR, the forests subject to the fire regimes change constantly, explaining the fluctuations in the behavior of the index values for the pixel 6 of 9 under study during the period 2000-2019. The time series of the NDFI index values respect to the twister that occurred in 2017 is shown in Figure 3. Certain types of disturbances may shift an ecosystem from a stable to one of instability, but provided these disturbances do not occur frequently, the forest will tend to return to its original stable state through succession [20].
The time series of the NDFI index values with respect to illicit crops established up until the first semester of 2019 is shown in Figure 4. altered. The first significant decrease in the values is observed during the month of October 2017, with a minimum value of -1. Subsequently, a further decline in the index values was detected to a minimum of -0.76 at the end of November 2018. The satellite corroboration allows this behavior to be attributed to the phenology of Cannabis sp. illicit crops, that is shown in Figure 5.  The forest fires that occurred in 2005 report the maximum magnitude value among the 3 disturbances analyzed, with a value of 132 ( Figure 6). The areas with reddish coloration are classified as critical and with a high magnitude of impact, therefore they have a lower resilience potential. These could coincide with those points from which the fires spread or areas where the fires were concentrated for a longer time. Other factors also influence, such as floristic composition, atmospheric humidity, temperature, wind direction and the slope of the terrain.
On the other hand, the intermediate to low magnitudes are registered in areas in which the fires did not have a serious impact, so their resilience potential is higher. They are located around the areas classified as critical. Due to the location of the fires, it can be said that they are closely related to the burning of weedy areas for the establishment of crops in the Reserve buffer zone, because these forests are bordering large extensions of areas destined for agricultural use.
Next, the second most significant disturbance in terms of registered values is the illicit crops observed until the first semester of 2019, which reported a maximum value of 124. The areas with reddish coloration are classified as critical and with a high magnitude of affectation, in which the implementation of these activities causes serious incidence. These are mainly observed in fallows and areas recently abandoned after being harvested, therefore they are devoid of vegetation. The yellowish-greenish coloration is observed in areas that presents less magnitude of affectation. Such behavior is due to the existence of vegetation cover, and depending on the phenological state of the crop, a high response in photosynthetic activity can occur, resembling forests in this aspect, but not in terms of composition and structure.
Lastly, the disturbance that registered the least significant value of magnitude is the twister that occurred in 2017, which reports a maximum value of 72. Areas with higher magnitudes are classified as critical, because they were the most affected and the climatic phenomenon could have manifested with higher intensity. These zones coincide with the central areas in which the tornado passed through and they have a lower resilience potential. On the other hand, areas that register intermediate to low magnitudes are classified as moderately critical, due to the pass of the storm did not directly affect these areas. They are located around critical areas, and they have the highest potential for resilience.
According to data from the Station of the Meteorology and Hydrology Direction (2017), located at the Air Base "Cnel. DEM José María Argaña" from the District of Capitán Meza, Itapúa, the maximum wind direction on 04/26/2017, in which the disturbance was registered at night, corresponds to 149 ° in the northwest direction [21]. This can be considered an important factor because it coincides with the trajectory followed by the climatic phenomenon.

Conclusions
The multitemporal analysis of the NDFI index with respect to the disturbances analyzed in the study area shows that the post-disturbance regeneration process tends to occur immediately in terms of remote sensing. The CODED algorithm detected degradation events after the disturbance of interest, in this case forest fires, for which the response in terms of degradation-regeneration is highly variable. On the other hand, in the case of the twister and illicit crops, the behavior is related to the factors that depends of the nature of each disturbance, in the first case, to the severity of the climatic phenomenon, and in the second, to the dynamism of the affected areas by illicit crops that responds both to the phenology of the crops and to the constantly migratory nature of these activities.
In terms of magnitude, forest fires reported the maximum values, indicating a high impact on the pixel selected for the study, taking into account that this phenomenon is of high occurrence within the Reserve. Next, the pixel selected to analyze illicit crops reported the maximum magnitude values during the thinning/burning stages of the forested areas and in the harvest stage.