Calculating NDWI

Calculating area changes using NDWI

In 2020, Lake Baringo, in Kenya, burst into the international limelight. Like bees let loose from a beehive, a multitude of stories were made about the lake’s phenomenal rise in water levels. All manner of theories deciphering the cause was laid on the table. They ranged from climate change and geologic faults to sedimentation. Pictures and documentaries of the lake’s devastation to the surrounding schools, hotels and houses were portrayed to all and sundry. It seemed Lake Bogoria was swallowing everything in its swath, and from locals, its expansion had superseded all historical limits. Generally, most Rift Valley lakes, including those in Ethiopia, had experienced an increase in their extent.

using NDWI
Source: Al Jazeera

Of course, it is terrifying to wake up and find water levels knocking at your doorstep unannounced and unwelcomed. More so if they show no signs of receding any time soon. According to the government report, Rising water levels in Kenya’s Rift Valley Lakes: Turkwel Gorge Dam and Lake Victoria 379, 935 persons have been affected across 13 counties. Surprisingly, even though Lake Baringo featured prominently in most news reports, the largest rise in water extent was from Turkwel Dam. It recorded a 158% change. Lake Victoria takes the crown with the largest increase in terms of area –542km^2. Lake Baringo ranks second in area extent increase in terms of percentage and third in terms of km^2.

Three satellite images were used in the analysis. These were for 2010, 2014 and 2020.

print.data.frame(rift_lakes)
##    S/N     Lake/Dam Toposheet\r\n(km2) 2010 Extent\r\n(km2)
## 1    1      Turkana            7485.48              7947.85
## 2    2       Logipi             129.88               183.46
## 3    3      Baringo 128.08000000000001               159.49
## 4    4      Bogoria              34.35                34.79
## 5    5        Solai               5.74                 9.54
## 6    6       Nakuru              40.04                45.01
## 7    7 Ol’ Bolossat               18.2                39.81
## 8    8   Elmenteita              21.24                21.54
## 9    9     Naivasha             135.32               154.16
## 10  10       Magadi             141.09               148.09
## 11  11     Victoria            3971.37              4030.22
## 12  12  Turkwel Dam              18.97                22.97
## 13  NA         <NA>              Total             12129.76
##    2014 Extent\r\n(km2) 2020 Extent\r\n(km2) Change\r\n(2020–2010) % Change
## 1               8064.09              8265.07                317.22     3.99
## 2                178.43               190.93                  7.47     4.07
## 3                221.77               268.06                108.57    68.07
## 4                 40.50                43.25                  8.46    24.32
## 5                 13.00                13.61                  4.07    42.66
## 6                 56.82                68.18                 23.17    51.48
## 7                 50.39                52.16                 12.35    31.02
## 8                 22.78                22.97                  1.43     6.64
## 9                164.86               193.48                 39.32    25.51
## 10               150.90               153.30                  5.21     3.52
## 11              4421.32              4572.70                542.48    13.46
## 12                31.68                59.27                 36.30   158.03
## 13             12796.93             13416.54              13902.98  1106.05
##    Change (2020–Normal) % Change from Normal
## 1                779.59                10.41
## 2                 61.05                47.00
## 3                139.98               109.29
## 4                  8.90                25.91
## 5                  7.87               137.11
## 6                 28.14                70.28
## 7                 33.96               186.59
## 8                  1.73                 8.15
## 9                 58.16                42.98
## 10                12.21                 8.65
## 11               601.33                15.14
## 12                40.30               212.44
## 13                 8.64              1773.22

So how can you measure the lake area using satellite images? One way is to draw a line around the lake edges as a vector polygon and automatically get the area using inbuilt tools in Qgis or ArcGIS. However, that is prone to error. You may be biased on edges that appear clear-cut compared to those interwoven with lacustrine vegetation, mud or indiscernible at all! A more sophisticated way is to use the Normalized Difference Water Index (NDWI). First proposed by Gao (1996) for measuring vegetation liquid water from space, NDWI is a remote-sensing derived index for monitoring changes in water content in leaves and with some tweaking, it also monitors water content in water bodies.

In this article, NDWI was used to derive the lake extent of Lake Baringo. The methodology for using the NDWI formula was largely borrowed from the Medium blog by Fransisco J. Penal. The NDWI formula is as follows:

NDWI = (Green – NIR)/(Green + NIR)

