Coastal Flooding in Norfolk, Virginia

Introduction

As a result of rising sea levels due to climate change, coastal communities are increasingly dealing with the effects of flooding. Coastal flooding occurs when storm surges or high tides cause coastal land that is typically dry to become inundated. Coastal flooding can increase erosion, harm wetland ecosystems, and threaten infrastructure such as roads, wastewater facilities, and drinking plants. In Norfolk, located on the Virginia Coastline, rising sea levels and sinking land endanger livelihood, property, and vitality in local communities. As Norfolk is facing the highest rate of sea level rise on the East Coast, it is crucial to understand the spatial variation of flooding in this area to understand its impacts on local communities and ecosystems.

Data Preparation

## Load Libraries 
library(tidyverse)
library(sf)
library(tigris)
library(tmap)
library(ggthemes)
library(wordcountaddin)
library(tmaptools)
library(kableExtra)
library(plotly)
library(spdep)
## Read in
nri_tct_va <- read_csv("../raw_data/NRI_Table_CensusTracts_Virginia.csv")

## Select using column names and dplyr style code
nri_tct_va <- nri_tct_va |>
  select(TRACTFIPS, AREA, CFLD_EXP_AREA)

## Select using column names and dplyr style code
nri_tct_va <- nri_tct_va |>
  mutate(TRACTFIPS = as.character(TRACTFIPS))

## Rename columns
nri_tct_va <- nri_tct_va |>
  rename(GEOID = TRACTFIPS,
         IMPCTAREA = CFLD_EXP_AREA)

## Add column with the 1st through 5th chars in GEOID column
nri_tct_va <- nri_tct_va |>
  mutate(FIPS = str_sub(GEOID, 1, 5))

## Extract a portion of a string
nri_tct_norfolk <- nri_tct_va |>
  filter(FIPS == "51710")

## Read in Norfolk census tracts 
norfolk_tracts <- read_sf("../modified_data/norfolk_tracts.shp")

## Keep important columns
norfolk_tracts <- norfolk_tracts |>
  select(GEOID, geometry)

## Perform table join using merge
norfolk_nri_sp <- left_join(norfolk_tracts,         
                       nri_tct_norfolk,               
                       by = c("GEOID" = "GEOID"))      

## calculate proportion of area impacted 
norfolk_nri_sp <- norfolk_nri_sp |>
  mutate(AREAPCT = 100 * IMPCTAREA / AREA)

## kable data table with scrolling 
kable(norfolk_nri_sp) |> 
  kable_styling(bootstrap_options = c("striped", 
                                      "hover", 
                                      "condensed", 
                                      "responsive"), 
                full_width = T) %>% 
  scroll_box(width = "60%", 
             height = "300px", 
             fixed_thead = TRUE)
