Distance Metric-based Forest Cover Change Detection using MODIS ...

Distance Metric-based Forest Cover Change Detection using MODIS Time Series Xiaoman Huang∗, Mark A. Friedl Department of Earth and Environment, Boston University, 675 Commonwealth Avenue, Boston, MA 02215, USA

Abstract More than 12 years of global observations are now available from NASA’s Moderate Resolution Imaging Spectroradiometer (MODIS). As this time series grows, the MODIS archive provides new opportunities for identification and characterization of land cover at regional to global spatial scales and interannual to decadal temporal scales. In particular, the high temporal frequency of MODIS provides a rich basis for monitoring land cover dynamics. At the same time, the relatively coarse spatial resolution of MODIS (250-500 m) presents significant challenges for land cover change studies. In this paper we present a distance metric-based change detection method for identifying changed pixels at annual time steps using 500 m MODIS time series data. The approach we describe uses distance metrics to measure 1) the similarity between a pixel’s annual time series to annual time series for pixels of the same land cover class, and 2) the similarity between annual time series from different years at the same pixel. Pre-processing, including gap-filling, smoothing and temporal subsetting of MODIS 500 m Nadir BRDF-adjusted Reflectance (NBAR) time series is essential to the success of our method. We evaluated our approach using three case studies. We first explored the ability of our method to detect change in temperate and boreal forest training sites in North America and Eurasia. We applied our method to map regional forest change in the Pacific Northwest region of the United States, and in tropical forests of the Xingu River Basin in Mato Grosso, Brazil. Results from these case studies show that the ∗

Corresponding author. Tel.: +1 617 233 0491. Email address: [email protected] (Xiaoman Huang)

Published in International Journal of Applied Earth Observation and Geoinformation (2014), Vol.29, pp. 78–92. http://dx.doi.org/10.1016/j.jag.2014.01.004

method successfully identified pixels affected by logging and fire disturbance in temperate and boreal forest sites. Change detection results in the Pacific Northwest compared well with a Landsat-based disturbance map, yielding a producer’s accuracy of 85%. Assessment of change detection results for the Xingu River Basin demonstrated that detection accuracy improves as the fraction of deforestation within a MODIS pixel increases, but that relatively small changes in forest cover were still detectable from MODIS. Annually, over 80% of pixels with >20% deforested area were correctly identified and the timing of change showed good agreement with reference data. Errors of commission were largely associated with pixels located at the edges of disturbance events and inadequate characterization of land cover changes unrelated to deforestation in the reference data. Although our case studies focused on forests, this method is not specific to detection of forest cover change and has the potential to be applied to other types of land cover change including urban and agricultural expansion and intensification. Keywords: MODIS, Change detection, Distance metrics, Time series, Forest disturbance

1

1. Introduction

2

Optical remote sensing is an important source of information for characterizing the state

3

and dynamics of the vegetated land surface. Coarse spatial resolution sensors with high

4

revisit capability, such as the Moderate Resolution Imaging Spectroradiometer (MODIS),

5

provide synoptic views of regional and global changes in vegetation cover that complement

6

lower temporal frequency, higher spatial resolution sensors such as the Landsat Thematic

7

Mapper (TM) and Enhanced Thematic Mapper Plus (ETM+). Such comprehensive charac-

8

terization, including the location, nature, and intensity of vegetation change, is essential to

9

process studies focused on human modification of terrestrial ecosystems, carbon budgets, and

10

biosphere-atmosphere interactions (Dickinson and Kennedy, 1992; Vitousek, 1994; Copeland

11

et al., 1996; Pielke et al., 2002; Achard et al., 2004; Foley et al., 2005; Baccini et al., 2012).

12

Over the past 4 decades, observations from sensors onboard the Landsat series of satellites

13

have been widely used to map land cover changes at local scales related to deforestation,

14

forest mortality, urbanization, and wetland loss (e.g., Skole and Tucker, 1993; Collins and

15

Woodcock, 1994; Jensen et al., 1995; Ji et al., 2001). To monitor land cover changes at 2

16

regional and global scales, early studies utilized coarse spatial resolution data from the

17

Advanced Very High Resolution Radiometer (AVHRR) (e.g., Ehrlich et al., 1994; Lambin,

18

1997; Lambin et al., 2003; Mucher et al., 2000; Young and Wang, 2001; Lepers et al., 2005).

19

Since 2000, the availability of MODIS data with vastly superior spatial, spectral, geometric,

20

and radiometric qualities relative to AVHRR, has provided an improved basis for regional

21

and global land cover change mapping.

22

Previous efforts focusing on mapping large scale land cover change with MODIS data have

23

used two main approaches: hybrid methods that use Landsat type data with MODIS and

24

wall-to-wall mapping. In hybrid studies (e.g., Hansen et al., 2008a, 2009, 2010; Hayes et al.,

25

2008; Potapov et al., 2008; Broich et al., 2011), land cover change such as net forest cover loss

26

over a period of ∼5 years is estimated using regression models calibrated with training data

27

from higher spatial resolution sensors (e.g., Landsat TM/ETM+). This approach requires

28

a careful sample design and reliable interpretation of Landsat data to generate high quality

29

training data. Wall-to-wall mapping, on the other hand, is designed to detect change at

30

the pixel level using either change indices (e.g., Linderman et al., 2005; Mildrexler et al.,

31

2009; Coops et al., 2009), temporal trajectory-based change detection algorithms (e.g., Sulla-

32

Menashe et al., 2013), or classification (e.g., Carroll et al., 2011). In contrast to sampling-

33

based approaches, these methods provide binary change information at native MODIS spatial

34

resolution, but often do not fully utilize the seasonal information in the time series.

35

Despite these efforts, improved methods to monitor regional to global land cover change at

36

moderate spatial resolution are urgently needed. The majority of change detection methods

37

use two-date or multi-date methods that were originally developed using Landsat data (e.g.,

38

Angelici and Bryant, 1977; Singh, 1989; Jensen et al., 1995; Collins and Woodcock, 1996;

39

Hayes and Sader, 2001; Kennedy et al., 2007; Masek et al., 2008; Huang et al., 2010). An

40

alternative strategy, which we pursue here, is to develop methods that better utilize temporal

41

information in MODIS data. In particular, better use of multi-temporal information has the

42

potential to provide more efficient and accurate identification of change locations.

43

Recognizing this potential, a variety of studies have begun to emphasize the use of high

44

temporal frequency information in MODIS data (e.g., Verbesselt et al., 2010; Lhermitte

45

et al., 2011). However, such methods frequently require gap-free input data, and also tend 3

46

to be sensitive to noise that can lead to spurious detection of change. This poses substantial

47

challenges in many areas of the world, since the availability and quality of MODIS data is

48

affected by numerous factors including clouds, aerosols, the presence of snow, large view

49

zenith angles, and large solar zenith angles at high latitudes (Tarpley et al., 1984; Holben,

50

1986; Moody and Strahler, 1994). In the tropics, persistent cloud cover creates prolonged

51

data gaps, and biomass burning during the dry season leads to aerosol contamination that

52

can significantly reduce data quality. To address this, compositing techniques are widely used

53

to reduce the amount of missing data by aggregating data over longer time periods (Holben,

54

1986; Viovy et al., 1992; Huete et al., 2002). Other studies have used Fourier analysis,

55

least-squares distribution fitting, or moving window filters to gap-fill and smooth the time

56

series (Sellers et al., 1994; Moody and Johnson, 2001; Jonsson and Eklundh, 2002, 2004; Chen

57

et al., 2004; Beck et al., 2006; Hocke and Kaempfer, 2009). More recently, hybrid approaches

58

that gap-fill missing data using both temporal trajectories and spatial information based on

59

land surface types have also been developed (Gao et al., 2008; Fang et al., 2008). The results

60

from these studies reinforce the need for pre-processing to reduce gaps and noise in the

61

MODIS time series.

62

With these issues in mind, the main objective of this paper is to describe and assess

63

a new method for continuous monitoring of forest change from 500 m MODIS data. Our

64

approach includes 2 key elements: 1) pre-processing of 500 m MODIS time series data to

65

gap-fill and reduce noise, and 2) a distance metric-based change detection method designed

66

to detect annual forest change. To evaluate our approach we present three case studies for

67

the period of 2003-2010. First, we applied our method to MODIS time series for temperate

