Synopsis

Weather events cause damages for communities and municipalities; in this essay we will explore the NOAA (National Oceanic and Atmospheric Administration) database with the aim of finding out what type of events are most harmful in terms of population health (e.g. fatalties/injuries) and economic damage. The database reports data about severe events happened across the US and spanning from 1950 to 2011; due to the 61 years long span, we will take into account the depreciation of the USD in our study: due to inflation, one dollar in 1961 is equivalent to 8.43 dollars in 2018 (see [1]), moreover, we will also explain the strategy used in order to take care of data quality issues present in the dataset. In the following, we will assume that the dataset is stored as stormdata.csv in the same directory where the R markdown is hosted; Libraries used are shown in the snippet below:

library(knitr)
library(kableExtra)
library(plyr)
library(ggplot2)
library(gridExtra)
library(grid)

Data Analisys and preparation

Let’s start by loading the dataset:

data <- read.csv("stormdata.csv")

A quick glance to the column of the data set gives:

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"

Being interested in the damage caused across the all US, and since we want to take into account the time variable in order to keep into accunt inflation, we shall adapt the dataset to our scope, we start with removing the columns that are not of interest:

columnsToBeRemoved <- c("STATE__","BGN_TIME","TIME_ZONE","COUNTY","COUNTYNAME","STATE","BGN_RANGE","BGN_AZI","BGN_LOCATI","END_DATE","END_TIME","COUNTY_END","COUNTYENDN","END_RANGE","END_AZI","END_LOCATI","LENGTH","WIDTH","F","MAG","WFO","STATEOFFIC","ZONENAMES","LATITUDE","LONGITUDE","LATITUDE_E","LONGITUDE_","REMARKS","REFNUM")
data<-data[ , !(names(data) %in% columnsToBeRemoved)]

A quick look to the dataset after the filtering gives:

head(data) %>%
  kable() %>%
  kable_styling(full_width=TRUE)
BGN_DATE EVTYPE FATALITIES INJURIES PROPDMG PROPDMGEXP CROPDMG CROPDMGEXP
4/18/1950 0:00:00 TORNADO 0 15 25.0 K 0
4/18/1950 0:00:00 TORNADO 0 0 2.5 K 0
2/20/1951 0:00:00 TORNADO 0 2 25.0 K 0
6/8/1951 0:00:00 TORNADO 0 2 2.5 K 0
11/15/1951 0:00:00 TORNADO 0 2 2.5 K 0
11/15/1951 0:00:00 TORNADO 0 6 2.5 K 0

Types of columns can be retrieved as follows:

sapply(data, class)
##   BGN_DATE     EVTYPE FATALITIES   INJURIES    PROPDMG PROPDMGEXP 
##   "factor"   "factor"  "numeric"  "numeric"  "numeric"   "factor" 
##    CROPDMG CROPDMGEXP 
##  "numeric"   "factor"

As we can see the field BGN_DATE is classifed as factor; Since we want to take into account inflation, we shall covert this column into an integer taking into account only the year, we will also change the name of the column accordingly:

data$BGN_DATE <-  as.numeric(format(as.Date(data$BGN_DAT,'%m/%d/%Y'),'%Y'))
colnames(data)[1] <- "YEAR"

A quick look to the dataset shall confirm that the data type is now numeric

##       YEAR     EVTYPE FATALITIES   INJURIES    PROPDMG PROPDMGEXP 
##  "numeric"   "factor"  "numeric"  "numeric"  "numeric"   "factor" 
##    CROPDMG CROPDMGEXP 
##  "numeric"   "factor"

Data Cleansing

Before proceding with our analisys, we shall check various aspects about data quality, we start from checking the presence of NAs in the dataset.

sapply( data, function(x) sum(is.na(x)))%>%
  kable(col.names = c("Number of NAs"), caption="Number of NAs per column") %>%
  kable_styling()
Number of NAs per column
Number of NAs
YEAR 0
EVTYPE 0
FATALITIES 0
INJURIES 0
PROPDMG 0
PROPDMGEXP 0
CROPDMG 0
CROPDMGEXP 0

No NAs are present in the dataset, let’s now have a look to the unique values per column, we start with the PROPDMGEXP:

