knitr::opts_chunk$set(echo = TRUE)
In this work we analyze the U.S. National Oceanic and Atmospheric Administration’s (NOAA) storm database which tracks characteristics of major storms and weather events in the United States. The data includes when and where these events occur, as well as estimates of any fatalities, injuries, property and corn damage. The goal of the work is to find out the most harmful weather events with respect to population health and economic consequences.
library(R.utils) #can't read origin file from bz2 to csv
## Warning: package 'R.utils' was built under R version 3.4.2
## Loading required package: R.oo
## Loading required package: R.methodsS3
## R.methodsS3 v1.7.1 (2016-02-15) successfully loaded. See ?R.methodsS3 for help.
## R.oo v1.21.0 (2016-10-30) successfully loaded. See ?R.oo for help.
##
## Attaching package: 'R.oo'
## The following objects are masked from 'package:methods':
##
## getClasses, getMethods
## The following objects are masked from 'package:base':
##
## attach, detach, gc, load, save
## R.utils v2.5.0 (2016-11-07) successfully loaded. See ?R.utils for help.
##
## Attaching package: 'R.utils'
## The following object is masked from 'package:utils':
##
## timestamp
## The following objects are masked from 'package:base':
##
## cat, commandArgs, getOption, inherits, isOpen, parse, warnings
download.file("https://d396qusza40orc.cloudfront.net/repdata%2Fdata%2FStormData.csv.bz2", "StormData.csv.bz2")
bunzip2("D:\\DATA\\coursera\\Data Scientist\\datasciencecoursera\\Reproducible Research\\Assignment 2\\StormData.csv.bz2")
stormdata <- read.csv("D:\\DATA\\coursera\\Data Scientist\\datasciencecoursera\\Reproducible Research\\Assignment 2\\StormData.csv")
head(stormdata)
## STATE__ BGN_DATE BGN_TIME TIME_ZONE COUNTY COUNTYNAME STATE
## 1 1 4/18/1950 0:00:00 0130 CST 97 MOBILE AL
## 2 1 4/18/1950 0:00:00 0145 CST 3 BALDWIN AL
## 3 1 2/20/1951 0:00:00 1600 CST 57 FAYETTE AL
## 4 1 6/8/1951 0:00:00 0900 CST 89 MADISON AL
## 5 1 11/15/1951 0:00:00 1500 CST 43 CULLMAN AL
## 6 1 11/15/1951 0:00:00 2000 CST 77 LAUDERDALE AL
## EVTYPE BGN_RANGE BGN_AZI BGN_LOCATI END_DATE END_TIME COUNTY_END
## 1 TORNADO 0 0
## 2 TORNADO 0 0
## 3 TORNADO 0 0
## 4 TORNADO 0 0
## 5 TORNADO 0 0
## 6 TORNADO 0 0
## COUNTYENDN END_RANGE END_AZI END_LOCATI LENGTH WIDTH F MAG FATALITIES
## 1 NA 0 14.0 100 3 0 0
## 2 NA 0 2.0 150 2 0 0
## 3 NA 0 0.1 123 2 0 0
## 4 NA 0 0.0 100 2 0 0
## 5 NA 0 0.0 150 2 0 0
## 6 NA 0 1.5 177 2 0 0
## INJURIES PROPDMG PROPDMGEXP CROPDMG CROPDMGEXP WFO STATEOFFIC ZONENAMES
## 1 15 25.0 K 0
## 2 0 2.5 K 0
## 3 2 25.0 K 0
## 4 2 2.5 K 0
## 5 2 2.5 K 0
## 6 6 2.5 K 0
## LATITUDE LONGITUDE LATITUDE_E LONGITUDE_ REMARKS REFNUM
## 1 3040 8812 3051 8806 1
## 2 3042 8755 0 0 2
## 3 3340 8742 0 0 3
## 4 3458 8626 0 0 4
## 5 3412 8642 0 0 5
## 6 3450 8748 0 0 6
stormdamages <- stormdata[c("STATE", "EVTYPE", "FATALITIES", "INJURIES", "PROPDMG", "PROPDMGEXP", "CROPDMG", "CROPDMGEXP")]
head(stormdamages)
## STATE EVTYPE FATALITIES INJURIES PROPDMG PROPDMGEXP CROPDMG CROPDMGEXP
## 1 AL TORNADO 0 15 25.0 K 0
## 2 AL TORNADO 0 0 2.5 K 0
## 3 AL TORNADO 0 2 25.0 K 0
## 4 AL TORNADO 0 2 2.5 K 0
## 5 AL TORNADO 0 2 2.5 K 0
## 6 AL TORNADO 0 6 2.5 K 0
unique(stormdamages$PROPDMGEXP)
## [1] K M B m + 0 5 6 ? 4 2 3 h 7 H - 1 8
## Levels: - ? + 0 1 2 3 4 5 6 7 8 B h H K m M
unique(stormdamages$CROPDMGEXP)
## [1] M K m B ? 0 k 2
## Levels: ? 0 2 B k K m M
stormdamages[stormdamages$PROPDMGEXP == "K", ]$PROPDMG <- stormdamages[stormdamages$PROPDMGEXP == "K", ]$PROPDMG * 1000
stormdamages[stormdamages$PROPDMGEXP == "M", ]$PROPDMG <- stormdamages[stormdamages$PROPDMGEXP == "M", ]$PROPDMG * 1000000
stormdamages[stormdamages$PROPDMGEXP == "m", ]$PROPDMG <- stormdamages[stormdamages$PROPDMGEXP == "m", ]$PROPDMG * 1000000
stormdamages[stormdamages$PROPDMGEXP == "B", ]$PROPDMG <- stormdamages[stormdamages$PROPDMGEXP == "B", ]$PROPDMG * 1000000000
stormdamages[stormdamages$CROPDMGEXP == "K", ]$CROPDMG <- stormdamages[stormdamages$CROPDMGEXP == "K", ]$CROPDMG * 1000
stormdamages[stormdamages$CROPDMGEXP == "k", ]$CROPDMG <- stormdamages[stormdamages$CROPDMGEXP == "k", ]$CROPDMG * 1000
stormdamages[stormdamages$CROPDMGEXP == "M", ]$CROPDMG <- stormdamages[stormdamages$CROPDMGEXP == "M", ]$CROPDMG * 1000000
stormdamages[stormdamages$CROPDMGEXP == "m", ]$CROPDMG <- stormdamages[stormdamages$CROPDMGEXP == "m", ]$CROPDMG * 1000000
stormdamages[stormdamages$CROPDMGEXP == "B", ]$CROPDMG <- stormdamages[stormdamages$CROPDMGEXP == "B", ]$CROPDMG * 1000000000
# change unit
head(stormdamages)
## STATE EVTYPE FATALITIES INJURIES PROPDMG PROPDMGEXP CROPDMG CROPDMGEXP
## 1 AL TORNADO 0 15 25000 K 0
## 2 AL TORNADO 0 0 2500 K 0
## 3 AL TORNADO 0 2 25000 K 0
## 4 AL TORNADO 0 2 2500 K 0
## 5 AL TORNADO 0 2 2500 K 0
## 6 AL TORNADO 0 6 2500 K 0
#create subset of fatality and injuries related to event type
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
injuries <- stormdamages %>% filter(stormdamages$INJURIES >0) %>% group_by(EVTYPE) %>% summarise(totalInjuries = sum(INJURIES)) %>% arrange(desc(totalInjuries))
fatalities <- stormdamages %>% filter(stormdamages$FATALITIES >0) %>% group_by(EVTYPE) %>% summarise(totalFatalities = sum(FATALITIES)) %>% arrange(desc(totalFatalities))
head(fatalities)
## # A tibble: 6 x 2
## EVTYPE totalFatalities
## <fctr> <dbl>
## 1 TORNADO 5633
## 2 EXCESSIVE HEAT 1903
## 3 FLASH FLOOD 978
## 4 HEAT 937
## 5 LIGHTNING 816
## 6 TSTM WIND 504
head(injuries)
## # A tibble: 6 x 2
## EVTYPE totalInjuries
## <fctr> <dbl>
## 1 TORNADO 91346
## 2 TSTM WIND 6957
## 3 FLOOD 6789
## 4 EXCESSIVE HEAT 6525
## 5 LIGHTNING 5230
## 6 HEAT 2100
Plot
par(mfrow = c(1,2), mar = c(10,5,3,2))
topFatalities <- fatalities[1:10,]
barplot(topFatalities$totalFatalities, las = 2, cex.axis = .7, names.arg = topFatalities$EVTYPE, main = "Top 10 Storms by Death Count", ylab = "number of fatalities", col = "sky blue")
topInjuries <- injuries[1:10,]
barplot(topInjuries$totalInjuries, las = 2, cex.axis = .7, names.arg = topInjuries$EVTYPE, main = "Top 10 Storms by Death Count", ylab = "number of fatalities", col = "sky blue")
crops <- stormdamages %>% filter(stormdamages$CROPDMG >0) %>% group_by(EVTYPE) %>% summarise(totalcrops = sum(CROPDMG)) %>% arrange(desc(totalcrops))
property <- stormdamages %>% filter(stormdamages$PROPDMG >0) %>% group_by(EVTYPE) %>% summarise(totalproperty = sum(PROPDMG)) %>% arrange(desc(totalproperty))
head(crops)
## # A tibble: 6 x 2
## EVTYPE totalcrops
## <fctr> <dbl>
## 1 DROUGHT 13972566000
## 2 FLOOD 5661968450
## 3 RIVER FLOOD 5029459000
## 4 ICE STORM 5022113500
## 5 HAIL 3025954473
## 6 HURRICANE 2741910000
head(property)
## # A tibble: 6 x 2
## EVTYPE totalproperty
## <fctr> <dbl>
## 1 FLOOD 144657709807
## 2 HURRICANE/TYPHOON 69305840000
## 3 TORNADO 56937160779
## 4 STORM SURGE 43323536000
## 5 FLASH FLOOD 16140812067
## 6 HAIL 15732267048
Merge the data
total <- merge(property, crops, by = "EVTYPE")
total$totaldamage = total$totalcrops + total$totalproperty
total <- arrange(total, desc(totaldamage))
Plot
top <- total[1:15,]
barplot(top$totaldamage/1000000000, las = 2, cex.axis = .7, names.arg = top$EVTYPE, main = "Top 15 Storms by Economics consequences", ylab = "Damage Caused (billion dollars)", col = "sky blue")