It should be noted that other bands, such as SWIR have also been used to calculate NDWI. According to McFeeters (1996) the advantages of using the Green and NIR bands above are as follows:

  1. maximize the typical reflectance of water features by using green light wavelengths
  2. minimize the low reflectance of NIR by water features; and
  3. take advantage of the high reflectance of NIR by terrestrial vegetation and soil features.

When the above equation is used, water features have positive values. Soil and terrestrial values on the other end of the spectrum have zero to negative values owing to their typically higher reflectance of NIR than Green light.

Having taken enough of all this, it is time to wade into some exercise of NDWI.

Calculating the area extent of Lake Baringo in 2010, 2014 and 2020

As always, an area of interest (AOI) has to be defined. Here it is:

var baringoArea = 
    /* color: #d63000 */
    /* shown: false */
    /* displayProperties: [
      {
        "type": "rectangle"
      }
    ] */
    ee.Geometry.Polygon(
        [[[36.00121112109377, 0.7351363293842054],
          [36.00121112109377, 0.5010039974835908],
          [36.174589111816424, 0.5010039974835908],
          [36.174589111816424, 0.7351363293842054]]], null, false);

Thereafter, we get the coordinates of Lake Baringo. Actually getting the coordinates comes first because you have to know the location of the lake before drawing the AOI. However, because geometries appear at the top of the GEE code editor, it is also logical that our coding steps will proceed in the same manner they appear in the code editor.

// get the coordinates of lake baringo

Map.setCenter( 36.080862,  0.611892, 12);

// years to study from scoping report - 2010, 2014, 2020

First, we will load the satellite images for Lake Baringo. Thereafter we will dive deeper by performing NDWI operations on them. For all three years, Surface Reflectance (SR) images were used. This is because SR images have effects of atmospheric scattering already removed, otherwise, they would affect the NDWI values.

We will start with 2010, in tandem with the first year of analysis in the government scoping report.

Step 1: Acquire satellite images

// landsat image for 2010 - clearest from december to february for L5
var baringo2010 = ee.ImageCollection('LANDSAT/LT05/C02/T1_L2')
    .filterDate('2010-01-01', '2010-02-28') // driest months for 2010
    .filterBounds(baringoArea)   // the region we want our sat images to intersect
    .sort('CLOUD_COVER') 
   
print("Metdata for L. Baringo 2010: ", baringo2010);
var baringo2010Image = ee.Image('LANDSAT/LT05/C02/T1_L2/LT05_169060_20100130') .clip(baringoArea);  // clip to our aoi(baringoArea)
Map.addLayer(baringo2010Image, {bands: ['SR_B3', 'SR_B2', 'SR_B1'], min: 8000, max: 15000},
 'Lake Baringo 2010');
print('Metadata for selected Image for L. Baringo 2010: ', baringo2010Image);

We have used an image from Landsat 5. Landsat 5 was one of the Landsat series in operation in 2010.

Why not Landsat 7? You may ask.

Good question.

It should be known that Landsat 7 suffered a Scan Line Error in 2003. This resulted in the satellite images having line gaps that result in a 22% image loss. For this reason, Landsat 5 was preferred.

Another thing that may hardly escape a hawk-eyed GEE user was, why select an image using ee.Image rather than choose the least cloudy image from the ee.ImageCollection method? Two criteria had to be met before using any image over Lake Baringo.

  1. It must be the least cloudy image
  2. It must cover the entire lake for proper analysis to be done.

The least cloudy images for 2010 and 2014 only covered 3/4 of the lake. Even though criterion 1 had been met, this necessitated choosing the second-least cloudy images which met both criteria with substantial clarity over Lake Baringo.

Having fished out the best image for analysis, the above code script for 2010 image acquisition could be replicated for the remaining two years. The only change needed was changing the satellite type. For the 2014 image, a Sentinel-2 image was used.

Below is the code script for acquiring images for 2014 and 2020.

// landsat image for 2014 - clearest from december to february for L8
var baringo2014 = ee.ImageCollection('LANDSAT/LC08/C02/T1_L2')
    .filterDate('2014-01-01', '2014-02-28') // driest months for 2014
    .filterBounds(baringoArea)   // the region we want our sat images to intersect
    .sort('CLOUD_COVER')
   
print("Metdata for L. Baringo 2014: ", baringo2014);

var baringo2014Image = ee.Image('LANDSAT/LC08/C02/T1_L2/LC08_169060_20140125')
    .clip(baringoArea);
  