68

and boreal forest training sites used to create the MODIS collection 5 land cover type prod-

69

uct (MCD12Q1, Friedl et al., 2010). Second, we tested our method for coniferous forests in

70

the Pacific Northwest region of the contiguous United States. Third, we applied our method

71

to detect and map forest cover change in the Xingu River Basin located in Mato Grosso,

72

Brazil.

4

73

2. Study Areas and Data

74

2.1. Temperate and Boreal Forest Training Sites

75

In the first case study, we use our change detection method (Section 3) to identify changes

76

occurring between 2003-2010 in a set of temperate and boreal forest sites. These sites were

77

selected from the training site database used to produce the MODIS land cover type prod-

78

uct (MCD12Q1), which is known as the System for Terrestrial Ecosystem Parameteriza-

79

tion (STEP) (Muchoney et al., 1999; Friedl et al., 2002, 2010). This database is designed

80

to capture the geographic variability of global land surface types using the International

81

Geosphere-Biosphere Programme (IGBP) land cover classification scheme. Each site in the

82

database consists of a polygon of homogeneous land cover delineated from Landsat or higher

83

spatial resolution imagery and is assumed to be stable over time. To ensure the accuracy of

84

land cover representation provided by STEP, each site needs to be monitored for changes.

85

For this work we used 269 STEP sites located in temperate and boreal forests of North

86

America and Eurasia by selecting evergreen needleleaf forest (ENF), deciduous needleleaf

87

forest (DNF), deciduous broadleaf forest (DBF), and mixed forest (MXF) sites located within

88

the biogeographic realms ’Nearctic’ and ’Palearctic’ defined in Olson et al. (2001). The size

89

of sites ranged from 2 MODIS 500 m pixels (∼0.4 km2 ) to 66 pixels (∼14.2 km2 ), with a

90

median size of 12 pixels (∼2.6 km2 ). We further defined subgroups within the IGBP forest

91

classes based on Olson biogeographic realms and biomes, and excluded subgroups that were

92

under-represented in STEP database ( 1250 m according to the NED.

132

The final mask included an area of ∼93000 km2 . Disturbance information from a study

133

by Kennedy et al. (2012) based on Landsat from 1984-2008 were aggregated to MODIS 500

134

m spatial resolution and used as reference data.

135

2.3. Tropical Forests of the Xingu River Basin

136

In the third case study, we applied our change detection method to map changes in the

137

Xingu River Basin (Mato Grosso, Brazil) for an area covered by 9 Landsat scenes (WRS 2

138

path 224-226, row 67-69, ∼260000 km2 ). The study area has a May-September dry season

139

and an October-April wet season, with several distinct types of natural vegetation including

140

moist tropical rainforest, cerrado, and deciduous forest (Furley et al., 1988). The study area

141

also includes substantial indigenous land and the headwaters of the Xingu River. As a result

142

of ongoing efforts to balance conservation with economic incentives for soybean production

143

and cattle ranching, this area has become a focus for many land cover mapping and change

144

studies (Morton et al., 2005; Hansen et al., 2008b; Arvor et al., 2010; Walker et al., 2010).

145

Here we focused on changes in evergreen broadleaf forest (EBF) using STEP EBF sites

146

in South America located in the Olson ’Tropical and Subtropical Moist Broadleaf Forests’

147

biome between 0 and 20◦ S. MODIS data used for our analysis included NBAR band 7 and

148

QA time series for the Xingu River Basin and the selected STEP EBF sites for 2001-2011.

149

To assess our results, we used the PRODES (Monitoring the Brazilian Amazon Gross De-

150

forestation) dataset produced by Brazil’s National Institute for Space Research (INPE) (INPE,

151

2012). This annual product provides polygon-based maps of deforestation delineated from

152

Landsat data using a combination of methods including linear spectral unmixing, image seg-

153

mentation, unsupervised classification, and manual interpretation (Valeriano et al., 2004).

154

This dataset was used for two purposes. First, we derived annual sub-pixel fractions of de-

155

forestation for each MODIS pixel by resampling the 30 m data to 500 m spatial resolution

156

corresponding to MODIS pixels (Figure 1). This allowed us to assess the detectability and

157

timing of change with respect to the proportion of sub-pixel deforestation. The PRODES

158

information was also used to generate a mask of evergreen forests in the study area, which

7

159

we define as 500 m pixels with ≥60% forest cover in 2003 according to PRODES. Second,

160

we derived binary change maps at 500 m to compare with change results from MODIS. The

161

identified year of change at 500 m from PRODES depends on the percentage deforestation

162

within each pixel. We tested two thresholds in our analysis (0% and 20%), which correspond

163

to the year when deforestation in any 500 m pixel exceeds 0% (any change), or 20% (small

164

change). [Figure 1 about here.]

165

166

3. Methods

167

3.1. Pre-processing

168

To pre-process NBAR time series extracted for STEP sites, we used a combined spatial

169

and temporal gap-filling approach, utilizing knowledge of land cover type for the sites. For

170

each missing value, we first searched for ’candidate’ fill values on the same date from nearby

171

pixels of the same class. If suitable pixels were not found, the search window was expanded

172

until suitable fill values were found. Specifically, we successively searched within the training

173

site, within the MODIS tile, in a 3×3 tile window, in a 3-tile wide swath restricted to the

174

same Olson biogeographic realm, and finally, within all training site pixels of the same class

175

globally. The median of high quality candidate values was then used to fill the data gap.

176

To assess the quality of the resulting gap-filled values, we evaluated their representative-

177

ness based on current and previous time series at the pixel. Fill values were determined

178

to be unrealistic if they differ substantially from observations of the same time from past

179

years, and were replaced using interpolation. An additional ’despiking’ procedure was then

180

applied to the time series to eliminate sudden increases or decreases in all values at each

181

pixel. Finally, a 3-point median filter was applied to all time series values to further reduce

182

noise in the data.

183

Although the STEP forest sites encompass a wide range of geographic and climate

184

regimes, the most prominent periods with persistent data gaps occur in the Northern Hemi-

185

sphere winter. During these periods, missing values were often filled initially with candidate

186

values from a very large spatial window (e.g., 3×3 tile window) with lower quality. To reduce 8

187

noise from clouds and large solar zenith angles in the Northern Hemisphere winter we used

188

only data from May through September (19 consecutive NBAR observations).

189

For NBAR time series in the NWFP area and Xingu River Basin, where unlike STEP,

190

we do not necessarily know the land cover at each pixel, we employed a gap-filling strategy

191

that only uses temporal information. To do this, we first extracted subsets for May through

192

September, which correspond to the snow/cloud-free season in each region. Missing values

193

were replaced with linear interpolation if good data were available within a 5-point time

194

window. Otherwise, missing values were filled using a quadratic polynomial fitted to the

195

nearest available values. The same despiking and median filter smoothing that we used for

196

STEP data was also applied following gap-filling.

197

3.2. Change Detection

198

Distance Metrics

199

We combined two distance metrics calculated from pre-processed, gap-free MODIS time

200

series to identify pixels that were 1) distinct from their assumed land cover class in each year,

201

and 2) showed change in their annual time series across years. Both distances make explicit

202

use of temporal information at each pixel. For a given year, we define the within-class dis-

203

