index.Rmd

---
title: "Montgomery County, PA - 911 Calls and COVID-19 Exploratory Data Analysis"
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```


```{r, include=FALSE}
library(tidyverse)
library(lubridate)
library(data.table) # need this for the "%like%"function
source("./src/funcs.R")
```

```{r, echo=FALSE}
tz_latest <- read.csv("https://storage.googleapis.com/montco-stats/tz.csv")
jhu_latest <- read.csv("https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_time_series/time_series_covid19_confirmed_US.csv")
jhu_deaths <- read.csv("https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_time_series/time_series_covid19_deaths_US.csv")
```


```{r, echo=FALSE}
# data cleanup
tz_latest$datetime <- ymd_hms(tz_latest$timeStamp, tz=NULL)
tz_latest$desc <- as.character(tz_latest$desc)
tz_latest$zip <- as.character(tz_latest$zip)

```
**This is a work in progress. Please check in regularly for more updated plots. I try to update this page at least twice a week or more.**

## Historical Trends In 911 Calls From 2017 To Present For Select Complaints

### Fever
Very notable spike in calls at the beginning of 2020 compared to previous years.  Scroll down to see a more focused plot specific for 2020 calls.
*The ongoing drop in calls is concerning because it may imply people are now reluctant to call 911 and go to hospital for fear of getting exposed to COVID-19. We see a slight pick-up beginning of May, will be interesting to see if it starts to trend up again.*


```{r, echo=FALSE, warning=FALSE, message=FALSE}

plot_by_complaint(by_complaint(tz_latest, "fever", "EMS: FEVER"), "Fever")

```

### Respiratory Emergencies
A very slight increase here which surprised me. I was expecting to see a much higher spike considering the initial understanding of SARS-CoV2 to be a respiratory virus.  

```{r, echo=FALSE, warning=FALSE, message=FALSE}
# Calls for Respiratory Emergencies
plot_by_complaint(by_complaint(tz_latest, "respiratory", "EMS: RESPIRATORY EMERGENCY"), "Respiratory Emergencies")

```

### Cardiac Emergencies
No noticeable spike here around COVID-19. However, the ongoing drop in calls is very concerning because it may imply people are now reluctant to call 911 and go to hospital for fear of getting exposed to COVID-19.  This trend is now being reported in the media as well (see references section below)


```{r, echo=FALSE, warning=FALSE, message=FALSE}
# Calls for Cardiac Emergencies
plot_by_complaint(by_complaint(tz_latest, "cardiac", "EMS: CARDIAC EMERGENCY"), "Cardiac Emergencies")

```

### Cardiac Arrests
Fairly consistent. 

```{r, echo=FALSE, warning=FALSE, message=FALSE}
# Calls for Cardiac Arrests
plot_by_complaint(by_complaint(tz_latest, "cardiac_arrests", "EMS: CARDIAC ARREST"), "Cardiac Arrests")

```

### Medical Alert Alarms
A noticeable increase in medic alert alarms starting late 2019, may be completely unrelated but would be interesting to triangulate and see if most of these are from nursing home locations?

```{r, echo=FALSE, warning=FALSE, message=FALSE}
# Calls for Medic Alert Alarms
plot_by_complaint(by_complaint(tz_latest, "medic_alert_alarms", "EMS: MEDICAL ALERT ALARM"), "Medical Alert Alarms")

```

### All EMS related Calls
A drop off in all EMS calls, likely both a strong reflection of the quarantine measures and of people not wanting to go to hospital due to fear of contracting Sars-CoV-2? (See linked media article in the references section.)

```{r, echo=FALSE, warning=FALSE, message=FALSE}
# Calls for all EMS complaints
plot_by_complaint(by_complaint_like(tz_latest, "ems", "EMS"), "EMS Services (all)")
```


### Traffic-related incidents
A total drop off in 911 calls for traffic-related, likely a strong reflection of the quarantine measures, fewer people out and about.

```{r, echo=FALSE, warning=FALSE, message=FALSE}
# Calls for traffic incidents
plot_by_complaint(by_complaint_like(tz_latest, "traffic", "Traffic"), "Traffic-Related Incidents")

```


## SINCE PRESUMED START OF COVID-19 PANDEMIC
Narrowing down the focus on the data to just 2020, from around the presumed start of COVID-19 spread in the US.


```{r, echo=FALSE, warning=FALSE, message=FALSE}

# fever
plot_by_complaint_2020(by_complaint_2020(tz_latest, "fever", "EMS: FEVER"), "Fever")

# Calls for Respiratory Emergencies
plot_by_complaint_2020(by_complaint_2020(tz_latest, "respiratory", "EMS: RESPIRATORY EMERGENCY"), "Respiratory Emergencies")

# Calls for Cardiac Emergencies
plot_by_complaint_2020(by_complaint_2020(tz_latest, "cardiac", "EMS: CARDIAC EMERGENCY"), "Cardiac Emergencies")

# Calls for Cardiac Arrests
plot_by_complaint_2020(by_complaint_2020(tz_latest, "cardiac_arrests", "EMS: CARDIAC ARREST"), "Cardiac Arrests")

# All EMS
plot_by_complaint_2020(by_complaint_like_2020(tz_latest, "ems", "EMS"), "All EMS")

```

```{r, echo=FALSE, warning=FALSE,message=FALSE}
# Extract Covid-19 data for Montco PA only
  
