Remote Sensing Evaluation of vegetation cover from IVDN vegetation indexes: possible effects of climate change at Catimbau National Park

The Catimbau National Park (Parna Catimbau) is an environmental conservation unit (CU) created by Decree Law No. 4,340, dated 08/22/2002, in accordance with Federal Law No. 9,985 (SNUC Law). This CU is characterized by being a National Archaeological Heritage by IPHAN and presenting a floristic diversity, which makes it an area of extreme biological importance. However, despite being a region defined by the Ministry of Environment as an integral protection unit, we can see in this area the presence of families and consequent anthropic activities such as agriculture, irrigated production, pastures and other agents that have been contributing to the modification of this Parna. This article aimed to map and evaluate the spectrum-temporal dynamics of the vegetation cover in the Catimbau National Park (PE) area and its possible relationships with climate change, in the period of 2003 and 2016. Remote Sensing techniques such as Normalized Difference Vegetation Index IVDN, whose values vary in the range of -1 to 1, were applied to the Vegetation study assisted by Landsat 5, TM sensor and Landsat 8, OLI sensor satellite images. The software SPRING was manipulated to create a geographic database and for the development of programming aided by LEGAL in order to correct atmospheric influences, among other parameters. The responses of these analyzes and their respective problematic of the potential contributing factors indicated changes in the vegetation cover with the highest increase for the anthropic area (0.107-0.207) and Soil / Outcrop (0.207-0.307) for the year 2016. These expansions are characteristics of physical, environmental and anthropogenic factors, where they all contribute for hazardous situations in the area. These results indicate the need for social and environmental studies to support environmental management of this conservation unit of the Caatinga biome.