tance for a pixel assumed to belong to class c using the Mahalanobis distance (Mahalanobis,

204

1936): 2 ¯ r,b,c )0 Σ−1 (Xi − X ¯ r,b,c ) Dmah = (Xi − X

(2)

205

where Xi is an s×1 vector of annual MODIS observations at pixel i in realm r and biome

206

¯ r,b,c is the mean seasonal curve of the population (i.e., a ’baseline’ characterization of b, X

207

the annual time series for a stable class c pixel in the same realm, biome, and year), and Σ

208

is the s×s covariance matrix of the population, based on a sample for the class of interest.

209

This distance measures the dissimilarity between pixel i and pixels belonging to the same

210

class that are located in the same realm and biome.

211

212

We define the between-year distance as the magnitude of a change vector between two annual time series from the same pixel using the Euclidean distance: 2 Deuc = (Xyr+,i − Xyr−,i )0 (Xyr+,i − Xyr−,i )

9

(3)

213

¯ yr−,i and X ¯ yr+,i are the annual time series at pixel i preceding and following the where X

214

year of interest. Thus the between-year distance measures the dissimilarity between a pixel’s

215

temporal profile across time.

216

Our change detection approach is different from previous land cover change studies using

217

distance metrics and change vector analysis (e.g., Malila, 1980; Lambin and Strahler, 1994;

218

Johnson and Kasischke, 1998; Ridd and Liu, 1998; Zhan et al., 1998; Bontemps et al., 2008;

219

Bucha and Stibig, 2008; Xian et al., 2009) in three main aspects. First, annual observations

220

for a pixel, Xi , are not restricted to a single band, spectral feature, or a set of features

221

at a single time. Rather, it is composed of annual time series of multiple bands/features

222

selected for the change detection. Second, using land cover information, we formulate the

223

2 within-class distance Dmah to characterize how different a pixel is from the class population

224

for each year. Third, by explicitly combining within-class and between-year distances, we

225

reduce the probability of false detection by requiring each pixel’s spectral-temporal behavior

226

to be both distinct from pixels of the same land cover class and to exhibit changes in time.

227

Figure 2 illustrates how the within-class and between-year distances are calculated for

228

an EBF pixel located in the Xingu River Basin. Visual inspection of high resolution images

229

for this pixel shows that it was converted to agriculture in 2007 (Figure 2(a)). For the

230

within-class distance, the pixel’s annual time series is compared with the EBF population

231

mean scaled by the covariance matrix of the population (Figure 2(b)). For the between-

232

year distance, we compare the pixel’s temporal profile one year before and after the year of

233

interest (Figure 2(c)). [Figure 2 about here.]

234

235

Distance Thresholding

236

Thresholds to identify changes based on within-class and between-year distances were

237

estimated for each forest class-realm-biome combination by pooling distances calculated for

238

multiple years (2003-2010). For the temperate and boreal forest case study, we estimated

239

empirical distributions of within-class and between-year distances for all subgroups of ENF,

240

DNF, DBF, and MXF classes. For the NWFP, we estimated the distance distribution for

241

a subgroup of the ENF class in the Pacific Northwest. For the STEP sites, we identified 10

242

changed pixels using thresholds based on 90% and 95% quantiles of empirical distributions

243

for both within-class and between-year distances. We note that the selection of a specific

244

threshold involves a trade-off. By applying a less conservative threshold, we lower the errors

245

of omission but increase commission errors (and vice versa).

246

For the Xingu River Basin, we estimated the distance distribution using a subgroup of

247

the EBF class in the Amazon (Figure 3). A range of distance thresholds was tested in

248

order to assess the balance between errors of omission and commission. Since the within-

249

class distance is formulated as a Mahalanobis distance, which assumes that the data is

250

multivariate normal, the distribution of the within-class distance should closely follow a χ2

251

distribution. We tested how selection of threshold affected omission and commission errors

252

by varying the within-class and between-year distance thresholds. We tested within-class

253

distance thresholds between 90% quantile of the empirical distribution and χ2p=0.75,df =19 =

254

14.562, and between-year distance thresholds between 90% quantile and the maximum value

255

of the empirical distribution.

256

[Figure 3 about here.]

257

With the distance thresholds established, change detection was performed on a per pixel

258

basis, in all case studies. For each year, pixels with both within-class and between-year

259

distances exceeding prescribed thresholds were flagged as ’potential change’. To minimize

260

error of omission caused by transient conditions or noise, changed pixels were required to

261

exceed the selected thresholds in two consecutive years, and the timing of change was defined

262

as the first year when the pixel was flagged.

263

3.3. Assessment of Change Detection Results

264

Assessment of change detection results from remote sensing relies on independent refer-

265

ence datasets. For the STEP temperate and boreal forest sites, we did not have reference

266

information related to change. Our assessment was therefore focused on commission errors

267

by visually examining each forest pixel detected as changed by our method. To do this, we

268

combined manually interpreted information from MODIS time series, Landsat TM/ETM+

269

image stacks, and Google EarthTM high spatial resolution images, where available. For the 11

270

NWFP area, we used reference information from a Landsat-based disturbance map for 1985-

271

2008 (Kennedy et al., 2012; Sulla-Menashe et al., 2013). We assessed errors of omission and

272

commission at MODIS pixel scale using repeated random samples. To do this, we drew

273

multiple (1000) random samples with change and no change (n=500) within the ENF mask

274

area. Due to differences in the time period associated with the MODIS and Landsat NWFP

275

maps, we excluded areas with change after year 2008 according to our map, as well as areas

276

with changes prior to 2002 according to the Landsat NWFP map. The random samples were

277

drawn from about 90% of the ENF mask areas. In addition, we assessed overall agreement

278

between our change map and Landsat-derived change map by aggregating them to ∼10 km

279

spatial scale.

280

For the Xingu case study, we compared change detection results from MODIS with the

281

PRODES deforestation map. For each year, we examined the proportion of change pixels that

282

were identified as changed, using 10% deforestation intervals between 0 and 100%. We then

283

assessed the accuracy of the final change map, focusing on new deforestation that occured

284

between 2003 and 2010, by looking at pixels that were 100% forest in 2002. We assessed both

285

spatial (location of change) and temporal (timing of change) agreement between our MODIS

286

results and the PRODES change maps and calculated the producer’s and user’s accuracy

287

for spatial agreement. For pixels that were detected as changed from MODIS which agreed

288

with data from PRODES, we calculated the temporal accuracy. Specifically, if the year of

289

change determined from MODIS was within ±1 year from PRODES, we considered it to be

290

in agreement.

291

While we used the PRODES dataset as the benchmark for assessing the performance of

292

our method for Xingu River Basin, it is important to note a few limitations of this dataset

293

that may affect our comparison. First, the types of change characterized in the PRODES

294

dataset are restricted to forest clear cuts. Other types of change such as burning, selective

295

logging, and degradation in forest fragments are not included (Souza et al., 2005). Also the

296

minimum mapping unit for PRODES (6.25 ha) may have minor effects on the estimated

297

proportion of deforestation for some pixels. Finally, since the year of change in PRODES is

298

defined as the first year it is detected, the presence of clouds and other sources of noise and

299

missing data can cause pixels to be flagged as changed later in PRODES than the actual 12

300

year of deforestation.

301

4. Results

302

4.1. Pre-processing

303

To assess the effectiveness our gap-filling approach for the STEP time series data, we

304

used a sample of 100 forest sites to simulate missing data. For each site, data for 5 dates

305

were removed at random for each year between 2002-2011, where the random selections were

306

weighted by the proportion of missing values in NBAR data for each date and land cover

307

class of interest within the corresponding Olson biogeographic realm. Thus, dates with higher

308

proportions of missing data were more likely to be selected, and vice versa. By removing

309

values of select dates for entire sites, we forced the candidate fill values to be selected from

310

other sites.

311

For each year, we evaluated the gap-filling strategy using the root-mean-square error

312

(RMSE) of gap-filled values relative to the actual observations. For comparison, the natural

313

variability within each class was quantified using the standard deviation of un-filled data.

314

Figure 4 presents barplots showing the standard deviation within each land cover class, the

315

RMSE with just gap-filling, and the RMSE with gap-filling and filtering. The error bars

316

indicate one standard deviation for Stdevc and RM SEc over 2002-2011. These results show

317

that the pre-processing strategy we used introduced errors that were smaller than the natural

318

variablity within each land cover class. This is observed for all forest classes, especially the

319

DBF and MXF classes. We therefore conclude that our gap-filling procedure provides a

320

good basis for imputing missing values and does not introduce additional noise that would

321

negatively affect our change detection results.

322

323

[Figure 4 about here.] 4.2. Change Detection in Temperate and Boreal Forest Training Sites

324

Using the distance thresholding approach described in Section 3.2, a total of 16 temperate

325

and boreal forest sites were identified as containing pixels affected by change events between

326

2003-2010 that were also confirmed by manual intepretation (Table 1). Change sites were

327

found in all 4 IGBP forest classes included in our analysis, with the majority of change 13

328

occurring in temperate coniferous forests. Specifically the ENF class had the highest number

329

of detected change sites (7), followed by the MXF class (4). [Table 1 about here.]

330

331

Detected changes were associated with two major types of disturbance: fire and logging.

332

Fire disturbance was detected in needleleaf forest sites including an ENF site in Alaska

333

and two DNF sites in Eastern Siberia. Logging affected all classes with the exception of

334

DNF, with change sites found across the US, Canada, and Russia. Figure 5 shows MODIS

335

EVI2 and band 7 time series for selected change pixels affected by change events. Both

336

EVI2 and band 7 exhibit signatures related to disturbance events that are complementary

337

to each other. For example, for a DNF pixel disturbed by fire in 2003 (Figure 5(b)), EVI2

338

values recovered to pre-disturbance levels over the course of 7 years, but peak reflectance in

339

band 7 remained higher than those for stable forest, suggesting development of ground-cover

340

vegetation post-disturbance while the fraction of exposed soil remained higher than the pre-

341

disturbance state. In some cases, changes in band 7 are more evident than in EVI2, since

342

even small fractions of soil visible in the field of view can greatly increase MODIS band 7

343

reflectance (Figure 5(a), 5(c)). [Figure 5 about here.]

344

345

In addition to the successful detection of changes described above, our method also

346

identified some spurious changes, which were mainly caused by high levels of noise in MODIS

347

time series. This was especially common for pixels with very high proportions of missing

348

data (e.g., >75%) due to persistent clouds or snow cover. For these pixels, high frequencies

349

of missing and low quality data reduced our ability to effectively filter and smooth MODIS

350

time series.

351

4.3. Change Detection in Coniferous Forests of the NWFP Area

352

In the NWFP area, our change detection results compare well with the Landsat-derived

353

disturbance map. In particular, major fire events during the study period, such as the 2002

354

Biscuit Fire in Siskiyou Mountains in southwest Oregon (Figure 6a), are evident. Table 2 14

355

shows the confusion matrix averaged across 1000 random samples at MODIS pixel scale.

356

The producer’s accuracy for the change class is 96% indicating low levels of omission errors,

357

but the user’s accuracy (62%) indicates considerable amount of commission errors compared

358

to the Landsat change map.

359

[Figure 6 about here.]

360

[Table 2 about here.]

361

Figures 6b and 6c show the locations of commission and omission errors in the MODIS

362

change map. One source of commission errors is related to high elevation areas in the North-

363

ern Cascades and Olympic Mountains, where our training sites do not adequately represent

364

the ENF vegetation signature. In addition, visual examination showed some disturbance

365

events that are not documented in the Landsat change map in the Western Cascades and

366

Coastal Range as another source of commission. We drew additional repeated random sam-

367

ples from MODIS tile h09v04 for bins at 10% cumulative proportions of sub-pixel disturbance

368

to assess the nature of omission errors. As expected, pixels containing small proportions of

369

disturbance in the NWFP area are challenging to detect, with levels of omission errors de-

370

creasing fastest between 20-70% of sub-pixel disturbance.

371

Given the heterogeneous landscape and complex disturbance history of this area, we

372

also assessed the MODIS change map (produced using 95% distance thresholds) in spatially

373

aggregated form. Figure 7 shows a comparison between Landsat-derived disturbance and

374

MODIS-derived disturbance at ∼10 km spatial resolution, for all 10×10 km blocks that

375

contain ENF forests according to our forest mask. The proportions summarized here reflect

376

the amount of disturbed pixels identified between 2002-2008 in each block. During the 2002-

377

2008 period, most of the blocks in the area contain 20% change threshold (Table 4, Figure 9). However, the user’s accuracy for this scenario is

440

lower, because pixels with small sub-pixel deforestation detected from MODIS are considered

441

commission errors.

442

5. Discussion

443

Pre-processing

444

The change detection method we describe in this paper relies on pre-processing to provide

445

gap-free and smoothed MODIS time series data. While various methods are available for gap-

446

filling remotely sensed time series, we devised our strategy for pre-processing based on data

447

availability, data quality, and other case study-specific information. For the STEP training

448

sites, utilizing land cover class information to fill data gaps using high quality observations

449

from nearby sites of the same class proved to be feasible and effective. In the NWFP study

450

area and Xingu River Basin, excluding periods of persistent clouds/snow which affect NBAR

451

retrievals helped to remove a large proportion of noise in annual time series that are unrelated

452

to vegetation conditions and land cover changes.

453

The distance metrics we used were sensitive to noise in the time series data, especially

454

sudden spikes. We used two measures to minimize the effects of such noise. First, we removed

455

spikes based on temporal context in the time series. Any sudden large increase or decrease

456

in consecutive values (greater than one standard deviation of the pixel’s time series) were

457

considered unrealistic and were replaced with interpolated values. This ensures that a single

458

anomalous data point will not dramatically increase the distance metric and lead to spurious

459

detection of change. Second, additional smoothing was applied using a 3-point median filter.

460

This approach performed well in all study areas, with the exception of a few STEP sites

461

with very high proportions of missing data.

462

Change Detection

463

Our analysis demonstrated that the distance metric-based change detection method was

464

effective for identifying forest change across a wide range of forest ecosystems at 500 m

465

spatial resolution. For STEP forest sites, MODIS EVI2 and MODIS band 7 data provide 18

466

complementary information that support detection of changes such as fire disturbance fol-

467

lowed by gradual recovery and deforestation by clear cutting and thinning. This provided a

468

proof-of-concept for broader application of the method to a wide range of forest ecosystems.

469

It does not, however, provide a comprehensive assessment of the method for all forest types,

470

different landscape patterns, and different types of change with varying intensity. For exam-

471

ple, we did not observe forest conversion to land use types such as urban and agriculture in

472

the STEP sites, nor did we specifically explore insect infestation which is a significant source

473

of disturbance in temperate and boreal forests.

474

Results from the NWFP study area show good agreement with Landsat-based results,

475

but also identified challenges associated with small area land cover changes and complex

476

management histories (Sulla-Menashe et al., 2013). Commission errors were relatively high.

477

While this can be partially explained by gridding artifacts in MODIS data, where distur-

478

bances in adjacent pixels lead to false detection of change, altitudinal gradient is the most

479

substantial source of commission error. First, the relatively limited set of training examples

480

we used for this work did not fully represent the diversity of the ENF spectral and temporal

481

signatures within the NWFP region, especially the ENF at higher elevations. Second, the

482

choice of a May-September temporal subset applying to these high elevation locations could

483

result in inclusion of snow cover and problematic gap-filling values. These suggests that fur-

484

ther stratification of the ENF class, varying length of the temporal subset, as well as better

485

characterization of the ENF class population, is needed to yield the best results from our

486

method.

487

In the Xingu River Basin case study, where the the disturbance regime is dominated by

488

forest clear cuts, we performed a more detailed analysis using varying thresholds in our dis-

489

tance metrics to identify change. The results confirmed that using a more generous threshold

490

(allowing more pixels to be flagged as potential change) results in smaller omission but higher

491

commission errors. This suggests that it may be necessary to adjust the choice of thresholds

492

used to detect change, and perhaps modify the definition of land cover change according to

493

the intended use of the change map. For instance, if we want to use change results from

494

MODIS to guide acquisition of higher spatial resolution data, we would choose thresholds

495

that yield fewer omission errors at the expense of allowing higher rates of commission errors. 19

496

It is important to note that several limitations of the PRODES dataset (described in

497

Section 3.3) affected our results for Xingu River Basin. In particular, inadequate represen-

498

tation of fire disturbance and forest degradation in PRODES produced apparent, but false,

499

errors of commission. While a large proportion of the ’false positives’ were located next to

500

the edges of deforestation patches, visual examination of higher spatial resolution images

501

indicated that considerable areas of fire disturbance, especially in 2010, were incorrectly

502

counted as commission errors since they were not accounted for in the PRODES dataset.

503

Additionally, regrowth in some historically deforested areas in the Xingu River Basin leads

504

to apparent omission errors. Finally, since the PRODES methodology uses two single-date

505

images, it identifies the timing of change incorrectly if clouds were present in the actual

506

year of change or if the change occurred after the image acquisition date, which can also

507

incorrectly introduce apparent errors of commission.

508

6. Conclusions

509

We developed and evaluated a distance metric-based change detection method that iden-

510

tifies annual forest cover change using MODIS 500 m NBAR time series. Our method exploits

511

two aspects of change that are reflected in each pixel’s annual time series: 1) deviation from

