-
Notifications
You must be signed in to change notification settings - Fork 1
/
wetlands_gee.js
419 lines (346 loc) · 13.9 KB
/
wetlands_gee.js
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
// ------------- IMAGE ACQUISITION ------------- //
// ------------- LANDSAT 7 2002 --------------- //
/**
* Function to mask clouds using the Landsat-7 QA band
* @param {ee.Image} image Landsat-7 image
* @return {ee.Image} cloud masked Landsat-7 image
*/
var cloudMaskL457 = function(image) {
var qa = image.select('QA_PIXEL');
// If the cloud bit (5) is set and the cloud confidence (7) is high
// or the cloud shadow bit is set (3), then it's a bad pixel.
var cloud = qa.bitwiseAnd(1 << 5)
.and(qa.bitwiseAnd(1 << 7))
.or(qa.bitwiseAnd(1 << 3));
// Remove edge pixels that don't occur in all bands
var mask2 = image.mask().reduce(ee.Reducer.min());
return image.updateMask(cloud.not()).updateMask(mask2);
};
var data = ee.ImageCollection("LANDSAT/LE07/C02/T1_TOA")
.filterDate('1999-01-01', '2002-12-31')
.filterBounds(aoi)
.filterMetadata("CLOUD_COVER", "less_than", 10)
.map(cloudMaskL457);
var landsat7_2002 = data.median().clip(aoi);
var Vis = {
bands: ['B3', 'B2', 'B1'],
min: 0.0,
max: 0.4,
gamma: 1.2,
};
Map.centerObject(aoi, 8);
// Map.addLayer(landsat7_2002, Vis, 'Landsat 7 2002');
// Map.addLayer(image02, {min: 0, max: 3, palette: ['pink', 'green', 'red', 'blue']}, 'image02');
// Map.addLayer(image12, {min: 0, max: 3, palette: ['pink', 'green', 'red', 'blue']}, 'image12');
// Map.addLayer(image22, {min: 0, max: 3, palette: ['pink', 'green', 'red', 'blue']}, 'image22');
// ------------- LANDSAT 5 2012 --------------- //
/**
* Function to mask clouds using the Landsat-5 QA band
* @param {ee.Image} image Landsat-5 image
* @return {ee.Image} cloud masked Landsat-5 image
*/
var cloudMaskL457 = function(image) {
var qa = image.select('QA_PIXEL');
// If the cloud bit (5) is set and the cloud confidence (7) is high
// or the cloud shadow bit is set (3), then it's a bad pixel.
var cloud = qa.bitwiseAnd(1 << 5)
.and(qa.bitwiseAnd(1 << 7))
.or(qa.bitwiseAnd(1 << 3));
// Remove edge pixels that don't occur in all bands
var mask2 = image.mask().reduce(ee.Reducer.min());
return image.updateMask(cloud.not()).updateMask(mask2);
};
var dataset = ee.ImageCollection("LANDSAT/LT05/C02/T1_TOA")
.filterDate('2009-01-01', '2012-05-31')
.filterBounds(aoi)
.filterMetadata("CLOUD_COVER", "less_than", 10)
.map(cloudMaskL457);
var landsat5_2012 = dataset.median().clip(aoi);
var Vis = {
bands: ['B3', 'B2', 'B1'],
min: 0.0,
max: 0.4,
gamma: 1.2,
};
Map.centerObject(aoi, 8);
// Map.addLayer(landsat5_2012, Vis, 'Landsat 5 2012');
//------------------ LANDSAT 8 2022 --------------------- //
/**
* Function to mask clouds using the Landsat-8 PIXEL_QA band
* @param {ee.Image} image Landsat-8 image
* @return {ee.Image} cloud masked Landsat-8 image
*/
var getQABits = function(image, start, end, newName) {
// Compute the bits we need to extract.
var pattern = 0;
for (var i = start; i <= end; i++) {
pattern += Math.pow(2, i);
}
// Return a single band image of the extracted QA bits, giving the band
// a new name.
return image.select([0], [newName])
.bitwiseAnd(pattern)
.rightShift(start);
};
// A function to mask out cloudy pixels.
var cloud_shadows = function(image) {
// Select the QA band.
var QA = image.select(['QA_PIXEL']);
// Get the internal_cloud_algorithm_flag bit.
return getQABits(QA, 4,4, 'cloud_shadows').eq(0);
// Return an image masking out cloudy areas.
};
// A function to mask out cloudy pixels.
var clouds = function(image) {
// Select the QA band.
var QA = image.select(['QA_PIXEL']);
// Get the internal_cloud_algorithm_flag bit.
return getQABits(QA, 3,3, 'Cloud').eq(0);
// Return an image masking out cloudy areas.
};
var maskClouds = function(image) {
var cs = cloud_shadows(image);
var c = clouds(image);
image = image.updateMask(cs);
return image.updateMask(c);
};
var visParams = {
bands: ['B4', 'B3', 'B2'],
min:0,
max: 0.3
};
var collection = ee.ImageCollection("LANDSAT/LC08/C02/T1_TOA")
.filterDate('2020-01-01', '2021-05-31')
.filterBounds(aoi)
.map(maskClouds);
var landsat8_2022 = collection.median().clip(aoi);
var visualization = {
bands: ['B4', 'B3', 'B2'],
min: 0.0,
max: 0.3,
};
Map.centerObject(aoi, 8);
// Map.addLayer(landsat8_2022, visualization, 'Landsat 8 2022');
// ---------- Feature extraction ------------ //
// lets extract mndwi, ndvi and evi indices to use in classification
// mndwi = (Green - SWIR)/(Green + SWIR)
var mndwi7_02 = landsat7_2002.normalizedDifference(['B2', 'B7']).rename('NDWI');
var landsat7_2002 = landsat7_2002.addBands(mndwi7_02);
var mndwi5_12 = landsat5_2012.normalizedDifference(['B2', 'B7']).rename('NDWI');
var landsat5_2012 = landsat5_2012.addBands(mndwi5_12);
var mndwi8_22 = landsat8_2022.normalizedDifference(['B3', 'B7']).rename('NDWI');
var landsat8_2022 = landsat8_2022.addBands(mndwi8_22);
// ndvi = (NIR - Red)/(NIR + Red)
var ndvi7_02 = landsat7_2002.normalizedDifference(['B3', 'B4']).rename('NDVI');
var landsat7_2002 = landsat7_2002.addBands(ndvi7_02);
var ndvi5_12 = landsat5_2012.normalizedDifference(['B3', 'B4']).rename('NDVI');
var landsat5_2012 = landsat5_2012.addBands(ndvi5_12);
var ndvi8_22 = landsat8_2022.normalizedDifference(['B4', 'B5']).rename('NDVI');
var landsat8_2022 = landsat8_2022.addBands(ndvi8_22);
// EVI = G * {(NIR - Red)/(NIR + C1 * Red - C2 * Blue + L)}
// evi landsat 7 2002
var EVI02 = landsat7_2002.expression(
'2.5 * ((NIR - RED) / (NIR + 6 * RED - 7.5 * BLUE + 1))', {
'NIR': landsat7_2002.select('B4'),
'RED': landsat7_2002.select('B3'),
'BLUE': landsat7_2002.select('B1')
}).rename("EVI")
var landsat7_2002 = landsat7_2002.addBands(EVI02)
// print('Landsat 7 2002', landsat7_2002)
// evi landsat 5 2012
var EVI12 = landsat5_2012.expression(
'2.5 * ((NIR - RED) / (NIR + 6 * RED - 7.5 * BLUE + 1))', {
'NIR': landsat5_2012.select('B4'),
'RED': landsat5_2012.select('B3'),
'BLUE': landsat5_2012.select('B1')
}).rename("EVI")
var landsat5_2012 = landsat5_2012.addBands(EVI12)
// print('Landsat 5 2012', landsat5_2012)
// evi landsat 8 2022
var EVI22 = landsat8_2022.expression(
'2.5 * ((NIR - RED) / (NIR + 6 * RED - 7.5 * BLUE + 1))', {
'NIR': landsat8_2022.select('B5'),
'RED': landsat8_2022.select('B4'),
'BLUE': landsat8_2022.select('B2')
}).rename("EVI")
var landsat8_2022 = landsat8_2022.addBands(EVI22)
// print('Landsat 8 2022', landsat8_2022)
// ------------- Deriving training points ---------------//
// merge classes to feature collection
var FCmerged = wetlands.merge(vegetation).merge(urban).merge(water);
//Stratified Random Sampling:
var FCimage = ee.Image().byte().paint(FCmerged, "class").rename("LC")
// print(FCimage, 'image');
var stratifiedsample = FCimage.stratifiedSample({
numPoints:10000,
classBand:'LC',
region:aoi,
scale:30,
classValues:[0,1,2,3],
classPoints:[1500,1500,1500,1500],
geometries:true
});
print('Stratified samples', stratifiedsample);
print(stratifiedsample.reduceColumns(ee.Reducer.frequencyHistogram(),['LC']).get('histogram','No of points'));
Map.addLayer(stratifiedsample, {}, 'Stratified Samples', false);
// ------------------ Extract pixel values ---------------------//
var stratifiedTraining02 = landsat7_2002.sampleRegions({
tileScale: 3,
collection: stratifiedsample,
properties:['LC'],
geometries: true,
scale:30,
});
var stratifiedTraining12 = landsat5_2012.sampleRegions({
tileScale: 3,
collection: stratifiedsample,
properties:['LC'],
geometries: true,
scale:30,
});
var stratifiedTraining22 = landsat8_2022.sampleRegions({
tileScale: 3,
collection: stratifiedsample,
properties:['LC'],
geometries: true,
scale:30,
});
// ------------- Split Training and Validation data --------------- //
// add a random column
var stratifiedTraining22 = stratifiedTraining22.randomColumn();
// split your **full** sample into training and validation points to keep them independent of each other
var trainingSample = stratifiedTraining22.filter('random <= 0.8');
var validationSample = stratifiedTraining22.filter('random > 0.8');
// --------------- Train Random Forest Classifier ----------------//
// select bands for classification
var bands57 = ['B1', 'B2', 'B3', 'B4', 'NDWI', 'NDVI', 'EVI']
var bands8 = ['B2', 'B3', 'B4', 'B5', 'NDWI', 'NDVI', 'EVI']
// RF classifier
var RFclassifier02 = ee.Classifier.smileRandomForest(100).train({
features: trainingSample,
classProperty: 'LC',
inputProperties: bands57,
});
var RFclassifier12 = ee.Classifier.smileRandomForest(100).train({
features: trainingSample,
classProperty: 'LC',
inputProperties: bands57,
});
var RFclassifier22 = ee.Classifier.smileRandomForest(100).train({
features: trainingSample,
classProperty: 'LC',
inputProperties: bands8,
});
// --------------------- Classify the images -----------------------//
// landsat 7 2002 classification
var classifiedImage2002 = landsat5_2002.classify(RFclassifier02);
//confusion matrix about the resubstitution accuracy
var trainAccuracy_2002 = RFclassifier02.confusionMatrix()
print('2002 Resubstitution error matrix: ', trainAccuracy_2002)
print('2002 Training overall accuracy: ', trainAccuracy_2002.accuracy())
print('2002 Training Kappa index:', trainAccuracy_2002.kappa())
// landsat 5 2012 classification
var classifiedImage2012 = landsat5_2012.classify(RFclassifier12);
// confusion matrix about the resubstitution accuracy
var trainAccuracy_2012 = RFclassifier12.confusionMatrix()
print('2012 Resubstitution error matrix: ', trainAccuracy_2012)
print('2012 Training overall accuracy: ', trainAccuracy_2012.accuracy())
print('2012 Training Kappa index:', trainAccuracy_2012.kappa())
// landsat 8 2022 classification
var classifiedImage2022 = landsat8_2022.classify(RFclassifier22);
var trainAccuracy_2022 = RFclassifier22.confusionMatrix()
print('2022 Resubstitution error matrix: ', trainAccuracy_2022)
print('2022 Training overall accuracy: ', trainAccuracy_2022.accuracy())
print('2022 Training Kappa index:', trainAccuracy_2022.kappa())
Map.addLayer(classifiedImage2002, {min: 0, max: 3, palette: ['pink', 'green', 'red', 'blue']}, 'Landsat 7 2002 Classified');
Map.addLayer(classifiedImage2012, {min: 0, max: 3, palette: ['pink', 'green', 'red', 'blue']}, 'Landsat 5 2012 Classified');
Map.addLayer(classifiedImage2022, {min: 0, max: 3, palette: ['pink', 'green', 'red', 'blue']}, 'Landsat 8 2022 Classified');
// -------------- Classify validation set ------------------//
// classify the validation sample
var validation = validationSample.classify(RFclassifier22)
// Get a confusion matrix representing expected accuracy of **test data (validation sample)**
var testAccuracy = validation.errorMatrix('LC', 'classification')
// print(testAccuracy)
// print('2022 Testing accuracy', testAccuracy.accuracy())
// --------------- EXPORT TO ASSET ------------------ //
// Export to asset
Export.image.toAsset({
image: classifiedImage2022.clip(aoi),
description: 'classified2022',
assetId: 'projects/ee-kimeu11/assets/classified2022', // <> modify these
region: aoi,
scale: 30,
crs: 'EPSG:4326',
});
// ------------------------- LEGEND --------------------------//
// Class color and label info.
var classInfo = [
{name: 'Wetlands', color: 'pink'},
{name: 'Vegetation', color: 'green'},
{name: 'Urban/Bare', color: 'red'},
{name: 'Water', color: 'blue'}
];
// Makes a legend entry: color and label side-by-side in a panel.
function legendEntry(info) {
var color = ui.Panel({style: {
width: '20px',
height: '20px',
backgroundColor: info.color,
margin: '6px 0px 0px 0px'
}});
var label = ui.Label({
value: info.name,
});
return ui.Panel({
widgets: [color, label],
layout: ui.Panel.Layout.flow('horizontal'),
style: {
stretch: 'horizontal',
margin: '-6px 0px 0px 0px'
}});
}
// Define a panel to hold all legend entries.
var legend = ui.Panel({
style: {
position: 'top-left',
padding: '8px 8px 0px 8px'
}
});
// Legend title
// Create legend title
var legendTitle = ui.Label({
value: '2012',
style: {
fontWeight: 'bold',
fontSize: '18px',
margin: '0 0 4px 0',
padding: '0'
}});
// Add the title to the panel
legend.add(legendTitle);
// Loop through the map classes, add each entry to the legend panel.
for (var i = 0; i < classInfo.length; i++) {
legend.add(legendEntry(classInfo[i]));
}
// Show legend on the map.
Map.add(legend);
// ---------------------- CHANGE DETECTION-------------------- //
// ---------------------- Image Difference -------------------- //
// compute image difference between 2002 and 2012
var difference02_12 = image12.subtract(image02);
// print(image02)
// print(difference02_12)
var reclassified = difference02_12.expression('b(0) < 0 ? -1 : b(0) == 0 ? 0 : 1');
Map.addLayer(reclassified.clip(aoi), {min: -1, max: 1, palette: ['#0096a0', '#ffbb22', '#fa0000']}, 'image difference');
// compute image difference between 2012 and 2022
var difference12_22 = image22.subtract(image12);
// print(image02)
print(difference12_22)
var reclassified = difference12_22.expression('b(0) < 0 ? -1 : b(0) == 0 ? 0 : 1');
Map.addLayer(reclassified.clip(aoi), {min: -1, max: 1, palette: ['#0096a0', '#ffbb22', '#fa0000']}, 'image difference');
// compute image difference between 2012 and 2022
var difference02_22 = image22.subtract(image02);
// print(image02)
print(difference02_22)
var reclassified = difference02_22.expression('b(0) < 0 ? -1 : b(0) == 0 ? 0 : 1');
Map.addLayer(reclassified.clip(aoi), {min: -1, max: 1, palette: ['#0096a0', '#ffbb22', '#fa0000']}, 'image difference');