Introduction
The Catimbau National Park (Parna Catimbau) was created by Decree Law No. 4,340, dated 08/22/2002, accordingly with Federal Law Nº. 9,985 (SNUC Law), which defines, in its Article 11 th , that national parks "have a basic objective, the preservation of natural ecosystems of great ecological relevance and scenic beauty, enabling scientific research and development of activities of education and environmental interpretation, recreation in contact with nature and ecological tourism" (Santos, 2003).The area has several archaeological sites (dating to about 6,000 years old BP) and brings occurrences of inscriptions and rock paintings, being for this reason, transformed into national archaeological patrimony by IPHAN (Institute of National Historical and Architectural Heritage) (Martin, 1991;Silva and Maia, 2008;SIGEP, 2010;Mendes, 2012).In addition, as a result of what the area provides, as it has been mentioned, such as a unique geological and geomorphological scenic beauties, the Parna Catimbau was mentioned as an area to become a Geopark in the Environmental, Geomorphological and Archaeological categories by Unesco (SIGEP, 2010).
The Parna Catimbau covers part of the municipalities of Buíque, Tupanatinga and Ibimirim (Agreste, State of Pernambuco), being considered an area of extreme biological importance and priority for the conservation of the Caatinga biome due to presenting relevant numbers of endemisms and rare species (Sampaio et al., 2002).Figueiredo et al. (2000) and Andrade et al. (2004) describe vegetation physiognomy of this Conservation Unit (CU) varies according to the type of soil, therefore it can be found vegetation from arboreal caatinga, of non-thorny species, to shrubs with predominance of thorny species, which is a peculiar characteristic related to these xeric environments.
Although this CU presents a legal mechanism for its conservation and is defined by the SNUC (National System of Conservation Units) as an integral protection unit, the presence of families that live in the region is perceived (Azevedo and Bernard, 2015), as well as the presence of indigenous groups, agricultural production, livestock and ecotourism (Freire, 2013).A recently published research, Freire et al. (2017), reveals that from 14 conservation units (CU's) of integral protection of the Caatinga biome examined, Parna Catimbau is one of the CU's with the highest environmental degradation index.The authors concluded that, in addition to the socioenvironmental conflicts between farmers, squatters, indigenous populations and managers of ICMBIo (Chico Mendes Institute for Biodiversity Conservation), there were significant losses of caatinga physiognomies, especially in areas with rock vegetation in areas with natural wall of sandstone -typical from this CU and with difficult access.
The detection of spatio-temporal changes in vegetation cover can be facilitated by the use of Remote Sensing (RS) techniques.Wang (2012) describes that RS is a tool that contributes to the monitoring of protection areas, since it provides information that makes possible understanding and to realize the relationships between the past and present state.The use of RS techniques directed to the study of vegetation with the application of Vegetation Indexes (VI's), which are related, among other factors, to biophysical vegetation parameters such as biomass and leaf area index (Formigoni et al., 2011).These VI's allow obtaining different data about the biophysical conditions of the vegetation, allowing spatial and temporal comparisons of terrestrial photosynthetic activities that are influenced by the consequent actions of the anthropical actions and, as a consequence, promote desertification and loss of floristic variability in the ecosystem (Jensen, 2009).
Among all the indexes the most well-known and applied is the Normalized Difference Vegetation Index(IVDN).IVDN is used for improving the discrimination of vegetation cover, for measuring plant vigour and biomass production.Jensen (2009) describes some appropriate factors for the applicability of this index, these are the following: seasonal and yearly changes in the development and vegetation activity that can be monitored, and the reason that reduces many forms of multiplicative noise, for example differences of solar irradiance, cloud shadows, some atmospheric attenuations, some topographic variations, that are found in multiple bands of images of multiple data.
Large-scale and long-term ecological studies have been one of the main points in multidisciplinary research on global climate change.The variation of vegetation cover not only influences energy balance, climate, hydrological processes and biogeochemical cycles, but also indirectly illustrates the interaction between climate, environment and anthropogenic activities (Heute and Warrick, 1990).Thus, the main objective of this article is to map and evaluate the dynamics of vegetation cover in the Parna Catimbau area and its possible relationships with climatic changes using temporal indexes of IVDN, considering the interval between two dates of the years 2003 and 2016, obtained in the rainy season.
This article was elaborated by virtue of results obtained in the inter-institutional research Joaquim Nabuco Foundation and Federal University of Campina Grande entitled "Mapping and Spectro-Temporal Analysis of Conservation Units of Integral Protection of the Caatinga Biome of the Federal Administration" (Freire et al., 2017), carried out between 2013 and 2016, as well as the contents of the course "Remote Sensing and Environmental Monitoring".Ministered at the Integrated Center for Georeferenced Studies for Social Research Mário Lacerda de Melo (CIEG/ FUNDAJ, 2017).

Study Area
The Parna Catimbau is located in the northeast region of Brazil, in the central region of the State of Pernambuco and it represented by the geographical coordinates 8 ° 24 '00 "and 8 ° 36'35" S and 37 ° 09'30' and 37 ° 14'40 'W , with an area of 60,787ha, of which 40.81% belongs to Ibimirim; 38.73% to Tupanatinga and 20.46% to of the municipality of Buíque (Siqueira, 2006), as observed in Figure 1.Aragão and Peraro (1994) describe that this CU is characterized by presenting a landscape formed on sandstone rocks that is related to the Devonian period, presenting unique forms and reshaped by weathering.Parna Catimbau presents a semi-arid tropical climate with an average annual temperature of 23ºC, average rainfall of 300-500 mm annually and altitude between 600-1000m (Brazil, 2005).Due to the variations of relief and climate, besides the typical species of Caatinga, it also occurs in the area elements of the Cerrado, Rocky Fields, Mata Atlântica and Restinga (Brazil, 2008).The vegetation cover is associated with the Caatinga domain and presents five types of physiognomies, such as: arboreal-shrub caatinga, shrub caatinga with cerrado elements, shrub caatinga with elements of rocky fields, perennial forest vegetation and shrub caatinga subperennial.

Pluviometric Analyzes
In order to evaluate variations in vegetation from the Normalized Difference Vegetation Index (IVDN), it is necessary to analyze the rainfall of the municipalities that comprise the park (table 1).Thusly, the images were used in the period of August 25, 2003 andAugust 12, 2016, since this period corresponds to the end of the rainy season and cloud cover <20%.The period of greatest rainfall (60% to 75% of rainfall) occurs in the months of March / April and June / July, between September and January the driest period occurs (Siqueira, 2006).Kasmierczak (1996) describes that the IVDN equation shows a different behavior according to the variation of the spectral bands, that is, when the band 3 of the TM sensor is used, the green, dense and uniform vegetation has a darker color.While in the use of band 4 of this same sensor, the water bodies will have a darker color and the green vegetation, dense and uniform lighter.With the division between these original spectral bands, the result will be a new spectral band, which represents the ratio between the original bands, indicating in the light areas the presence of vegetation and in the darker tones the absence for the same date (Freire and Pacheco, 2005, p. 6).According to these authors, the resulting spectral band should be inserted in the range of -1 to 1, indicating for values closer to 1 the presence of areas with strong chlorophyll activity, whereas for results close to -1 express areas with low chlorophyll activity.

Image Processing
The process to calculate IVDN for TM and OLI sensor images requires the adoption of measures aimed at correcting the influences of external and internal factors that make it difficult to achieve satisfactory results.In this way, the TM sensor was initially performed to convert the raw digital values (digital numbers) to physical values (radiance and reflectance).Radiance represents the solar energy reflected by each pixel, per unit area, time, solid angle and wavelength through the sensors alongside the satellite (Bastianseen et al., 2000), which in the case of LandSAT 5 satellite TM sensors correspond to the measurements performed in channels 1, 2, 3, 4, 5, 6 and 7 as observed in table 1.
With the use of the Georeferenced Information Processing System (SPRING), which is a software provided by the National Institute of Space Research (INPE), a programming language was developed in LEGAL.For the process of transformation of digital number (ND) of each pixel of the image in monochromatic spectral radiance, the equation proposed by Markham and Baker (1987) was used (equation 1).And, with the reflectance resulting from the ratio between the emerging flux of the atmosphere and the incident flux at its top, the equation defined by Bastiaanssen (1995) (equation 2) was applied. (1) Where: L λ = spectral radiance of each band (Wm -2 sr -1 µm -1 ); "a" e "b" = minimum and maximum spectral radiances (Wm -2 sr -1 µm -1 ); ND = digital number of each pixel and Dnmáx = 255, since the radiometric resolution for the LandSAT 5 satellite images is 8 bits.Considering the values shown in the table above, the following calibration coefficients were used, as shown in Table 3.
Table 3 -Parameters considered for the calculation of the calibration coefficient for the LandSAT 5 satellite, TM sensor.

-
The reflectance enables the monitoring of the behavior of certain vegetation cover, since the changes in the spectral response are related to the structural or spectral changes of the vegetation cover and not the changes in the sunlight angle or variation of the sensors (Pereira et al, 1996).In this context, for the calculation of the monochromatic reflectance, the equation proposed by Bastiaanssen (1995) (2). (2) Where: ρ ʎi = banda monochromatic reflectance of each band (Wm -2 sr -1 µm -1 ); L λ = monochromatic spectral radiance (Wm -2 sr -1 µm -1 ); k ʎi= spectral solar irradiance of each band at the top of the atmosphere (Wm -2 sr -1 µm -1 ); Z = solar zenith angle at the moment of acquisition and dr = ratio between the average Earth-Sun distance (r0) and the Earth-Sun distance on a given day of the year squared.
In order to determine the Earth-Sun distance, the equation proposed by Duffie and Beckman (1980) was applied, as observed below.
(3) Where: DOY is the Julian day of the year and the angle is given in radians.
The DOY of the orbital image under study corresponds to 237 because the orbital image used dated 25 August 2003.Substituting it in the equation, the resulting value for dr is 1.032916377.This value is within the range established for dr, that is, between 0.97 and 1.03 being dimensionless.
For the conversion of the monochromatic radiance to reflectance, in other words, the transformation of the digital numbers from the image to physical values using the SPRING software, then it becomes necessary to establish the relationship between the equation proposed by Markham and Baker (1987) and Bastiaanssen's equation (1995).
(4) Where: ρʎ is the monochromatic reflectance to be obtained; dr is the inverse of the relative distance Earth-Sun squared in astronomical units; z is solar zenith angle (degrees) at the moment of acquisition; Eʎ is the mean of solar irradiance at the top of the atmosphere; L is the monochromatic spectral radiance where the minimum and maximum values are obtained for each band.Lastly, there is DNmax which is the maximum digital number for each band (255) and DN represents the maximum and the minimum values achieved in each band.
The result obtained by this equation will correspond to a new reflectance image which brings physical values.And, this will be replaced in the equation proposed by Rouse et al. (1973) for the calculation of the IVDN (equation 5).
(5) Where: ρ IV : reflectance value in the near infrared band and ρ V : reflectance value in the red range of the visible spectral region.
The observed reason for IVDN calculation is characterized by the intense absorption the vegetation presents in the red spectral regions and by an intense energy reflected in the near infrared region caused by the cellular structure of the leaves (Lourenço and Landim, 2003).Figure 2 shows the routine in LEGAL developed for the calculation of the radiance, reflectance and NVDI for the year 2003.The IVDN calculated for the OLI sensor considered the red and the near infrared spectral regions, that is, bands 4 and 5, respectively.The orbital images provided by the OLI sensor have 16-bit radiometric resolution and can be rescaled to Top of Atmosphere (TOA) reflectance and / or radiance.This conversion of digital number values to Top of Atmosphere (TOA) reflectance is performed from the reflectance rescaling coefficients.The values of these coefficients are provided in the metadata of the products.Equation 6shows how the TOA of the planetary reflectance is calculated without correction for solar angle.However, in order to determine the TOA reflectance with correction for solar angle it is necessary to establish the relation of TOA planetary reflectance without correction for solar angle with the solar elevation angle and the solar zenith angle, as observed in equation 7.

Detection of vegetation loss
For the detection of vegetation loss, Crósta (1993) and Moreira (2003) propose a technique to map these areas based on IVDN data, that is, with IVDN images resulting from the years 2003 and 2016, the mean and the variance values were obtained using the statistic grid (Table 4).

Table 4: Statistical data of IVDN images
The mean and variance values for the years 2003 and 2016 demonstrate some variation, in this case it is necessary that the reflectance levels of the most recent image (2016) coincide with the mean and variance of the oldest reflectance levels (2003).In this way, it will be determined the gain and the offset values.The gain is obtained from the difference between the variances as observed in equation 8.This result will be used to calculate the offset, applied in equation 9.With the values of Ganho and Offset the recent image (2016) will be adjusted so that its average and variance correspond to the same values of the oldest image ( 2003), according to equation 10.The result generated will correspond to a new image reflectance (NDVI_2016_Igual ) coinciding with the old one.From this resulting image it will be possible to detect the losses of vegetation.
(10) Loss detection will occur from the difference between the 2003 IVDN images with the resulting image (NVDI_2016_Igual), considering the value of Ganho wqual to 100 and offset equal to 150.This enhancement occurs because the pixels that present similar values for the two dates, that is, there were no changes in the vegetation, will present low values of gray levels in the resulting information plane, in contrast, in the regions where changes occurred, the resultant pixels values presented high gray levels (Freire and Pacheco, 2005).

Analysing
the Normalized Difference Vegetation Index (IVDN) applied in Parna Catimbau for the years 2003 and 2016 it was possible to identify five classes that were classified in the following ranges: a) 0,107-0,207 for anthropical areas, b) 0,207-0.307representing exposed soil and / or rocky outcrop, c) 0.307 to 0.407 herbaceous caatinga, d) 0.407-0.507shrub caatinga and e) greater than 0.507 arboreal caatinga.
Comparing years 2003 and 2016, it was identified a significant modification was identified in the anthropical areas (0,107 -0,207) and in the Soil / Outcrop (0.207 -0.307),where consequent vegetation reductions of the Caatinga biome in this CU for the year 2016 were found, according to Figures 4 and 5.These foment the deforestation of the caatinga vegetation for the appropriation of new extensive pasture areas, even when the soil is predominantly shallow, stony and with low organic nutrients, hence it has low aptitude for cattle pasture (Freire and Pacheco, 2011); 2) area "B" presents mostly a rocky vegetation loss in the sandstone walls.Table 5 and Graph  The considerable increase of the exposed soil / rock outcrop was identified for the year 2016 with an increase of 65.243 km².This growth brings considerable material losses and concerns to both the rural populations and the public managers involved with the areas affected by the anthropical phenomenon (Freire and Pacheco, 2011).
As a result of the physical interferences, climatic changes and anthropical conflicts in this region, it demonstrates the fragility in the conservation of the biodiversity of the caatinga biome in this conservation unit.For the year 2016, there is a variation in the presence of herbaceous caatinga vegetation (327.36km²),shrub caatinga (92.62km²) and arboreal caatinga (5.02km²) when compared to 2003, as observed in the graph above.This reduction is mainly due to the anthropogenic actions promoted by residents and smallholders in order to increase cultivated areas and / or pastures in order to increase the value of their lands with the purpose of indemnifying the land regularization in progress by ICMBio in the park area.Another decisive factor is the rainfall data.These show that the precipitation values for the years 2003 and 2016 present differences, 51.8 mm and 0 mm, respectively for 10 days before the scene image; however, when considering 60 days prior to the scene, the most current scene had a greater influence of the rainy season (54.7 mm), while the oldest scene was smaller (25.2 mm).In this way, we can affirm that there was loss of vegetation essentially due to the action of environmental degradation promoted by the intra park dwellers (area identified by nº 1 in Figure 6).
This identification of the areas of greatest loss of vegetation are determining points for the study, since they permit analyzing possible causes of these problematic modifications faced by this CU in the studied period.
The areas in the darker shades (blue) highlight the regions of environmental recovery, while those in green identify denser vegetations (green) and the red ones highlight the vegetation losses of caatinga biome that occurred in years 2003 and 2016.Area 1 indicates the growth of degradation spots from socioenvironmental conflicts due to anthropical interference in the cultivation and agricultural production both predominant in the region (Figure 6), and this dynamic scenario extends itself throughout the Park.
However, the areas identified in the same Figure 6 by nº 2 indicate the sandstone walls (Figure 8) -areas of hard access and improper for farming and grazing activities activities.Vegetation in these areas is rock type and its losses are, a priori, related to the possible changes in recorded in the sandstone walls (Figure 7) of the Parna Catimbau.In addition, there is also the prolonged dry and droughts in the region that have probably been intensified by global climate changes, which registered a considerable increase in temperature, with a maximum increase of 0.6 °C and minimum one of 0.5 ° C (Salati et al., 2007, p. 25).