512

the pixel’s own past temporal behavior (measured by the between-year distance), and 2)

513

deviation in annual time series from neighboring pixels of the same class (measured by the

514

within-class distance). We established region- and class-specific populations and their dis-

515

tance distributions using high quality land cover training sites. By applying the distance

516

thresholds to both within-class and between-year distances at annual steps, we successfully

517

identified changes within temperate and boreal forest sites. For the NWFP area, our method

518

showed good agreement with Landsat-based disturbance areas, but also had relatively high

519

levels of commission errors. For the Xingu River Basin, our change map showed 80-90%

520

spatial agreement and close to 80% temporal agreement with the PRODES deforestation

521

dataset.

522

A key advantage of our approach is that instead of solely relying on the time series at

523

each pixel to detect change, it incorporates land cover related information from nearby pixels.

524

Using Olson geographic realms and biomes to define subgroups within IGBP forest classes, 20

525

we established regional class populations that characterize the range of variability in the dis-

526

tance metrics for undisturbed pixels. Our analysis of MODIS forest training sites worldwide

527

showed that we were able to identify changes due to logging and fire in geographically diverse

528

temperate and boreal forest ecosystems. Regional applications for the NWFP area and the

529

Xingu River Basin showed overall good agreement with Landsat-derived results. Further,

