02-defid2_gee_code.RmdThe overall DEFID2 dataset consists of a core database containing disturbance geometries and their attributes, and an additional extension database containing derived metrics (temporal segments from LandTrendr algorithm, NBR time-series and co-occurrence with fire and windthrow datasets). These metrics were extracted for each ‘exact’ point or polygon geometry using the Google Earth Engine (GEE) computing platform. Code to reproduce these steps is explained below.
By visiting the URL: https://code.earthengine.google.com/?accept_repo=users/agataelia1991/DEFID2 the public users/agataelia1991/DEFID2 repository will be added to any GEE account. This repository contains two folders:
Codes folder contains the individual codes used to produce the DEFID2 extension outputs and statistics. These are briefly explained in the sections below.
Modules folder contains the modified LandTrendr JavaScript (js) module necessary to build LandTrendr input collections, run LandTrendr, and deal with the outputs. This module is briefly explained in the section below.
LandTrendrModified.js is the JavaScript (js) module used to generate the Landsat based yearly medoid composites, to run LandTrendr (LT) given a set of parameters and to retrieve the LandTrendr (LT) spectral-temporal segments information array.
This module is based on the LandTrendr.js module retrieved from https://github.com/eMapR/LT-GEE with in depth documentation at https://emapr.github.io/LT-GEE/index.html.
Some additional functions were added to the original LandTrendr.js module in order to generate yearly per pixel medoid composites exclusively based on the three consecutive greenness months per pixel centered on the average highest NDVI month identified with the 00_MODIS_max_NDVI_month code. In order to be used, the LandTrendrModified.js module must be imported into the scripts requiring the module by using the following line:
This code identifies the average highest NDVI month per pixel based on two inputs:
var modisCollection = ee.ImageCollection("MODIS/006/MOD13Q1"),
europeShp = ee.FeatureCollection("users/agataelia1991/aoi/europe_filtered");
Map.centerObject(europeShp, 7);
var europe = europeShp.geometry().bounds();
/////////////////////////////////
//// DEFINE USEFUL FUNCTIONS ////
/////////////////////////////////
// Define a function that rescales images
var rescale = function(image){
var rescaled = image.reproject({
crs: 'EPSG:4326',
scale: 250
});
return rescaled;
};
// Define a function that adds a band with indication of month of reference per pixel
var addMonthBand = function(image, month){
var monthImage = ee.Image.constant(month).rename('MONTH').cast({'MONTH': 'uint16'});
var newImage = ee.Image.cat(image, monthImage);
return newImage
}
//////////////////////////////////////////////////////////////////////////////////////////
//// CREATE COLLECTION OF AVERAGE MONTHLY MODIS NDVI AND QUALITY MOSAIC OVER MAX NDVI ////
//////////////////////////////////////////////////////////////////////////////////////////
// Create collection of average monthly moving average (3 months) of NDVI between 2000 and 2020
var modisMonthlyAverageCollection = ee.ImageCollection.fromImages(
ee.List.sequence(1, 12).map(function(month) {
var m = ee.Number(month);
var mp1;
var mm1;
if(m==1){
mp1 = m.add(1);
mm1 = 12;
} else if (m==12) {
mp1 = 1;
mm1 = m.subtract(1);
} else {
mp1 = m.add(1);
mm1 = m.subtract(1);
}
var modisMonthlyAverage = modisCollection
.filterBounds(europe)
.select('NDVI')
.filter(ee.Filter.calendarRange(2000, 2020, 'year'))
.filter(ee.Filter.calendarRange(mm1, mp1, 'month'))
.reduce(ee.Reducer.mean())
.rename('NDVI');
var modisMonthlyAverageWithMonth = addMonthBand(modisMonthlyAverage, m);
var modisMonthlyAverageRescaled = rescale(modisMonthlyAverageWithMonth);
return modisMonthlyAverageRescaled;
}));
// Inspect results
print(modisMonthlyAverageCollection);
// Quality mosaic where maximum NDVI to identify month with maximum average NDVI over 20 years
var maxNDVI = modisMonthlyAverageCollection.qualityMosaic('NDVI');
print(maxNDVI);
// Visualize results
var palette = ['aec3d4', '152106', '225129', '369b47', '30eb5b', '387242',
'6a2325', 'c3aa69', 'b76031', 'd9903d', '91af40', '30fb5c' ];
Map.addLayer(maxNDVI.select('MONTH'), {min: 1, max: 12, palette: palette}, 'maxNDVImonth');
////////////////
//// EXPORT ////
////////////////
Export.image.toAsset({
image: maxNDVI.select('MONTH').uint16(),
description: 'maxNDVImonthMAuint16Europe',
assetId: 'users/agataelia1991/MODIS/maxNDVImonthMAuint16Europe',
region: europe,
scale: 250,
crs: 'EPSG:4326',
maxPixels: 1e13});
Export.image.toDrive({
image: maxNDVI.select('MONTH').clip(europe).uint16(),
description: 'maxNDVImonthMAuint16Europe',
scale: 250,
crs: 'EPSG:4326',
region: europe,
fileFormat: 'GeoTIFF',
maxPixels: 1e13
});This code leverages the LandTrendrModified.js module to generate yearly maps of LandTrendr identified segments intersecting temporally and spatially the DEFID2 records. This code takes as inputs:
var hansen = ee.Image("users/agataelia1991/Hansen/hansenTreecover2000Mask");
///////////////////////////////////////////////
//// DEFINE COUNTRY, BOUNDING BOX AND DATA ////
///////////////////////////////////////////////
// Define country code to run this code by country to avoid memory issues
var country = 'RO';
// Retrieve the correspoding data bounding box from the assets
var bbFeature = ee.FeatureCollection('users/agataelia1991/DEFID2_Loic/AOI_by_country/' + country + '_aoi');
var bb = bbFeature.geometry();
// Retrieve the correspoding DEFID2 data from the assets
var defid = ee.FeatureCollection('users/agataelia1991/DEFID2_Loic/DEFID2_by_country/' + country + '_defid_exact_polygons_v1');
///////////////////////////////////////////////////////
//// IDENFITY DISTINCT SURVEY YEARS IN DEFID2 DATA ////
///////////////////////////////////////////////////////
var surveyYears = defid.aggregate_array('YEAR').distinct().getInfo();
/////////////////////////////////
//// DEFINE USEFUL FUNCTIONS ////
/////////////////////////////////
// Define a function that rescale images
var rescale = function(image, meter){
var rescaled = image.reproject({
crs: 'EPSG:4326',
scale: meter
});
return rescaled;
};
////////////////////////
//// RUN LANDTRENDR ////
////////////////////////
// Load the LandTrendr.js module
var ltgee = require('users/agataelia1991/DEFID2:Modules/LandTrendrModified.js');
// Define collection parameters
var startYear = 1995;
var endYear = 2021;
var startDay = '12-01';
var endDay = '01-31';
var aoi = bb;
var index = 'NBR';
var maskThese = ['cloud', 'shadow', 'snow', 'water'];
// Define landtrendr parameters
var runParams = {
maxSegments: 6,
spikeThreshold: 0.9,
vertexCountOvershoot: 3,
preventOneYearRecovery: true,
recoveryThreshold: 0.25,
pvalThreshold: 0.05,
bestModelProportion: 0.75,
minObservationsNeeded: 6
};
// Run landtrendr
var lt = ltgee.runLT(startYear, endYear, startDay, endDay, aoi, index, [], runParams, maskThese);
// Get segments informations
var segInfo = ltgee.getSegmentData(lt, index, 'all', true);
// Reproject LT output to 30 meters resolution
var segInfo30m = rescale(segInfo, 30);
////////////////////////////////////////////
//// START LOOPING THROUGH SURVEY YEARS ////
////////////////////////////////////////////
for(var i in surveyYears) {
var year = surveyYears[i];
print(year);
// Select only DEFID2 data with current survey year
var yearlyDefid = defid.filter(ee.Filter.eq('YEAR', year));
/////////////////////////////////////////////////////////////////////////
//// RASTERIZE AND REDUCE DEFID2 DATASET BY YEAR OF SURVEY PER PIXEL ////
/////////////////////////////////////////////////////////////////////////
// Rasterize DEFID2, reduce to year of survey per pixel and rescale to 30 meters resolution
var defidRaster = (yearlyDefid.reduceToImage({ properties: ['YEAR'], reducer: ee.Reducer.min()})).rename(['YEAR']);
var defidRaster30m = rescale(defidRaster, 30);
// Mask lt segments with DEFID2 layer
var segMasked = segInfo30m.updateMask(defidRaster30m.select(['YEAR']).mask());
// Create image collection of image stacks for each segment, each band is a segment attribute
var segments = [0, 1, 2, 3, 4, 5];
var segCol = ee.ImageCollection.fromImages(
segments.map(function(seg) {
var segment = segMasked.arraySlice(1, seg, seg+1);
var segImg = ee.Image.cat(segment.arraySlice(0,0,1).arrayProject([1]).arrayFlatten([['startYr']]),
segment.arraySlice(0,1,2).arrayProject([1]).arrayFlatten([['endYr']]),
segment.arraySlice(0,2,3).arrayProject([1]).arrayFlatten([['startVal']]),
segment.arraySlice(0,4,5).arrayProject([1]).arrayFlatten([['mag']]),
segment.arraySlice(0,5,6).arrayProject([1]).arrayFlatten([['dur']]),
segment.arraySlice(0,6,7).arrayProject([1]).arrayFlatten([['rate']]));
// Filter the segments and cast to int16
var segImgNotZero = segImg.select(['startYr']).gt(0);
var segImgInt = segImg.updateMask(segImgNotZero).int16();
return segImgInt;
}));
/////////////////////////////////////////////////////////////
//// GET COLLECTION OF LT INFORMATION ENRICHED DEFID2 DATA //
/////////////////////////////////////////////////////////////
// Map through the LT segments collection and select only the segments interesecting the year of DEFID2 survey
var yearlySegments = segCol.map(function(image){
var y = ee.Number(year);
// Mask segments in order to keep the one intersecting the surveyYr and if the surveyYr corresponds to two segments'
// vertex, keep the previous one
var maskedByYear = image.updateMask(image.select(['startYr']).lt(y)
.and(image.select(['endYr']).gte(y)));
return maskedByYear;
})
// Reduce the collection of segments with any reducer since they are univocal per pixel
.qualityMosaic('dur');
var defidLT = yearlySegments;
// Inspect results
print(defidLT);
//////////////////////////////////////////
//// EXPORT YEARLY MAP AS ASSETS DATA ////
//////////////////////////////////////////
// Select desired LT properties and mask with hansen 2000 forest mask
var defidLTExport = defidLT.select(['startYr', 'startVal', 'mag', 'dur', 'rate'])
.updateMask(hansen.select(['treecover2000']));
Export.image.toAsset({
image: defidLTExport,
description: country + '_' + year + '_DEFID2_LT_NBR',
assetId: 'users/agataelia1991/DEFID2_output/LT/' + country + '_' + year + '_DEFID2_LT_NBR',
region: bb,
scale: 30,
crs: 'EPSG:4326',
maxPixels: 1e13});
Export.image.toDrive({
image: defidLTExport.unmask(9999),
description: country + '_' + year + '_DEFID2_LT_NBR_drive',
scale: 30,
crs: 'EPSG:4326',
region: bb,
fileFormat: 'GeoTIFF',
maxPixels: 1e13});
}This code leverages the LandTrendrModified.js module to generate a band stack image representing yearly unmasked pixels count that are used to generate the yearly medoid composites used as input by LandTrendr. This code takes as inputs:
var hansen = ee.Image("users/agataelia1991/Hansen/hansenTreecover2000Mask");
///////////////////////////////////////////////
//// DEFINE COUNTRY, BOUNDING BOX AND DATA ////
///////////////////////////////////////////////
// Define country code
var country = 'RO';
// Retrieve the correspoding bb from the assets
var bbFeature = ee.FeatureCollection('users/agataelia1991/DEFID2_Loic/AOI_by_country/' + country + '_aoi');
var bb = bbFeature.geometry();
// Retrieve the correspoding DEFID2 data from the assets
var defid = ee.FeatureCollection('users/agataelia1991/DEFID2_Loic/DEFID2_by_country/' + country + '_defid_exact_polygons_v1');
////////////////////////////
//// LOAD THE LT MODULE ////
////////////////////////////
// Load the LandTrendr.js module
var ltgee = require('users/agataelia1991/DEFID2:Modules/LandTrendrModified.js');
// Define collection parameters
var startYear = 1995;
var endYear = 2021;
var startDay = '12-01';
var endDay = '01-31';
var aoi = bb;
var maskThese = ['cloud', 'shadow', 'snow', 'water'];
var yearList = ee.List.sequence(startYear, endYear).getInfo().map(function(year){return year.toString()});
//////////////////////////////////////////
//// GET UNMASKED PIXEL COUNT IMAGE //////
//////////////////////////////////////////
// Get unmkased pixels count image
var nClearCollection = ltgee.buildClearPixelCountCollection(startYear, endYear, startDay, endDay, aoi, maskThese);
// Reduce unmasked pixels collection to bands
var nClear = nClearCollection.toBands().rename(yearList);
// Mask unmasked pixels with the DEFID2 dataset, mask with hansen 2000 forest mask and assign 9999 to masked pixels
var nClearMasked = nClear.updateMask(defid.reduceToImage({properties: ['YEAR'], reducer: ee.Reducer.min()}).mask())
.updateMask(hansen.select(['treecover2000']))
.uint16().unmask(9999);
// Export to drive
Export.image.toDrive({
image: nClearMasked,
description: country + '_yearly_clear_pixels_count',
scale: 30,
crs: 'EPSG:4326',
region: bb,
fileFormat: 'GeoTIFF',
maxPixels: 1e13});This code identifies yearly spatial overlappings between DEFID2 records and wind throws and forest fires events, selecting the closest in time prior and following the survey date of each DEFID2 record. Is also summarizes statistics of the overlapping events at the record level. This code takes as inputs:
var modis = ee.ImageCollection("ESA/CCI/FireCCI/5_1"),
forwind = ee.FeatureCollection("users/agataelia1991/FORWIND/FORWIND_v2_clean"),
hansen = ee.Image("users/agataelia1991/Hansen/hansenTreecover2000Mask");
///////////////////////////////////////////////
//// DEFINE COUNTRY, BOUNDING BOX AND DATA ////
///////////////////////////////////////////////
// Define country code
var country = 'RO';
// Retrieve the correspoding bb from the assets
var bbFeature = ee.FeatureCollection('users/agataelia1991/DEFID2_Loic/AOI_by_country/' + country + '_aoi');
var bb = bbFeature.geometry();
// Retrieve the correspoding DEFID2 data from the assets
var defid = ee.FeatureCollection('users/agataelia1991/DEFID2_Loic/DEFID2_by_country/' + country + '_defid_exact_polygons_v1');
///////////////////////////////////////////////////////
//// IDENFITY DISTINCT SURVEY YEARS IN DEFID2 DATA ////
///////////////////////////////////////////////////////
var surveyYears = defid.aggregate_array('YEAR').distinct().getInfo();
/////////////////////////////////////////////////
//// CREATE EMPTY FEATURE COLLECTION TO FILL ////
/////////////////////////////////////////////////
var reducedMedianDifferenceAll = ee.FeatureCollection([]);
/////////////////////////////////
//// DEFINE USEFUL FUNCTIONS ////
/////////////////////////////////
// Function to select only pixels with confidence level higher that 80%
var maskConfidence = function(image){
var confidenceMask = image.select(['ConfidenceLevel']).gt(80);
var maskedImage = image.updateMask(confidenceMask);
return maskedImage;
};
// Define a function that rescale images to EPSG:4326
var rescale30 = function(image){
var rescaled = image.reproject({
crs: 'EPSG:4326',
scale: 30
});
return rescaled;
};
// Define a function that rescale images to desired scale
var rescale = function(image, meter){
var rescaled = image.reproject({
crs: 'EPSG:4326',
scale: meter
});
return rescaled;
};
////////////////////////////////
//// MANAGE MODIS FIRE DATA ////
////////////////////////////////
// Create collection of MODIS yearly products
var modisYearlyCollection = ee.ImageCollection.fromImages(
ee.List.sequence(2001, 2019).map(function(year) {
var modisYearly = modis
.filterBounds(bb)
.map(maskConfidence)
.select('BurnDate')
.filter(ee.Filter.calendarRange(year, year, 'year'))
.reduce(ee.Reducer.count());
var modisYearlyMask = modisYearly.mask();
var modisYearlyFinal = ee.Image.constant(year).int16().updateMask(modisYearlyMask).rename(['YEAR'])
.set('YEAR', year)
.clip(bb);
var modisYearlyFinalRescaled = rescale(modisYearlyFinal, 500);
return modisYearlyFinalRescaled;
}));
////////////////////////////////////////////
//// START LOOPING THROUGH SURVEY YEARS ////
////////////////////////////////////////////
for(var i in surveyYears) {
var year = surveyYears[i];
print(year);
// Select only DEFID2 data with current survey year
var yearlyDefid = defid.filter(ee.Filter.eq('YEAR', year));
// Rasterize DEFID2, reduce to year of survey per pixel
var defidRaster = (yearlyDefid.reduceToImage({ properties: ['YEAR'], reducer: ee.Reducer.min()})).int16().rename(['YEAR']);
// Rescale selected DEFID2 to 30 meters resolution, optional
var defidRaster = rescale(defidRaster, 30);
// Create image collection of raster difference between each yearly layer of EFFIS and FORWIND and the
// survey year of the LandTrendr results, then reduce by minimum (gt or lt 0 according to pre or post difference,
// if difference between years is 0 then it falls under pre) to identify the distance in year (pre and post)
// of the closest fire or windstorm event to the survey year for each LT pixel
var ForwindDifference = (ee.ImageCollection(
ee.List.sequence(2000, 2018).map(function(year){
var yearlyForwind = forwind.filter(ee.Filter.eq('YEAR', year));
var rasterizedYearlyForwind = (yearlyForwind.reduceToImage({properties: ['YEAR'], reducer: ee.Reducer.min()})).int16().rename(['YEAR']);
var rasterizedYearlyForwind500 = rescale(rasterizedYearlyForwind, 500);
var maskedRasterizedYearlyForwind = rasterizedYearlyForwind500.updateMask(defidRaster.select(['YEAR']).mask());
var difference = (defidRaster.select(['YEAR']).updateMask(maskedRasterizedYearlyForwind.select(['YEAR']).mask())).subtract(maskedRasterizedYearlyForwind);
var differencePre = difference.updateMask(difference.gte(0)).abs().rename(['forwind_same_pre']);
var differencePost = difference.updateMask(difference.lt(0)).abs().rename(['forwind_post']);
return ee.Image.cat(differencePre, differencePost);
}))).reduce(ee.Reducer.min());
print(ForwindDifference);
var ModisDifference = (ee.ImageCollection(
ee.List.sequence(2001, 2019).map(function(year){
var yearlyModis = modisYearlyCollection.filter(ee.Filter.eq('YEAR', year)).toBands().rename(['YEAR']);
var maskedYearlyModis = yearlyModis.updateMask(defidRaster.select(['YEAR']).mask());
var difference = (defidRaster.select(['YEAR']).updateMask(maskedYearlyModis.select(['YEAR']).mask())).subtract(maskedYearlyModis);
var differencePre = difference.updateMask(difference.gte(0)).abs().rename(['modis_same_pre']);
var differencePost = difference.updateMask(difference.lt(0)).abs().rename(['modis_post']);
return ee.Image.cat(differencePre, differencePost);
}))).reduce(ee.Reducer.min());
print(ModisDifference);
// Create final layer representing distance in year (pre and post) of the closest fire or windstorm event for each LT pixel
// and mask with hansen 2000 forest mask
var finalDifferenceLayer = ee.Image.cat(ModisDifference, ForwindDifference).updateMask(hansen.select(['treecover2000']))
.uint16();
print (finalDifferenceLayer);
// Calculate weighted mode of distance in years within DEFID2 polygons
var reducedMedianDifference = finalDifferenceLayer.reduceRegions({
collection: yearlyDefid,
reducer: ee.Reducer.median(),
scale: 30
//tileScale: 16
});
// Append created feature collections
var reducedMedianDifferenceAll = reducedMedianDifferenceAll.merge(reducedMedianDifference);
Export.image.toDrive({
image: finalDifferenceLayer.unmask(9999),
description: country + '_' + year + '_DEFID2_fire_wind',
scale: 30,
crs: 'EPSG:4326',
region: bb,
fileFormat: 'GeoTIFF',
maxPixels: 1e13});
}
// Export
Export.table.toDrive({
collection: reducedMedianDifferenceAll,
description: country + '_DEFID2_fire_wind',
selectors: ['PG_ID', 'COUNTRYN', 'YEAR', 'modis_same_pre_min', 'modis_post_min', 'forwind_same_pre_min', 'forwind_post_min']
});This code leverages the LandTrendrModified.js module to generate a band stack image of yearly masked medoid composites of NBR used by LandTrendr and save it as an asset. In a second run it summarizes NBR statistics at DEFID2 record level. This code takes as inputs:
var hansen = ee.Image("users/agataelia1991/Hansen/hansenTreecover2000Mask");
///////////////////////////////////////////////
//// DEFINE COUNTRY, BOUNDING BOX AND DATA ////
///////////////////////////////////////////////
// Define country code
var country = 'RO';
// Retrieve the correspoding bb from the assets
var bbFeature = ee.FeatureCollection('users/agataelia1991/DEFID2_Loic/AOI_by_country/' + country + '_aoi');
var bb = bbFeature.geometry();
// Retrieve the correspoding DEFID2 data from the assets
var defid = ee.FeatureCollection('users/agataelia1991/DEFID2_Loic/DEFID2_by_country/' + country + '_defid_exact_polygons_v1');
// Retrieve the corresponding yearly NBR band stack from the assets
var nbr = ee.Image('users/agataelia1991/DEFID2_output/NBR/' + country + '_yearly_medoid_NBR_stack');
/////////////////////////////////
//// DEFINE USEFUL FUNCTIONS ////
/////////////////////////////////
// Define a function to create new bands names
var newBandsNames = function(startYear, endYear, indexList){
var nIndex = ee.List(indexList).length().getInfo();
var list = ee.List([]);
for (var i = startYear; i <= endYear; i++){
for (var j = 0; j < nIndex; j++){
var newName = indexList[j] + '_' + i.toString();
list = list.add(newName);
}}
return list;
};
///////////////////////////
//// BUILT COLLECTIONS ////
///////////////////////////
// Load the LandTrendr.js module
var ltgee = require('users/agataelia1991/DEFID2:Modules/LandTrendrModified.js');
// Define collection parameters
var startYear = 1995;
var endYear = 2021;
var startDay = '12-01';
var endDay = '01-31';
var aoi = bb;
var index = 'NBR';
var maskThese = ['cloud', 'shadow', 'snow', 'water'];
var bandList = ['NBR'];
// Create list of new bands names
var bandNames = newBandsNames(startYear, endYear, bandList);
//print(bandNames);
// Build annual cloud and cloud shadow masked medoid composites of Landsat surface reflectance TM-equivalent
// bands 1,2,3,4,5,7 over the three months of growing season centered on max NDVI month
var annualSRcollection = ltgee.buildSRcollection(startYear, endYear, startDay, endDay, aoi, maskThese);
print(annualSRcollection);
// Transform the annual surface reflectance bands into a selected indices collection
var indexCollection = ltgee.transformSRcollection(annualSRcollection, bandList);
print(indexCollection);
// Rescale index values to -1-1
var indexCollectionRescaled = indexCollection.map(function(image){
return image.divide(1000)
.set('system:time_start', image.get('system:time_start'));
});
// Transform image collection of yearly indices to a image band stack and mask with hansen 2000 forest mask
var collectionBandStack = indexCollectionRescaled.toBands().rename(bandNames)
.updateMask(hansen.select(['treecover2000']));
print(collectionBandStack);
// Export image to asset before extracting statistics
Export.image.toAsset({
image: collectionBandStack,
description: country + '_yearly_medoid_NBR_stack',
assetId: 'users/agataelia1991/DEFID2_output/NBR/' + country + '_yearly_medoid_NBR_stack',
region: bb,
scale: 30,
crs: 'EPSG:4326',
maxPixels: 1e13});
/////////////////////////////////////////////////////
//// REDUCE AND EXPORT MULTI-TEMPORAL NBR AS CSV ////
/////////////////////////////////////////////////////
// Reduce the images band stack over polygons
var reducedMeanVarianceImage = nbr.reduceRegions({
collection: defid,
reducer: ee.Reducer.mean().combine({reducer2: ee.Reducer.variance(), sharedInputs: true}),
scale: 30,
//tileScale: 16
});
// Export
Export.table.toDrive({
collection: reducedMeanVarianceImage,
description: country + '_DEFID2_index_mean_variance_image',
selectors: ['PG_ID', 'COUNTRYN', 'YEAR', 'NBR_1995_mean', 'NBR_1996_mean', 'NBR_1997_mean', 'NBR_1998_mean', 'NBR_1999_mean', 'NBR_2000_mean',
'NBR_2001_mean', 'NBR_2002_mean', 'NBR_2003_mean', 'NBR_2004_mean', 'NBR_2005_mean', 'NBR_2006_mean', 'NBR_2007_mean',
'NBR_2008_mean', 'NBR_2009_mean', 'NBR_2010_mean', 'NBR_2011_mean', 'NBR_2012_mean', 'NBR_2013_mean', 'NBR_2014_mean',
'NBR_2015_mean', 'NBR_2016_mean', 'NBR_2017_mean', 'NBR_2018_mean', 'NBR_2019_mean', 'NBR_2020_mean', 'NBR_2021_mean',
'NBR_1995_variance', 'NBR_1996_variance', 'NBR_1997_variance', 'NBR_1998_variance', 'NBR_1999_variance', 'NBR_2000_variance',
'NBR_2001_variance', 'NBR_2002_variance', 'NBR_2003_variance', 'NBR_2004_variance', 'NBR_2005_variance', 'NBR_2006_variance', 'NBR_2007_variance',
'NBR_2008_variance', 'NBR_2009_variance', 'NBR_2010_variance', 'NBR_2011_variance', 'NBR_2012_variance', 'NBR_2013_variance', 'NBR_2014_variance',
'NBR_2015_variance', 'NBR_2016_variance', 'NBR_2017_variance', 'NBR_2018_variance', 'NBR_2019_variance', 'NBR_2020_variance', 'NBR_2021_variance']
});This code summarizes yearly the LandTrendr metrics from the outputs saved in the assets at DEFID2 record level. This code takes as inputs:
///////////////////////////////////////////////
//// DEFINE COUNTRY, BOUNDING BOX AND DATA ////
///////////////////////////////////////////////
// Define country code
var country = 'RO';
// Retrieve the correspoding DEFID2 data from the assets
var defid = ee.FeatureCollection('users/agataelia1991/DEFID2_Loic/DEFID2_by_country/' + country + '_defid_exact_polygons_v1');
///////////////////////////////////////////////////////
//// IDENFITY DISTINCT SURVEY YEARS IN DEFID2 DATA ////
///////////////////////////////////////////////////////
var surveyYears = defid.aggregate_array('YEAR').distinct().getInfo();
//////////////////////////////////////////////////
//// CREATE EMPTY FEATURE COLLECTIONS TO FILL ////
//////////////////////////////////////////////////
var reducedModeStartYrAll = ee.FeatureCollection([]);
var reducedMedianOtherLTPropsAll = ee.FeatureCollection([]);
var reducedCountLTAll = ee.FeatureCollection([]);
var reducedCountNegLTAll = ee.FeatureCollection([]);
////////////////////////////////////////////
//// START LOOPING THROUGH SURVEY YEARS ////
////////////////////////////////////////////
for(var i in surveyYears) {
var year = surveyYears[i];
print(year);
// Select only DEFID2 data with current survey year
var yearlyDefid = defid.filter(ee.Filter.eq('YEAR', year));
// Retrieve the correspoding LT output from the assets
var LT = ee.Image('users/agataelia1991/DEFID2_output/LT/' + country + '_' + year + '_DEFID2_LT_NBR');
// Mask the yearly LT outputs to keep only the negative segments
var LTNeg = LT.updateMask(LT.select(['mag']).lt(0));
// Select the different bands of the LT outputs
var startYr = LTNeg.select(['startYr']);
var otherLTProps = LTNeg.select(['startVal', 'mag', 'dur', 'rate']);
////////////////////////////////////////
//// REDUCE AND EXPORT STATS AS CSV ////
////////////////////////////////////////
// Calculate weighted mode of startYr of negative segments within DEFID2 polygons
var reducedModeStartYr = startYr.reduceRegions({
collection: yearlyDefid,
reducer: ee.Reducer.mode(), //(null, null ,100000),
scale: 30
//tileScale: 16
});
// Calculate weighted median of other properties of negative segments within DEFID2 polygons
// setting an high maxRaw in order not to take median from histogram in large polygons
var reducedMedianOtherLTProps = otherLTProps.reduceRegions({
collection: yearlyDefid,
reducer: ee.Reducer.median(), //(null, null ,100000),
scale: 30
//tileScale: 16
});
// Calculate weighted count of pixels within DEFID2 polygons
// by creating a mask and summing number of pixels
var reducedCountLT = LT.select(['startYr']).reduceRegions({
collection: yearlyDefid,
reducer: ee.Reducer.sum(),
scale: 30
//tileScale: 16
});
// Calculate weighted count of negative pixels within DEFID2 polygons
// by creating a mask and summing number of pixels
var reducedCountNegLT = LTNeg.select(['startYr']).reduceRegions({
collection: yearlyDefid,
reducer: ee.Reducer.sum(),
scale: 30
//tileScale: 16
});
// Append created feature collections
var reducedModeStartYrAll = reducedModeStartYrAll.merge(reducedModeStartYr);
var reducedMedianOtherLTPropsAll = reducedMedianOtherLTPropsAll.merge(reducedMedianOtherLTProps);
var reducedCountLTAll = reducedCountLTAll.merge(reducedCountLT);
var reducedCountNegLTAll = reducedCountNegLTAll.merge(reducedCountNegLT);
}
// Export
Export.table.toDrive({
collection: reducedModeStartYrAll,
description: country + '_DEFID2_LT_startYr_mode_negatives',
selectors: ['PG_ID', 'COUNTRYN', 'YEAR', 'mode']
});
// Export
Export.table.toDrive({
collection: reducedMedianOtherLTPropsAll,
description: country + '_DEFID2_LT_other_median_negatives',
selectors: ['PG_ID', 'COUNTRYN', 'YEAR', 'startVal', 'mag', 'dur', 'rate']
});
// Export
Export.table.toDrive({
collection: reducedCountLTAll,
description: country + '_DEFID2_LT_sum',
selectors: ['PG_ID', 'COUNTRYN', 'YEAR', 'sum']
});
// Export
Export.table.toDrive({
collection: reducedCountNegLTAll,
description: country + '_DEFID2_LT_sum_negatives',
selectors: ['PG_ID', 'COUNTRYN', 'YEAR', 'sum']
});DEFID2 data with survey year from 2000 onwards were taken into account for the LandTrendr analysis, in order to guarantee at least five years (from 1995) of stable Landsat imagery prior to the survey year fed to the temporal segmentation algorithm.
LandTrendr segments maps are calculated per survey year of DEFID2 dataset, in order to avoid issues generated by existing spatially overlapping DEFID2 records with different survey years.
The medoid compositing technique of Flood, 2013 https://www.mdpi.com/2072-4292/5/12/6481 was used to generate yearly composites. Since medoid was proved to not be robust against a single outlier in case of less than three observations (Flood, 2013), we provide a code to retrieve a band stack representing for each year the number of clear observations per pixel (meaning per pixel observations that were not masked due to cloud, shadow, snow, water or for not belonging to the desired per pixel greenness months data range). We then suggest to mask out pixels with less than three clear observations per year for a more reliable outcome.
When the start year of the identified segment corresponds with the beginning of the time series of imagery used for running the segmentation, it means that no breakpoint was identified in the time series by the LandTrendr algorithm. In this case, caution with data interpretation is advised.
If the survey year reported in the disturbance record attributes is delayed with respect to forest management interventions (e.g. salvage logging) this can cause the extracted LandTrendr spectral segment to be representing the regrowth phase. When segments with positive magnitude are encountered, it is advised to look into the index full time series of the whole disturbance record.
In general, salvage logging management is a recurrent practice in Europe. Temporal segmentation may not be sufficient in separating the spectral response between the insect disturbance and the salvage logging intervention into two separate segments, especially with a large scale (parameters wise) approach.
FORWIND wind disturbance records represent areas disturbed by wind in the period 2000-2018 (Forzieri et al., 2020 https://essd.copernicus.org/articles/12/257/2020).
FireCCI51: MODIS data represents burned area information in the period 2001-2019 (https://developers.google.com/earth-engine/datasets/catalog/ESA_CCI_FireCCI_5_1#description).
For running LandTrendr a standard set of values for the segmentation parameters was used, referring to Kennedy et al., 2018 https://www.mdpi.com/2072-4292/10/5/691.
No filtering on mmu (minimung mapping unit) or filtering on magnitude, duration and rate values was applied on the resulting LandTrendr segments (Rodman at al., 2021 https://www.sciencedirect.com/science/article/abs/pii/S0034425720306179 Senf & Seidl, 2020 https://www.nature.com/articles/s41893-020-00609-y?proof=t#Sec7)