Conclusions
The applicability of remote sensing in the study of vegetation cover using IVDN indexes contribute to environmental monitoring and to the study of the behavioral dynamics of vegetation distribution in areas of integral protection, due to the natural complexity of the CU's, to perceive and comprehend the conflicts involved by the farmers, squatters and indigenous people and to understand the historical of human occupation and land use in the Catimbau National Park.
The result of this study indicated that this CU has undergone alterations and consequent losses of caatinga biome vegetation.Thusly, a systematic monitoring of this Parna constitutes an important instrument of environmental management by allowing the understanding of the biological phenomena and the dynamics of human occupation and land use, that correspond to the anthropic actions in this area, inducing to multidisciplinary studies based on the relations of the conviviality social and environmental characteristics of the region and its implications on the environment and the populations that inhabit its interior and surroundings.

Figure 1 -
Figure 1 -Location of the Catimbau National Park.
Table 1 -Total monthly precipitation of the municipalities bordering the Catimbau National Park.Source: APAC (2017)2.3Normalized Difference Vegetation Index -IVDNFor this study, it was considered the use of the LandSAT 5, TM sensor and LandSAT 8, OLI sensor satellites.The orbital images acquired corresponded to the years 2003 and 2016, obtained for orbit 215 or 66 point.For obtaining IVDN, it was used spectral regions of the red (V) and near infrared (IV) bands, corresponding to the spectral bands B3 and B4 of the sensor TM and B4 and B5 of the OLI sensor, respectively.