530

as new observations become available, this type of analysis can be repeated to ensure the

531

representativeness of our estimated class populations.

532

While the results presented in this paper show that our method was effective for the case

533

studies we examined, it has not been extensively tested over large areas in temperate, bo-

534

real, and tropical forests. Further, we have assumed knowledge of land cover in our analysis,

535

which means that in order to apply this approach over large areas, a good characteriza-

536

tion of the initial state of land cover is required. Patterns of missing data are important

537

for determining the optimal pre-processing strategy, which may require flexible gap-filling,

538

smoothing, subsetting or compositing methods to maximize information content related to

539

land cover change. Additionally, our change detection results within MODIS forest train-

540

ing sites suggests that the best feature (i.e., spectral bands and indices) varies for different

541

land cover classes and geographic regions. Thus, feature selection may significantly influence

542

change detection performance, especially with specific types of target disturbance (e.g., pest

543

infestation in western U.S.).

544

Moving forward, we plan to evaluate our approach more extensively in temperate, boreal,

545

and tropical forests, and to extend the method beyond forest ecosystems to map changes in

546

agricultural and urban land use. We will also continue to improve methods to characterize

547

class populations, especially for under-represented sub-class groups. Finally, due to the

548

coarse spatial resolution of MODIS observations, our ability to detect change is affected

549

by the size and intensity of change events. A logical next step is to adapt our method

550

for use with higher spatial, but lower temporal resolution time series data such as Landsat

551

TM/ETM+ imagery.

21

552

553

Acknowledgements The research described in this paper was supported by NASA grant numbers NNX11AE75G

554

and NNX11AG40G. The authors also gratefully acknowledge the National Institute for Space

555

Research in Brazil who generate and distribute the PRODES data, as well as Dr. Robert

556

Kennedy, who provided us with the Landsat-based Northwest Forest Plan forest change data.

557

References

558

Achard, F., Eva, H. D., Mayaux, P., Stibig, H.-J., Belward, A., 2004. Improved estimates of

559

net carbon emissions from land cover change in the tropics for the 1990s. Global Biogeo-

560

chemical Cycles 18 (2), GB2008.

561

562

563

564

Angelici, G., Bryant, N., 1977. A Land Use Change Monitoring System Based on LANDSAT. LARS Symposia. Arvor, D., Meirelles, M., Vargas, R., 2010. Monitoring land use changes around the indigenous lands of the Xingu Basin in Mato Grosso, Brazil. Proceedings of IGARSS ’10.

565

Baccini, A., Goetz, S. J., Walker, W. S., Laporte, N. T., Sun, M., Sulla-Menashe, D.,

566

Hackler, J., Beck, P. S. A., Dubayah, R., Friedl, M. A., Samanta, S., Houghton, R. A.,

567

2012. Estimated carbon dioxide emissions from tropical deforestation improved by carbon-

568

density maps. Nature Climate Change 2 (3), 182–185.

569

Beck, P., Atzberger, C., Hogda, K., Johansen, B., Skidmore, A., 2006. Improved monitoring

570

of vegetation dynamics at very high latitudes: A new method using MODIS NDVI. Remote

571

Sensing of Environment 100 (3), 321–334.

572

Bontemps, S., Bogaert, P., Titeux, N., Defourny, P., 2008. An object-based change detection

573

method accounting for temporal dependences in time series with medium to coarse spatial

574

resolution. Remote Sensing of Environment 112 (6), 3181–3191.

575

Broich, M., Hansen, M. C., Potapov, P., Adusei, B., Lindquist, E., Stehman, S. V., 2011.

576

Time-series analysis of multi-resolution optical imagery for quantifying forest cover loss in

22

577

Sumatra and Kalimantan, Indonesia. International Journal Of Applied Earth Observation

578

And Geoinformation 13 (2), 277–291.

579

580

Bucha, T., Stibig, H.-J., 2008. Analysis of MODIS imagery for detection of clear cuts in the boreal forest in north-west Russia. Remote Sensing of Environment 112 (5), 2416–2429.

581

Carroll, M., Townshend, J. R. G., Hansen, M., DiMiceli, C., Sohlberg, R., Wurster, K.,

582

2011. MODIS Vegetative Cover Conversion and Vegetation Continuous Fields. In: Land

583

Remote Sensing and Global Environmental Change: NASA’s Earth Observing System and

584

the Science of ASTER and MODIS. Springer New York, New York, NY, pp. 725–745.

585

Chen, J., Jonsson, P., Tamura, M., Gu, Z., Matsushita, B., Eklundh, L., 2004. A simple

586

method for reconstructing a high-quality NDVI time-series data set based on the Savitzky-

587

Golay filter. Remote Sensing of Environment 91 (3-4), 332–344.

588

Cohen, W. B., Fiorella, M., Gray, J., Helmer, E., Anderson, K., 1998. An efficient and

589

accurate method for mapping forest clearcuts in the Pacific Northwest using Landsat

590

imagery. Photogrammetric Engineering And Remote Sensing 64 (4), 293–299.

591

Collins, J. B., Woodcock, C. E., 1994. Change detection using the Gramm-Schmidt trans-

592

formation applied to mapping forest mortality. Remote Sensing of Environment 50 (3),

593

267–279.

594

Collins, J. B., Woodcock, C. E., 1996. An assessment of several linear change detection

595

techniques for mapping forest mortality using multitemporal landsat TM data. Remote

596

Sensing of Environment 56 (1), 66–77.

597

Coops, N. C., Wulder, M. A., Iwanicka, D., Jun. 2009. Large area monitoring with a MODIS-

598

based Disturbance Index (DI) sensitive to annual and seasonal variations. Remote Sensing

599

of Environment 113 (6), 1250–1261.

600

Copeland, J., Pielke, R., Kittel, T., 1996. Potential climatic impacts of vegetation change: A

601

regional modeling study. Journal of Geophysical Research-Atmospheres 101 (D3), 7409–

602

7418.

23

603

604

605

606

607

608

de Beurs, K. M., Townsend, P. A., 2008. Estimating the effect of gypsy moth defoliation using MODIS. Remote Sensing of Environment 112 (10), 3983–3990. Dickinson, R. E., Kennedy, P., 1992. Impacts on regional climate of Amazon deforestation. Geophysical Research Letters 19 (19), 1947–1950. Ehrlich, D., Estes, J., Singh, A., 1994. Applications of NOAA-AVHRR 1 km data for environmental monitoring. International Journal of Remote Sensing 15 (1), 145–161.

609

Fang, H., Liang, S., Townshend, J. R. G., Dickinson, R. E., 2008. Spatially and temporally

610

continuous LAI data sets based on an integrated filtering method: Examples from North

611

America. Remote Sensing of Environment 112 (1), 75–93.

612

Foley, J. A., DeFries, R. S., Asner, G. P., Barford, C. C., Bonan, G. B., Carpenter, S. R.,

613

Chapin, III, F. S., Coe, M. T., Daily, G. C., Gibbs, H. K., Helkowski, J. H., Holloway, T.,

614

Howard, E. A., Kucharik, C. J., Monfreda, C., Patz, J. A., Prentice, I. C., Ramankutty,

615

