Libraries used for this project:

library(tidyverse)
library(sf)
library(raster)
library(getlandsat)
library(mapview)
library(leaflet)
library(stats)

Question 1: Identifying the AOI

This region defines the AOI for this analysis. We will be looking at Palo, Iowa and its surroundings in order to extract its flood extents. The bounding box is a 5000 meter buffer around the point.

Question 2: Loading in Images as Multiband Raster Object in R and Crop to Domain of AOI

## class      : RasterStack 
## dimensions : 7811, 7681, 59996291, 6  (nrow, ncol, ncell, nlayers)
## resolution : 30, 30  (x, y)
## extent     : 518085, 748515, 4506585, 4740915  (xmin, xmax, ymin, ymax)
## crs        : +proj=utm +zone=15 +datum=WGS84 +units=m +no_defs 
## names      : band1, band2, band3, band4, band5, band6 
## min values :     0,     0,     0,     0,     0,     0 
## max values : 65535, 65535, 65535, 65535, 65535, 65535

As shown above the dimensions of the stacked image are 7811, 7681, 59996291, 6. The CRS is +proj=utm +zone=15 +datum=WGS84 +units=m +no_defs. The cell resolution is 30, 30.

## class      : RasterBrick 
## dimensions : 340, 346, 117640, 6  (nrow, ncol, ncell, nlayers)
## resolution : 30, 30  (x, y)
## extent     : 594075, 604455, 4652535, 4662735  (xmin, xmax, ymin, ymax)
## crs        : +proj=utm +zone=15 +datum=WGS84 +units=m +no_defs 
## source     : /private/var/folders/pv/33xhzt056mv7mj5063lqq3d80000gn/T/RtmplmOhqa/raster/r_tmp_2020-09-07_002429_21370_15133.grd 
## names      : Coastal.Aerosol,  Blue, Green,   Red, Near.Infrared, SWIR.1 
## min values :            8566,  7721,  6670,  6057,          5141,   4966 
## max values :           18482, 19079, 19592, 21965,         24772,  23896

The dimensions of the cropped image stack are 340, 346, 117640, 6. The CRS is +proj=utm +zone=15 +datum=WGS84 +units=m +no_defs and the cell resolution is also 30, 30.

Question 3: Visualization Through RGB plots and Stretches

Landsat Images Without Color Stretch

## [1] "R-G-B (natural color)"

## [1] "NIR-R-G(fa) (color infrared)"

## [1] "NIR-SWIR1-R (false color water focus)"

## [1] "SWIR2-SWIR1-R"

Landsat Images With Color Stretch

## [1] "R-G-B (natural color)"

## [1] "NIR-R-G(fa) (color infrared)"

## [1] "NIR-SWIR1-R (false color water focus)"

## [1] "SWIR2-SWIR1-R"

The purpose of using a stretch is to put more contrast into the images by reducing the variability by stretching the image. Instead of using the maximum/minimum possible RGB values, it uses minimum and maximum values that exist in the image. This overall increases the contrast of the image effect and makes it look a little bit nicer.

Question 4: Utilizing Raster Algebra and Raster Thresholding

In terms of the similarity between all these images, they all split the image into two parts by color. The blue part is prominent in both the NDVI and SWI image while that same area of land is red in the NDWI, MNDWI, and WRI images. The NDWI, MNDWI, and WRI images show the water threshold at cells greater than 0 and 1 while the NDVI shows the water threshold at cells less than 0. NDVI is used to show the normalized difference vegetation index, NDWI shows the normalized difference water index, MNDWI is a modified version of the NDWI, WRI show the water ratio index, and the SWI show the simple water index.

Through raster thresholding I was able to create the graph above which is a binary flood raster showing the flooded areas in blue and the non-flooded areas in white for each of the different images.

Question 5: Using K-Means Algorithm to Group Raster Cells With Similar Spectral Properties

## [1] 117640      6

In terms of the dimensions of the extracted values, there are 117,640 rows and 6 columns, one column for each different band of the raster stack. This shows that the data was extracted cell-by-cell for each of the methods as there are 6,705,840 points (117,640 cells for each method).

In this, we used unsupervised classification to identify similar features in a continuous field. For this, I specified 12 different clusters and the algorithm created a categorization based on patterns in the data. Then I added this layer back into the raster stack to show it next to all the other methods of displayed the flood.

Question 6: Summary

Number of Flooded Cells per Image
Number
NDVI 6666
NDWI 7212
MNDWI 11939
WRI 8469
SWI 15207
K.Means 9074

Since we know that each cell has a resolution of 30 meters by 30 meters, each cell has an area of 900 meters squared. We then multiply 900 by the number of cells to get the total flooded area in each image.

Area of Flooded Cells (m^2)
Area
NDVI 5999400
NDWI 6490800
MNDWI 10745100
WRI 7622100
SWI 13686300
K.Means 8166600

We can visualize the uncertainty in our classifications by summing the entire stack using the calc function. The higher the count in each pixel, the more certain we can be about its flooded state. For example, if a cell has a value of 6, it indicates that every method identified the cell as flooded, if it has a value of 2 then we know that two of the methods identified the cell as flooded.

When looking at a mapview of raster, it is evident that in certain areas the pixel level is not a whole number. My explanation for this would be that this is caused by the overall resolution of the image. If we were to perform an aggregate function on the image in order to reduce the resolution so that each cell has solely one color, each cell would have a whole number for its pixel level. Because the image is higher quality, the colors can be a little mixed which gives us decimals for the pixel levels.

Question 7: Checking Flood Raster Against Real Drone Footage of Flood

Using context clues, I was able to find the exact location of the image above on the map with the bounding box of Palo, Iowa. Using the latitude and longitude of that point, I created an sfc object in order to extract the binary flood values at that specific location from the six layer flood map stack. As you can see below, in each of the layers, the value at this point is 1 which means that all of the layers capture flooding at this specific location.

## [1] "Lng = -91.78959, Lat = 42.06305"
##      NDVI NDWI MNDWI WRI SWI K.Means
## [1,]    1    1     1   1   1       1