Time Series Analysis of Land Cover Change : Developing Statistical Tools to Determine Significance of Land Cover Changes in Persistence Analyses

Despite the existence of long term remotely sensed datasets, change detection methods are limited and often remain an obstacle to the effective use of time series approaches in remote sensing applications to Land Change Science. This paper establishes some simple statistical tests to be applied to NDVI-derived time series of remotely sensed data products. Specifically, the methods determine the statistical significance of three separate metrics of the persistence of vegetation cover or changes within a landscape by comparison to various forms of ―benchmarks‖; directional persistence (changes in sign relative to some fixed reference value), relative directional persistence (changes in sign relative to the preceding value), and massive persistence (changes in magnitude relative to the preceding value). Null hypotheses are developed on the basis of serially independent, normally distributed random variables. Critical values are established theoretically through consideration of the numeric properties of those variables, application of extensive Monte Carlo simulations, and parallels to random walk processes. Monthly pixel-level NDVI values for the state of Florida are analyzed over 25 years, illustrating the techniques‘ abilities to identify areas and/or times of significant change, and facilitate a more detailed understanding of this landscape. The potential power and utility of such techniques is diverse within the area of remote sensing studies and Land Change Science, especially in the context of global change. OPEN ACCESS Remote Sens. 2014, 6 4474