montco_cases <- jhu_latest %>% filter(Combined_Key == "Montgomery, Pennsylvania, US") %>%
  select(starts_with("X")) 

montco_data_cases <- montco_cases %>% rbind(montco_cases, as.vector(unlist(str_split(colnames(montco_cases), " ")))) %>%
  t() %>%
  as_tibble(.name_repair = NULL) %>%
  mutate(date = str_sub(V3, 2)) %>%
  select(date, number_of_cases = V1) %>%
  mutate(date =  mdy(date), number_of_cases = as.numeric(number_of_cases)) %>%
  filter(wday(date) == 1) # to matchup with dates in the fevers table


ggplot(montco_data_cases, aes(x=as.Date(date), y=number_of_cases)) +
  geom_bar(stat = "identity") + scale_x_date(date_breaks = "2 weeks") + 
  geom_smooth(method = "gam", color = "red", se = FALSE) +
  theme_minimal() +
  labs(title = "COVID-19 Confirmed Cases (Cumulative)", subtitle = "Montgomery County, PA (tuesdays") +
  xlab("Week") +
  ylab("Total Number of Cases") +
  theme(axis.text.x=element_text(angle=60, hjust=1))
```

```{r, echo=FALSE, warning=FALSE, message=FALSE}
# Extract Covid-19 death data for Montco PA only

montco_deaths <- jhu_deaths %>% filter(Combined_Key == "Montgomery, Pennsylvania, US") %>%
  select(starts_with("X"))

montco_data_deaths <- montco_deaths %>% rbind(montco_deaths, as.vector(unlist(str_split(colnames(montco_deaths), " ")))) %>%
  t() %>%
  as_tibble(.name_repair = NULL) %>%
  mutate(date = str_sub(V3, 2)) %>%
  select(date, number_of_deaths = V1) %>%
  mutate(date =  mdy(date), number_of_deaths = as.numeric(number_of_deaths)) %>%
  filter(wday(date) == 1) # to matchup with dates in the fevers table

#montco_data_deaths
ggplot(montco_data_deaths, aes(x=as.Date(date), y=number_of_deaths)) +
  geom_bar(stat = "identity") + scale_x_date(date_breaks = "2 weeks") + 
  geom_smooth(method = "gam", color = "red", se = FALSE) +
  theme_minimal() +
  labs(title = "COVID-19 Deaths (Cumulative)", subtitle = "Montgomery County, PA (Tuesdays") +
  xlab("Week") +
  ylab("Total Number of Deaths") +
  theme(axis.text.x=element_text(angle=60, hjust=1))
```


### Comparing Calls for Fever with Confirmed Cases 
Exploring the correlation between number of calls for fever and number of confirmed cases we see a lag in the trends which may be indicative of the delay between when COVID-19 may have started affecting people in the area vs when confirmation of diagnosis (not an individual basis but from a general big-picture perspective)
```{r, echo=FALSE, warning=FALSE}
#merge the two
tz_latest_fever_2020 = by_complaint_2020(tz_latest, "fever", "EMS: FEVER")
tz_latest_fever_2020$week <- as.Date(tz_latest_fever_2020$week)
merged <- left_join(tz_latest_fever_2020, montco_data_cases, by = c("week" = "date"))
merged$number_of_cases[is.na(merged$number_of_cases)] <- 0

```

```{r, echo=FALSE, warning=FALSE}
coeff <- 50

# A few constants
callsColor <- "blue"
casesColor <- "magenta"

ggplot(merged, aes(x=week)) +
  
  geom_line( aes(y=number_of_calls), color=callsColor) +
  scale_x_date(date_breaks = "2 weeks") +
  theme_minimal() +
  geom_line( aes(y=number_of_cases / coeff), color=casesColor) +
  
  scale_y_continuous(
    
    # Features of the first axis
    name = "Number of Calls",
    
    # Add a second axis and specify its features
    sec.axis = sec_axis(~.*coeff, name="Total Confirmed Cases")
  ) + 

  theme(
    axis.title.y = element_text(color = callsColor, size=13),
    axis.title.y.right = element_text(color = casesColor, size=13),
    axis.text.x=element_text(angle=60, hjust=1)
  ) +

  labs(title = "Number of 911 Calls for Fever compared with Confirmed (cumulative) Cases", subtitle = "Montgomery County, PA")
```

### DATA SOURCES:
* COVID-19 Case Numbers - <a href="https://coronavirus.jhu.edu/map.html">John Hopkins University</a>
* Montgomery County EMS Data - Compiled and provided by <a href=https://www.kaggle.com/mchirico/montcoalert>  Mike Chirico</a> 

## REFERENCES:
* <a href="https://www.acc.org/latest-in-cardiology/articles/2020/04/22/13/36/the-new-pandemic-threat-people-may-die-because-theyre-not-calling-911-coronavirus-disease-2019-covid-19">The New Pandemic Threat: People May Die Because They’re Not Calling 911</a>
* <a href="https://www.npr.org/sections/health-shots/2020/05/06/850454989/eerie-emptiness-of-ers-worries-doctors-where-are-the-heart-attacks-and-strokes">Eerie Emptiness Of ERs Worries Doctors: Where Are The Heart Attacks And Strokes?</a>

<a href="https://github.com/marymlucas">My Github</a> 

Last Updated:
```{r echo=FALSE}
Sys.time()
```