N., Snyder, P. K., 2005. Global consequences of land use. Science 309 (5734), 570–574.

616

Franklin, J. F., Dyrness, C. T., 1988. Natural vegetation of Oregon and Washington, 2nd

617

Edition. Oregon State University Press, Corvallis, OR.

618

Friedl, M. A., McIver, D., Hodges, J., Zhang, X., Muchoney, D., Strahler, A., Woodcock,

619

C. E., Gopal, S., Schneider, A., Cooper, A., Baccini, A., Gao, F., Schaaf, C., 2002.

620

Global land cover mapping from MODIS: algorithms and early results. Remote Sensing of

621

Environment 83, 287–302.

622

Friedl, M. A., Sulla-Menashe, D., Tan, B., Schneider, A., Ramankutty, N., Sibley, A., Huang,

623

X., 2010. MODIS Collection 5 global land cover: Algorithm refinements and characteriza-

624

tion of new datasets. Remote Sensing of Environment 114 (1), 168–182.

625

Furley, P. A., Ratter, J. A., Gifford, D. R., 1988. Observations on the Vegetation of Eastern

626

Mato Grosso, Brazil. III. The Woody Vegetation and Soils of the Morro de Fumaca,

627

Torixoreu. Proceedings Of The Royal Society Of London Series B-Biological Sciences

628

235 (1280), 259–280.

24

629

Gao, F., Morisette, J., Wolfe, R., Ederer, G., Pedelty, J., Masuoka, E., Myneni, R., Tan,

630

B., Nightingale, J., 2008. An algorithm to produce temporally and spatially continuous

631

MODIS-LAI time series. IEEE Geoscience and Remote Sensing Letters 5 (1), 60–64.

632

Hansen, M. C., Roy, D. P., Lindquist, E., Adusei, B., Justice, C. O., Altstatt, A., 2008a.

633

A method for integrating MODIS and Landsat data for systematic monitoring of forest

634

cover and change in the Congo Basin. Remote Sensing of Environment 112, 2495–2513.

635

Hansen, M. C., Shimabukuro, Y., Potapov, P., Pittman, K., 2008b. Comparing annual

636

MODIS and PRODES forest cover change data for advancing monitoring of Brazilian

637

forest cover. Remote Sensing of Environment 112 (10), 3784–3793.

638

Hansen, M. C., Stehman, S., Potapov, P., 2010. Quantification of global gross forest cover

639

loss. Proceedings of the National Academy of Sciences of the United States Of America

640

107 (19), 8650–8655.

641

Hansen, M. C., Stehman, S., Potapov, P., Arunarwati, B., Stolle, F., Pittman, K., 2009.

642

Quantifying changes in the rates of forest clearing in Indonesia from 1990 to 2005 using

643

remotely sensed data sets. Environmental Research Letters 4 (3).

644

Hayes, D. J., Cohen, W., Sader, S. A., Irwin, D., 2008. Estimating proportional change in

645

forest cover as a continuous variable from multi-year MODIS data. Remote Sensing of

646

Environment 112 (3), 735–749.

647

Hayes, D. J., Sader, S. A., 2001. Comparison of change-detection techniques for monitor-

648

ing tropical forest clearing and vegetation regrowth in a time series. Photogrammetric

649

Engineering And Remote Sensing 67 (9), 1067–1075.

650

Hocke, K., Kaempfer, N., 2009. Gap filling and noise reduction of unevenly sampled data

651

by means of the Lomb-Scargle periodogram. Atmospheric Chemistry and Physics 9 (12),

652

4197–4206.

653

654

Holben, B. N., 1986. Characteristics of maximum-value composite images from temporal AVHRR data. International Journal of Remote Sensing.

25

655

Huang, C., Coward, S. N., Masek, J. G., Thomas, N., Zhu, Z., Vogelmann, J. E., 2010.

656

An automated approach for reconstructing recent forest disturbance history using dense

657

Landsat time series stacks. Remote Sensing of Environment 114 (1), 183–198.

658

Huete, A., Didan, K., Miura, T., Rodriguez, E., Gao, X., Ferreira, L., 2002. Overview of

659

the radiometric and biophysical performance of the MODIS vegetation indices. Remote

660

Sensing of Environment 83, 195–213.

661

662

663

664

Hunt, Jr, E. R., Rock, B. N., 1989. Detection of changes in leaf water content using Nearand Middle-Infrared reflectances. Remote Sensing of Environment 30 (1), 43–54. INPE, 2012. Project PRODES: Monitoring the Brazilian Amazon Forests by Satellite. Available at: http://www.dpi.inpe.br/prodesdigital/.

665

Jensen, J. R., Rutchey, K., Koch, M., Narumalani, S., 1995. Inland wetland change detection

666

in the Everglades Water Conservation Area 2A using a time-series of normalized remotely-

667

sensed data. Photogrammetric Engineering And Remote Sensing 61 (2), 199–209.

668

669

670

671

Ji, C., Liu, Q., Sun, D., Wang, S., Lin, P., Li, X., 2001. Monitoring urban expansion with remote sensing in China. International Journal of Remote Sensing 22 (8), 1441–1455. Jiang, Z., Huete, A., Didan, K., 2008. Development of a two-band enhanced vegetation index without a blue band. Remote Sensing of Environment.

672

Johnson, R., Kasischke, E., 1998. Change vector analysis: a technique for the multispectral

673

monitoring of land cover and condition. International Journal of Remote Sensing 19 (3),

674

411–426.

675

Jonsson, P., Eklundh, L., 2002. Seasonality extraction by function fitting to time-series of

676

satellite sensor data. IEEE Transactions on Geoscience and Remote Sensing 40 (8), 1824–

677

1832.

678

679

Jonsson, P., Eklundh, L., 2004. TIMESAT - a program for analyzing time-series of satellite sensor data. Computers & Geosciences 30 (8), 833–845.

26

680

Kennedy, R., Cohen, W., Schroeder, T., 2007. Trajectory-based change detection for auto-

681

mated characterization of forest disturbance dynamics. Remote Sensing of Environment

682

110 (3), 370–386.

683

Kennedy, R. E., Yang, Z., Cohen, W. B., Pfaff, E., Braaten, J., Nelson, P., 2012. Spatial and

684

temporal patterns of forest disturbance and regrowth within the area of the Northwest

685

Forest Plan. Remote Sensing of Environment 122, 117–133.

686

Key, C. H., Benson, N. C., 1999. Measuring and remote sensing of burn severity: the CBI

687

and NBR. In: Proceedings Joint Fire Science Conference and Workshop. Boise, ID, pp.

688

1–1.

689

690

691

692

Lambin, E. F., 1997. Modelling and monitoring land-cover change processes in tropical regions. Progress in Physical Geography 21 (3), 375–393. Lambin, E. F., Geist, H., Lepers, E., Nov. 2003. Dynamics of land-use and land-cover change in tropical regions. Annu. Rev. Environ. Resour. 28 (1), 205–241.

693

Lambin, E. F., Strahler, A., 1994. Change-Vector Analysis in Multitemporal Space: A Tool

694

to Detect and Categorize Land-Cover Change Processes Using High Temporal-Resolution

695

Satellite Data. Remote Sensing of Environment 48 (2), 231–244.

696

Lepers, E., Lambin, E. F., Janetos, A., DeFries, R. S., Achard, F., Ramankutty, N., Scholes,

697

R., 2005. A synthesis of information on rapid land-cover change for the period 1981-2000.

698

Bioscience 55 (2), 115–124.

699

Lhermitte, S., Verbesselt, J., Verstraeten, W. W., Coppin, P., 2011. A comparison of time

700

series similarity measures for classification and change detection of ecosystem dynamics.

701

Remote Sensing of Environment 115 (12), 3129–3152.

702

Linderman, M., Rowhani, P., Benz, D., Serneels, S., Lambin, E. F., 2005. Land-cover change

703

and vegetation dynamics across Africa. Journal of Geophysical Research-Atmospheres

704

110 (D12).

27

705

706

707

708

Mahalanobis, P. C., 1936. On the generalized distance in statistics. In: Proceedings National Institute of Science, India. pp. 49–55. Malila, W., 1980. Change vector analysis: an approach for detecting forest changes with Landsat. LARS Symposia.

709

Masek, J. G., Huang, C., Wolfe, R., Cohen, W., Hall, F., Kutler, J., Nelson, P., 2008. North