unique(data$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

Integer numbers are taken as the exponent to be given to 10 in order to obtain the actual amount (as an example, 3 would yield to 1000 e.g. \(10^3\)); according to [2] and [3] the following mapping shall be done for other values:

propDmgValues <- unique(data$PROPDMGEXP)
propDmgMeanings <- c(3,6,0,9,6,0,0,5,6,0,4,2,3,2,7,2,0,1,8)
data.frame(propDmgValues, propDmgMeanings) %>%
  kable(col.names = c("Value in the dataset","Exponent"), caption="Exponent associated to the PROPEXP column") %>%
  kable_styling()
Exponent associated to the PROPEXP column
Value in the dataset Exponent
K 3
M 6
0
B 9
m 6
0
0 0
5 5
6 6
? 0
4 4
2 2
3 3
h 2
7 7
H 2
0
1 1
8 8

We intentionallu set the exponent as 0 (e.g. a multiplier 1) for the ‘-’,‘+’ and ‘?’ values, same is done for the feature related to the agricoltural damanage:

propDmgValues <- unique(data$CROPDMGEXP)
propDmgMeanings <- c(0,6,3,6,9,0,0,3,2)
data.frame(propDmgValues, propDmgMeanings) %>%
  kable(col.names = c("Value in the dataset","Exponent"), caption="Exponent associated to the CROPEXP column") %>%
  kable_styling()
Exponent associated to the CROPEXP column
Value in the dataset Exponent
0
M 6
K 3
m 6
B 9
? 0
0 0
k 3
2 2

We can now create a simple function converting the values above in absolute values and adding a new column

expConverted <- function(x,y) {
    if (x=='M') return(y*1000000);
    if (x=='m') return(y*1000000);
    if (x=='B') return(y*1000000000);
    if (x=='K') return(y*1000);
    if (x=='0') return(y);
    if (x=='5') return(y*100000);
    if (x=='6') return(y*1000000);
    if (x=='4') return(y*10000);
    if (x=='2') return(y*100);
    if (x=='3') return(y*1000);
    if (x=='h') return(y*100);
    if (x=='7') return(y*10000000);
    if (x=='H') return(y*100);
    if (x=='1') return(y*10);
    if (x=='8') return(y*100000000);
    if (x=='k') return(y*1000);
  return(y);
}
data$CROPDMGVALUE <- mapply(expConverted, data$CROPDMGEXP, data$CROPDMG)
data$PROPDMGVALUE <- mapply(expConverted, data$PROPDMGEXP, data$PROPDMG)

Now that we have the absolute values in the two new columns we created (CROPDMGVALUE and PROPDMGVALUE) we can drop the others:

data<-data[ , !(names(data) %in% c("PROPDMG","PROPDMGEXP","CROPDMG","CROPDMGEXP"))]

We can also sum those columns (CROPDMGVALUE and PROPDMGVALUE) in a new column (ECONOMIC_DAMAGE)

data$ECONOMIC_DAMAGE_ABSOLUTE_VALUE <- data$CROPDMGVALUE + data$PROPDMGVALUE
data<-data[ , !(names(data) %in% c("CROPDMGVALUE","PROPDMGVALUE"))]

The dataset now looks like:

head(data) %>%
  kable() %>%
  kable_styling()
YEAR EVTYPE FATALITIES INJURIES ECONOMIC_DAMAGE_ABSOLUTE_VALUE
1950 TORNADO 0 15 25000
1950 TORNADO 0 0 2500
1951 TORNADO 0 2 25000
1951 TORNADO 0 2 2500
1951 TORNADO 0 2 2500
1951 TORNADO 0 6 2500

We now have to deal with the inflation, as said above, comparing USD of 1950 and USD of 2011 does not make any sense (a damage of 10.000 USD in 1950 would be equivalent to more than 90K in 2011); in order to take care of this, we will run another function calculating the value of the ECONOMIC_DAMAGE column depending on the YEAR one (see [1] for details), for the sake of this essay, we will convert the amounts in 2011 USDs.

adjustInflation <- function(x,y) {
    if (x==1950) return(y*9.33);if (x==1951) return(y*8.65);if (x==1952) return(y*8.49);if (x==1953) return(y*8.42);if (x==1954) return(y*8.36);if (x==1955) return(y*8.39);
    if (x==1956) return(y*8.27);if (x==1957) return(y*8.00);if (x==1958) return(y*7.78);if (x==1959) return(y*7.73);if (x==1960) return(y*7.60);if (x==1961) return(y*7.52);
    if (x==1962) return(y*7.45);if (x==1963) return(y*7.35);if (x==1964) return(y*7.26);if (x==1965) return(y*7.14);if (x==1966) return(y*6.94);if (x==1967) return(y*6.73);    
    if (x==1968) return(y*6.46);if (x==1969) return(y*6.13);if (x==1970) return(y*5.80);if (x==1971) return(y*5.55);if (x==1972) return(y*5.38);if (x==1973) return(y*5.07);
    if (x==1974) return(y*4.56);if (x==1975) return(y*4.18);if (x==1976) return(y*3.95);if (x==1976) return(y*4.43);if (x==1977) return(y*4.16);if (x==1978) return(y*3.87);
    if (x==1979) return(y*3.47);if (x==1980) return(y*3.06);if (x==1981) return(y*2.77);if (x==1982) return(y*2.61);if (x==1983) return(y*2.53);if (x==1984) return(y*2.43);    
    if (x==1985) return(y*2.34);if (x==1986) return(y*2.30);if (x==1987) return(y*2.22);if (x==1988) return(y*2.13);if (x==1989) return(y*2.03);if (x==1990) return(y*1.93);      
    if (x==1991) return(y*1.85);if (x==1992) return(y*1.80);if (x==1993) return(y*1.74);if (x==1994) return(y*1.70);if (x==1995) return(y*1.65);if (x==1996) return(y*1.61);
    if (x==1997) return(y*1.57);if (x==1998) return(y*1.55);if (x==1999) return(y*1.51);if (x==2000) return(y*1.46);if (x==2001) return(y*1.42);if (x==2002) return(y*1.40);    
    if (x==2003) return(y*1.37);if (x==2004) return(y*1.33);if (x==2005) return(y*1.29);if (x==2006) return(y*1.25);if (x==2007) return(y*1.22);if (x==2008) return(y*1.17);
    if (x==2009) return(y*1.18);if (x==2010) return(y*1.16);
    return(y);
}

data$ECONOMIC_DAMAGE <- mapply(adjustInflation, data$YEAR, data$ECONOMIC_DAMAGE_ABSOLUTE_VALUE)
data<-data[ , !(names(data) %in% c("ECONOMIC_DAMAGE_ABSOLUTE_VALUE"))]

As we can see above, we also added a new column and, for consistency, removed the buffer one previously created, the dataset now looks like:

head(data) %>%
  kable() %>%
  kable_styling()
YEAR EVTYPE FATALITIES INJURIES ECONOMIC_DAMAGE
1950 TORNADO 0 15 233250
1950 TORNADO 0 0 23325
1951 TORNADO 0 2 216250
1951 TORNADO 0 2 21625
1951 TORNADO 0 2 21625
1951 TORNADO 0 6 21625

We now have to deal with the factor variable referring to the event type, a quick glance of the column EVTTYPE (in terms of different values) gives the following:

length(unique(data$EVTYPE))
## [1] 985

There are 985 different values as event types, however, in [3] we can see that the possible values shall be in the 38 elements set given below

numbers <- c(1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38)
values <- c("Dust Storm","Excessive Heat","Extreme Cold/Wind Chill","Flash Flood","Flood","Frost/Freeze","Funnel Cloud","Freezing Fog","Hail","Heat","Heavy Rain","Heavy Snow","High Surf","High Wind","Hurricane (Typhoon)","Ice Storm","Lake-Effect Snow","Lakeshore Flood","Lightning","Marine Hail","Marine High Wind","Marine Strong Wind","Marine Thunderstorm Wind","Rip Current","Seiche","Sleet","Storm Surge/Tide","Strong Wind","Thunderstorm Wind","Tornado","Tropical Depression","Tropical Storm","Tsunami","Volcanic Ash","Waterspout","Wildfire","Winter Storm","Winter Weather")
data.frame(numbers, values) %>%
  kable(col.names = c("#","Event Type"), caption="Possible Event Types according to the NOAAs (see [3] for details)") %>%
  kable_styling()
Possible Event Types according to the NOAAs (see [3] for details)
# Event Type
1 Dust Storm
2 Excessive Heat
3 Extreme Cold/Wind Chill
4 Flash Flood
5 Flood
6 Frost/Freeze
7 Funnel Cloud
8 Freezing Fog
9 Hail
10 Heat
11 Heavy Rain
12 Heavy Snow
13 High Surf
14 High Wind
15 Hurricane (Typhoon)
16 Ice Storm
17 Lake-Effect Snow
18 Lakeshore Flood
19 Lightning
20 Marine Hail
21 Marine High Wind
22 Marine Strong Wind
23 Marine Thunderstorm Wind
24 Rip Current
25 Seiche
26 Sleet
27 Storm Surge/Tide
28 Strong Wind
29 Thunderstorm Wind
30 Tornado
31 Tropical Depression
32 Tropical Storm
33 Tsunami
34 Volcanic Ash
35 Waterspout
36 Wildfire
37 Winter Storm
38 Winter Weather

So we will have to map each of the 985 possible values that are in the dataset with one of the values shown above, let’s start by finding out the frequency of event types in the dataset intentionally limiting the list to 38 events:

frequencyEvts<-count(data, 'EVTYPE')
frequencyEvts<-frequencyEvts[with(frequencyEvts, order(-freq)), ]
head(frequencyEvts,38) %>%
  kable() %>%
  kable_styling()
EVTYPE freq
238 HAIL 288661
854 TSTM WIND 219940
759 THUNDERSTORM WIND 82563
830 TORNADO 60652
147 FLASH FLOOD 54277
164 FLOOD 25326
783 THUNDERSTORM WINDS 20843
354 HIGH WIND 20212
452 LIGHTNING 15754
304 HEAVY SNOW 15708
284 HEAVY RAIN 11723
972 WINTER STORM 11433
978 WINTER WEATHER 7026
210 FUNNEL CLOUD 6839
486 MARINE TSTM WIND 6175
485 MARINE THUNDERSTORM WIND 5812
933 WATERSPOUT 3796
672 STRONG WIND 3566
917 URBAN/SML STREAM FLD 3392
953 WILDFIRE 2761
22 BLIZZARD 2719
88 DROUGHT 2488
424 ICE STORM 2006
123 EXCESSIVE HEAT 1678
368 HIGH WINDS 1533
956 WILD/FOREST FIRE 1457
206 FROST/FREEZE 1342
82 DENSE FOG 1293
980 WINTER WEATHER/MIX 1104
872 TSTM WIND/HAIL 1028
134 EXTREME COLD/WIND CHILL 1002
269 HEAT 767
342 HIGH SURF 725
844 TROPICAL STORM 690
154 FLASH FLOODING 682
133 EXTREME COLD 655
47 COASTAL FLOOD 650
434 LAKE-EFFECT SNOW 636

The strategy we want to follow to address this issue (e.g. to map the 985 possible values that are in the dataset to the 38 provided by NOAA) is to map them into one of the allowed one until the ones that are not mapped have a frequency low enough to make the error we will make neglicible in order to find out the most harmful ones (hence the choice to limit the table to the first 38); to do so, we want to add a new column to the table generated above whose aim is to state whether the event (actually its description) falls within the ones provided by NOAA.

isInNoaa <- function(x) {
    noaaValues <- c("Dust Storm","Excessive Heat","Extreme Cold/Wind Chill","Flash Flood","Flood","Frost/Freeze","Funnel Cloud","Freezing Fog","Hail","Heat","Heavy Rain","Heavy Snow","High Surf","High Wind","Hurricane (Typhoon)","Ice Storm","Lake-Effect Snow","Lakeshore Flood","Lightning","Marine Hail","Marine High Wind","Marine Strong Wind","Marine Thunderstorm Wind","Rip Current","Seiche","Sleet","Storm Surge/Tide","Strong Wind","Thunderstorm Wind","Tornado","Tropical Depression","Tropical Storm","Tsunami","Volcanic Ash","Waterspout","Wildfire","Winter Storm","Winter Weather");
    noaaValues <- sapply(noaaValues, toupper);
     if (toupper(x) %in% noaaValues) {
        return('YES');   
     } 
     return('NO');   
}
frequencyEvts$IN_NOAA <- mapply(isInNoaa, frequencyEvts$EVTYPE)

The frequency dataset now looks like:

head(frequencyEvts,38) %>%
  kable() %>%
  kable_styling()
EVTYPE freq IN_NOAA
238 HAIL 288661 YES
854 TSTM WIND 219940 NO
759 THUNDERSTORM WIND 82563 YES
830 TORNADO 60652 YES
147 FLASH FLOOD 54277 YES
164 FLOOD 25326 YES
783 THUNDERSTORM WINDS 20843 NO
354 HIGH WIND 20212 YES
452 LIGHTNING 15754 YES
304 HEAVY SNOW 15708 YES
284 HEAVY RAIN 11723 YES
972 WINTER STORM 11433 YES
978 WINTER WEATHER 7026 YES
210 FUNNEL CLOUD 6839 YES
486 MARINE TSTM WIND 6175 NO
485 MARINE THUNDERSTORM WIND 5812 YES
933 WATERSPOUT 3796 YES
672 STRONG WIND 3566 YES
917 URBAN/SML STREAM FLD 3392 NO
953 WILDFIRE 2761 YES
22 BLIZZARD 2719 NO
88 DROUGHT 2488 NO
424 ICE STORM 2006 YES
123 EXCESSIVE HEAT 1678 YES
368 HIGH WINDS 1533 NO
956 WILD/FOREST FIRE 1457 NO
206 FROST/FREEZE 1342 YES
82 DENSE FOG 1293 NO
980 WINTER WEATHER/MIX 1104 NO
872 TSTM WIND/HAIL 1028 NO
134 EXTREME COLD/WIND CHILL 1002 YES
269 HEAT 767 YES
342 HIGH SURF 725 YES
844 TROPICAL STORM 690 YES
154 FLASH FLOODING 682 NO
133 EXTREME COLD 655 NO
47 COASTAL FLOOD 650 NO
434 LAKE-EFFECT SNOW 636 YES

We have a few values marked as ‘NO’ (e.g. they are not among the values provided by NOAA), namely:

firstMapping <- function(x) {
    if (toupper(x)=='MARINE TSTM WIND') {
      return(toupper('Marine Thunderstorm Wind'));   
    }
    if (toupper(x)=='URBAN/SML STREAM FLD') {
      return(toupper('Flood'));   
    }
    if (toupper(x)=='BLIZZARD') {
      return(toupper('Heavy Snow'));   
    }
    if (toupper(x)=='DROUGHT') {
      return(toupper('Heat'));   
    }
    if (toupper(x)=='HIGH WINDS') {
      return(toupper('High Wind'));   
    }
    if (toupper(x)=='WILD/FOREST FIRE') {
      return(toupper('Wildfire'));   
    }
    if (toupper(x)=='DENSE FOG') {
      return(toupper('Freezing Fog'));   
    }
    if (toupper(x)=='WINTER WEATHER/MIX') {
      return(toupper('Winter Weather'));   
    }
    return(toupper(x));
}
data$EVTYPE <-  mapply(firstMapping, data$EVTYPE)


frequencyEvts<-count(data, 'EVTYPE')
frequencyEvts<-frequencyEvts[with(frequencyEvts, order(-freq)), ]
frequencyEvts$IN_NOAA <- mapply(isInNoaa, frequencyEvts$EVTYPE)
head(frequencyEvts,38) %>%
  kable() %>%
  kable_styling()
EVTYPE freq IN_NOAA
203 HAIL 288661 YES
772 TSTM WIND 219942 NO
679 THUNDERSTORM WIND 82564 YES
749 TORNADO 60652 YES
129 FLASH FLOOD 54277 YES
145 FLOOD 28719 YES
312 HIGH WIND 21747 YES
703 THUNDERSTORM WINDS 20843 NO
265 HEAVY SNOW 18427 YES
404 LIGHTNING 15754 YES
435 MARINE THUNDERSTORM WIND 11987 YES
245 HEAVY RAIN 11742 YES
881 WINTER STORM 11433 YES
886 WINTER WEATHER 8149 YES
181 FUNNEL CLOUD 6844 YES
865 WILDFIRE 4218 YES
846 WATERSPOUT 3796 YES
595 STRONG WIND 3569 YES
234 HEAT 3255 YES
379 ICE STORM 2006 YES
106 EXCESSIVE HEAT 1678 YES
178 FROST/FREEZE 1343 YES
168 FREEZING FOG 1339 YES
790 TSTM WIND/HAIL 1028 NO
116 EXTREME COLD/WIND CHILL 1002 YES
301 HIGH SURF 734 YES
763 TROPICAL STORM 690 YES
136 FLASH FLOODING 682 NO
115 EXTREME COLD 657 NO
41 COASTAL FLOOD 656 NO
387 LAKE-EFFECT SNOW 636 YES
149 FLOOD/FLASH FLOOD 625 NO
554 SNOW 617 NO
390 LANDSLIDE 600 NO
60 COLD/WIND CHILL 539 NO
162 FOG 538 NO
515 RIP CURRENT 470 YES
431 MARINE HAIL 442 YES

A further mapping is now advisable to get rid of the values marked as ‘NO’ in the above table:

secondMapping <- function(x) {
    if (toupper(x)=='TSTM WIND') {
      return(toupper('Thunderstorm Wind'));   
    }
    if (toupper(x)=='THUNDERSTORM WINDS') {
      return(toupper('Thunderstorm Wind'));   
    }
    return(toupper(x));
}
data$EVTYPE <-  mapply(secondMapping, data$EVTYPE)


frequencyEvts<-count(data, 'EVTYPE')
frequencyEvts<-frequencyEvts[with(frequencyEvts, order(-freq)), ]
frequencyEvts$IN_NOAA <- mapply(isInNoaa, frequencyEvts$EVTYPE)
head(frequencyEvts,38) %>%
  kable() %>%
  kable_styling()
EVTYPE freq IN_NOAA
679 THUNDERSTORM WIND 323349 YES
203 HAIL 288661 YES
748 TORNADO 60652 YES
129 FLASH FLOOD 54277 YES
145 FLOOD 28719 YES
312 HIGH WIND 21747 YES
265 HEAVY SNOW 18427 YES
404 LIGHTNING 15754 YES
435 MARINE THUNDERSTORM WIND 11987 YES
245 HEAVY RAIN 11742 YES
879 WINTER STORM 11433 YES
884 WINTER WEATHER 8149 YES
181 FUNNEL CLOUD 6844 YES
863 WILDFIRE 4218 YES
844 WATERSPOUT 3796 YES
595 STRONG WIND 3569 YES
234 HEAT 3255 YES
379 ICE STORM 2006 YES
106 EXCESSIVE HEAT 1678 YES
178 FROST/FREEZE 1343 YES
168 FREEZING FOG 1339 YES
788 TSTM WIND/HAIL 1028 NO
116 EXTREME COLD/WIND CHILL 1002 YES
301 HIGH SURF 734 YES
762 TROPICAL STORM 690 YES
136 FLASH FLOODING 682 NO
115 EXTREME COLD 657 NO
41 COASTAL FLOOD 656 NO
387 LAKE-EFFECT SNOW 636 YES
149 FLOOD/FLASH FLOOD 625 NO
554 SNOW 617 NO
390 LANDSLIDE 600 NO
60 COLD/WIND CHILL 539 NO
162 FOG 538 NO
515 RIP CURRENT 470 YES
431 MARINE HAIL 442 YES
96 DUST STORM 427 YES
11 AVALANCHE 386 NO

Now all the 38 events that are supplied by NOAA are on top as frequency (except for ‘OTHER’), for the sake of correctness and to quantify how much all others accounts in terms of damage and to have an idea of how frequent they are, we can now introdude a new event type (OTHER) to take care of all the event types that are not accounted for:

secondMapping <- function(x) {
    noaaValues <- c("Dust Storm","Excessive Heat","Extreme Cold/Wind Chill","Flash Flood","Flood","Frost/Freeze","Funnel Cloud","Freezing Fog","Hail","Heat","Heavy Rain","Heavy Snow","High Surf","High Wind","Hurricane (Typhoon)","Ice Storm","Lake-Effect Snow","Lakeshore Flood","Lightning","Marine Hail","Marine High Wind","Marine Strong Wind","Marine Thunderstorm Wind","Rip Current","Seiche","Sleet","Storm Surge/Tide","Strong Wind","Thunderstorm Wind","Tornado","Tropical Depression","Tropical Storm","Tsunami","Volcanic Ash","Waterspout","Wildfire","Winter Storm","Winter Weather");
     noaaValues <- sapply(noaaValues, toupper);
     if (toupper(x) %in% noaaValues) {
        return(toupper(x));   
     } 
     return('OTHER');   
}
data$EVTYPE <-  mapply(secondMapping, data$EVTYPE)

One last glance to the frequency:

frequencyEvts<-count(data, 'EVTYPE')
frequencyEvts<-frequencyEvts[with(frequencyEvts, order(-freq)), ]
frequencyEvts$IN_NOAA <- mapply(isInNoaa, frequencyEvts$EVTYPE)
head(frequencyEvts,38) %>%
  kable() %>%
  kable_styling()
EVTYPE freq IN_NOAA
29 THUNDERSTORM WIND 323349 YES
9 HAIL 288661 YES
30 TORNADO 60652 YES
4 FLASH FLOOD 54277 YES
5 FLOOD 28719 YES
14 HIGH WIND 21747 YES
12 HEAVY SNOW 18427 YES
18 LIGHTNING 15754 YES
23 OTHER 14414 NO
22 MARINE THUNDERSTORM WIND 11987 YES
11 HEAVY RAIN 11742 YES
37 WINTER STORM 11433 YES
38 WINTER WEATHER 8149 YES
8 FUNNEL CLOUD 6844 YES
36 WILDFIRE 4218 YES
35 WATERSPOUT 3796 YES
28 STRONG WIND 3569 YES
10 HEAT 3255 YES
15 ICE STORM 2006 YES
2 EXCESSIVE HEAT 1678 YES
7 FROST/FREEZE 1343 YES
6 FREEZING FOG 1339 YES
3 EXTREME COLD/WIND CHILL 1002 YES
13 HIGH SURF 734 YES
32 TROPICAL STORM 690 YES
16 LAKE-EFFECT SNOW 636 YES
24 RIP CURRENT 470 YES
19 MARINE HAIL 442 YES
1 DUST STORM 427 YES
27 STORM SURGE/TIDE 148 YES
20 MARINE HIGH WIND 135 YES
31 TROPICAL DEPRESSION 60 YES
26 SLEET 59 YES
21 MARINE STRONG WIND 48 YES
17 LAKESHORE FLOOD 23 YES
34 VOLCANIC ASH 23 YES
25 SEICHE 21 YES
33 TSUNAMI 20 YES

We can now proceed to generate the results of our essay.

Results

First of all, we have to group by event type injuries, fatalities and economic damage, in this regards, we will create three new dataset:

injuriesData <- aggregate(INJURIES ~ EVTYPE, data, sum);
fatalitiesData <- aggregate(FATALITIES ~ EVTYPE, data, sum);
economicData <- aggregate(ECONOMIC_DAMAGE ~ EVTYPE, data, sum);

And then we can proceed with plotting the results:

injuriesData<-injuriesData[order(-injuriesData$INJURIES),]
injuriesPlot <- ggplot(data=head(injuriesData,10), aes(x=reorder(EVTYPE, INJURIES), y=INJURIES)) +
   geom_bar(fill="red",stat="identity")  + coord_flip() + 
    ylab("Number of Injuries") + xlab("") +
    ggtitle("Health Impacts of weather events") +
    theme(legend.position="none")

fatalitiesData<-fatalitiesData[order(-fatalitiesData$FATALITIES),]
fatalitiesPlot <- ggplot(data=head(fatalitiesData,10), aes(x=reorder(EVTYPE, FATALITIES), y=FATALITIES)) +
    geom_bar(fill="yellow",stat="identity") + coord_flip() +
    ylab("Number of Fatalities") + xlab("") +
    theme(legend.position="none")

grid.arrange(injuriesPlot, fatalitiesPlot, nrow =2)

We can conclude that, as far as health impacts are concerned, tornadoes are by far the most harmful.

As per the economic damage, we can proceed in a similar way:

economicData<-economicData[order(-economicData$ECONOMIC_DAMAGE),]
economicsPlot <- ggplot(data=head(economicData,10), aes(x=reorder(EVTYPE, ECONOMIC_DAMAGE), y=ECONOMIC_DAMAGE)) +
   geom_bar(fill="blue",stat="identity")  + coord_flip() + 
    ylab("Damage in USD (2011 inflation adjusted values)") + xlab("") +
    ggtitle("Financial Impacts of weather events") +
    theme(legend.position="none")
grid.arrange(economicsPlot, nrow =1)

We can say that, among the NOAA values, flood is the most harmful by a financial point of view.

References

[1] https://www.usinflationcalculator.com/
[2] http://www5.csudh.edu/oliver/chemdata/prefix.htm
[3] https://d396qusza40orc.cloudfront.net/repdata%2Fpeer2_doc%2Fpd01016005curr.pdf