Analysis of Grassland Degradation in Zona da Mata, MG, Brazil, Based on NDVI Time Series Data with the Integration of Phenological Metrics

There is currently a lot of interest in determining the state of Brazilian grasslands. Governmental actions and programs have recently been implemented for grassland recovery in Brazilian states, with the aim of improving production systems and socioeconomic indicators. The aim of this study is to evaluate the vegetative growth, temporal vigor, and long-term scenarios for the grasslands in Zona da Mata, Minas Gerais State, Brazil, by integrating phenological metrics. We used metrics derived from the normalized difference vegetation index (NDVI) time series from moderate resolution imaging spectroradiometer (MODIS) data, which were analyzed in a geographic information system (GIS), using multicriteria analysis, the analytical hierarchy process, and a simplified expert system (ESS). These temporal metrics, i.e., the growth index (GI) for 16-day periods during the growing season; the slope; and the maximum, minimum, and mean for the time series, were integrated to investigate the grassland vegetation conditions and degradation level. The temporal vegetative vigor was successfully described using the rescaled range (R/S statistic) and the Hurst exponent, which, together with the metrics estimated for the full time series, imagery, and field observations, indicated areas undergoing degradation or areas that were inadequately managed (approximately 61.5%). Time series analysis revealed that most grasslands showed low or moderate vegetative vigor over time with long-term persistence due to farming practices associated with burning and overgrazing. A small part of the grasslands showed high and sustainable plant densities (approximately 8.5%). A map legend for grassland management guidelines was developed using the proposed method with remote sensing data, which were applied using GIS software and a field campaign.


Introduction
Minas Gerais State is the most important milk-producing area in Brazil.It produces approximately nine billion liters of milk per year, and the Zona da Mata region is responsible for 10% of the total production [1].Zona da Mata consists of many small farms in the traditional milk regions, with typical intensive management.Grassland degradation in Zona da Mata is a current research topic due to its importance to the state's dairy chain, resulting in surveys, monitoring, and planning of public policies for the sector.Several Brazilian government programs such as the Low Carbon Agriculture Program (ABC) of the Ministry of Agriculture, which aims to reduce greenhouse gas emissions, and the rural credit-financing program of the joint Federation States have designed actions to reduce and recover degraded grasslands (http://www.agricultura.gov.br/assuntos/sustentabilidade/plano-abc/).However, these government programs need technical support and methodologies that are compatible with decision-making mechanisms for large pasture areas, which are present across most of the country.The Zona da Mata region, located in the state of Minas Gerais, has a total area of 3,574,800 ha, and due to its status as one of the most traditional regions for extensive livestock farming, this region presents very important demands from the agricultural sector in relation to estimates of spatiotemporal conditions in grasslands.Monitoring grassland conditions through remote sensing demands temporal regularity and quality to generate cartographic results regarding the vegetative development scenario because of phenological dynamics and long-term effects, in addition to management and climatic factors that influence the plant conditions and survey complexities.It has been shown that caution should be taken when using image data series in relation to grasslands in the peak milk production period due to the coexistence among vegetal biomass quantity, herd consumption, rainfall periods, high photosynthetic activity, and the resulting broad vegetal soil cover in the summer season.Phenological metrics by remote sensing aim to characterize general or specific transformations of plants over time by using the spectral responses of plant tissues caused by their interactions with incident energy.Generally, exploratory phenological metrics are better suited for periods of higher plant mass production, while fine metrics aim to diagnose irregular pasture biomass production.Vegetation indices adequately respond to demands that involve vegetation condition mapping, and some types of indices can be selected in relation to variations in vegetation density due to climatic conditions and for the species studied.Vegetation indices generally indicate saturation in the Brazilian summer, and the use of soil-adjusted indices that describe the canopy structure is recommended, especially in cases of very dense vegetation, such as in forests.However, for herbaceous and sparse vegetation, the energy absorption levels of chlorophyll, specifically in the red band, are more effective for biomass evaluations and biophysical characteristic indicators when combined with the availability of high temporal resolution images [2,3].Normalized difference vegetation index (NDVI) data obtained from the Moderate Resolution Imaging Spectroradiometer (MODIS) sensor data from the Terra satellite (available from the MODIS database) may yield synoptic information regarding forage photosynthetic activity.Vegetation and photosynthetically active radiation indices are not linearly related to the leaf area index (LAI), with an LAI saturation of approximately 3.0; however, they are nearly linearly related to the canopy photosynthetic activity [4][5][6].Grasslands have been studied and monitored using satellite data in several different ecosystems and over seasonal and interannual periods to verify variability, production gains, and resilience.Researchers have correlated grassland net primary productivity, as measured in the field with MODIS-derived productivity, and have obtained better results for the growing season [7][8][9], indicating that extreme weather conditions make it difficult to analyze biophysical correlations, requiring time series for estimates of comparative metrics.Several vegetation indices have been derived from hyperspectral reflectance data to measure biomass due to their ability to provide more information on vegetation variability.However, there is also an increasing need for advanced techniques to analyze large datasets and spectral signatures.Thus, further research is needed to improve biophysical estimates from multispectral data [10].Using multispectral sensors and hypertemporal data, the accuracy of detecting invasive plants and other targets increases, and the results from the use of classification methods and vegetative growth analysis metrics thus improve.However, invasive species in grasslands have also been identified using NDVI/MODIS data, such as in semiarid regions [11].With the sparse vegetation due to the rainfall regime, climatic conditions are decisive in defining the methods and data to be used for identification and phenological analyses, since the influence of climate on grasslands, and their degradation aspects and recovery, are long-term.Thus, the interannual variation of grassland growth trends has been analyzed using NDVI/MODIS MOD13Q1 data, and a direct relationship with accumulated rainfall has been shown in some studies.Researchers have found a relationship between seasonal accumulations in the NDVI/MODIS data (cool or hot seasons) and rainfall for natural pastures in the USA and Australia [3,12,13].Phenological metrics have been used to detect differences between crops and differences in the vigor of different vegetation types, e.g., crops, pastures, and forest vegetation.Therefore, temporal metrics, particularly the start of the season (SOS), the peak of the season (POS), the end of the season (EOS), and the slope of the curve or the series, as well as differences between the maximum, minimum, and mean, have provided satisfactory results [14][15][16][17].The grassland vegetation growth index (GI) from an NDVI/MODIS time series has been estimated to map grassland degradation in the spring in China [18].Additionally, ordinary least squares regression has been applied to evaluate the maximum annual NDVI linear trend during the study period and to estimate the greenness rate of change (GRC) for each pixel.The GRC is estimated by using the slope of the regression line for each pixel, where a slope > 0 indicates that the vegetation is recovering, and a slope < 0 indicates ongoing degradation [19,20].Trend measurements are fundamental in describing geographic phenomena that are related to vegetation phenology and accumulate over time.Grassland degradation processes often begin to appear as subtle variations in vegetation and depend on plant species, grazing, herd stocking rates, and phytosanitary controls.Thus, long-term memory evaluation based on time series analysis can enable the description of grassland trends of degradation or maintain vegetative growth conditions.The Hurst exponents, estimated from rescaled time series, can reveal the persistence levels of vegetation patterns, which can be very useful for decision-making regarding other biophysical phenomena identified by the phenological metrics of vegetation.Researchers have analyzed vegetation trends in North China using the slope of the NDVI/SPOT vegetation [21] and the Hurst exponent [22,23], and have found less degradation, lower sustainability, and lower stability for natural pastures than for forests.The estimate of the slope metric from the NDVI data, integrated annually, has been combined with the Hurst exponent using rescaled range (R/S) analysis and is used to create classes of vegetation degradation based on vegetation changes and long-term trends [24].R/S analysis and the Hurst exponent have been used for trend analysis in different fields of study, such as economics and computation [25][26][27][28][29]. Hurst exponents have also been estimated for a region in Tibet, China, and showed future trends in vegetation dynamics with persistence for a large portion of the forests and lower persistence (nonpersistence) for pastures, both of which are affected by topography and altitude [30].The authors also observed a linear trend for the NDVI annual mean values, confirming the previously observed pasture degradation.Traditionally, current phenological conditions are identified to evaluate pasture degradation, for which MODIS data have been shown to be well suited [31].Weighted linear combination multicriteria analysis (MCA) has been used to combine different levels of information using appropriate weights for the phenomena under investigation and to generate eigenvalues to determine the degrees of suitability and to make decisions based on the factors or criteria [32][33][34].However, when layers of complex data are inserted or when the data require specific restrictions, hybrid methods may be applied.Computational blocks in an Expert System apply human analysis and knowledge based on defined conditions to available databases and constitute a powerful tool to introduce information and make decisions in several areas of knowledge [35][36][37][38].The Expert System inspires the creation of algorithms to aid decision-making in complex and simple contexts with multiple alternatives or a few proposed solutions.We introduce a simplified expert system (ES S ), from which few solutions may resolve the problem posed where a layer involves complex classes and transcends weight combinations, such as the occurrence with Hurst exponents.The use of weighted MCA and ESs could fill the gap between analysis automation, eigenvalue generation for spatial analysis, and appropriate control over the factors related to statistical indices with values that do not represent the phenomenon linearly.The preparation of hypertemporal databases involves the elimination of invalid values while maintaining the time series behavior.This principle was adopted for the analysis of NDVI data series in pastures of the studied region and considered the specifics related to herbaceous vegetation and decumbent, tussock-forming foliage, and it can reflect photosynthetic activity for all soil covers, despite the potential for overgrazing, which is typically caused by local livestock.Time series data treatment methods can be used to eliminate spurious data points; to define the pattern of temporal and spatial profiles using wavelet functions, regression windows, and the Savitzky-Golay filter; they are also methods to replace low-quality pixels and to smooth and filter time series [39][40][41].The NDVI/MODIS database for grasslands was filtered via the cubic spline wavelet method with second-level decomposition, which was found to exhibit the equivalence with other methods to detect outliers, and a quadratic mean approach was employed to smooth the NDVI profiles.The aim of the present study was to analyze grassland development in Zona da Mata based on seasonal metrics and statistics of time series, integrating the weighted MCA and ES S .This methodology was employed in a raster-based geographic information system (GIS) and generated management guidelines from the temporal vigor classes, which may enable decision-making and support public policies in the agricultural sector.The decision-making considers the vegetation condition classes mapped using GIS methodology, as well as the edaphoclimatic characteristics of the region, pasture species, stocking rates, and socioeconomic aspects related to management practices.The pastures of Brachiaria spp.(Brachiaria brizantha and Brachiaria decumbens) are the most widespread species in the region and are perennial due to their higher biomass production, efficient bud break, and rapid tiller growth (https://www.infoteca.cnptia.embrapa.br/handle/doc/783242).Areas classified as persistently degraded pastures require greater attention, the adoption of herd stocking rate-monitoring procedures, promotion of pasture renovation, burn disuse and mowing, liming and fertilization due to aluminum presence, and low availability of phosphorus in the soils.For the areas classified with a high vegetative density and persistence over time, continuous stocking can be adopted, considering the rapid leaf growth of tillers and saturation of leaf biomass production.

