This paper presents some insights regarding effects of storms and other severe weather events on both public health and economic problems for communities and municipalities. Many severe events can result in fatalities, injuries, and property damage, and preventing such outcomes to the extent possible is a key concern.
For this research we will explore U.S. National Oceanic and Atmospheric Administration’s (NOAA) storm database.The database covers the time period between 1950 and November 2011.
The research aims to address the following questions:
Across the United States, which types of events are most harmful with respect to population health?
Across the United States, which types of events have the greatest economic consequences?
Storms and other severe weather events can cause both public health and economic problems for communities and municipalities. Many severe events can result in fatalities, injuries, and property damage, and preventing such outcomes to the extent possible is a key concern.
This project involves exploring the U.S. National Oceanic and Atmospheric Administration’s (NOAA) storm database. This database tracks characteristics of major storms and weather events in the United States, including when and where they occur, as well as estimates of any fatalities, injuries, and property damage.
The data for this assignment come in the form of a comma-separated-value file compressed via the bzip2 algorithm to reduce its size.
You can download the file from the course web site:
Let’s start with data download. We copy link used in project description and assign it to fileurl.
# We copy link used in project description and assign it to
# fileurl
filename <- "repdata_data_StormData.csv"
# Controll for already existing files If folder doesn't exist
# proceed with download
if (!file.exists(filename)) {
fileURL <- "https://d396qusza40orc.cloudfront.net/repdata%2Fdata%2FStormData.csv.bz2"
download.file(fileURL, filename, method = "curl")
}
# If file exists proceed with unzip
if (!file.exists("repdata_data_StormData.csv")) {
unzip(filename)
}
Now we load our downloaded data:
# Read repdata_data_StormData.csv
data <- read.csv("repdata_data_StormData.csv", sep = ",", header = TRUE,
stringsAsFactors = FALSE)
names(data)
## [1] "STATE__" "BGN_DATE" "BGN_TIME" "TIME_ZONE" "COUNTY"
## [6] "COUNTYNAME" "STATE" "EVTYPE" "BGN_RANGE" "BGN_AZI"
## [11] "BGN_LOCATI" "END_DATE" "END_TIME" "COUNTY_END" "COUNTYENDN"
## [16] "END_RANGE" "END_AZI" "END_LOCATI" "LENGTH" "WIDTH"
## [21] "F" "MAG" "FATALITIES" "INJURIES" "PROPDMG"
## [26] "PROPDMGEXP" "CROPDMG" "CROPDMGEXP" "WFO" "STATEOFFIC"
## [31] "ZONENAMES" "LATITUDE" "LONGITUDE" "LATITUDE_E" "LONGITUDE_"
## [36] "REMARKS" "REFNUM"
Let’s see how are data look like:
# A first look to data
head(data)
## STATE__ BGN_DATE BGN_TIME TIME_ZONE COUNTY COUNTYNAME STATE EVTYPE
## 1 1 4/18/1950 0:00:00 0130 CST 97 MOBILE AL TORNADO
## 2 1 4/18/1950 0:00:00 0145 CST 3 BALDWIN AL TORNADO
## 3 1 2/20/1951 0:00:00 1600 CST 57 FAYETTE AL TORNADO
## 4 1 6/8/1951 0:00:00 0900 CST 89 MADISON AL TORNADO
## 5 1 11/15/1951 0:00:00 1500 CST 43 CULLMAN AL TORNADO
## 6 1 11/15/1951 0:00:00 2000 CST 77 LAUDERDALE AL TORNADO
## BGN_RANGE BGN_AZI BGN_LOCATI END_DATE END_TIME COUNTY_END COUNTYENDN
## 1 0 0 NA
## 2 0 0 NA
## 3 0 0 NA
## 4 0 0 NA
## 5 0 0 NA
## 6 0 0 NA
## END_RANGE END_AZI END_LOCATI LENGTH WIDTH F MAG FATALITIES INJURIES PROPDMG
## 1 0 14.0 100 3 0 0 15 25.0
## 2 0 2.0 150 2 0 0 0 2.5
## 3 0 0.1 123 2 0 0 2 25.0
## 4 0 0.0 100 2 0 0 2 2.5
## 5 0 0.0 150 2 0 0 2 2.5
## 6 0 1.5 177 2 0 0 6 2.5
## PROPDMGEXP CROPDMG CROPDMGEXP WFO STATEOFFIC ZONENAMES LATITUDE LONGITUDE
## 1 K 0 3040 8812
## 2 K 0 3042 8755
## 3 K 0 3340 8742
## 4 K 0 3458 8626
## 5 K 0 3412 8642
## 6 K 0 3450 8748
## LATITUDE_E LONGITUDE_ REMARKS REFNUM
## 1 3051 8806 1
## 2 0 0 2
## 3 0 0 3
## 4 0 0 4
## 5 0 0 5
## 6 0 0 6
# Let's check structure of dataframe
glimpse(data)
## Observations: 902,297
## Variables: 37
## $ STATE__ <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,...
## $ BGN_DATE <chr> "4/18/1950 0:00:00", "4/18/1950 0:00:00", "2/20/1951 0:0...
## $ BGN_TIME <chr> "0130", "0145", "1600", "0900", "1500", "2000", "0100", ...
## $ TIME_ZONE <chr> "CST", "CST", "CST", "CST", "CST", "CST", "CST", "CST", ...
## $ COUNTY <dbl> 97, 3, 57, 89, 43, 77, 9, 123, 125, 57, 43, 9, 73, 49, 1...
## $ COUNTYNAME <chr> "MOBILE", "BALDWIN", "FAYETTE", "MADISON", "CULLMAN", "L...
## $ STATE <chr> "AL", "AL", "AL", "AL", "AL", "AL", "AL", "AL", "AL", "A...
## $ EVTYPE <chr> "TORNADO", "TORNADO", "TORNADO", "TORNADO", "TORNADO", "...
## $ BGN_RANGE <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,...
## $ BGN_AZI <chr> "", "", "", "", "", "", "", "", "", "", "", "", "", "", ...
## $ BGN_LOCATI <chr> "", "", "", "", "", "", "", "", "", "", "", "", "", "", ...
## $ END_DATE <chr> "", "", "", "", "", "", "", "", "", "", "", "", "", "", ...
## $ END_TIME <chr> "", "", "", "", "", "", "", "", "", "", "", "", "", "", ...
## $ COUNTY_END <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,...
## $ COUNTYENDN <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, ...
## $ END_RANGE <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,...
## $ END_AZI <chr> "", "", "", "", "", "", "", "", "", "", "", "", "", "", ...
## $ END_LOCATI <chr> "", "", "", "", "", "", "", "", "", "", "", "", "", "", ...
## $ LENGTH <dbl> 14.0, 2.0, 0.1, 0.0, 0.0, 1.5, 1.5, 0.0, 3.3, 2.3, 1.3, ...
## $ WIDTH <dbl> 100, 150, 123, 100, 150, 177, 33, 33, 100, 100, 400, 400...
## $ F <int> 3, 2, 2, 2, 2, 2, 2, 1, 3, 3, 1, 1, 3, 3, 3, 4, 1, 1, 1,...
## $ MAG <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,...
## $ FATALITIES <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 4, 0, 0, 0,...
## $ INJURIES <dbl> 15, 0, 2, 2, 2, 6, 1, 0, 14, 0, 3, 3, 26, 12, 6, 50, 2, ...
## $ PROPDMG <dbl> 25.0, 2.5, 25.0, 2.5, 2.5, 2.5, 2.5, 2.5, 25.0, 25.0, 2....
## $ PROPDMGEXP <chr> "K", "K", "K", "K", "K", "K", "K", "K", "K", "K", "M", "...
## $ CROPDMG <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,...
## $ CROPDMGEXP <chr> "", "", "", "", "", "", "", "", "", "", "", "", "", "", ...
## $ WFO <chr> "", "", "", "", "", "", "", "", "", "", "", "", "", "", ...
## $ STATEOFFIC <chr> "", "", "", "", "", "", "", "", "", "", "", "", "", "", ...
## $ ZONENAMES <chr> "", "", "", "", "", "", "", "", "", "", "", "", "", "", ...
## $ LATITUDE <dbl> 3040, 3042, 3340, 3458, 3412, 3450, 3405, 3255, 3334, 33...
## $ LONGITUDE <dbl> 8812, 8755, 8742, 8626, 8642, 8748, 8631, 8558, 8740, 87...
## $ LATITUDE_E <dbl> 3051, 0, 0, 0, 0, 0, 0, 0, 3336, 3337, 3402, 3404, 0, 34...
## $ LONGITUDE_ <dbl> 8806, 0, 0, 0, 0, 0, 0, 0, 8738, 8737, 8644, 8640, 0, 85...
## $ REMARKS <chr> "", "", "", "", "", "", "", "", "", "", "", "", "", "", ...
## $ REFNUM <dbl> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 1...
We will focus to thosse observations related to injuries, fatalities, property damage and crop damage. Also we will turn variable names to lower case
# Create df with focus data Name of variables in the lower
# case
df <- select(data, EVTYPE, FATALITIES, INJURIES, PROPDMG, PROPDMGEXP,
CROPDMG, CROPDMGEXP) %>% rename_all(tolower)
names(df)
## [1] "evtype" "fatalities" "injuries" "propdmg" "propdmgexp"
## [6] "cropdmg" "cropdmgexp"
Next step is to understand variable type and check for missing values:
# Let's check for the number of missing values
sum(is.na(df))
## [1] 0
# We check classes of every column inside data
lapply(df, class)
## $evtype
## [1] "character"
##
## $fatalities
## [1] "numeric"
##
## $injuries
## [1] "numeric"
##
## $propdmg
## [1] "numeric"
##
## $propdmgexp
## [1] "character"
##
## $cropdmg
## [1] "numeric"
##
## $cropdmgexp
## [1] "character"
There are no missing values in the variables inside df. We have 3 variables of type character (evtype, propdmgexp, cropdmgexp).
For evtype we will clean cases from leading and trailing whitespaces and dublications to avoid problems with classification:
# Clean evtype
df <- df %>% mutate(evtype = gsub(" +", " ", trimws(evtype)))
Let’s check if there are other symbols other than K, M, B, H inside variables propdmgexp, cropdmgexp
# For variable propdmgexp
unique(df$propdmgexp)
## [1] "K" "M" "" "B" "+" "0" "5" "6" "?" "4" "2" "3" "H" "7" "-" "1" "8"
# For variable crodmgexp
unique(df$cropdmgexp)
## [1] "" "M" "K" "B" "?" "0" "2"
Based in National Weather Service Storm Data Documentation we have to recode all these symbols to their appropriate values. The result of recoding we will save in 2 new variables propdmgcost, cropdmgcost, to represent damage created in dollars for properties and damage created in dollars for crop.
We will also create a new variable propcroptotal to aggregrate total damage.
Let’s start process step by step.
# Create 3 new variables and initiate them
df$propdmgcost = 0
df$cropdmgcost = 0
df$economy_damage_total = 0
# Let's fill new variable propdmgcost with data from recoding
# of propdmgexp
df[df$propdmgexp == "H", ]$propdmgcost = df[df$propdmgexp ==
"H", ]$propdmg * 100
df[df$propdmgexp == "K", ]$propdmgcost = df[df$propdmgexp ==
"K", ]$propdmg * 1000
df[df$propdmgexp == "M", ]$propdmgcost = df[df$propdmgexp ==
"M", ]$propdmg * 1e+06
df[df$propdmgexp == "B", ]$propdmgcost = df[df$propdmgexp ==
"B", ]$propdmg * 1e+09
# Let's fill new variable cropdmgcost with data from recoding
# of cropdmgexp
df[df$cropdmgexp == "H", ]$cropdmgcost = df[df$cropdmgexp ==
"H", ]$cropdmg * 100
df[df$cropdmgexp == "K", ]$cropdmgcost = df[df$cropdmgexp ==
"K", ]$cropdmg * 1000
df[df$cropdmgexp == "M", ]$cropdmgcost = df[df$cropdmgexp ==
"M", ]$cropdmg * 1e+06
df[df$cropdmgexp == "B", ]$cropdmgcost = df[df$cropdmgexp ==
"B", ]$cropdmg * 1e+09
# Total cost
df$totalcostdmg = df$propdmgcost + df$cropdmgcost
## Aggregate cost of damages
economic_effect <- df %>% select(evtype, totalcostdmg) %>% group_by(evtype) %>%
summarise(costdmg = sum(totalcostdmg)) %>% arrange(desc(costdmg))
head(economic_effect, 10)
## # A tibble: 10 x 2
## evtype costdmg
## <chr> <dbl>
## 1 FLOOD 150319678250
## 2 HURRICANE/TYPHOON 71913712800
## 3 TORNADO 57340613590
## 4 STORM SURGE 43323541000
## 5 HAIL 18752904670
## 6 FLASH FLOOD 17562178610
## 7 DROUGHT 15018672000
## 8 HURRICANE 14610229010
## 9 RIVER FLOOD 10148404500
## 10 ICE STORM 8967041310
Now let’s see if we can answer to question about types of events (as indicated in the Evtype} variable) are most harmful with respect to population health across the United States ?
For this we aggregrate fatalities by evtype to calculate sums and see top values
## Aggregate fatalities
agg_fatalities <- select(df, evtype, fatalities) %>% group_by(evtype) %>%
summarise(sum_fatalities = sum(fatalities)) %>% arrange(desc(sum_fatalities))
head(agg_fatalities, 10)
## # A tibble: 10 x 2
## evtype sum_fatalities
## <chr> <dbl>
## 1 TORNADO 5633
## 2 EXCESSIVE HEAT 1903
## 3 FLASH FLOOD 978
## 4 HEAT 937
## 5 LIGHTNING 816
## 6 TSTM WIND 504
## 7 FLOOD 470
## 8 RIP CURRENT 368
## 9 HIGH WIND 248
## 10 AVALANCHE 224
Let’s visualize our numeric results regarding top fatalities:
## Plot top fatalities
top_fatalities <- ggplot(head(agg_fatalities, 10), aes(x = reorder(evtype,
sum_fatalities), y = sum_fatalities)) + geom_bar(fill = "green",
stat = "identity") + coord_flip() + ylab("Fatalities") +
xlab("Event") + ggtitle("Most harmful weather events by fatality") +
theme(legend.position = "none")
top_fatalities
Let’s repeat our previous steps for identifying top injuries
## Aggregate injuries
agg_injuries <- select(df, evtype, injuries) %>% group_by(evtype) %>%
summarise(sum_injuries = sum(injuries)) %>% arrange(desc(sum_injuries))
head(agg_injuries, 10)
## # A tibble: 10 x 2
## evtype sum_injuries
## <chr> <dbl>
## 1 TORNADO 91346
## 2 TSTM WIND 6957
## 3 FLOOD 6789
## 4 EXCESSIVE HEAT 6525
## 5 LIGHTNING 5230
## 6 HEAT 2100
## 7 ICE STORM 1975
## 8 FLASH FLOOD 1777
## 9 THUNDERSTORM WIND 1488
## 10 HAIL 1361
And the plot:
## Plot top injuries
top_injuries <- ggplot(head(agg_injuries, 10), aes(x = reorder(evtype,
sum_injuries), y = sum_injuries)) + geom_bar(fill = "green",
stat = "identity") + coord_flip() + ylab("Injuries") + xlab("Event") +
ggtitle("Most harmful weather events by injuries") + theme(legend.position = "none")
top_injuries
When it comes to economic consequences, We will vizualize top 10 factors judging by damages in dollars.
## Plot top events which damages have highest cost
top_costdmg <- ggplot(head(economic_effect, 10), aes(x = reorder(evtype,
costdmg), y = costdmg)) + geom_bar(fill = "green", stat = "identity") +
coord_flip() + ylab("Cost in dollars") + xlab("Event") +
ggtitle("Top events which damages have highest cost") + theme(legend.position = "none")
top_costdmg
Now let’s replot the above plot using logarithmic scale. We will put both plots in the same grid to compare if the used scale helps in inproving understanding.
## Plot top events which damages have highest cost
top_costdmg <- ggplot(head(economic_effect, 10), aes(x = reorder(evtype,
costdmg), y = costdmg)) + geom_bar(fill = "green", stat = "identity") +
coord_flip() + ylab("Cost in dollars") + xlab("Event") +
ggtitle("Top events which damages have highest cost") + theme(legend.position = "none")
top_costdmg
## Plot the same plot using logarithmic scale
top_costdmg_log <- ggplot(head(economic_effect, 10), aes(x = reorder(evtype,
costdmg), y = log10(costdmg))) + geom_bar(fill = "yellow",
stat = "identity") + coord_flip() + ylab("Cost in dollars(log scale)") +
xlab("Event") + ggtitle("Top events which damages have highest cost") +
theme(legend.position = "none")
top_costdmg_log
## Arrange in grid
grid.arrange(top_costdmg, top_costdmg_log, ncol = 1, nrow = 2)
Analysis shows that Tornados and Excessive Heat are the top weather events when it comes to fatalities.
Tornados are also the top weather event with respect to injuries. Suprising fact is that no matter how fatal can be Excessive Heat, this event is only in the fourth place when it comes to injuries, and is precedded by Thunderstorms and Flood.
Events that mostly affect economics and largest property damages are flood and hurricanes.