Map.addLayer(baringo2014Image, {bands: ['SR_B4', 'SR_B3', 'SR_B2'], min: 7000.0, max: 13000},
'Lake Baringo 2014');

print('Metadata for selected Image for L. Baringo 2014: ', baringo2014Image);

// sentinel2 image for 2020 - clearest from december to february for S2
var baringo2020 = ee.ImageCollection('COPERNICUS/S2_SR')
    .filterDate('2020-01-01', '2020-02-28') // driest months for 2020
    .filterBounds(baringoArea)   // the region we want our sat images to intersect
    .sort('CLOUDY_PIXEL_PERCENTAGE')
    .first();
   
var baringo2020Image = baringo2020.clip(baringoArea);  // no need to select an image here from list since it cover our aoi perfectly
   
Map.addLayer(baringo2020Image, {bands: ['B4', 'B3', 'B2'], min: 200, max: 1900},
'Lake Baringo 2020');

print('Metadata for selected Image for L. Baringo 2020: ', baringo2020Image);

Having done all the above, Phase I is closed.

Now Phase II is where the rubber meets the road. Phase II involves calculating the NDWI and subsequently calculating the lake area from this index. For each year both procedures of NDWI and area calculation will be done back-to-back. A friendly note. Once both operations of NDWI calculation and area estimation are done for 2010, replicating the code for the remaining two years is as easy as ABC. Follow on to find out.

Calculate the NDWI

Below is the code script for calculating the NDWI for the Landsat-5 2010 image.

// ndwi for 2010 L5 image

// ndwi = (Green - NIR)/(Green + NIR)
var green2010 = 'SR_B2';
var nir2010 = 'SR_B4';

var ndwi2010 = baringo2010Image.normalizedDifference([green2010, nir2010])
    .rename('ndwi2010');

print("NDWI 2010", ndwi2010);