Study Area
The Zona da Mata mesoregion has a grassland area of 1,213,506 ha, approximately 2.5 million inhabitants, a predominance of strong, hilly relief, granitic, and gneissic soils, and a humid subtropical climate (Cwa).Historically used for livestock, the main economic activity in the region is milk production, which plays an important part in milk production in the state of Minas Gerais.The state of Minas Gerais produces 9 billion liters of milk annually [42], constituting 27% of the national milk production.Zona da Mata is responsible for 10% of the total milk production from the state and is an important and traditional area for milk production.Milk production in the state of Minas Gerais has increased markedly since the first survey was performed in 1974, which estimated 2 billion liters of production.Most farms that involve extensive animal farming are largely covered with perennial grasslands of different forage species and have confined, and to a lesser extent, semiconfined production systems.These features were observed in a field survey and through an analysis of pixel time series from NDVI/MODIS images (Figure 1).A field campaign was carried out at specific points, with visits designed according to random sampling criteria and logistical convenience.Some pure pixels were selected in relation to a typical pasture temporal pattern for estimating correlations and comparisons with metrics and qualitative field information.MODIS images from 2012 were used for mapping the grassland areas in the region of interest.

Materials
Relevant data were extracted from the MODIS/Terra database (https://modis.gsfc.nasa.gov/data/dataprod/mod13.php) to obtain a series of images of Zona da Mata, map the grasslands, and estimate phenological metrics.Figure 2a shows the articulation of the MOD13Q1 tiles (H13V10, H13V11, H14V10, and H14v11) that were used to build mosaics for the dates and periods presented in Figure 2b.

Materials
Relevant data were extracted from the MODIS/Terra database (https://modis.gsfc.nasa.gov/data/dataprod/mod13.php)to obtain a series of images of Zona da Mata, map the grasslands, and estimate phenological metrics.Figure 2a shows the articulation of the MOD13Q1 tiles (H13V10, H13V11, H14V10, and H14v11) that were used to build mosaics for the dates and periods presented in Figure 2b.
RapidEye images (5 m spatial resolution) were used to assess the accuracy of the classifications, to check for changes in land use, and to support the field campaign as a source of high-resolution data.These images were obtained from Embrapa due to an agreement with the Brazilian Ministry of the Environment, which holds the image database.

Materials
Relevant data were extracted from the MODIS/Terra database (https://modis.gsfc.nasa.gov/data/dataprod/mod13.php) to obtain a series of images of Zona da Mata, map the grasslands, and estimate phenological metrics.Figure 2a shows the articulation of the MOD13Q1 tiles (H13V10, H13V11, H14V10, and H14v11) that were used to build mosaics for the dates and periods presented in Figure 2b.RapidEye images (5 m spatial resolution) were used to assess the accuracy of the classifications, to check for changes in land use, and to support the field campaign as a source of high-resolution data.These images were obtained from Embrapa due to an agreement with the Brazilian Ministry of the Environment, which holds the image database.

Methods
The full NDVI/MODIS series was first filtered using the cubic spline wavelet algorithm and second-level decomposition, which smoothed the time series and thus eliminated outliers and spurious data.This action corrected the time series for anomalous atmospheric effects while retaining the typical pasture frequencies.To analyze grassland development and degradation in the study area, several information layers were produced from the relevant seasonal metrics.Seasonal GIs were estimated for 16-day periods in spring 2012 and reflected the evolution of grassland conditions relative to the same period in previous years [18], with the SOS, POS, and season of vegetation growth taken into account.The slope was calculated based on the annual maximum NDVI data from 2000 to 2013, which could have indicated a linear growth or senescence trend [19,20].The exploratory metrics, i.e., the mean, minimum, and maximum NDVIs [14][15][16][17], were estimated due to the importance of comparing absolute values to estimated indices and trends.In addition, the temporal statistic Hurst exponent was calculated by R/S analysis [23] to rank the weighted criteria in terms of bias or spatiotemporal trends.For the spatial evaluation of grasslands, decision criteria were adopted based on weighted information layers, and conditions were defined based on technical knowledge of the temporal trends obtained through estimations of the Hurst exponent.The degree of grassland development was estimated by assigning weights to the geographical layers by the use of the weighted multicriteria analytical hierarchy process (AHP) [43] and by integrating ESS applied to the GIS called SIS MULT 4, with a single knowledge block and solution.The purpose of this procedure was to define the nonlinear trends through the Hurst exponents and apply them to the MCA results, expressed as eigenvalues.A normalization procedure based on a function available in GIS that applies histogram equalization to an 8-bit linear stretch, i.e., assigning values from 0 to 255, was used to standardize the floating point map database produced from phenological metrics.The flowchart for the analysis process and the application developed in the GIS for the multicriteria and detailed expert knowledge layers are presented in Figure 3a,b.

Methods
The full NDVI/MODIS series was first filtered using the cubic spline wavelet algorithm and second-level decomposition, which smoothed the time series and thus eliminated outliers and spurious data.This action corrected the time series for anomalous atmospheric effects while retaining the typical pasture frequencies.To analyze grassland development and degradation in the study area, several information layers were produced from the relevant seasonal metrics.Seasonal Gis were estimated for 16-day periods in spring 2012 and reflected the evolution of grassland conditions relative to the same period in previous years [18], with the SOS, POS, and season of vegetation growth taken into account.The slope was calculated based on the annual maximum NDVI data from 2000 to 2013, which could have indicated a linear growth or senescence trend [19,20].The exploratory metrics, i.e., the mean, minimum, and maximum NDVIs [14][15][16][17], were estimated due to the importance of comparing absolute values to estimated indices and trends.In addition, the temporal statistic Hurst exponent was calculated by R/S analysis [23] to rank the weighted criteria in terms of bias or spatiotemporal trends.For the spatial evaluation of grasslands, decision criteria were adopted based on weighted information layers, and conditions were defined based on technical knowledge of the temporal trends obtained through estimations of the Hurst exponent.The degree of grassland development was estimated by assigning weights to the geographical layers by the use of the weighted multicriteria analytical hierarchy process (AHP) [43] and by integrating ES S applied to the GIS called SIS MULT 4, with a single knowledge block and solution.The purpose of this procedure was to define the nonlinear trends through the Hurst exponents and apply them to the MCA results, expressed as eigenvalues.A normalization procedure based on a function available in GIS that applies histogram equalization to an 8-bit linear stretch, i.e., assigning values from 0 to 255, was used to standardize the floating point map database produced from phenological metrics.The flowchart for the analysis process and the application developed in the GIS for the multicriteria and detailed expert knowledge layers are presented in Figure 3a

Growth Index
The GI of the grasslands was calculated as the ratio between the NDVI/MODIS 16-day images during the spring of 2012 and the means for spring in 2000-2011.Each 16-day period was treated as one fortnight of the month, approximately, the first (1st) and second of the month (2nd), generating a map or factor, totaling eight factors between September 2012 and January 2013.
The beginning of spring, in September, corresponds to the SOS metric, and the beginning of January corresponds to the POS metric.These metrics are important for pasture phenology and describe the vegetation recovery.The GI is calculated using the following equation: where GI is the pasture vegetative growth rate, NDVIm is the vegetation index for day m, and NDVIn is the vegetation index mean for the same period in previous years.Therefore, NDVIm is the vegetation index for each 16-day period between September 2012 and January 2013, and NDVIn is the mean NDVI for the periods between September and January from 2000-2011.The GIs were classified based on the pasture growth observed in the field, temporal profiles, and high-resolution images.
The resulting indices were compared to the sampled pixels and to the knowledge of the visited regions.Several points were monitored in the field between September 2014 and April 2015, data were collected, and the resulting maps and legends were evaluated.The field points were continuously monitored.Twenty sites were visited, and the pasture conditions were recorded in all four cardinal directions for a total of 100 field observations that were located using a global positioning system (GPS).These sites were visited using state and federal roads during the dry and rainy seasons.The GIs were classified into the following categories: Very low (<-0.05),low (-0.05-0.00),balanced (0.00-0.05), high (0.05-0.10), and very high (>0.10).This classification was based on the distribution of the calculated indices, reference materials, field observations, histograms, and statistical agreement, according to the chi-square confidence interval.In total, eight information layers were generated and showed the evolution of pasture growth from the beginning to its peak.The degree of importance of each layer was evaluated using multicriteria evaluation and validation.

Linear Trend: Slope
The linear trends of the intra and inter-annual metrics, at either annual or seasonal scales, were used to indicate the degree of vegetation changes.Phenological metrics, combined with long or shortterm trends, were used to define management or planning classes to continuously monitor pasture conditions.Linear trends were estimated using the slope metric.The maximum annual NDVIs were

Growth Index
The GI of the grasslands was calculated as the ratio between the NDVI/MODIS 16-day images during the spring of 2012 and the means for spring in 2000-2011.Each 16-day period was treated as one fortnight of the month, approximately, the first (1 st ) and second of the month (2 nd ), generating a map or factor, totaling eight factors between September 2012 and January 2013.
The beginning of spring, in September, corresponds to the SOS metric, and the beginning of January corresponds to the POS metric.These metrics are important for pasture phenology and describe the vegetation recovery.The GI is calculated using the following equation: where GI is the pasture vegetative growth rate, NDVI m is the vegetation index for day m, and NDVI n is the vegetation index mean for the same period in previous years.Therefore, NDVI m is the vegetation index for each 16-day period between September 2012 and January 2013, and NDVI n is the mean NDVI for the periods between September and January from 2000-2011.The Gis were classified based on the pasture growth observed in the field, temporal profiles, and high-resolution images.The resulting indices were compared to the sampled pixels and to the knowledge of the visited regions.Several points were monitored in the field between September 2014 and April 2015, data were collected, and the resulting maps and legends were evaluated.The field points were continuously monitored.Twenty sites were visited, and the pasture conditions were recorded in all four cardinal directions for a total of 100 field observations that were located using a global positioning system (GPS).These sites were visited using state and federal roads during the dry and rainy seasons.The Gis were classified into the following categories: Very low (<-0.05),low (-0.05-0.00),balanced (0.00-0.05), high (0.05-0.10), and very high (>0.10).This classification was based on the distribution of the calculated indices, reference materials, field observations, histograms, and statistical agreement, according to the chi-square confidence interval.In total, eight information layers were generated and showed the evolution of pasture growth from the beginning to its peak.The degree of importance of each layer was evaluated using multicriteria evaluation and validation.

Linear Trend: Slope
The linear trends of the intra and inter-annual metrics, at either annual or seasonal scales, were used to indicate the degree of vegetation changes.Phenological metrics, combined with long or short-term trends, were used to define management or planning classes to continuously monitor pasture conditions.Linear trends were estimated using the slope metric.The maximum annual NDVIs were calculated for the time series from 2000 to 2012.The angular coefficient of the linear regression line was estimated using the following equation: where Slope is the slope obtained by the linear regression of the maximum NDVI (NDVI max ) for year I with n = 13.

Exploratory Metrics: Maximum, Minimum, and Mean
The analyzed vegetation development indices and trends along the hypertemporal series provide information regarding pasture development conditions.However, different typologies presenting higher, lower, or even different photosynthetic activity levels may exhibit similar behaviors.Exploratory temporal statistics may therefore supply important data regarding pasture stages and physiognomy that depend on different edaphoclimatic conditions, with the overall NDVI statistics being an indicator of the vegetation patterns.The maximum, minimum, and mean NDVIs for the entire series show the vegetation index amplitudes and levels and may support the classification of pasture conditions along with the linear and nonlinear trends.Exploratory metrics were obtained using algorithms in the GIS, and layers for the NDVI maximum , NDVI minimum , and NDVI mean were generated for the hypertemporal series, producing a total of three layers of information, which were then weighted and included in the MCA.
where n is the number of years, and NDVIi is the vegetation index for year i.

Hurst Exponent
The evaluation of long-term memory is useful for insertion into future scenarios and is a result of integrating seasonal metrics into the grassland development legend.R/S analysis (variable amplitude divided by standard deviation) is commonly used in financial and economic time series and describes temporal behavior.R/S analysis was introduced by Harold Hurst (1951) to determine the long-term water storage capacity of reservoirs, and its usage expanded to predict time series behavior.It was first applied by Hurst to the construction of a reservoir in the Nile River and was later consolidated into a method for the investigation of long and short-term memory.R/S analysis and the Hurst exponent have been structured to determine memory in time series, where H < 0.5 indicates anti-persistence or noncorrelation, with the phenomena tending to revert in the future; H = 0.5 indicates random motion, also called Brownian motion, with values not reaching or defining a trend; and H > 0.5 indicates persistence or memory, with past values affecting future results [21].In economics, when applied to price series, investments, and capital markets, the Hurst exponent is used to detect memory effects and trends [26,[44][45][46].It has also been used in academic, environmental, and remote sensing applications [47].Power-law relations for spatiotemporal transects in two types of pasture in China have been estimated [48].Ashutosh et al. [27] and Soterroni et al. [49] concluded that R/S analysis is equivalent to the wavelet method for estimating the Hurst exponent for both long and short time series.Other studies [24,30] have used series of SPOT-VGT and AVHRR NDVI images to evaluate the future trends of pasture and overall vegetation degradation.NDVI/MODIS images have been measured and reported from February 2000 for images with a 250 m spatial resolution constituting medium and short time series; however, these images do not meet the minimum recommendation for R/S analysis of 500-1000 data points.Recent studies of short series [50] have developed a new method for databases with 100 data points.Additionally, short biotic data series have been used for R/S analysis by adopting lags or arithmetic scales instead of the classic geometric or binary scales [29].R/S analysis is biased if filtering is not adequately applied and is less affected by periodic data because of the calculations of the cumulative deviations and means [51].In addition, all methods for estimating the Hurst exponent may result in over or underestimation, depending on annual variations; R/S analysis results in underestimation with increasing annual variances or amplitudes.The normal distribution of the estimated Hurst exponents, which occurs for artificial or real series, has also been studied [28].Quasinormality has been reported for adimensional NDVI data, and logarithmic returns were used in the analysis of medium and short series of seasonal data with clear periodicity due to the Gaussian distribution of the deviations calculated by the R/S analysis.
The time series were initially filtered by eliminating outliers and using the quadratic mean to replace spurious values.According to the classic methodology, to determine the measuring scale or division of time series of image cells or pixels, the total series length (τ) must first be divided by two until rounded to the smallest subdivision.Images from the 305 NDVI/MODIS (MOD13Q1) product were collected and analyzed using the 16-day means from 18 February 2000 to 9 May 2013 (LP-DAAC, 2013).The algorithm was developed using a tool and an algorithm contained in the GIS, which allowed us to review the grassland map for consistency with the MODIS and RapidEye images.
Hurst exponents were estimated for the matrices arranged into time series using the principles established in the equations below and were numbered from 6 to 10 [26].The calculations for each set of pixels in time were performed by binary division.In principle, the mean of the returns, or observed r values for each segment τ block or size, and σ (Si; standard deviation) within the segment were calculated.The algorithm considers the maximum or minimum deviation around the mean along a block, and the maximum amplitude divided or normalized by the standard deviation is obtained.After calculating the R/S for each block, the mean is obtained and is related to the H exponent according to Equation (9).For τ blocks or geometric lags with base 2n, where n = 1, 2, 3,... I, R/S analysis assumes values that result in random motion with H = 0.5 = β in the linear regression equation.
H is therefore the angular coefficient of the line obtained by a linear regression of Log(R/S) x Log (τ) and may assume values close to 0.5 (± 0.02) for a random series.H < 0.5 indicates the reversibility of the observed trend through the different metrics, whereas H > 0.5 indicates continuance of the current trend of geographical phenomena.For the case of grasslands, to better use the results of the temporal analysis of the NDVI/MODIS images and to estimate the Hurst exponents, the seasonal, intra-annual, and interannual phenological parameters and vegetation growth metrics were estimated.For the series of 305 NDVI/MODIS images, because they presented a medium-sized series, to better fit the lags, binary division was performed from the image total (τ, as the first block), followed by τ/2n blocks, with n = 1, 2, 4, 8, and 16, and rounding of the smallest block.

Multicriteria Analysis and Expert System
As there are several layers of information regarding pasture conditions, the criteria and methods are needed to extract the maximum weightings for the importance, integration, and spatial relationships between the factors.This process results in a single map with scores or eigenvalues that define the phenomena with a legend that enables management decisions.Often used in spatial analyses, MCA is based on the mapping of factors in information layers and defining the relevance level of the layers for the construction of a final map with criteria that help decision-makers [52,53].
The weighted multicriteria method was applied using the Gis (16-day periods from September to December 2012, with a total of eight information layers), slopes, maxima, means, and minima via the following equation: where M is the eigenvalue or score produced for the pixel in column I and row j from the sum of F factors multiplied by the respective weight W, n varies between 1 and 12 in the present study, and c is the Boolean factor (0 or 1).The pairwise hierarchization method [43] was applied along with AHP analysis to the degree of favoring grassland development using the numeric scale presented in Table 1.Intermediaries between the criteria Weighted gradient between factors The numerical scale was chosen according to the decision criteria, wherein the factors are close to or overlap other factors in terms of importance for the analysis of pasture conditions.While up to nine factors are typically used, 12 factors were used in the present study, with a correction for the random consistency index for the number of factors added [54] (Table 2).After tabulating the cross-comparison between the factors, using each classification criterion and its reciprocal value, the data were normalized by dividing them by the sum of each real and reciprocal criterion.The weights were then calculated as the means of the normalized values.To validate the weights, the products of the matrices of pairwise values and weights were generated, and the final mean of the proportionality ratio (product divided by the weight) was used to estimate the consistency index (CI).The consistency ratio (CR) was obtained by dividing CI by RI (random index or random consistency index [43,54]).CR < 0.1 indicates coherence between the values attributed to the pairwise criteria and the estimated weights for the factors.
Weight estimation and evaluation were performed using the following matrices and equations: where F is the factor, n is the number of factors, v is the value associated with the pairwise evaluation criteria, and r is the reciprocal value, arranged in an n x n matrix.
where M is the criterion value associated with factor n of matrix A and is normalized by summing each column.Matrix C where W is the final weight or mean of the sum of row I of matrix B, and n is the number of factors.
(AC) n, 1 n r=1 a nr c r1 = a n1 c 11 + . . .+ a nn c n1 (15) where AC is the product matrix of matrices A and C, with n rows and one column.
The mean of the ratio between the product matrix and its respective weight yields the estimated λ max , from which the CI is calculated as follows: CR = CI/RI (17) After validating the weights to be applied to the factors represented by the information layers in the GIS, eigenvalues were estimated using the weighted multicriteria method on the normalized layers by the binarization of 8-bit images, using values from 0 to 255, via a built-in function of the GIS [55].The resulting map represents the degree of pasture development that depends on the seasonal factors or metrics listed above.However, because we needed to understand future scenarios using temporal statistics, a nonlinear trend factor, represented by the Hurst exponent (H), was inserted into the evaluation process.To accomplish this, an artificial, or simplified, ES S was created to insert the conditions of the applied H into the map generated by MCA in a single block, with only one solution, but containing technical knowledge of H and the factors used in the MCA.An ES S is a variation of an artificial intelligence system, where the system analysis rules insert conditions that are inherent to the studied phenomena.This system is characterized by the data range, information, and ability to simplify the expected conclusion and explain decisions based on the conditions or "thought" [56].The mechanism cited is similar to MCA, except with rules that are necessary to understand and interpret the complex information layer and to allow for several solutions and explanations.In MCA, the rules originating from the knowledge base are inserted into the mathematical expressions through a careful evaluation of the weights and Boolean constants.Inferences are drawn from a knowledge base translated into rules through processing the decision-making rules for one or more solutions.In the context of a GIS, the ES S is represented by the flow chart presented in Figure 4.An application was implemented in the ArcInfo 9.0 interface (the GIS product of Esri) to integrate map M resulting from the statistically weighted MCA.The system block is based on rules, and the legend items corresponding to the grassland development gradation are generated according to combinations of eigenvalues with H.
M eigenvalues were classified based on the natural distribution of the data in a histogram, field campaign data (Figure 5), and images with high spatial resolution according to temporal vegetative conditions; growth capacity and vigor were classified as low, medium, or high.An application was implemented in the ArcInfo 9.0 interface (the GIS product of Esri) to integrate map M resulting from the statistically weighted MCA.The system block is based on rules, and the legend items corresponding to the grassland development gradation are generated according to combinations of eigenvalues with H.
M eigenvalues were classified based on the natural distribution of the data in a histogram, field campaign data (Figure 5), and images with high spatial resolution according to temporal vegetative conditions; growth capacity and vigor were classified as low, medium, or high.The low class corresponds to grasslands with little vegetation development, a pattern of degradation over time, sparse vegetation, compacted soil, and overgrazing.The medium class corresponds to grasslands with reasonably balanced development, intermediate vigor, and sparse to semidense vegetation.The high class corresponds to forage with good vegetative growth, good tiller density, and high photosynthetic activity.These classes therefore reflect growth over time with vegetation size.The Gis and slopes, in addition to the exploratory metrics, which indicate vertical stratification at the herbaceous level and density through the NDVI, play an important part in class definition.An application was implemented in the ArcInfo 9.0 interface (the GIS product of Esri) to integrate map M resulting from the statistically weighted MCA.The system block is based on rules, and the legend items corresponding to the grassland development gradation are generated according to combinations of eigenvalues with H.
M eigenvalues were classified based on the natural distribution of the data in a histogram, field campaign data (Figure 5), and images with high spatial resolution according to temporal vegetative conditions; growth capacity and vigor were classified as low, medium, or high.The low class corresponds to grasslands with little vegetation development, a pattern of degradation over time, sparse vegetation, compacted soil, and overgrazing.The medium class corresponds to grasslands with reasonably balanced development, intermediate vigor, and sparse to semidense vegetation.The high class corresponds to forage with good vegetative growth, good tiller density, and high photosynthetic activity.These classes therefore reflect growth over time with vegetation size.The GIs and slopes, in addition to the exploratory metrics, which indicate vertical stratification at the herbaceous level and density through the NDVI, play an important part in class definition.The H and M combinations were included in the final legend, which was derived from inferences of the rules of integration of the H statistic and M classes obtained through the multicriteria crossing of metrics, as shown in Table 3.The definition of the classes was based on the vegetative conditions observed in the field, temporal profiles of pure and mixture pixels, and high-resolution imagery; the metrics between known and visited regions were compared.Approximately 900 km were covered in the field campaign, with 80 ground points observed, and data were collected between July and December 2014; 20 additional points were visited between March and April 2015.In this research, qualitative data, such as soil cover levels, vegetative conditions, tiller formation, and occurrences of weeds and termites, were tabulated and analyzed.Field annotations, reference iconography, and pixel values were analyzed using the GIS, and coordinates were determined using GPS for the definition of the GI and M classes as well as to associate the grassland evolution with the phenological metrics and Hurst exponents.The suitability of the estimations and field-to-map associations were tested and established based on further observations of histograms and statistical agreements through intervals of confidence using a chi-square distribution table.
The following 12 classes were defined, with names and explanations derived from the field survey, pixel samplings, and analysis of imagery time series:

•
Degraded and reversible: The vegetation index decreases over time, with low plant cover; the future trend is an increase in plant cover.

•
Degraded and unstable: A decreasing vegetation index with low plant cover and a small or unobservable trend.

•
Degraded and slightly persistent: The decreasing vegetation index indicates low plant cover, a moderate trend for maintenance of the degradation process, sparse grasslands, compacted soils, and overgrazing.

•
Degraded and persistent: The decrease in the vegetation index shows low plant cover and a strong trend for maintenance of the grassland degradation.

•
Semidense and reversible: This class has a moderate to high vegetation index and sparse to semidense plant cover, and it shows a trend toward plant density instability from fallow, burnings, or crop rotation.

•
Semidense and unstable: This class is characterized by moderate plant cover and a low or no temporal trend due to complex farming practices such as variations in the cattle stocking rate.

•
Semidense and slightly persistent: Areas in this class have moderate plant cover and show a reasonable trend for sustained grassland vigor.

•
Semidense and persistent: This class has reasonable plant cover and exhibits a strong trend for sustaining this vegetation condition.

•
Dense and reversible: There is a high level of plant cover, but with a trend toward reversal of the plant density over time, from the adjacent management or degradation.

•
Dense and unstable: There is a high level of plant cover and a low or no temporal trend.

• Dense and slightly persistent:
There is a high level of plant cover and a moderate trend for sustained grassland vigor.

•
Dense and persistent: This class has a high level of plant cover and a strong indication for sustained density or a high vegetation index.
The analyzed indices and trends of vegetation development along the hypertemporal series provide information regarding grassland development conditions.However, different typologies with different photosynthetic activities may present similar behavior.The aim of determining a higher number of metrics, combinations, and classes is to differentiate the temporal profiles of the NDVI and enable decision-making using vegetation vigor and development scenarios.The decision-making map is intended to guide field activities, while the temporal vigor classes aim to lay the groundwork for general pasture management and public policy planning.According to the identified vegetative condition classes and considering the long-term memory aspects, management actions may be taken by stakeholders, producers, public technical assistance agents, and decision-makers in general.Processes with persistent character in areas that present degradation or even in grasslands classified as dense vegetation, but with instability or reversibility bias, can be adopted for incisive management measures that will also allow for long-term conservation and observations of local pest infestation conditions.These processes can also address the occurrences of overgrazing, burning, or degradation caused by soil compaction or low fertility.Grasslands classified as semidense also require grazing control and corrective actions after local verification.Therefore, for those classes that indicate degradation or semidense vegetation, stocking rates, fertilization, and renovation control actions are important due to the eventual need for soil protection or vegetation recovery.These grasslands may be suffering from long-term damage due to scale effects in the time series, indicating the persistence aspects of the degradation phenomenon.Classes that are reversible point to antipersistent scale phenomena that are related to ongoing practices that may reverse degradation, such as controlling stocking rates, grassland renovation, avoidance of traditional fire use, mowing, and fertilization.There may also be indications of improper practices, such as disproportionate increases in herd and grazing levels, which cause the decrease in vegetation density and necessitate the adoption of grassland monitoring, canopy height measurements, light interception, and adequate herd sizing.The decumbent, tussock-forming canopy that is predominant in the grasslands in the region allows for good evaluations from MODIS sensor images due to good soil cover and plant volume.

Results
The grassland vegetative conditions observed in the field campaign were compared with the results obtained by the GI with a confidence level of 95%.Strong qualitative and quantitative correlations were estimated from the data and chi-square tests.The accuracies of the producer and user for the grassland/nongrassland binary map were approximately 74% and 100%, respectively, with low errors of omission and commission and no errors related to the masked areas of nongrasslands.We obtained a Kappa index of 0.90 and an overall accuracy of 91.5%, meaning that the resulting grassland map exhibited very good quality.
Coherent distribution patterns were observed for the values between the GI and slope layers, enabling an equivalent distribution of weights between maps.These, together with a hierarchical analysis of the remaining information layers, yielded good consistency for the factor weights by AHP analysis.In the factor normalization, maps were generated with scales that were compatible with MCA.Normalization was applied to all maps, except for the H map due to the characteristics of nonlinear analysis and the heterogeneity of the results, as explained in the Methods Section.The ES S was used for factor treatment and integration with map M.

GI for Spring
An analysis of the GI maps for spring 2012 using NDVI/MODIS for 16-day periods showed that 60% of the area showed symptoms of grassland degradation based on field data, high-resolution images, and pure pixels.Areas with good growth patterns occupied approximately 4% of the area.Balanced or stable growth during the analyzed period was observed in 36% of the total area.A qualitative-quantitative analysis using the chi-square test at a 95% confidence level (calculated p-value > tabulated p-value (χ 2 0.95 )) was used to classify the GI as follows: Very low (<-0.05),low (-0.05-0.00),balanced (0.00-0.05), high (0.05-0.10), and very high (>0.10)(results of the estimates, as percentages, are shown in Figure 6a).However, the results were normalized to compare the studied factors or metrics, and the growth levels were redistributed by a growing linear function between 0 and 255.The three main maps (out of eight) produced better-visualized growth evolution from the initial 16-day period in September to the end of October and December 2012 and are presented in Figure 6b.These maps show the regions with different levels of degradation or decreased grassland growth capacity as darker colors.The lighter colors show those areas with higher growth maintenance or recovery during spring 2012 when compared to the NDVI data between 2000 and 2011.

Slope Metrics
The slopes indicated that approximately 44% of the area presented low or very low linear trends, indicating low growth.This finding reveals important degradation aspects that are related to grassland development over time.Approximately 27% of the area remained balanced, with stable growth, and 29% presented a tendency toward high growth (according to the chi-square test with 95% confidence) compared with the field data and remote sensing images for the same samples used for the GI evaluation (Figure 6c).The slopes for the maximum annual NDVI were classified according to the same classes used for the GI to evaluate grassland growth: Very low (<-0.004),low (-0.004--0.001),balanced (-0.001-0.001),high (0.001-0.004), and very high (>0.004).
The normalized map of the slopes, between 0 and 255, showed a degraded, stable growth and high development areal distribution similar to that estimated for the GI at the end of spring (Figure 6d).Overall, the percentage distributions between the different classes were similar.However, the percentage for the balanced class was slightly lower, and for the high and very high classes, the slopes were much higher than the GI.

Maximum, Minimum, and Mean NDVI
The Gis and slopes are complemented by the maximum, minimum, and mean NDVIs along the entire time series because these statistics reflect aspects of the vegetation structure.After smoothing the series to eliminate spurious data and restructuring the temporal pattern or signature, these metrics provide additional information about eventual anomalous atmospheric effects.These metrics, therefore, indicate the vegetation state through the NDVI amplitude, even if it presents vegetative growth or degradation.NDVI maximum was between 0.40 and 0.93, NDVI minimum was between -0.06 and 0.70, and NDVI mean was between 0.28 and 0.81 for the entire study area.The distribution maps of the normalized metrics are presented in Figure 6e.NDVI minimum and NDVI mean demonstrated very similar patterns, with similar distributions and regions of occurrence in the temporal curves.The combination of these three layers of information may, for example, show a predominance of sparse herbaceous vegetation.If the Gis or slopes indicate higher vegetative growth, the amplitudes between the maximum and minimum NDVIs and the NDVI means are low.The opposite may also occur with dense forage, i.e., high NDVI amplitudes with lower growth trends.Thousands of different combinations of the pixels of the analyzed layers can therefore be produced, indicating different degrees of development, which can be revealed as persistent or volatile via R/S analysis through the integration of MCA and the ES S .(g) (h)

Temporal Statistic (H)
The nonlinear analysis of the dynamic processes based on R/S analysis and scale phenomena through H show future scenarios with nonperiodic cycles for the trends of sustainability, persistence, and long-term memory, as well as reversibility, randomness, and unpredictability.The estimates for the pixels in the hypertemporal series presented a high degree of fitting (Figure 6f), and the determination coefficient (R 2 ) indicated high accuracy between the independent and dependent variables via the log-log representation of τ versus R/S.The areas and percentages occupied by H classes showed a predominance of areas with a slightly persistent memory effect (approximately 69%) and are presented in Table 4. Areas with strong persistence (higher than 0.65) are highlighted and represent 20% of the total area.The antipersistent and random behavior classes constituted slightly more than 11% of the total area.According to the literature, these areas are usually associated with more complex land management, involving crop rotation, fallow, burning, degradation, or crop replacement.In other scientific topics, such as economics, the environment, computation, or geophysics, the Hurst exponent is used to detect memory effects as well as trends in time series, spatial sets, and geometry.H may be employed for the interpretation of long and short time series from satellite image databases.To estimate, filtering and stabilization were adopted for the NDVI MODIS time series, in addition to the necessary R/S estimates.
A map with the classes of the estimated H is presented in Figure 6g.These results compose the system of rules for the implemented ES S in the GIS for the treatment of map M and result from the weighted crossing of the 12 layers, or factors, in the MCA.The results of the integration of the MCA and ES S are described in the next section.

Multicriteria and Expert Systems
Both MCA and ES S can yield similar results, considering that conditions can be included in the multicriteria through the determination of weights by the AHP and factor knowledge.Moreover, layer crossing based on information regarding the studied phenomena can be incorporated in an algorithm or application in an information system such as a GIS.The weights were generated based on the relationship to the importance defined through the AHP, which is essential for the application of MCA, and a management legend with the H exponent in the ES S was created using GIS.

MCA
The first stage after the normalization of the information layers was to populate a table with the numerical relationships of importance using the AHP.Technical criteria for each parameter and pairwise comparisons were considered and used to form a matrix, followed by the calculations necessary to define and evaluate the consistency of the estimated weights.The results from the pairwise analyses between factors and their relative importance were estimated, and the factors considered very close in terms of importance were assigned a value of one and were real and reciprocal.The indices and weights estimated through the matrices are presented in Table 5.The estimated CR was 0.05, showing coherence between the weights because it is recommended that the CR be less than 0.1.Saaty [54] recalculated the Ris for several factors greater than nine, which is valid for the present study, where 12 factors were used.The distribution of the vegetation vigor classes over time, based on layers linearly weighted by M, which was normalized and classified based on histogram fitting to natural breaks, is presented in Figure 6h.Areas presenting lower development, especially for the Gis in the early growing season and lower slope values indicated by the exploratory metrics, were classified as low vigor, with lower eigenvalues than the average.The group around the mean was classified as medium vigor and included the intermediate values of the temporal metrics or compensation values.Grasslands with higher Gis, slopes, and NDVI maximum eigenvalues and lower compensatory effects of the metrics for the end of spring or NDVI minimum were classified as having high temporal vegetative vigor.The percentages of these classes are presented in Table 6.The MCA and ES s results with the NDVI time series were integrated to map the trends of grassland development, degradation, and a long-term scenario.The map summarizing the vegetative vigor over time for the region of interest is presented in Figure 7 and shows a predominance of sparse vegetation and sparse to semidense areas with low persistence.The data in the figure show areas of magnification in the North and South of Zona da Mata, with a prevalence of grassland degradation and higher leaf density.The predominance of degraded to semidense development over time indicates low pasture vigor, soil compaction, and overgrazing, as was observed in the field.The MCA and ESs results with the NDVI time series were integrated to map the trends of grassland development, degradation, and a long-term scenario.The map summarizing the vegetative vigor over time for the region of interest is presented in Figure 7 and shows a predominance of sparse vegetation and sparse to semidense areas with low persistence.The data in the figure show areas of magnification in the North and South of Zona da Mata, with a prevalence of grassland degradation and higher leaf density.The predominance of degraded to semidense development over time indicates low pasture vigor, soil compaction, and overgrazing, as was observed in the field.The grassland development classes obtained from the crossing of metrics treated by MCA and the H statistic, using the algorithm in the GIS based on the rules of inference imposed on the ESs, i.e., SIS MULT 4, are presented in Table 7. Dominance of a slight to strong memory effect was observed in both vegetation undergoing degradation and vegetation with good development.These results indicate that degradation or conservation processes tend to persist over time and are reversing in low proportions, based on observations of the obtained scenario.
Most of the area is covered by the semidense class, which was designated as the class that generally contains grasslands with intermediate soil cover and degradation aspects, indicating that they are potentially in the process of loss of relevant grazing biomass.There is also a predominance of the semidense slightly persistent class, indicating that the conditions of sparse to dense vegetation are accompanied by moderate or slight long-term persistence, which imposes difficulties for the restoration of these grasslands in the short term.This condition requires adequate management The grassland development classes obtained from the crossing of metrics treated by MCA and the H statistic, using the algorithm in the GIS based on the rules of inference imposed on the ESs, i.e., SIS MULT 4, are presented in Table 7. Dominance of a slight to strong memory effect was observed in both vegetation undergoing degradation and vegetation with good development.These results indicate that degradation or conservation processes tend to persist over time and are reversing in low proportions, based on observations of the obtained scenario.Most of the area is covered by the semidense class, which was designated as the class that generally contains grasslands with intermediate soil cover and degradation aspects, indicating that they are potentially in the process of loss of relevant grazing biomass.There is also a predominance of the semidense slightly persistent class, indicating that the conditions of sparse to dense vegetation are accompanied by moderate or slight long-term persistence, which imposes difficulties for the restoration of these grasslands in the short term.This condition requires adequate management practices such as soil fertility improvement, pasture renewal, and herd stocking rate control, which are traditional difficulties for dairy farming in the region.

Discussion
NDVI/MODIS data yield important information for the study of vegetation.The analysis of herbaceous vegetation demands a spatial and temporal approach because of the changes occurring at the vegetation scale and the dynamics.The temporal, seasonal, annual, and phenological metrics obtained through vegetation indices from remote sensing images fill a gap in the larger-scale characteristics of vegetation conditions and high territorial extension in grassland management areas.Research carried out in China has shown that the GI can effectively identify areas and levels of vegetation temporal evolution during the growth season [18].Similarly, the slope factor has been observed to be effective in investigations of grassland linear growth trends in Brazil and worldwide [21,58].However, the combination of seasonal and annual metrics with descriptive metrics reveals information regarding vegetation growth, status, and structural typology [17].MCA, viewed as effective in several studies of decision-making, shows or combines vegetation, photosynthetic conditions and profiles over time and defines temporal vigor classes.Some studies of grassland degradation or management have reported that the typologies, biomass, and conditions may be estimated through remote sensing [2,3,7,15,16].Simultaneously, processes may persist over long periods because of traditional production systems and long-term cultivation and practices with low efficiency on farms.The persistence of temporal phenomena was estimated by the R/S analysis and Hurst exponents-the memory effect that has been extensively studied in several fields of study [26,29,30].An ES S was built to obtain solutions for action planning and to diagnose several processes involving decision-making from a knowledge base of the study area.The ESs enables the use of, and inferences to be made from, algorithms that attempt to mimic human thought [35,38,56,57].In some situations, only one deterministic solution may be sufficient to integrate the classified information layers, which cannot be efficiently normalized and weighted by MCA.
Areas with low vegetative vigor were mapped by MCA over approximately 24% of Zona da Mata and were mainly distributed in the Northwest, central area, and Southeast, where signs of degradation were also observed in high-resolution images and in the field.Approximately 46% of the area was evaluated as having reasonable growth conditions with sparse leaf density and was classified as semidense or sparse vegetation; these areas were mostly distributed in the central area, West, and South.Approximately 30% of the region, in the West and South, was classified as having dense herbaceous vegetation, with more recurrent tillering and more efficient clump development over time, regardless of grazing, indicating adequate management of the production system, together with favorable edaphoclimatic conditions.The raster algorithm built in the ES s in the GIS revealed the predominance of low long-term persistence, as indicated by the H exponent for the identified vegetation vigor stages.However, temporal analysis revealed no anticorrelations or notable reversibility or instability trends and low antipersistence and showed a degree of balance or sustainment of the ongoing degradation or recovery processes.Areas with persistent and unstable degraded grasslands were concentrated in the North and Southeast and constituted approximately 21.5% and 2% of the region, respectively.These areas were near coffee plantations with hilly relief and a tropical climate, which may partly explain the temporal instability and low growth.Approximately 40% of the area was mapped as presenting sparse or semidense vegetation, which is persistent over time, especially in the West, South, and center of the study area.Of this area, 6.5% was included in the unstable class and could therefore possibly change its status in the long-term.The higher percentage observed for the semidense and slightly persistent class (approximately 33%) arose from the region's agricultural practices and edaphoclimatic conditions.The tradition of raising dairy cattle in subtropical climates may consolidate production systems with grasslands with a lower propensity for degradation stricto sensu, although these grasslands present some characteristics or aspects of loss of productivity or capacity for vegetation self-recovery.Approximately 27% of the pastures presented a higher vegetation density with long-term persistence in the West-central area, South, and Southeast, and 3% of the area was included in the unstable dense grassland class.Additionally, for this case, the favorable edaphoclimatic conditions and efficient production practices in extensive systems explained the distribution of the classes dense and slightly dersistent, dense and persistent, dense and unstable, and dense and reversible.
The map of grassland development for Zona da Mata may guide decision-making in management planning, land use policies, grassland recovery, and monitoring.Using MCA and raster algorithms based on an ES S , the spatiotemporal phenological metrics from NDVI are combined and degraded areas or areas undergoing degradation are identified.From this information stored in the GIS and presented in maps, it is possible to constantly monitor the pastures and support management decisions in the field.In general, control of the stocking rate for cattle in extensive livestock areas aims to correct consumption of pasture plant biomass in terms of the number of leaves originating from tillers, which regrow and continuously generate new material for grazing.Therefore, a decline in leaf biomass occurs after approximately 35 days, resulting in increases in the number of internodes and plant structures with lower feed conversion.Thus, the occurrences of excessive rest periods can be as harmful as overgrazing, as dead material may increase on grasslands.Weeds and even the natural regeneration of forest species, formation of forest fragments, and changes in land cover may also occur in extreme cases of grassland abandonment.The study region is in the Atlantic Forest domain region with a humid subtropical climate, which favors rapid forest regeneration.However, due to the aluminous soils with low fertility, inadequate stocking rates have led to soil compaction and losses of pasture productivity.These factors promote long-term memory processes according to the information obtained from field data and that observed in a small series of high-resolution satellite images.

Conclusions
The present study used the NDVI/MODIS Terra data series to estimate pasture growth-related metrics and integrated them through multicriteria and categorization by Hurst exponents using a simplified expert system employed to demonstrate long-term memory processes in conjunction with vegetation phenological parameters expressed in time series.Several other metrics and parameters related to NDVI variations may be evaluated through Hurst statistics, e.g., aboveground biomass [47].
Thus, we believe that this methodology may be employed for constant monitoring of grassland vigor using the MODIS data or could be adapted to other time series from other sensors.
The effects of pasture degradation detected in this study region may influence a reduction in livestock areas, as pointed out in a recent study, despite a relative increase of 17% in pasture areas in the Brazilian territory during the 21 st century [59].A prior study also suggested that social aspects, such as low levels of human development presented by municipalities, may influence the environmental and grassland degradation associated with decreased rainfall and the adoption of improper soil conservation practices [60].
Based on the NDVI/MODIS data and the resulting metrics, the following conclusions are drawn: • Most grasslands in Zona da Mata (approximately 61.5%) are degraded or are undergoing degradation with long-term persistence; • Approximately 27% of the grasslands presented adequately sustainable density conditions, according to the temporal analysis methods; and The phenological metrics based on remote sensing and temporal statistics were used effectively to build a management legend to provide guidelines for decision-making in relation to grassland development derived from NDVI/MODIS data and field evaluations.
Ultimately, knowledge of the dimensions and tendencies of pasture degradation makes it possible to outline actions and policies for regional development, given that dairy cattle farming is of paramount importance for the Zona da Mata region.In addition to monitoring the grassland dynamics for the studied region, this methodology may be applied to several other regions in the country and to other parts of the world to indicate the application of specific restoration strategies for each location with varying resource intensities, as exhibited by temporal vigor.

Figure 1 .
Figure 1.Location of the study area and field observations.

Figure 1 .
Figure 1.Location of the study area and field observations.

25 Figure 1 .
Figure 1.Location of the study area and field observations.

Figure 3 .
Figure 3. (a) Overall flowchart of the proposed method of operation; (b) application developed in the geographic information system (GIS) to cross inference layers.

Figure 3 .
Figure 3. (a) Overall flowchart of the proposed method of operation; (b) application developed in the geographic information system (GIS) to cross inference layers.

Figure 5 .
Figure 5. Overall conditions observed in the field: (a) Degraded grasslands with sparse herbaceous vegetation; (b) grasslands in good condition, with semidense or medium-sparse vegetation; and (c) grasslands with dense vegetation and adequate herd size.

Figure 4 .
Figure 4. ES S block (source: Adapted from [57]), procedures, and GIS operations used in the methodology.

Figure 5 .
Figure 5. Overall conditions observed in the field: (a) Degraded grasslands with sparse herbaceous vegetation; (b) grasslands in good condition, with semidense or medium-sparse vegetation; and (c) grasslands with dense vegetation and adequate herd size.

Figure 5 .
Figure 5. Overall conditions observed in the field: (a) Degraded grasslands with sparse herbaceous vegetation; (b) grasslands in good condition, with semidense or medium-sparse vegetation; and (c) grasslands with dense vegetation and adequate herd size.

Figure 6 .
Figure 6.(a) GI class percentages for Zona da Mata; (b) normalized maps of the GI for selected 2012 periods compared with the hypertemporal series; (c) slope metric percentages; (d) slope normalized map for the hypertemporal series; (e) maximum, minimum, and mean NDVI normalized maps for the hypertemporal series; (f) typical example of fitting for a pixel from the hypertemporal series for the region of interest, with H showing a memory effect; (g) map of the H classification distribution; and (h) map of the spatial distribution of vegetative vigor resulting from MCA.

Figure 6 .
Figure 6.(a) GI class percentages for Zona da Mata; (b) normalized maps of the GI for selected 2012 periods compared with the hypertemporal series; (c) slope metric percentages; (d) slope normalized map for the hypertemporal series; (e) maximum, minimum, and mean NDVI normalized maps for the hypertemporal series; (f) typical example of fitting for a pixel from the hypertemporal series for the region of interest, with H showing a memory effect; (g) map of the H classification distribution; and (h) map of the spatial distribution of vegetative vigor resulting from MCA.

Figure 7 .
Figure 7. Grassland development as indicated by vegetative vigor and a long-term scenario.

Figure 7 .
Figure 7. Grassland development as indicated by vegetative vigor and a long-term scenario.

Table 3 .
Decision process rules and inferences.

Table 4 .
Areas and percentages of H classes.

Table 5 .
Estimated indices and weights.

Table 6 .
Results of the MCA of grassland conditions over time.

Table 4 .
Areas and percentages of H classes.

Table 5 .
Estimated indices and weights.

Table 6 .
Results of the MCA of grassland conditions over time.

Table 7 .
Results of the ES S analysis of the long-term scenario of grassland development.