Show the code
library(sf)
library(lidR)
library(RCSF)
library(terra)
library(raster)
library(tmap)
library(RColorBrewer)In this report, I will be focusing on working with LiDAR data to analyze whether the salt marsh in Jamaica Bay, Queens increased or decreased in elevation between 2014 and 2017. The LiDAR data will be obtained from the National Oceanic and Atmospheric Administration (NOAA). This will be done by processing point cloud data into DEM maps and comparing the differences in elevation data from both years. This will help me determine whether the marsh is healthy and growing or slowly eroding over time.
library(sf)
library(lidR)
library(RCSF)
library(terra)
library(raster)
library(tmap)
library(RColorBrewer)I read in the data files and reduced the number of points to improve processing speed. Afterwards, I classified the LiDAR point cloud into ground and non-ground points, isolating ground level elevation from vegetation and other above ground features.
# read las file
las_2014 <- readLAS("/Users/Ally/Desktop/LIDAR_Project/25170_2014.las")Warning: There are 3874701 points flagged 'withheld'.
las_2017 <- readLAS("/Users/Ally/Desktop/LIDAR_Project/25170_2017.las")# decreasing number of points
las_2014 <- decimate_points(las_2014, homogenize(1, res = 1))
las_2017 <- decimate_points(las_2017, homogenize(1, res = 1))
# classify LiDAR point clouds into ground and non-ground points
las_2014 <- classify_ground(las_2014, algorithm = csf(time_step = 5))Original dataset already contains 1132774 ground points. These points were reclassified as 'unclassified' before performing a new ground classification.
las_2017 <- classify_ground(las_2017, algorithm = csf(time_step = 5))Original dataset already contains 162414 ground points. These points were reclassified as 'unclassified' before performing a new ground classification.
# filtering ground point classification
ground_2014 <- filter_poi(las_2014, Classification == 2)
ground_2017 <- filter_poi(las_2017, Classification == 2) Using the classified ground points, I generated a DEM for each year displaying the salt marsh elevation across Jamaica Bay. I then subtracted the two grids to produce a difference map showing elevation change between 2014 and 2017.
# clean and convert scattered LiDAR points into an elevation grid
dem_2014 <- rasterize_terrain(ground_2014, res = 1, algorithm = tin())
# clean and convert scattered LiDAR points into an elevation grid
dem_2017 <- rasterize_terrain(ground_2017, res = 1, algorithm = tin())
crs(dem_2017) <- "EPSG:32618"
# reproject 2017 to match 2014 coordinate system
dem_2017 <- project(dem_2017, crs(dem_2014))
dem_2017 <- resample(dem_2017, dem_2014)plot(dem_2014,
main = "2014 Salt Marsh Elevation",
cex.main = 0.95,
xlab = "Longitude",
ylab = "Latitude",
plg = list(title = "meters", title.cex = 0.8))
contour(dem_2014, add = TRUE, nlevels = 5, col = "black", lwd = 0.5)In the 2014 Salt Marsh Elevation map, the marsh platform is clearly visible in the light green and teal areas, with elevations between 5–10 meters, indicating a well established marsh structure at this time.
plot(dem_2017,
main = "2017 Salt Marsh Elevation",
cex.main = 0.95,
xlab = "Longitude",
ylab = "Latitude",
plg = list(title = "meters", title.cex = 0.8))
contour(dem_2017, add = TRUE, nlevels = 5, col = "black", lwd = 0.5)By 2017, the marsh platform has become less distinct, with ground level dropping closer to 0 meters elevation which are similar to surrounding water levels. This suggests the marsh may be losing elevation relative to sea level over time.
dem_difference <- dem_2017 - dem_2014
plot(dem_difference,
main = "Elevation Change\n(2014-2017)",
cex.main = 0.9,
col = brewer.pal(11, "RdYlGn"),
xlab = "Longitude",
ylab = "Latitude",
plg = list(title = "meters", title.cex = 0.8))
contour(dem_difference, add = TRUE, nlevels = 5, col = "black", lwd = 0.5)The difference grid (2017 - 2014) reveals mixed patterns of elevation change across Jamaica Bay. The red areas centered around the marsh island indicate elevation loss indicating erosion, while the green areas suggests water gains.
Organic matter inventory was estimated by applying the linear regression equation reported in Wang et al. (2023), relating marsh accretion rate to organic matter inventory. I will first covert the difference elevation change by 3 since its three years to find the rate of change per year.
# convert elevation change to annual accretion (cm/yr)
accretion <- (dem_difference * 100) / 3
# wang equation
organic_matter <- (accretion - 0.12) / 0.15
plot(organic_matter,
main = "Relating Marsh Accretion Rate\nto Organic Matter",
cex.main = 0.9,
col = brewer.pal(11, "RdYlGn"),
xlab = "Longitude",
ylab = "Latitude",
plg = list(title = "2 g cm⁻²", title.cex = 0.8))
contour(organic_matter, add = TRUE, nlevels = 5, col = "black", lwd = 0.5)This map shows similar pattern to the difference chart supporting negative values (0 to −2000) indicate organic matter loss, concentrated at the marsh island, while surrounding areas show modest gains.
This project analyzed LiDAR elevation data from 2014 and 2017 to assess salt marsh change in Jamaica Bay, finding that the marsh platform lost significant elevation over the three year period. The organic matter analysis further confirmed these losses, suggesting the Jamaica Bay marsh is increasingly vulnerable to sea level rise and long term coastal erosion.