var minNDWI2010 = ndwi2010.reduceRegion({
  reducer: ee.Reducer.min(),
  scale: 30,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Min NDWI of 2010: ', minNDWI2010);

var maxNDWI2010 = ndwi2010.reduceRegion({
  reducer: ee.Reducer.max(),
  scale: 30,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Max NDWI of 2010: ', maxNDWI2010);

Map.addLayer(ndwi2010, {palette: ['red', 'yellow', 'green', 'blue'], min: -0.46, max: 0.08}, 'NDWI 2010 image');

The purpose of calculating the NDWI min and max values was two-fold. One, it established the NDWI range. Two, the NDWI min-max values were used in creating a well-balanced NDWI map with the water features enhanced. See Map.addLayer for this.

NDWI Image for 2010
NDWI image for 2010

Calculating the area of the lake

Calculating the area is not so straightforward. First, you have to derive the water values, which according to McFrees(1996) range from 0 to 1. Therefore, in this case, and also borrowing from Fransisca’s blog, values above 0 were masked out. The unmasked features left behind belong to the water feature –the lake. It is from this unmasked water feature that the lake’s area will be calculated.

The script to perform the above tasks is shown in the code script below.

Mask-out water features

Here is the script for masking out non-water features.

// create ndwi water mask for 2010
var ndwi2010Threshold = ndwi2010.gte(0.0);
var ndwi2010Mask = ndwi2010Threshold.updateMask(ndwi2010Threshold);
//Map.addLayer(ndwi2010Threshold, {palette: ['black', 'white']}, 'NDWI 2010 binary mask');
Map.addLayer(ndwi2010Mask, {palette: ['blue']}, 'NDWI 2010 Mask');

The meaning of the method updateMask confuses the hoi polloi. Let’s unmask it. Think of updateMask as a method of removing masked-out pixels, with masked-out referring to those pixels that do not meet the zero-to-one criterion. In other words, updateMask removes all pixels ranging from 0 to -1 while values above the 0 thresholds are saved as ndwi2010Mask. The updateMask method does the same job as any other variable.

NDWI
L5 image 2010
NDWI
Mask overlay over Lake Baringo in the L5 RGB image

Calculate the area of the water mask image

If the masked area fits over the lake perfectly, the area will be a close estimate, if not the true estimate, of the lake size. At this point, it will be a good protocol to check and uncheck the 2010 L5 image and NDWI 2010 Mask image to check if they line up perfectly. In the case of this exercise, the mask fit the underlying lake of the RGB image perfectly. Yay!

Below is the script for calculating the area of the water mask – ndwi2010Mask.

// get area of ndwi 2010 in km2 which will be size of lake baringo 2010

var ndwi2010Area = ndwi2010Mask.multiply(ee.Image.pixelArea());

var ndwi2010Stats = ndwi2010Area.reduceRegion({
  reducer: ee.Reducer.sum(),
  scale: 30,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Area of L. Baringo in 2010: ', ndwi2010Stats);

var ndwi2010Final = ee.Number(ndwi2010Stats.get('ndwi2010')).divide(1e6);
print('Area of L. Baringo in Km2 in 2010: ', ndwi2010Final);

The area of the mask comes to 129.67 km^2. The government report puts the area of the same lake at 159.49km^2.

Once the area of Lake Baringo from the 2010 image has been successfully estimated, the same process can be replicated for the two images of 2014 and 2020. It will simply involve replacing the Green and NIR bands with the names of the two respective satellite images. For example, Green and NIR bands for Landsat 7 are SR_B3 and SR_B5 while those for Sentinel-2 are B3 and B8. Don’t confuse them! Additionally, ensure for each different satellite image, the dates for the NDWI thresholds and masks are updated.

The code script below shows how the lake area extent for the 2014 and 2020 images was calculated.

// ndwi for 2014 image L8 - baringo2014Image

// ndwi = (Green - NIR)/(Green + NIR)
// note that green band in L8 is 'SR_B3' while NIR is 'SR_B5'

var green2014 = 'SR_B3';
var nir2014 = 'SR_B5';

var ndwi2014 = baringo2014Image.normalizedDifference([green2014, nir2014])
    .rename('ndwi2014');

print("NDWI 2014", ndwi2014);

var minNDWI2014 = ndwi2014.reduceRegion({
  reducer: ee.Reducer.min(),
  scale: 30,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Min NDWI of 2014: ', minNDWI2014);

var maxNDWI2014 = ndwi2014.reduceRegion({
  reducer: ee.Reducer.max(),
  scale: 30,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Max NDWI of 2014: ', maxNDWI2014);

Map.addLayer(ndwi2014, {palette: ['red', 'yellow', 'green', 'blue'], min: -0.44, max: 0.16}, 'NDWI 2014 image');

// get area of L. Baringo in 2014

// create ndwi water mask for 2014
var ndwi2014Threshold = ndwi2014.gte(0.0);
var ndwi2014Mask = ndwi2014Threshold.updateMask(ndwi2014Threshold);
//Map.addLayer(ndwi2014Threshold, {palette: ['black', 'white']}, 'NDWI 2014 binary mask');
Map.addLayer(ndwi2014Mask, {palette: ['blue']}, 'NDWI 2014 Mask');

// get area of ndwi 2014 in km2 which will be size of lake baringo 2014

var ndwi2014Area = ndwi2014Mask.multiply(ee.Image.pixelArea());

var ndwi2014Stats = ndwi2014Area.reduceRegion({
  reducer: ee.Reducer.sum(),
  scale: 30,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Area of L. Baringo in 2014: ', ndwi2014Stats);

var ndwi2014Final = ee.Number(ndwi2014Stats.get('ndwi2014')).divide(1e6);
print('Area of L. Baringo in Km2 in 2014: ', ndwi2014Final);

///////////////////////////////////////////////
///////////////////////////////////////////////

// ndwi for 2020 image S2 - baringo2020Image

// ndwi = (Green - NIR)/(Green + NIR)
// note that green band in S2 is 'B3' while NIR is 'B8'

var green2020 = 'B3';
var nir2020 = 'B8';

var ndwi2020 = baringo2020Image.normalizedDifference([green2020, nir2020])
    .rename('ndwi2020');

print("NDWI 2020", ndwi2020);

var minNDWI2020 = ndwi2020.reduceRegion({
  reducer: ee.Reducer.min(),
  scale: 10,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Min NDWI of 2020: ', minNDWI2020);

var maxNDWI2020 = ndwi2020.reduceRegion({
  reducer: ee.Reducer.max(),
  scale: 10,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Max NDWI of 2020: ', maxNDWI2020);

Map.addLayer(ndwi2020, {palette: ['red', 'yellow', 'green', 'blue'], min: -0.81, max: 0.66}, 'NDWI 2020 image');

// get area of L. Baringo in 2020

// create ndwi water mask for 2020
var ndwi2020Threshold = ndwi2020.gte(0.0);
var ndwi2020Mask = ndwi2020Threshold.updateMask(ndwi2020Threshold);
//Map.addLayer(ndwi2014Threshold, {palette: ['black', 'white']}, 'NDWI 2020 binary mask');
Map.addLayer(ndwi2020Mask, {palette: ['blue']}, 'NDWI 2020 Mask');

// get area of ndwi 2020 in km2 which will be size of lake baringo 2020

var ndwi2020Area = ndwi2020Mask.multiply(ee.Image.pixelArea());

var ndwi2020Stats = ndwi2020Area.reduceRegion({
  reducer: ee.Reducer.sum(),
  scale: 10,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Area of L. Baringo in 2020: ', ndwi2020Stats);

var ndwi2020Final = ee.Number(ndwi2020Stats.get('ndwi2020')).divide(1e6);
print('Area of L. Baringo in Km2 in 2020: ', ndwi2020Final);

The estimated lake area for Lake Baringo in 2014 was 189.22km^2 while in 2020 it was 189.23km^2. The government report puts the area extent for the two years at 221.77km^2 and 268.06km^2 respectively.

Related Article: GIS Applications: Forest Mapping in Kenya

Create a chart visualizing the change in the lake area

Though not that important, one can create a GEE chart to visually show the changes in the derived lake area extents.

// create a datatable for the area of Lake Baringo in 2010, 2014 and 2020 in km2
//

var dataTable = [
  [
    {label: 'Year', role: 'domain', type: 'string'},
    {label: 'Area_km2', role: 'data', type: 'number'},
    {label: 'Satellite Used', role: 'domain', type: 'string'}
  ],
  ['2010', 129.67, 'Landsat 5'],
  ['2014', 189.22, 'Landsat 8'],
  ['2020', 189.23, 'Sentinel 2'],
];

// Define the chart and print it to the console.
var chart = ui.Chart(dataTable).setChartType('ColumnChart').setOptions({
  title: 'Variation of area of Lake Baringo in km2 2010, 2014, 2020',
  legend: {position: 'none'},
  hAxis: {title: 'Year', titleTextStyle: {italic: false, bold: true}},
  vAxis: {title: 'Area_km2', titleTextStyle: {italic: false, bold: true}},
  colors: [ '#800080']
});


print(chart);
NDWI Charts
Area changes of Lake Baringo in 2010, 2014 and 2020

Conclusion

The advantage of GEE can never be overemphasized. Where expensive field methods would have been used in estimating the change in lake area extent, GEE does so with a few repetitive lines of code. Furthermore, the derived indices can be subjected to further assessment, or the same images used for a different analysis altogether. For example, the satellite images used can serve as raw data for studying the algae content of the lake over the three periods. GEE continues to be a useful tool to scientists and ecologists alike, allowing the user to perform tasks just as well as those contracted by large organizations, including governments. Finally, the data collected from GEE exercises during the preliminary study can serve as a good determinant of whether in-situ fieldwork is necessary or not. In our case, it is necessary since there was almost a double increase in lake area size from 2000 – 2014.

Full code script for lake area extent estimation – 2010, 2014 and 2020

// get the coordinates of lake baringo

Map.setCenter( 36.080862,  0.611892, 12);

// years to study from scoping report - 2010, 2014, 2020

// let's begin with year 2010
// driest months are december to february
// take clearest image from between these three months

// landsat image for 2010 - clearest from december to february for L5
var baringo2010 = ee.ImageCollection('LANDSAT/LT05/C02/T1_L2')
    .filterDate('2010-01-01', '2010-02-28') // driest months for 2010
    .filterBounds(baringoArea)   // the region we want our sat images to intersect
    .sort('CLOUD_COVER')
    
    
print("Metdata for L. Baringo 2010: ", baringo2010);

var baringo2010Image = ee.Image('LANDSAT/LT05/C02/T1_L2/LT05_169060_20100130')
    .clip(baringoArea);  // clip to our aoi(baringoArea)

Map.addLayer(baringo2010Image, {bands: ['SR_B3', 'SR_B2', 'SR_B1'], min: 8000, max: 15000}, 
 'Lake Baringo 2010');
 
print('Metadata for selected Image for L. Baringo 2010: ', baringo2010Image);

// landsat image for 2014 - clearest from december to february for L8
var baringo2014 = ee.ImageCollection('LANDSAT/LC08/C02/T1_L2')
    .filterDate('2014-01-01', '2014-02-28') // driest months for 2014
    .filterBounds(baringoArea)   // the region we want our sat images to intersect
    .sort('CLOUD_COVER')
    
print("Metdata for L. Baringo 2014: ", baringo2014);

var baringo2014Image = ee.Image('LANDSAT/LC08/C02/T1_L2/LC08_169060_20140125')
    .clip(baringoArea);
   
Map.addLayer(baringo2014Image, {bands: ['SR_B4', 'SR_B3', 'SR_B2'], min: 7000.0, max: 13000}, 
'Lake Baringo 2014');

print('Metadata for selected Image for L. Baringo 2014: ', baringo2014Image);

// sentinel2 image for 2020 - clearest from december to february for S2
var baringo2020 = ee.ImageCollection('COPERNICUS/S2_SR')
    .filterDate('2020-01-01', '2020-02-28') // driest months for 2020
    .filterBounds(baringoArea)   // the region we want our sat images to intersect
    .sort('CLOUDY_PIXEL_PERCENTAGE')
    .first();
    
var baringo2020Image = baringo2020.clip(baringoArea);  // no need to select an image here from list since it cover our aoi perfectly
    
Map.addLayer(baringo2020Image, {bands: ['B4', 'B3', 'B2'], min: 200, max: 1900}, 
'Lake Baringo 2020');

print('Metadata for selected Image for L. Baringo 2020: ', baringo2020Image);


// Phase 2 - Calculate ndwi for each of the three years.

// ndwi for 2010 L5 image

// ndwi = (Green - NIR)/(Green + NIR)
var green2010 = 'SR_B2';
var nir2010 = 'SR_B4';

var ndwi2010 = baringo2010Image.normalizedDifference([green2010, nir2010])
    .rename('ndwi2010');

print("NDWI 2010", ndwi2010);

var minNDWI2010 = ndwi2010.reduceRegion({
  reducer: ee.Reducer.min(),
  scale: 30,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Min NDWI of 2010: ', minNDWI2010);

var maxNDWI2010 = ndwi2010.reduceRegion({
  reducer: ee.Reducer.max(),
  scale: 30,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Max NDWI of 2010: ', maxNDWI2010);

Map.addLayer(ndwi2010, {palette: ['red', 'yellow', 'green', 'blue'], min: -0.46, max: 0.08}, 'NDWI 2010 image');

// create ndwi water mask for 2010
var ndwi2010Threshold = ndwi2010.gte(0.0);
var ndwi2010Mask = ndwi2010Threshold.updateMask(ndwi2010Threshold);
//Map.addLayer(ndwi2010Threshold, {palette: ['black', 'white']}, 'NDWI 2010 binary mask');
Map.addLayer(ndwi2010Mask, {palette: ['blue']}, 'NDWI 2010 Mask');

// get area of ndwi 2010 in km2 which will be size of lake baringo 2010

var ndwi2010Area = ndwi2010Mask.multiply(ee.Image.pixelArea());

var ndwi2010Stats = ndwi2010Area.reduceRegion({
  reducer: ee.Reducer.sum(),
  scale: 30,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Area of L. Baringo in 2010: ', ndwi2010Stats);

var ndwi2010Final = ee.Number(ndwi2010Stats.get('ndwi2010')).divide(1e6);
print('Area of L. Baringo in Km2 in 2010: ', ndwi2010Final);


// ndwi for 2014 image L8 - baringo2014Image

// ndwi = (Green - NIR)/(Green + NIR)
// note that green band in L8 is 'SR_B3' while NIR is 'SR_B5'

var green2014 = 'SR_B3';
var nir2014 = 'SR_B5';

var ndwi2014 = baringo2014Image.normalizedDifference([green2014, nir2014])
    .rename('ndwi2014');

print("NDWI 2014", ndwi2014);

var minNDWI2014 = ndwi2014.reduceRegion({
  reducer: ee.Reducer.min(),
  scale: 30,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Min NDWI of 2014: ', minNDWI2014);

var maxNDWI2014 = ndwi2014.reduceRegion({
  reducer: ee.Reducer.max(),
  scale: 30,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Max NDWI of 2014: ', maxNDWI2014);

Map.addLayer(ndwi2014, {palette: ['red', 'yellow', 'green', 'blue'], min: -0.44, max: 0.16}, 'NDWI 2014 image');

// get area of L. Baringo in 2014

// create ndwi water mask for 2014
var ndwi2014Threshold = ndwi2014.gte(0.0);
var ndwi2014Mask = ndwi2014Threshold.updateMask(ndwi2014Threshold);
//Map.addLayer(ndwi2014Threshold, {palette: ['black', 'white']}, 'NDWI 2014 binary mask');
Map.addLayer(ndwi2014Mask, {palette: ['blue']}, 'NDWI 2014 Mask');

// get area of ndwi 2014 in km2 which will be size of lake baringo 2014

var ndwi2014Area = ndwi2014Mask.multiply(ee.Image.pixelArea());

var ndwi2014Stats = ndwi2014Area.reduceRegion({
  reducer: ee.Reducer.sum(),
  scale: 30,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Area of L. Baringo in 2014: ', ndwi2014Stats);

var ndwi2014Final = ee.Number(ndwi2014Stats.get('ndwi2014')).divide(1e6);
print('Area of L. Baringo in Km2 in 2014: ', ndwi2014Final);

///////////////////////////////////////////////
///////////////////////////////////////////////

// ndwi for 2020 image S2 - baringo2020Image

// ndwi = (Green - NIR)/(Green + NIR)
// note that green band in S2 is 'B3' while NIR is 'B8'

var green2020 = 'B3';
var nir2020 = 'B8';

var ndwi2020 = baringo2020Image.normalizedDifference([green2020, nir2020])
    .rename('ndwi2020');

print("NDWI 2020", ndwi2020);

var minNDWI2020 = ndwi2020.reduceRegion({
  reducer: ee.Reducer.min(),
  scale: 10,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Min NDWI of 2020: ', minNDWI2020);

var maxNDWI2020 = ndwi2020.reduceRegion({
  reducer: ee.Reducer.max(),
  scale: 10,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Max NDWI of 2020: ', maxNDWI2020);

Map.addLayer(ndwi2020, {palette: ['red', 'yellow', 'green', 'blue'], min: -0.81, max: 0.66}, 'NDWI 2020 image');

// get area of L. Baringo in 2020

// create ndwi water mask for 2020
var ndwi2020Threshold = ndwi2020.gte(0.0);
var ndwi2020Mask = ndwi2020Threshold.updateMask(ndwi2020Threshold);
//Map.addLayer(ndwi2014Threshold, {palette: ['black', 'white']}, 'NDWI 2020 binary mask');
Map.addLayer(ndwi2020Mask, {palette: ['blue']}, 'NDWI 2020 Mask');

// get area of ndwi 2020 in km2 which will be size of lake baringo 2020

var ndwi2020Area = ndwi2020Mask.multiply(ee.Image.pixelArea());

var ndwi2020Stats = ndwi2020Area.reduceRegion({
  reducer: ee.Reducer.sum(),
  scale: 10,
  bestEffort: true,
  maxPixels: 1e7,
});

print('Area of L. Baringo in 2020: ', ndwi2020Stats);

var ndwi2020Final = ee.Number(ndwi2020Stats.get('ndwi2020')).divide(1e6);
print('Area of L. Baringo in Km2 in 2020: ', ndwi2020Final);

////////////////////////////////////
////////////////////////////////////

// create a datatable for the area of Lake Baringo in 2010, 2014 and 2020 in km2
//

var dataTable = [
  [
    {label: 'Year', role: 'domain', type: 'string'},
    {label: 'Area_km2', role: 'data', type: 'number'},
    {label: 'Satellite Used', role: 'domain', type: 'string'}
  ],
  ['2010', 129.67, 'Landsat 5'],
  ['2014', 189.22, 'Landsat 8'],
  ['2020', 189.23, 'Sentinel 2'],
];

// Define the chart and print it to the console.
var chart = ui.Chart(dataTable).setChartType('ColumnChart').setOptions({
  title: 'Variation of area of Lake Baringo in km2 2010, 2014, 2020',
  legend: {position: 'none'},
  hAxis: {title: 'Year', titleTextStyle: {italic: false, bold: true}},
  vAxis: {title: 'Area_km2', titleTextStyle: {italic: false, bold: true}},
  colors: [ '#800080']
});

print(chart);

Leave a Reply

Your email address will not be published. Required fields are marked *

More Reading

Post navigation