710

American forest disturbance mapped from a decadal Landsat record. Remote Sensing of

711

Environment 112 (6), 2914–2926.

712

713

714

715

716

717

Mildrexler, D., Zhao, M., Running, S., 2009. Testing a MODIS global disturbance index across North America. Remote Sensing of Environment 113 (10), 2103–2117. Moody, A., Johnson, D., 2001. Land-surface phenologies from AVHRR using the discrete fourier transform. Remote Sensing of Environment 75 (3), 305–323. Moody, A., Strahler, A., 1994. Characteristics of composited AVHRR data and problems in their classification. International Journal of Remote Sensing 15 (17), 3473–3491.

718

Morton, D., DeFries, R. S., Shimabukuro, Y., Anderson, L., Espirito-Santo, F., Hansen,

719

M. C., Carroll, M., 2005. Rapid assessment of annual deforestaion in the Brazilian Amazon

720

using MODIS data. Earth Interactions 9, 22.

721

Mucher, C., Steinnocher, K., Kressler, F., Heunks, C., 2000. Land cover characterization

722

and change detection for environmental monitoring of pan-Europe. International Journal

723

of Remote Sensing 21 (6-7), 1159–1181.

724

Muchoney, D., Strahler, A., Hodges, J., LoCastro, J., 1999. The IGBP DISCover confidence

725

sites and the system for terrestrial ecosystem parameterization: Tools for validating global

726

land-cover data. Photogrammetric Engineering And Remote Sensing 65 (9), 1061–1067.

727

Olson, D., Dinerstein, E., Wikramanayake, E., Burgess, N., Powell, G., Underwood, E.,

728

D’Amico, J., Itoua, I., Strand, H., Morrison, J., Loucks, C., Allnutt, T., Ricketts, T., Kura,

729

Y., Lamoreux, J., Wettengel, W., Hedao, P., Kassem, K., 2001. Terrestrial ecoregions of

730

the worlds: A new map of life on Earth. Bioscience 51 (11), 933–938. 28

731

Pielke, R., Marland, G., Betts, R., Chase, T., Eastman, J., Niles, J., Niyogi, D., Running,

732

S., 2002. The influence of land-use change and landscape dynamics on the climate system:

733

relevance to climate-change policy beyond the radiative effect of greenhouse gases. Philo-

734

sophical Transactions of the Royal Society of London Series A: Mathematical Physical and

735

Engineering Sciences 360 (1797), 1705–1719.

736

Potapov, P., Hansen, M. C., Stehman, S., Loveland, T., Pittman, K., 2008. Combining

737

MODIS and Landsat imagery to estimate and map boreal forest cover loss. Remote Sensing

738

of Environment 112 (9), 3708–3719.

739

740

Ridd, M. K., Liu, J., 1998. A comparison of four algorithms for change detection in an urban environment. Remote Sensing of Environment 63 (2), 95–100.

741

Schaaf, C., Gao, F., Strahler, A., Lucht, W., Li, X., Tsang, T., Strugnell, N., Zhang, X., Jin,

742

Y., Muller, J., Lewis, P., Barnsley, M., Hobson, P., Disney, M., Roberts, G., Dunderdale,

743

M., Doll, C., d’Entremont, R., Hu, B., Liang, S., Privette, J., Roy, D. P., 2002. First

744

operational BRDF, albedo nadir reflectance products from MODIS. Remote Sensing of

745

Environment 83 (1-2), 135–148.

746

Sellers, P. J., Tucker, C. J., COLLATZ, G., Los, S., Justice, C. O., Dazlich, D., Randall,

747

D., 1994. A global 1 ◦ by 1 ◦ NDVI data set for climate studies. Part 2: The generation of

748

global fields of terrestrial biophysical parameters from the NDVI. International Journal of

749

Remote Sensing 15 (17), 3519–3545.

750

751

752

753

Singh, A., 1989. Digital change detection techniques using remotely-sensed data. International Journal of Remote Sensing. Skole, D., Tucker, C. J., 1993. Tropical Deforestation and Habitat Fragmentation in the Amazon: Satellite Data from 1978 to 1988. Science 260 (5116), 1905–1910.

754

Souza, C., Roberts, D., Cochrane, M., 2005. Combining spectral and spatial information to

755

map canopy damage from selective logging and forest fires. Remote Sensing of Environment

756

98 (2-3), 329–343.

29

757

Sulla-Menashe, D., Kennedy, R. E., Yang, Z., Braaten, J., Krankina, O. N., Friedl, M. A.,

758

2013. Detecting forest disturbance in the Pacific Northwest from MODIS time series using

759

temporal segmentation. Remote Sensing of Environment (in press).

760

Tan, B., Woodcock, C. E., Hu, J., Zhang, P., Ozdogan, M., Huang, D., Yang, W.,

761

Knyazikhin, Y., Myneni, R. B., 2006. The impact of gridding artifacts on the local spatial

762

properties of MODIS data: Implications for validation, compositing, and band-to-band

763

registration across resolutions. Remote Sensing of Environment 105 (2), 98–114.

764

765

Tarpley, J., Schneider, S., Money, R., 1984. Global vegetation indices from the NOAA-7 meteorological satellite. Journal of Climate and Applied Meteorology 23, 491–494.

766

Thomas, J. W., Franklin, J. F., Gorden, J., Johnson, K. N., Apr. 2006. The Northwest Forest

767

Plan: Origins, Components, Implementation Experience, and Suggestions for Change.

768

Conservation Biology 20 (2), 277–287.

769

Valeriano, D. M., Mello, E. M. K., Moreira, J. C., Shimabukuro, Y. E., Duarte, V., e Souza,

770

I. M., dos Santos, J. R., Barbosa, C. C. F., de Souza, R. C. M., 2004. Monitoring tropical

771

forest from space: the Prodes Digital project. Proceedings of ISPRS ’04.

772

Verbesselt, J., Hyndman, R., Newnham, G., Culvenor, D., 2010. Detecting trend and sea-

773

sonal changes in satellite image time series. Remote Sensing of Environment 114 (1),

774

106–115.

775

Vermote, E. F., Kotchenova, S., 2011. MODIS Directional Surface Reflectance Product:

776

Method, Error Estimates and Validation. In: Land Remote Sensing and Global Environ-

777

mental Change: NASA’s Earth Observing System and the Science of ASTER and MODIS.

778

Springer New York, New York, NY, pp. 533–547.

779

Viovy, N., Arino, O., Belward, A., 1992. The Best Index Slope Extraction (BISE): A method

780

for reducing noise in NDVI time-series. International Journal of Remote Sensing 13 (8),

781

1585–1590.

782

783

Vitousek, P. M., 1994. Beyond global warming: ecology and global change. Ecology 75 (7), 1861–1876. 30

784

Walker, W., Stickler, C., Kellndorfer, J., Kirsch, K., Nepstad, D., 2010. Large-Area Classifi-

785

cation and Mapping of Forest and Land Cover in the Brazilian Amazon: A Comparative

786

Analysis of ALOS/PALSAR and Landsat Data Sources. IEEE Journal of Selected Topics

787

in Applied Earth Observations and Remote Sensing 3 (4), 594–604.

788

789

Waring, R. H., Franklin, J. F., 1979. Evergreen coniferous forests of the Pacific Northwest. Science 204 (4400), 1380–1386.

790

Xian, G., Homer, C. G., Fry, J. A., Jun. 2009. Updating the 2001 National Land Cover

791

Database land cover classification to 2006 by using Landsat imagery change detection

792

methods. Remote Sensing of Environment 113 (6), 1133–1147.

793

Young, S., Wang, C., 2001. Land-cover change analysis of China using global-scale Pathfinder

794

AVHRR Landcover (PAL) data, 1982-92. International Journal of Remote Sensing 22 (8),

795

1457–1477.

796

Zhan, X., Huang, C., Townshend, J. R. G., DeFries, R. S., Hansen, M. C., Dimiceli, C.,

797

Sohlberg, R., Hewson-Scardelletti, J., Tompkins, A., 1998. Land cover change detection

798

with change vector in the red and near-infrared reflectance space. Proceedings of IGARSS

799

’98 2, 859–861.

31

800

801

List of Figures 1

802 803 804 805

2

806 807 808 809 810 811 812

3