GEOID geometry AREA IMPCTAREA FIPS AREAPCT
51710004400 MULTIPOLYGON (((-76.25925 3… 0.2763212 0.0000683 51710 0.0247101
51710006400 MULTIPOLYGON (((-76.24552 3… 1.1473991 0.3256889 51710 28.3849719
51710002200 MULTIPOLYGON (((-76.29326 3… 0.9241701 0.2877624 51710 31.1373829
51710004300 MULTIPOLYGON (((-76.27218 3… 0.4628282 0.0094488 51710 2.0415246
51710002100 MULTIPOLYGON (((-76.27796 3… 0.7102228 0.1804203 51710 25.4033453
51710006902 MULTIPOLYGON (((-76.21013 3… 1.1262393 0.3285644 51710 29.1735834
51710005100 MULTIPOLYGON (((-76.2735 36… 1.4174726 0.0817573 51710 5.7678211
51710007001 MULTIPOLYGON (((-76.22602 3… 0.6642082 0.0385023 51710 5.7967168
51710005702 MULTIPOLYGON (((-76.25835 3… 0.6936517 0.0272363 51710 3.9265060
51710004100 MULTIPOLYGON (((-76.28708 3… 0.1597471 0.0041904 51710 2.6231634
51710003700 MULTIPOLYGON (((-76.30169 3… 0.2431408 0.0610361 51710 25.1031810
51710006603 MULTIPOLYGON (((-76.21943 3… 0.6804853 0.0844220 51710 12.4061413
51710002700 MULTIPOLYGON (((-76.30205 3… 0.3874264 0.0446409 51710 11.5224237
51710002300 MULTIPOLYGON (((-76.30491 3… 0.6102213 0.2672848 51710 43.8012827
51710000300 MULTIPOLYGON (((-76.26356 3… 0.7157623 0.0115168 51710 1.6090280
51710000500 MULTIPOLYGON (((-76.27132 3… 0.4287847 0.0308805 51710 7.2018683
51710000202 MULTIPOLYGON (((-76.24556 3… 0.5091204 0.0060690 51710 1.1920469
51710006502 MULTIPOLYGON (((-76.20103 3… 0.6754866 0.2038408 51710 30.1768791
51710000800 MULTIPOLYGON (((-76.27574 3… 0.4469770 0.0357428 51710 7.9965535
51710003600 MULTIPOLYGON (((-76.29115 3… 0.3330928 0.1035538 51710 31.0885852
51710000600 MULTIPOLYGON (((-76.26077 3… 0.4869879 0.0086898 51710 1.7844053
51710004001 MULTIPOLYGON (((-76.29532 3… 0.1464680 0.0951905 51710 64.9906340
51710005602 MULTIPOLYGON (((-76.25322 3… 0.6100801 0.0101081 51710 1.6568524
51710000400 MULTIPOLYGON (((-76.30047 3… 2.9590014 0.3037776 51710 10.2662214
51710006604 MULTIPOLYGON (((-76.20341 3… 0.4097964 0.1304187 51710 31.8252504
51710980200 MULTIPOLYGON (((-76.24971 3… 0.8772063 0.1186578 51710 13.5267806
51710006100 MULTIPOLYGON (((-76.25426 3… 1.2287473 0.1581068 51710 12.8673156
51710005903 MULTIPOLYGON (((-76.23301 3… 0.5653327 0.0000000 51710 0.0000000
51710003000 MULTIPOLYGON (((-76.27476 3… 0.4750578 0.2173103 51710 45.7439645
51710000700 MULTIPOLYGON (((-76.26262 3… 0.5246623 0.0615018 51710 11.7221604
51710002500 MULTIPOLYGON (((-76.31676 3… 0.5586875 0.1796799 51710 32.1610704
51710006601 MULTIPOLYGON (((-76.23628 3… 0.3118729 0.0128240 51710 4.1119378
51710003200 MULTIPOLYGON (((-76.26631 3… 0.3461359 0.0583796 51710 16.8660934
51710003300 MULTIPOLYGON (((-76.26118 3… 0.5029339 0.0702680 51710 13.9716252
51710005000 MULTIPOLYGON (((-76.29586 3… 1.5290909 0.2541360 51710 16.6200686
51710001300 MULTIPOLYGON (((-76.28824 3… 0.2557990 0.0088394 51710 3.4556052
51710004600 MULTIPOLYGON (((-76.26425 3… 0.5734240 0.1603263 51710 27.9594655
51710003501 MULTIPOLYGON (((-76.28391 3… 0.4481478 0.0315097 51710 7.0311012
51710003100 MULTIPOLYGON (((-76.2612 36… 0.4480437 0.1292025 51710 28.8370277
51710001700 MULTIPOLYGON (((-76.28805 3… 0.5142883 0.0721866 51710 14.0362168
51710990000 MULTIPOLYGON (((-76.29906 3… NA NA NA NA
51710004500 MULTIPOLYGON (((-76.25908 3… 0.3016264 0.0327613 51710 10.8615464
51710005701 MULTIPOLYGON (((-76.26911 3… 0.8354736 0.0133135 51710 1.5935244
51710004200 MULTIPOLYGON (((-76.27972 3… 0.1799140 0.0766840 51710 42.6225893
51710005500 MULTIPOLYGON (((-76.27185 3… 0.8639806 0.0155387 51710 1.7985034
51710006607 MULTIPOLYGON (((-76.22815 3… 0.9264286 0.0349386 51710 3.7713195
51710002400 MULTIPOLYGON (((-76.31604 3… 1.4809528 0.4657953 51710 31.4524060
51710001100 MULTIPOLYGON (((-76.31746 3… 0.2574374 0.1121736 51710 43.5731604
51710004800 MULTIPOLYGON (((-76.28535 3… 0.4349904 0.1263917 51710 29.0562124
51710002000 MULTIPOLYGON (((-76.26816 3… 0.4627884 0.1273999 51710 27.5287617
51710005601 MULTIPOLYGON (((-76.25321 3… 0.8069619 0.0884761 51710 10.9640979
51710005800 MULTIPOLYGON (((-76.25042 3… 0.6124271 0.0000000 51710 0.0000000
51710006901 MULTIPOLYGON (((-76.21094 3… 1.6822884 0.0241012 51710 1.4326464
51710006602 MULTIPOLYGON (((-76.22764 3… 0.6053732 0.0199575 51710 3.2967192
51710006200 MULTIPOLYGON (((-76.23161 3… 0.8958245 0.0968652 51710 10.8129712
51710005902 MULTIPOLYGON (((-76.23852 3… 0.7937471 0.0083974 51710 1.0579413
51710000901 MULTIPOLYGON (((-76.31738 3… 1.3254450 0.1674073 51710 12.6302730
51710000902 MULTIPOLYGON (((-76.33316 3… 10.0050167 0.5839221 51710 5.8362929
51710003400 MULTIPOLYGON (((-76.27916 3… 0.4687913 0.0432488 51710 9.2255969
51710002800 MULTIPOLYGON (((-76.29409 3… 0.6541583 0.2794145 51710 42.7135844
51710004700 MULTIPOLYGON (((-76.27588 3… 0.6047765 0.0785308 51710 12.9851012
51710006605 MULTIPOLYGON (((-76.22231 3… 0.8930169 0.0370087 51710 4.1442338
51710006000 MULTIPOLYGON (((-76.25359 3… 0.6546599 0.1350215 51710 20.6246876
51710002900 MULTIPOLYGON (((-76.28921 3… 0.7495909 0.1976408 51710 26.3664835
51710000100 MULTIPOLYGON (((-76.22989 3… 0.3834527 0.0926094 51710 24.1514473
51710007002 MULTIPOLYGON (((-76.22808 3… 1.2987815 0.4273198 51710 32.9015920
51710006501 MULTIPOLYGON (((-76.21624 3… 0.3911918 0.2554775 51710 65.3074622
51710004002 MULTIPOLYGON (((-76.30743 3… 0.6665272 0.2586562 51710 38.8065464
51710980100 MULTIPOLYGON (((-76.3336 36… 1.9157681 0.0970578 51710 5.0662604
51710006606 MULTIPOLYGON (((-76.23849 3… 0.8668136 0.0027232 51710 0.3141583
51710002600 MULTIPOLYGON (((-76.30326 3… 0.3964619 0.1674730 51710 42.2419014
51710003800 MULTIPOLYGON (((-76.31721 3… 0.5511566 0.1991870 51710 36.1398179
51710005901 MULTIPOLYGON (((-76.24277 3… 0.5760109 0.0000000 51710 0.0000000
51710001200 MULTIPOLYGON (((-76.31876 3… 1.3328593 0.2977807 51710 22.3414983
51710001400 MULTIPOLYGON (((-76.27372 3… 0.5469617 0.0592775 51710 10.8375932
51710000201 MULTIPOLYGON (((-76.24569 3… 0.2748552 0.0038547 51710 1.4024330
51710980300 MULTIPOLYGON (((-76.21495 3… 2.3721508 0.0010261 51710 0.0432556
51710001500 MULTIPOLYGON (((-76.29227 3… 0.3717749 0.0410818 51710 11.0501876
51710006800 MULTIPOLYGON (((-76.20641 3… 1.1581729 0.0049174 51710 0.4245858
51710001600 MULTIPOLYGON (((-76.29546 3… 0.1989382 0.0197966 51710 9.9511381
51710004900 MULTIPOLYGON (((-76.29925 3… 0.7219774 0.3157660 51710 43.7362799

Coastal flooding risk data was taken from National Risk Index data at the census tract level. To prepare this data for analysis, the data was reduced to only necessary fields: FIPS, area, and area impacted by coastal flooding. These variables were renamed to increase clarity. The area and impacted area were measured in square miles. To create a GEOID field for joining purposes, the first five numbers of the TRACTFIPS field were taken and turned into character variables. Observations were then filtered to only census tracts from the observation area: Norfolk, Virginia. Census Tract data for Norfolk were from a cartographic boundary shapefile retrieved from the Census Bureau. This shapefile was subset to only GEOID and geometry. NRI data was table joined to census tract data using the GEOID columns so spatial analysis could be performed. The percentage of each census tract impact was calculated using the impacted area over the total area for normalization purposes.

Exploratory Data Science Analysis

Data Description and Summary

## row count
row_count <- nrow(norfolk_nri_sp)

## Identify complete rows (equal to TRUE)
complete_obs <- norfolk_nri_sp |> 
  na.omit()

## Identify rows with NAs (equal to FALSE in output)
na_obs <- row_count - nrow(complete_obs)

## Maximum Percentage Area Impacted
AREAPCT_MAX <- max(norfolk_nri_sp$AREAPCT, na.rm = TRUE)

## Minimum Percentage Area Impacted  
AREAPCT_MIN <- min(norfolk_nri_sp$AREAPCT, na.rm = TRUE)

## Average Percentage Area Impacted 
AREAPCT_AVE <- mean(norfolk_nri_sp$AREAPCT, na.rm = TRUE) 

## Histogram
norfolk_nri_sp |> 
  ggplot(mapping = aes(x = AREAPCT)) +
  geom_histogram(binwidth = 2,
                 na.rm = TRUE,
                  color = "black",
                 fill = "darkgray") +
  labs(title = "Percentage of Area Impacted by Coastal Flooding (by Census Tract)",
       subtitle = "Norfolk, VA",
       caption = "GEOG 215 (Tsehai Ricketts)",
       y = "Count",
       x = "Area Impacted (%)") +
  theme_tufte() +                                        ## Apply theme 
        theme(plot.title = element_text(size = 14,       ## Title font size
                                        face = "bold"),  ## Title font face
              plot.subtitle = element_text(size = 12))   ## Subitle font size

My variable of interest was the percentage of area impacted by coastal flooding by census tract in Norfolk, Virginia. Impacted area is defined by the NRI as the area in a census block that falls within minor, moderate, and major high tide flooding thresholds, which the NRI derives from NOAA high tide flooding data.

The dataset included 81 census tracts with 1 NA observation. The average proportion of area impacted was 17.06%, with a range from 0% to 65.31%. The data is skewed right, and the modal value is closer to zero. There seem to be around a break between the majority of the data and two outliers where over 60% of the tract was impacted by flooding, higher than the other observations.

Geographic Distribution and Clustering

## Load basemap
bg_basemap <- 'https://server.arcgisonline.com/ArcGIS/rest/services/World_Topo_Map/MapServer/tile/{z}/{y}/{x}'

## Put t-map in view mode
tmap_mode("view")

## Chloropleth map
areapct_tmap <- tm_shape(norfolk_nri_sp,
                     projection = "ESRI:102003") + 
  tm_basemap(bg_basemap) + 
  tm_polygons("AREAPCT",
              title = "Area Impacted (%)",
              style = "jenks",
              palette = "Reds",
              colorNA = "white",
              border.col = "black",
              border.alpha = 0.25,
              popup.vars = c("(%)" = "AREAPCT")) +          
  tm_scale_bar() +                                             ## Scale bar 
    tm_view(view.legend.position = c("right", "top"))         ## Location of legend

## Filter NA values
norfolk_nri_sp_com <- norfolk_nri_sp|>
  filter(!is.na(AREAPCT))

## Turn off scientific notation 
options(scipen = 999)

## Create Queen case neighbors
norfolk_nb_queen <- poly2nb(norfolk_nri_sp_com, 
                       queen = TRUE)

## Convert neighbor object to weight matrix
norfolk_wm_queen <-  nb2listw(norfolk_nb_queen,
                         style = "B",           # binary
                         zero.policy = TRUE)
## Moran's I
norfolk_moran <- moran.test(norfolk_nri_sp_com$AREAPCT,         # Variable we're analyzing
                           listw = norfolk_wm_queen,       # Sp weights matrix
                           alternative = "two.sided", # Clustering or Dispersion
                           randomisation = TRUE,      # Compare to randomized values
                           zero.policy = TRUE)        # Allow obs with 0 neighbors

## LISA -- Local Moran's I
norfolk_nri_lisa <- localmoran(norfolk_nri_sp_com$AREAPCT,         # Variable we're analyzing
                          listw = norfolk_wm_queen,       # Weights object
                          alternative = "two.sided", # Clustering or Dispersion
                          zero.policy = TRUE) |>     # Best to keep TRUE for LISA
  as_tibble() |>                                     # Better object type
  mutate(across(everything(), as.vector))            # Remove junk from localmoran output

## Add values required for LISA category
norfolk_nri_lisa <- norfolk_nri_lisa |>
  mutate(SCVAR =  scale(norfolk_nri_sp_com$AREAPCT) |> as.vector(),  # Original data column
         LAGVAR = lag.listw(norfolk_wm_queen, 
                            scale(norfolk_nri_sp_com$AREAPCT)),     # Lag of original data column
         LISACAT = case_when(SCVAR >= 0 & LAGVAR >= 0 & `Pr(z != E(Ii))` <= 0.05 ~ 1,
                             SCVAR <= 0 & LAGVAR <= 0 & `Pr(z != E(Ii))` <= 0.05 ~ 2,
                             SCVAR >= 0 & LAGVAR <= 0 & `Pr(z != E(Ii))` <= 0.05 ~ 3,
                             SCVAR <= 0 & LAGVAR >= 0 & `Pr(z != E(Ii))` <= 0.05 ~ 4,
                             `Pr(z != E(Ii))` > 0.05 ~ 5))

## This adds a label based on the values
norfolk_nri_lisa <- norfolk_nri_lisa |>
  mutate(CATNAME = case_when(LISACAT == 1 ~ "HH",
                             LISACAT == 2 ~ "LL",
                             LISACAT == 3 ~ "HL",
                             LISACAT == 4 ~ "LH",
                             LISACAT == 5 ~ "Not Significant"))

## Quick Table of LISA output
lisa_table <- table(norfolk_nri_lisa$CATNAME)

## Kable table of Lisa Summary
lisa_kable <- 
  lisa_table |>
  kable(col.names = c("Category", "Frequency")) |>
  kable_styling(bootstrap_options = c("striped", 
                                      "hover", 
                                      "condensed",
                                      "responsive"), 
                full_width = T)

## Add LISA category column to the spatial data for mapping
norfolk_nri_sp_com <- norfolk_nri_sp_com |>
  mutate(LISACAT = norfolk_nri_lisa$LISACAT,
         CATNAME = norfolk_nri_lisa$CATNAME)

## Plot maps together
### LISA map
lisa_tmap <- tm_shape(norfolk_nri_sp_com,
                      projection = "ESRI:102003") + 
  tm_basemap(bg_basemap) + 
  tm_polygons("LISACAT", 
              title = "LISA Category",
              breaks = c(1, 2, 3, 4, 5, 6),
              palette =  c("red", 
                           "blue", 
                           "lightpink", 
                           "skyblue", 
                           "grey90"),
              colorNA = "white",
              labels = c("High-High", 
                         "Low-Low",
                         "High-Low",
                         "Low-High", 
                         "Not significant"),
              border.col = "black", 
              border.alpha = 0.25,
              popup.vars = c("(%)" = "AREAPCT")) +
   tm_scale_bar() +   ## Scale bar 
  tm_view(view.legend.position = c("right", "top"))        

### Map both maps together
tmap_arrange(areapct_tmap, 
             lisa_tmap,
             ncol = 2,
             nrow = 1,
             sync = TRUE) 
## Print Lisa kable table 
lisa_kable
Category Frequency
HH 9
LH 2
LL 9
Not Significant 60

The southwest of the city, particularly near the Elizabeth and Lafayette rivers, seems to have the tracts that have the highest proportion of area impacted. There also seems to be another area of high values towards the northeastern edge. The interior portion of Norfolk seems largely unimpacted, which is logical as these areas are the furthest away from water.

The Moran’s I statistic is 0.3856931 and the p-value is 0. The Moran value tells us that the data is moderately clustered and the p-value tells us this result is statistically significant. There are 9 high-high clustered census tracts and 9 low low clustered census tracts. There are 2 low-high census tracts and 60 census tracts did not return a significant result. The LISA results support clustering in the southwest and northeast and suggests that the eastern interior portion of Norfolk has a cluster of largely unaffected census tracts surrounded by other unaffected census tracts.

Conclusions

This project wrangled and analyzed NRI Coastal flooding data to visualize and explore patterns and variations across Norfolk, Virginia. Mapping the distributions of the proportion of impacted area and spatial autocorrelation analysis suggests that southwestern and northeastern census tracts in Norfolk are more vulnerable to the impacts of coastal flooding. It also suggests that the interior section of Norfolk is associated with a lower risk of flooding. As a future extension, it would be interesting to research why some tracts directly along the coast have a high percentage of their area impacted by flooding while other tracts with the same proximity to water do not. While this study focused on one variable, future analysis to determine if there is an association between tracts highly impacted by high tide flooding and different social factors, such as poverty or race, could be beneficial in understanding the socio-economic consequences of the growing impact of climate change in the region.

Document Statistics

#### Word Count
```{r , echo = FALSE, message = FALSE}
wordcountaddin:::text_stats() %>%
  kable_styling(bootstrap_options = c("striped",
                                      "hover",
                                      "condensed",
                                      "responsive"),
                full_width = F)
```
Method koRpus stringi
Word count 745 759
Character count 4926 5094
Sentence count 39 Not available
Reading time 3.7 minutes 3.8 minutes