Introduction
Interest in human impact on the natural environment has increased dramatically over the last quarter century [1,2], as evidenced by such multi-disciplinary, multi-national programs as the International Geosphere Biosphere Programme [3].Over roughly the same period, repeat digital, synoptic measures of the Earth's surface have stimulated many Land Change Science (LCS) research questions designed to improve understanding of human-environment interactions.Satellite remote sensing is ideal for this task, providing consistent, repeatable measurements across a suite of spatial scales, and is well-suited to capture processes of land surface change, such as human activity, fire, floods, etc. [4].Detection and characterization of such changes is often the first step in understanding the mechanisms and identifying their drivers [5].Implicit in the use of remote sensing in LCS is the assumption that environmental parameters which define human-environment interactions can be detected in this way [6].Successful examination of environmental changes therefore depends on repetitive collection of data at varying temporal and spatial scales, [4,7].Identification of change within a system, and the type of change-seasonal, annual, step, monotonic, etc., is crucial to this process.Extended time series measured over short time increments are ideal [1], yet, despite the existence of longer term remotely sensed datasets, change detection methods are limited [5] and remain a major stumbling block in the effective use of time series approaches to remotely sensed data in research.
Continuous data products, most commonly vegetation or spectral indices, emerged in response to the inability of land cover classification analyses to explore within-class changes and the need for objective means by which emitted and reflected energy could be linked to surficial biophysical processes.Additionally, measures of biophysical properties such as vegetation cover can be more directly related to biogeochemical processes and therefore vegetation indices are increasingly used to measure change in land change analyses [8].The relationship between vegetation indices and biomass is well established [4,[9][10][11], thus indices are commonly used to examine inter-and intra-annual changes in the quantity and distribution of green biomass.The Normalized Difference Vegetation Index, or NDVI [12], based on the principal that vegetation is highly reflective in the near infrared and absorptive in the visible red [13], is widely accepted as an indicator of vegetation properties and ecological responses to environmental changes over large geographic regions [4,11,[14][15][16][17][18][19][20][21] although caution and local interpretation is also necessary given that this relationship is not always a simple, linear one.Three major advantages of continuous analyses of remotely sensed data are: (1) A closer reliance on the original data collected by the satellite that can be linked directly to land surface processes [22], so producing more realistic depictions of surface processes [12]; (2) Greater flexibility and incorporation of variation as a pixel no longer must belong to an individual class but can represent a gradient or mosaic; and (3) Greater utility as representation of initial boundary conditions [12] to many Earth system models.Time series approaches to continuous Earth Observation based variables are being adopted by many scholars [14,[23][24][25][26][27][28][29][30][31] to advance understanding of intra-and inter-annual variations in vegetation.The Advanced Very High Resolution Radiometer (AVHRR) on board the NOAA satellite series and the newer MODIS sensor onboard the Terra (EOS AM) satellite provide an excellent opportunity to employ these approaches to vegetation research because of their extensive spatial and temporal coverage.The AVHRR sensor aboard the NOAA polar orbiting satellite series (NOAA-7, 9, 11, 14, 16 and 17) was designed for weather studies, however, it is used extensively to monitor terrestrial vegetation type and condition.MODIS was designed to support land, atmospheric, and ocean imaging in a single instrument.MODIS data offer excellent land cover change products, often used in concert with AVHRR data, and yields a continuation of the AVHRR series.Therefore the large archive of AVHRR and the newer MODIS products provides a rich source of information for multi-temporal studies [32].
Some time series approaches to continuous indices have been developed, such as wavelet [33], Fourier [14,34], and PCA analyses [35,36], although some of these may transform the data under study.A suite of techniques which specifically assess and model seasonal and phenological changes within an image time-series [1,5,37] encompass those to detect changes of timing, magnitude and direction of change [5] over both seasonal and longer time periods.Washington-Allen et al. [38] argue for the concept of -benchmark conditions‖, and the interpretation of time series relative to some reference or benchmark.Such critical values may include an initial or baseline condition, change from a mean or other measure of central tendency, or some form of critical boundary condition (such as percentile, maximum or minimum, etc.).This approach provides a comparative basis for time series behaviors, such as -stable‖, -increasing‖ or -decreasing‖, relative to the selected benchmark.The approach permits an assessment of ecological resilience, or system response, [39] to some disturbance [38].However, within this framework there is still a need to determine the relevance of detected changes and ultimately their statistical significance [24].
This paper builds on the approach employed by Lanfredi and co-workers [24] to establish simple statistical tests to be applied to three separate -benchmark‖ metrics of the persistence of vegetation cover across a landscape extracted from the time series of remotely sensed data products: (1) Directional persistence addresses the direction of landscape change (increasing/decreasing NDVI) relative to a fixed benchmark condition, which could be related to a climatic event, policy change, etc.; (2) Relative directional persistence examines the direction of change of every observation relative to the preceding one; and (3) Massive persistence considers not only the direction of change from the previous observation, but also the magnitude of that change relative to the ‗control' area or boundary conditions (e.g., an area of unchanged forest cover over time, perhaps in a protected area).The latter two use the preceding value of NDVI as the -benchmark‖, investigating the sign and magnitudes of changes respectively.Together, some or all of these persistence metrics at a pixel level scale, applied to monthly time series of up to 30 years, yield a more detailed understanding of landscapes [24].Application of these techniques to remote sensing facilitates the identification of areas across the landscape in which statistically significant changes have occurred, and, within the concept of benchmark conditions [38], provides insight into the ecological resilience in the region.The potential power and utility of such techniques is diverse within the area of LCS and remote sensing studies, particularly in the context of global change studies.

Study Area
Florida displays quite diverse climate, resulting in varied ecosystems and natural resources [40][41][42] (Figure 1).Its vegetation cover is influenced by numerous factors not least of which is the 600 percent increase in population since 1940.Land use and land cover (LULC) change resulting from climate change/variability is a major concern for its direct and indirect impacts on vegetation cover.This has been studied intensively in several limited regions such as Everglades National Park [43][44][45], South Florida [42,46,47], or specific ecosystems [48][49][50].However a statewide, objective analysis of vegetation patterns and trends remains to be completed despite the fact that the natural landscape has been transformed extensively by agriculture, urbanization, and the diversion of surface water features [51].

Data Sets
AVHRR Global Inventory Modeling and Mapping Studies (GIMMS) NDVI dataset is available from July 1981 to December 2006 and has been used for numerous regional to global scale vegetation studies [30,[52][53][54][55][56][57][58].It has been corrected for a variety of potential errors and sources of noise, provides a spatial resolution of 8 km, and the data record, based on 15-day composites, facilitates the construction of cloud-free views of the Earth with the maximum NDVI at regularly spaced intervals.The Florida state boundary is subset from this intensive 25-year, monthly coverage, for the period with the larger NDVI value is chosen.
Florida land cover classification data are available from the Coastal Change Analysis Program (C-CAP) NOAA and are developed, primarily, from Landsat Thematic Mapper (TM) satellite imagery.The smallest feature size (spatial resolution) that can be mapped is 30 m pixels on the ground.This classification is available for Florida in 1996, 2001, and 2006 for 23 land use and land cover classes, which were hierarchically grouped and reduced into 8; agriculture, grassland, forest, scrub, estuarine wetlands, palustrine wetlands, barren and water.Additionally, to capture the complexity of climate in Florida (secondary winter rainfall peak diminishing towards the south), the area was subset based on six of NOAA's seven climate divisions-excluding the Keys (Figure 1).Given the absence of marked topographic features, climate tends to change slowly over the state and, for the sake of this study, is considered to be reasonably homogeneous within each division.Related research questions to be addressed include whether there are marked abrupt changes in NDVI resulting from land cover conversion, or longer term seasonal trends in particular land cover classes related to observed increases in fall and winter precipitation noted throughout the southeast [59].

Metrics of Persistence
To detect inter-annual changes of NDVI in a spatially explicit manner, three simple metrics are computed and predetermined critical levels of statistical significance applied to maps of the metrics derived over the period of record.

Directional Persistence
A measure, D, of the cumulative direction of change over the time series relative to the fixed benchmark observation of NDVI, during a sequence of observations of a particular month, j, is derived as: where V i,j is the monthly NDVI, in year i, month j, e.g., V crit is the critical benchmark NDVI value for observations made in October (j = 10).A value of +1 is assigned to t i,j when the pixel for that year and month records a value of NDVI greater than the critical observation, and −1 for observations less than.Theoretically, as NDVI is a continuous variable bounded by −1 and +1, identical values of NDVI are impossible.The benchmark will be selected by the investigator relative to the topic of research.It could be a particular value of theoretical interest, one determined by a specific empirically defined threshold, or, as in this case where persistence is of special interest, by the initial value of the time series.

Relative Directional Persistence
This is similar to the previous metric except directional comparisons are made to the monthly observation in the preceding year rather than the fixed benchmark:

Massive Persistence
The drawback of the other two measures is that they do not account for magnitude, only direction, of change, and are therefore fairly insensitive to abrupt -step‖ changes in series as opposed to monotonic changes.Massive persistence, M, incorporates changes in the magnitude of NDVI from the monthly value in the preceding year: (7)

Theoretical Development of a Time Series Persistence NDVI Analysis
Throughout this development, null hypotheses are defined by the assumption that, in the absence of any significant perturbations to the land use land cover (LULC) (e.g., climate variability/change, deforestation, urbanization, or fire) the values of NDVI are normally distributed and serially independent.That is to say, that they are independent and identically distributed normal.Theoretically, the ratio and bounded nature of NDVI contradicts this assumption and its distribution should become increasingly skewed as values approach the bounds.However, extensive empirical observations of monthly NDVI values drawn from pixels across Florida retaining a consistent land cover/use under an independent, a priori, classification over the time period of study (see for example, Figure 2; and Section 3.2 for more information on the specifics of the datasets used) indicate that, in reality, the assumption of normality is not unreasonable.Similarly, the fact that analyses are carried out on observations of NDVI for the same calendar month in successive years minimizes the problem of serial dependence.Values on the y-axis are empirical estimates of cumulative distribution calculated using the Weibul plotting position, rank/(sample size +1), where the 25 observations in the historic sample size are ranked from lowest (1) to highest (25).
It should be noted that the mean (the approximate central horizontal position of the data) and variance (their horizontal dispersion) of these observed NDVI values vary with month, geographic location and LULC.Thus, for instance, palustrine forest shows a lower mean NDVI and greater variance in winter months, than it does in summer months.Means and variances over developed land are both small.These temporal, spatial and land use changes in both mean and variance do not present a challenge in modeling the first two metrics which only consider the sign or direction of change, but they do necessitate a more complex set of considerations in dealing with the massive persistence.Derivation of the directional persistence variable matches the classic random walk process and its well-known associated statistical properties (see for example Wilson and Kirkby, 1980 [60]).A random walk process merely considers the accumulation of consecutive outcomes of a simple Bernoulli process in which the probabilities of a success (step = +1) and a failure (step = −1) are fixed and independent of the position of the random walk at that time.That is to say, the probability of a success is the same, regardless of the number, or combination, of preceding successes and failures.Under these circumstances the probability distribution of the number of successes (failures) after a specified number of trials is given by Pascal's triangle.In the case of changes relative to a fixed benchmark, the probability of success (failure) is determined by the comparative likelihood of that critical value under the null hypothesis.In the specific case of the examples employed, this critical value was selected as the first in the time series, and, under these circumstances the value of p is unknowable.Instead we resort to the notion that the value of any random variable drawn at random from the population under the conditions of the null hypothesis is the expectation, or mean, which has a probability of success/failure of 0.5.Figure 3 shows expected relative frequencies of this metric under the null conditions, p = q = 0.5, for varying time-steps (total years of record-1).Note also that when the total number of steps is even (odd), only even (odd) outcomes are possible.Cumulative frequency can be used to determine critical values of D, corresponding to the upper and lower desired significance level.As D is a discrete random variable, it is unlikely that its integer values and associated cumulative probabilities will correspond exactly to the levels of significance generally employed in statistical testing, e.g., 0.05 or 0.01.However, the smallest value of D, which equals, or is less than, the generally used levels of significance can be extracted from the frequencies of directional persistence generated under the conditions of the null hypothesis (Table 1).

Relative Directional Persistence
A success (failure) is redefined as occurring when a value of monthly NDVI is greater (less) than the one in the preceding year.Thus the probabilities of a success or failure change conditioned upon the value of the previous observation.The greater the absolute value of NDVI the preceding year, the less likely the series is to continue in that direction (increasing/decreasing). Figure 4 shows hypothetical normally distributed values of NDVI in the lower panel.The upper panel displays the probability of a success (positive step) or failure (negative step), conditioned upon the previous value of NDVI (horizontal axis).The lower (higher) the preceding NDVI value the less likely that the next value will be lower (higher), thereby restricting the range of likely scores in comparison to Pascal's triangle and the directional statistic, D. A simple Monte Carlo simulation of the hypothesized null conditions permits the computation of a sampling distribution.Figure 3 compares the frequency distributions (simulated n = 10,000) of anticipated values of R resulting after a number of transitions to those expected under conditions of the standard random walk (p = 0.5).Constriction of the range of anticipated sums of directional changes is apparent, and becomes more acute as the number of transitions increases.Note the continued correspondence between the odd/even nature of the simulated transitions and the possible outcomes.Table 2 contains critical values of R. As D and R only measure directional change, all applications, regardless of the original mean and variance of NDVI and LULC type, can be reduced to a standard normal distribution and handled in this fashion.

Massive Persistence
The fact that the massive persistence variable, M, attempts to capture some measure of the cumulative magnitude and direction of change presents a difficulty in establishing statistical significance bounds as this introduces a need to represent, or retain, the variance of the observed data under the null conditions.The expected mean should be zero after any number of transitions.Null conditions can be evaluated assuming unit variance, but, unlike the directional persistence variable, the original variance has to be reincorporated ultimately into the simulated values to establish statistical limits appropriate for the particular type of LULC, the time of year (wet season/dry season), and the particular climatic region within which the tests are to be completed.
When two independent normally distributed random variables X~N(μ x ,σ 2 x ) and Y~N(μ y ,σ y 2 ) are differenced, Z = Y − X, the resultant difference is itself normally distributed, such that μ z = μ y − μ x, and σ 2 z = σ 2 y + σ 2 x .Monthly NDVI values are differenced in successive years in this fashion.Assuming that all cases can be reduced to a standard normal distribution, the differences should be distributed: N (0 − 0, 1 + 1).In computing M, these differences incrementally constitute a partial sum of differences, i.e., The sums of two independent, normally distributed random variables, Q = R + S, are themselves normally distributed with mean μ q = μ r + μ s and variance σ 2 q = σ 2 r + σ 2 s .However, if the two random variables being summed are correlated, then the variance of the new variable becomes, σ 2 q = σ 2 r + σ 2 s + (2ρ r,s σ r + σ s ), where, ρ r,s is the correlation between the two variables being combined, and σ r and σ s are their respective standard deviations.The partial sums created in the massive persistence, are not serially independent.Numerically, the second random variable appears in the computation of both the first and second difference, Z 1 and Z 2 .Figure 3 also illustrates that if the previous difference, Z t − 1 , is negative then the next, Z t is more likely to be positive, and vice-versa.
Monte Carlo simulations (1,000,000 realizations), each of length 24 (the example time series consist of 25 observations) differenced variables, Z, are performed.They confirm that values of Z are negatively correlated with a value −0.5.Substitution of this property leads to an interesting and useful result.As both R and S are ~N (0, 2), the variance of the partial sum of two consecutive differences, σ 2 q , reduces to 2. The same process is repeated with each step of the partial sum, and its variance remains 2 regardless of the step in the summing process, or, in terms of the observed NDVI data, the length of the annual time series.Therefore, under the null hypothesis, the mean of the standardized massive persistence variable should be 0, and its variance, 2 (standard deviation, 2 ).Desired significance levels can be established easily as 2 times those shown in a regular z-table.Thus the equivalent critical (α = 0.05) standardized massive persistence value is calculated as 2 ± 1.645, or ± 2.326).
The problem remains as how to translate this critical standardized value of M back into the terms of the actual NDVI, which likely varies with season, geographic location (e.g., climate division) and LULC (see Figure 2).To this end, observations can be extracted by time of year or for each month; for the same climatic zone or region (if this varies across a study region) and from stable (unchanging) pixels e.g., those which remained in the same independent, a priori, land cover class classification during the study period, to constitute, -control data sets‖ (Figure 2).Although means and variances vary markedly, across months (m), climatic zones or divisions (d) and LULC class (l), the critical massive persistence value, M m,d,l,crit ,at some desired significance value, α, can be estimated as: where σ m,d,l is the observed standard deviation of the values of NDVI in the -control‖ pixels.

Results
Critical values of the persistence metrics provide the null circumstances against which observed values may be compared temporally and spatially, in an objective and repeatable fashion.These differences might arise from climatically induced trends in NDVI, or step functions resulting from rapid changes in LULC [39], but the methodology highlights those areas and times where and when values differ significantly from random, which are therefore worthy of further investigation.

Directional Persistence
Raw scores across this extensive area (Figure 5a) support the notion of an overall increasing trend in vegetation productivity in this quarter century, but also indicate some smaller sub-regions and months where an opposing declining trend is found, relative to the initial boundary conditions of the time-series.Areas of persistent decline are most common in late winter/early spring (February, March, April) and late summer/early fall (August and September).Most noticeable are the large negative values in northwest and north Florida during March, and August, and the negative values in the northeast Florida in September.All regions of the state record negative/declining scores relative to the boundary conditions at some point during the year, however, the overwhelming trend indicates that for the majority of the months (January, May to July, October to December) large expanses of the state experience increasing NDVI values.Statistically significant values (two-tailed test) of the metric (Figure 5b,c) highlight regions where patterns of change are most marked within the landscape.Spatial clustering of statistically significant declining NDVI evident again in the months of March, August and September, is less apparent in February and March than the raw scores (Figure 5a) indicate.Changing the significance level from 0.005 to 0.05 (Figure 5b vs. Figure 5c) further emphasizes the spatial contiguity of the increasing and declining NDVI patterns, as well as the dominance of the increasing trend throughout the majority of the months.Broadly speaking, winter precipitation (November-March) constitutes a significant proportion of annual totals (25%-35%) in the two northern climatic sub-divisions, declining southwards to less than 15% in the extreme south.Overall patterns of significant changes should be viewed within this broad context.A recent National Climate Assessment for the southeast region [61] notes a general increasing trend in winter precipitation in the southeastern United States.Increased values of D are found through much of the state, particularly northern areas during this period, however March and to a lesser extent April show significantly large negative scores, which in combination with a statewide set of positive scores in October suggest that patterns of -winter‖ precipitation may have shifted to earlier in the season in the last quarter of the century.August and September also show notable statewide patterns of significantly declining scores.These months correspond to the peak of the Atlantic -hurricane season‖, which underwent an extended period of lower than normal activity in association with several decades during which the Atlantic Multi-decadal Oscillation (AMO) was in its negative phase [62].Enfield et al. [62] showed a positive correlation (0.10 level) between annual rainfall totals in most climate divisions in Florida and the AMO.This is for annual, not seasonal, rainfall totals but does highlight that during warm phases of AMO (since about 2000) the positive association between Pacific sea surface temperatures has increased in the Southeast.ENSO impacts winter rainfalls here in Florida with warmer sea surface temperatures leading to more winter rains.At a finer spatial scale, the coastal conurbation in the southeastern of the state is generally an area of significant negative scores in most months suggesting continued urbanization.Directional persistence, based on the benchmark conditions of 1982, provides a potentially important analytical tool, especially in conjunction with the spatially explicit statistical significance of the results.

Relative Directional Persistence
Figure 6 maps values of R for NDVI throughout Florida, 1982-2006, for January (Figure 6a) and August (Figure 6b).This represents an overwhelming quantity of information to be assimilated in searching for areas and times, in and during which, environmental change may have occurred.Masking out those values which cannot be considered to be significantly different (α = 0.05) from random (maps on the right) quickly clarifies the spatial and temporal picture and allows the researcher to focus on them alone.Geographically, areas of north central and northeast Florida seem to show persistent patterns of negative relative steps (declines in NDVI).Areas in the northwest indicate opposing seasonal trends, with more negative scores in summer and more positive ones in winter.South and central Florida offer mixed patterns, which may be highly localized, but with more significant results in January than August, and the majority of those January changes, except in the extreme south, are positive.Potential biophysical explanations to be explored more fully are consistent with the patterns evinced by the Directional Persistence metric, D, include a noted increase in the percentage of annual precipitation falling in the winter, the relative dominance of cool and warm season rainfalls spatially, deforestation for commercial purposes, fires, and urban development.

Massive Persistence Changes in Palustrine Wetlands in Florida
The Massive Persistence variable (M) is specific to climatic zone, month, and land cover; however consideration of significance permits an objective comparison of the full landscape.We show results for a key land cover class in Florida, where the loss of wetlands is of considerable ecological significance [41,44].Figure 7 shows only results for pixels classified as palustrine wetlands by C-CAP 1996 (Figure 8), and for which the massive persistence scores exceeded the critical values (α = 0.05), for August-November.Many areas return significant persistent declines in NDVI during August, particularly in coastal regions.In September, the northeast and east central coastal areas of the state are dominated by declines, while much of the rest of the state reports significant increases.This trend of increasing monthly NDVI values almost completely dominates the statewide picture in October.The number of pixels indicating significantly large scores in November declines, but only one of them indicates a significant negative value.The major concentrations of significantly positive values of M are located in the east central and southeastern coastal regions.They appear to coincide with urban Remote Sens. 2014, 6 4489 fringes of heavily developed regions.The decline in summer (August) and increases in fall (October, November) NDVI in the wetlands match closely reported shifts in seasonal precipitation [59].

Discussion
This paper set out to establish statistical tests that can be applied to NDVI-derived time series of remotely sensed data products, to determine significant change as measured by three separate -benchmark‖ metrics; directional persistence, relative directional persistence, and massive persistence.Previous research [24] used only one such measure applied to only 5 years of data, and furnished no objective tests of significance.Through the use of our test site of Florida, and monthly NDVI values analyzed at a pixel level over 25 years, the tests identify areas and/or times in the landscape of significant change and facilitate a more detailed understanding of this landscape by expanding the change analysis beyond trend detection to include identification of changes that are different from what might be expected at random.
Lengths of some remote sensing records now permit the systematic quantification of vegetation change over sufficient temporal and spatial scales to begin operationalizing -state-and-transition‖ or -multiple stable state‖ models [63].Landsat data are of limited application in biologically complex systems [64], and present difficulties in obtaining series of images representative of the same phenological period from which to make objective statements about land cover change.However, monthly AVHRR-and MODIS-derived vegetation indices, available since the early 1980s, hold considerable promise for the large-scale quantification of complex vegetation-climate dynamics, and for regional analyses of landscape change as related to global environmental changes for example, [24,65].The research presented here highlights the utility of the time-series approach, with initial benchmark conditions [38] and the application of statistical significance at a pixel level to an entire study area.Considering the dramatic increase in research in global change and the evaluation of conditions based on some initial or boundary conditions, this type of methodology augurs well for further development and application of spatially explicit statistical testing to the intersection of remote sensing and LCS.
A limitation of AVHRR (and MODIS) products is their spatial resolution; especially those indices derived over 1-8 km pixel sizes.This coarse resolution is offset by a finely resolved (MODIS is every 8 days and AVHRR every 14 days, and most are extracted as monthly maximum values) continuous time series of up to 30 years in length, ideal for studies addressing vegetated surfaces, where phenological and seasonal impacts need to be understood and accounted for within longer term changes.The ability to unambiguously discern the location and timing of significant changes (both physical and statistical) is of paramount importance in developing an understanding of landscape-level processes.Unfortunately not all types of landscape changes can be tracked by means of coarse scale remote sensing data alone, and this does have to be considered when analyzing and interpreting the results.However, in this analysis, the incorporation of higher resolution static date land cover change analyses (Figure 8), in concert with monthly AVHRR NDVI coarse spatial scale (but high temporal fidelity) imagery was an ideal combination.In the study region a great deal of pixel level change was evaluated, but consideration of only statistically significant changes focused attention on landscape changes of import and requiring further study, management interventions or understanding.Such an objective and well-defined set of analytical tools allow managers, policy-makers, conservationists and institutions to better evaluate the changes and to discern their importance.
The proposed metrics are easily calculated and fairly easy to test, without recourse to particularly high-powered statistical methods like wavelets.In this regard they are somewhat reminiscent of the three methodologies proposed by De Beurs and Henebry [25], but are free of some of those restrictions.These authors proposed three methods to detect what they called discontinuities and trends in similar time series of NDVI.In the first of these, two statistical tests, the C-method and Fligner-Policello, are proposed to investigate potential step changes in following some a priori defined event, in the case of the paper the fall of the former Soviet Union.No such judgment is required for our methods, nor is there a requirement of -large‖ sample sizes.The directional persistence, D, could be used to measure changes relative to some fixed event in time.In the example shown in this paper the initial year of record was used as such a -benchmark‖, but it could equally easily be the observed values of NDVI at the times of other major events in the landscape, such as devastation from a hurricane, extended drought, insect infestation or fungal disease, all of which might constitute a priori -benchmarks‖.
De Beurs and Henebry [25] recommend the use of the seasonal Mann-Kendall test to detect trends.In essence the test is very similar to both the D and R statistics proposed in this paper as it constitutes some simple integer scaling according to comparative sizes of seasonally constant (10 day periods in the former paper, monthly values in this) values of NDVI.Determination of the significance level of the Mann-Kendall statistic becomes more complex as it requires the variances of each seasonal subset and the covariances of every combination of these subsets for each pixel analyzed.The simple tabulated values proposed in this paper make no such demands on the user.
The final method proposed by De Beurs and Henebry [25] involves the computation of quadratic relationships between NDVI and Accumulated Growing Degree Days.In order to compute the latter the authors assume the base for growing degree day calculations to be 0 °C for wheat, and require maximum and minimum temperatures over the period of record.The approach poses difficulties when mixed land uses/land covers are under consideration and requires arbitrary decisions by the user as to an appropriate base.Although maximum and minimum temperatures for the state of Florida could be derived at similar spatial resolutions from gridded data sets, similar, reliable sets are often not available in other regions of the world where this ability to examine large swathes of a fairly inaccessible landscape are most valuable.
The use of a readily accessible and commonly used measure for land cover changes, here NDVI, illustrates that the state of Florida when considered at the landscape scale has a dominant and statistically significant pattern of increasing NDVI.The spatially explicit nature of the analysis allows the identification of significant patterns of land cover change, indicating that there are regions with differing long term patterns of change in NDVI, the most apparent of which occurs in the northwest of the state.Additionally, the temporal extent of the dataset combined with the short time increments of observations (monthly) support a statewide analysis which highlights the temporal patterns of long term change in NDVI.This enables the identification of months with NDVI patterns that differ from the remainder of the year (e.g., March, August, and September).When considered in combination with static land cover classifications, such as the C-CAP, the methodology employed here enables the examination of the relationship between the long term spatial and temporal dynamics of NDVI and land cover classes of relevance/importance to the study site.For example, in Florida the palustrine wetlands, which are ecologically and socially important, and the patterns of change in NDVI relative to baseline observations can be further examined through the isolation of this land cover class.Similarly, the relationship between developed land along the coast (particularly in southeastern Florida) and long term NDVI trends can be isolated and assessed revealing information regarding not only the effect of development on NDVI but also, potentially, about changes in development pattern and landscape composition within this land use category as indicated by changes in NDVI over time.

Conclusions
The first two measures can be used in any location in which the assumption of approximately normally distributed values of NDVI is not severely violated.The use of the massive persistence is a little more restrictive as it requires secondary information from an independent classification to identify the statistical properties of -control‖ pixels that have remained unchanged in their LULC.These metrics and their accompanying significance levels provide a very valuable and broadly applicable means of objectively analyzing very large data sets and extracting the major patterns of spatial and temporal change.
Although the principle objective of this study is to illustrate the development of some relatively simply computed and statistically tested measures of persistence with reference to the state of Florida, the results also indicate the types of multi-scalar investigations that might be established as a result of this sort of analysis.At the largest space and time scales studied here, there are consistent indications of changes in NDVI which may well have their basis in independently observed long-term changes in the seasonality of precipitation.On smaller spatial scales it is possible that the consistent patterns of no change or negative changes may well be associated with localized drivers of LULC change, such as continued urban growth in the southeast conurbation, the Tampa-Orlando -I4-corridor‖ and the Jacksonville area.Similarly the temporal dimensions of change could be investigated at smaller scales, thereby investigating during which periods the various metrics became significant and in which regions.A signal of this type might be more apparent in northern parts of the state in association with the various phases of ENSO to which winter precipitation is known to be sensitive.
Through the use of such measures of NDVI directional persistence, relative directional persistence and massive persistence, we can now better evaluate and understand landscape level changes in vegetation amounts.This improves the characterization of the landscape upon which examinations of the change mechanisms and drivers are based.Rather than showing all changes occurring within a region, these methods allow us to identify and focus on just those areas of significant changes.This is of unparalleled use to managers and researchers alike and is a powerful tool within the remote sensing arsenal.Application of these newly developed statistical techniques to remote sensing facilitates the identification of areas across the landscape in which statistically significant changes have occurred.Evaluated within the concept of benchmark conditions [38], these significant changes provide a greater understanding of ecological resilience in the region.Benchmarks and associated trend-detection permit an assessment of ecological resilience, or system response, to some disturbances [38].However, within this framework, there is still a need to determine the relevance of trends and ultimately their statistical significance, as is developed here in this research.The potential power and utility of such techniques is diverse within the area of physical geography, especially in the context of global change studies.

Figure 1 .
Figure 1.Study area map of the state of Florida with the climate divisions labeled in bold.

Figure 2 .
Figure 2. Observed monthly distributions of NDVI values of pixels in (a) climatic division four, classified as Scrubland by C-CAP in 1996, 2002 and 2006; (b) division five, Agricultural Land; (c) division six, Developed Land; (d) division three, Palustrine Wetland.

4 N 6 N 3 N
S c r u b la n d s -D iv .u la ti v e D i s tr i b u ti o u la ti v e D i s tr i b u ti o u la ti v e D i s tr i b u u la ti v e D i s tr i b u ti o e lo p e d L a n d -D iv .u la ti v e D i s tr i b u ti o u la ti v e D i s tr i b u ti o s tr in e W e tla n d s -D iv .u la ti ve D i s tr i b u ti o u la ti ve D i s tr i b u ti o

Figure 3 .
Figure 3.Comparison of Monte Carlo simulations of the distribution of the relative directional persistence variable, R, following various numbers of transitions (a) 5 transitions, (b) 10 transitions, (c) 15 transitions, (d) 20 transitions, (e) 25 transitions, and (f) 28 transitions and the anticipated results from a random walk (p = 0.5) process.

Figure 4 .
Figure 4. Hypothetical normally distributed NDVI values (lower panel) and the changing probabilities of success and failure in the following step depending upon the value of the current value of NDVI, as used in defining the relative persistence variable, R.

Figure 5 .Figure 5
Figure 5. Monthly results from the application of the spatially explicit persistence analysis of NDVI from a 25-year time series for (a) Directional persistence relative to 1982, (Green represents increases in the variable, yellow areas of little to no change and red decreases); (b) regions of statistical significance at 0.005 significance level relating to >±14 years; and (c) those significant at 0.05 level relating to >±10 years

Figure 6 .
Figure 6.The NDVI relative directional persistence scores for the state of Florida for (a) January and (b) August showing all scores and only the statistically significant ones for January (c) and August (d).

Figure 7 .
Figure 7.The spatial distribution of statistically significant massive persistence scores across the state of Florida, for areas of palustrine wetland for (a) August, (b) September, (c) October and (d) November.

Figure 8 .
Figure 8. Land cover change across Florida as determined from the C-CAP analysis for 1996-2001-2006 only, highlighting some dominant land cover conversions across the state.

Table 1 .
Critical values of the directional persistence, D, derived from a simple Random Walk (p = 0.5) process.n is the number of observed transitions, or total record length minus 1.

Table 2 .
Critical values of the relative directional persistence statistic, R. The number of observed transitions is n, equal to total record length minus 1.