813 814

4

815 816

5

817 818 819 820 821 822

6

823 824

7

825 826

8

827 828 829 830 831 832 833

9

834 835 836 837 838

10 11

839 840 841 842

12

Cumulative percent deforestation in 2003 and 2010 within MODIS 500 m pixels located in the Xingu River Basin based on PRODES dataset (1(a), 1(b)). Distribution of sub-pixel deforestation size for new deforested pixels that appeared between 2003-2010 (1(c)). . . . . . . . . . . . . . . . . . . . . Calculation of within-class and between-year distances for an EBF pixel affected by clear cutting. Band 7 time series for this pixel is shown in 2(a), with the year of interest (2007) highlighted. The contrast between the pixel’s annual time series and those of the EBF population is captured by withinclass distance (2(b), annual mean EBF time series ± one standard deviation plotted in red). The contrast between annual time series before and after the year of interest (2(c)) is captured by between-year distance. . . . . . . . . . . Distribution of within-class distances calculated for band 7 for EBF population and 95% distance threshold. . . . . . . . . . . . . . . . . . . . . . . . . . . . Standard deviations of the forest populations versus RMSEs of gap-filled values in EVI2. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Examples of changed pixels detected in STEP temperate and boreal forest sites. The unprocessed NBAR EVI2 (green) and band 7 (blue) time series are plotted for 2001-2010, with 4 grey dashed lines indicating acquisition times of Landsat images shown below. Within each Landsat snapshot, the white polygon shows the boundaries of the STEP site, and the yellow polygon corresponds to the 500 m MODIS grid for the pixel plotted. . . . . . . . . . . . Change map for the NWFP and locations of the omission and commission errors compared to the Landsat-based disturbance map. . . . . . . . . . . . . Comparison between proportions of change pixels detected using Landsat and MODIS at 10×10 km block scale. The red line is the 1:1 line. . . . . . . . . Within-class distance, between-year distance, and potential change masks for Xingu River Basin in 2003 (left column) and 2010 (right column). The withinclass distance (top row) compares all pixels to forest baseline in 2003 and 2010. The between-year distance (middle row) compares each pixel’s own temporal profiles, between 2002 and 2004, and 2009 and 2011. The potential change masks (bottom row) show the results from combining within-class and between-year distances using 95% threshold. . . . . . . . . . . . . . . . . . . Producer’s accuracies for annual potential change masks for Xingu River Basin. Each point represents mean accuracy over 2003-2010 for each 10% deforestation interval. The error bars shown are one standard deviation from the mean accuracy. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Change map of Xingu River Basin using 95% distance thresholds. . . . . . . Zoom-in comparison between MODIS year of change (left column) and PRODES year of change aggregated to 500 m (right column) for two 30×30 km blocks. Producer’s and user’s accuracy for final change maps of Xingu River Basin, using varying thresholds for within-class distance and between-year distance. Thresholds for the between-year distance are displayed in log10 scale. . . . .

32

33

34 35 36

37 38 39

40

41 42 43

44

(a)

(c)

(b) Figure 1: Cumulative percent deforestation in 2003 and 2010 within MODIS 500 m pixels located in the Xingu River Basin based on PRODES dataset (1(a), 1(b)). Distribution of sub-pixel deforestation size for new deforested pixels that appeared between 2003-2010 (1(c)).

33

(a)

(b)

(c)

Figure 2: Calculation of within-class and between-year distances for an EBF pixel affected by clear cutting. Band 7 time series for this pixel is shown in 2(a), with the year of interest (2007) highlighted. The contrast between the pixel’s annual time series and those of the EBF population is captured by within-class distance (2(b), annual mean EBF time series ± one standard deviation plotted in red). The contrast between annual time series before and after the year of interest (2(c)) is captured by between-year distance.

34

Figure 3: Distribution of within-class distances calculated for band 7 for EBF population and 95% distance threshold.

35

Figure 4: Standard deviations of the forest populations versus RMSEs of gap-filled values in EVI2.

36

(a) Logging in an evergreen needleleaf forest pixel.

(b) Forest fire in a deciduous needleleaf forest pixel.

(c) Logging in a deciduous broadleaf forest pixel.

(d) Logging in a mixed forest pixel.

Figure 5: Examples of changed pixels detected in STEP temperate and boreal forest sites. The unprocessed NBAR EVI2 (green) and band 7 (blue) time series are plotted for 2001-2010, with 4 grey dashed lines indicating acquisition times of Landsat images shown below. Within each Landsat snapshot, the white polygon shows the boundaries of the STEP site, and the yellow polygon corresponds to the 500 m MODIS grid for the pixel plotted.

37

(a) ENF change in NWFP between 2002-2010.

(b) Locations of omission errors. (c) Locations of commission errors.

Figure 6: Change map for the NWFP and locations of the omission and commission errors compared to the Landsat-based disturbance map.

38

0.8

Proportion of change pixels from MODIS

0.7

0.6

0.5

0.4

0.3

0.2

0.1

0 0

0.1

0.2 0.3 0.4 0.5 0.6 Proportion of change pixels from Landsat

0.7

0.8

Figure 7: Comparison between proportions of change pixels detected using Landsat and MODIS at 10×10 km block scale. The red line is the 1:1 line.

39

Figure 8: Within-class distance, between-year distance, and potential change masks for Xingu River Basin in 2003 (left column) and 2010 (right column). The within-class distance (top row) compares all pixels to forest baseline in 2003 and 2010. The between-year distance (middle row) compares each pixel’s own temporal profiles, between 2002 and 2004, and 2009 and 2011. The potential change masks (bottom row) show the results from combining within-class and between-year distances using 95% threshold.

40

% pixels detected in annual masks

Proportion of total change pixels detected w.r.t. change size, 2003−2010 100 using 90% threshold 90

using 95% threshold

80 70 60 50 40 30

(0,10] (10,20] (20,30] (30,40] (40,50] (50,60] (60,70] (70,80] (80,90](90,100]

% deforestation in 500 m pixels Figure 9: Producer’s accuracies for annual potential change masks for Xingu River Basin. Each point represents mean accuracy over 2003-2010 for each 10% deforestation interval. The error bars shown are one standard deviation from the mean accuracy.

41

Figure 10: Change map of Xingu River Basin using 95% distance thresholds.

42

Figure 11: Zoom-in comparison between MODIS year of change (left column) and PRODES year of change aggregated to 500 m (right column) for two 30×30 km blocks.

43

Figure 12: Producer’s and user’s accuracy for final change maps of Xingu River Basin, using varying thresholds for within-class distance and between-year distance. Thresholds for the between-year distance are displayed in log10 scale.

44

843

844 845 846 847

List of Tables 1 2 3 4

Temperate and boreal disturbance sites detected. . . . . . . . . . . . . . . . Error matrix of sample counts of NWFP. . . . . . . . . . . . . . . . . . . . . Accuracy assessment for annual potential change masks of Xingu River Basin. Accuracy assessment for MODIS change map of Xingu River Basin. . . . . .

45

46 47 48 49

IGBP Class

Number of Sites

Change Type

7 2 3 4

logging, fire fire logging logging

Evergreen Needleleaf Forest (ENF) Deciduous Needleleaf Forest (DNF) Deciduous Broadleaf Forest (DBF) Mixed Forest (MXF)

Table 1: Temperate and boreal disturbance sites detected.

46

Landsat(reference)

Change No change

Change

No change

Total

309 13

191 487

500 500

Table 2: Error matrix of sample counts of NWFP.

47

Overall producer’s accuracy

95% threshold 90% threshold

Overall user’s accuracy

>0% >20% >40% >60% >80%

>0%

78.37 83.04

94.83 93.26

82.78 87.08

85.09 89.08

86.74 90.46

88.09 91.54

>20% >40% >60% >80% 90.13 88.00

85.74 83.32

Table 3: Accuracy assessment for annual potential change masks of Xingu River Basin.

48

80.42 77.84

72.98 70.39

Spatial agreement >0% P accuracy U accuracy 95% threshold 90% threshold

82.29 86.84

Temporal agreement

>20% P accuracy U accuracy

72.32 66.67

90.33 93.54

>0% >20% P accuracy P accuracy

63.27 57.23

Table 4: Accuracy assessment for MODIS change map of Xingu River Basin.

49

77.56 78.33

78.97 79.11