Table 2 -
Data related to the channels for the LandSAT 5 satellite, TM sensor with respective calibration coefficients (a, b e K).Source: Adapted from Chander et al. (2009)

Figure 2 -
Figure 2 -Legal programming fot the year 2003 = TOA planetary reflectance, without correction for solar angle; Mρ = Band-specific multiplicative rescaling factor from the metadata; Aρ = Band-specific additive rescaling factor from the metadata e Qcal = Quantized and calibrated standard product pixel values (DN).
= TOA planetary reflectance; θ SE: = Local sun elevation angle; and θ SZ: = Local solar zenith angle.The Figure3describes the algorithm used to calculate IVDN using SPRING software.

Figure 3 -
Figure 3 -Routine developed in the Legal program for calculating the IVDN fot the year 2016.
of the IVDN image for the year 2003; M2016 = average of the IVDN image for the year 2016.

Figure 4 :
Figure 4: Thematic Map of the IVDN Index for the year 2003.

Figure 5 -
Figure 5 -Thematic Map for the IVDN Index for the year 2016.
1 represent the changes obtained from the class measures in the studied periods.Table 5 -Changes in Class Measures for the year 2003 and 2016 in km².Graph 1 -Analysis of class measures for the study period in km².

Figure 6 -
Figure 6 -Dynamics of the vegetation in the conservation unit.