Doris Chen
2023-07-28
COVID-19 pandemic has swept across the world, affecting each of us in varying degrees over the past three years. We are looking into the data from Johns Hopkins University, trying to discover the trends and patterns regarding global COVID-19 cases and deaths. Data Source
library(tidyverse)
library(lubridate)
library(maps)
library(scales)
url_in <- "https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_time_series/"
file_names <- c("time_series_covid19_confirmed_global.csv", "time_series_covid19_deaths_global.csv")
urls <- str_c(url_in, file_names) ## string concatenate, element-wise joined
global_cases <- read_csv(urls[1])
global_deaths <- read_csv(urls[2])
str(global_cases)
str(global_deaths)
The results are too long, so I set results='hide' to
avoid a super long list in the report.
There are 289 observations of 1147 variables in both tables. Most of the variables are date. We are going to convert both tables from wide to long.
Pivot data from wide to long, and remove Lat, Long from the data.
global_cases <- global_cases %>%
pivot_longer(cols = -c('Province/State', 'Country/Region', Lat, Long),
names_to = "date", values_to = "cases" ) %>%
select(-c(Lat, Long))
global_deaths <- global_deaths %>%
pivot_longer(cols = -c('Province/State', 'Country/Region', Lat, Long),
names_to = "date", values_to = "deaths" ) %>%
select(-c(Lat, Long))
Join two tables into one and filter the data with cases greater than 0.
global <- global_cases %>%
full_join(global_deaths) %>%
dplyr::rename(Country_Region = "Country/Region", Province_State = "Province/State") %>%
mutate(date = mdy(date))
global <- global %>%
filter(cases > 0)
Add population data to global.
uid_lookup_url <- "https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/UID_ISO_FIPS_LookUp_Table.csv"
uid <- read_csv(uid_lookup_url, show_col_types = FALSE) %>%
select(-c(Lat, Long_, Combined_Key, code3, iso2, iso3, Admin2))
global <- global %>%
left_join(uid, by = c("Province_State", "Country_Region")) %>%
select(-c(UID, FIPS)) %>%
select(Province_State, Country_Region, date, cases, deaths, Population)
Create a new table “global_trend” to see the global trend of new cases and new deaths by time.
global_trend <- global %>%
group_by(date) %>%
summarize(cases = sum(cases), deaths = sum(deaths)) %>%
mutate(new_cases = cases - lag(cases), new_deaths = deaths - lag(deaths)) %>%
select(date, new_cases, new_deaths) %>%
ungroup()
Create a new table “country_total” to compare total cases and deaths in different countries.
country_total <- global %>%
filter(!is.na(Population)) %>%
group_by(Country_Region, date) %>%
summarize(cases = sum(cases), deaths = sum(deaths), Population = sum(Population)) %>%
group_by(Country_Region) %>%
summarize(cases = max(cases), deaths = max(deaths), Population = max(Population)) %>%
mutate(cases_per_million = cases * 1000000 / Population,
deaths_per_million = deaths * 1000000 / Population) %>%
select(Country_Region, cases, deaths, cases_per_million, deaths_per_million, Population) %>%
ungroup()
We also create a new table “country_map”, which we need to modify some country names to match the name in the world map data from “maps” package. So we can draw a heat map later.
world_map <- map_data("world")
world_map <- subset(world_map, region != "Antarctica")
global$Country_Region[global$Province_State == "Greenland"] <- "Greenland"
country_map <- global %>%
filter(!is.na(Population)) %>%
group_by(Country_Region, date) %>%
summarize(cases = sum(cases), deaths = sum(deaths), Population = sum(Population)) %>%
group_by(Country_Region) %>%
summarize(cases = max(cases), deaths = max(deaths), Population = max(Population)) %>%
select(Country_Region, cases, deaths) %>%
ungroup()
country_map$Country_Region[country_map$Country_Region == "US"] <- "USA"
country_map$Country_Region[country_map$Country_Region == "Congo (Kinshasa)"] <- "Democratic Republic of the Congo"
country_map$Country_Region[country_map$Country_Region == "Congo (Brazzaville)"] <- "Republic of Congo"
country_map$Country_Region[country_map$Country_Region == "Korea, North"] <- "North Korea"
country_map$Country_Region[country_map$Country_Region == "Korea, South"] <- "South Korea"
country_map$Country_Region[country_map$Country_Region == "Burma"] <- "Myanmar"
country_map$Country_Region[country_map$Country_Region == "Cote d'Ivoire"] <- "Ivory Coast"
country_map$Country_Region[country_map$Country_Region == "United Kingdom"] <- "UK"
We are going to discover the data to answer a few questions.
Question 1: What is the global trend of COVID-19 new cases and new deaths over the last 3 years?
global_trend %>%
ggplot(aes(x = date, y = new_cases)) +
geom_line(aes(color = "new_cases")) +
geom_point(aes(color = "new_cases")) +
geom_line(aes(y = new_deaths, color = "new_deaths")) +
geom_point(aes(y = new_deaths, color = "new_deaths")) +
scale_y_log10() +
theme(legend.position = "bottom", legend.title = element_blank(),
axis.text.x = element_text(angle = 90)) +
labs(title = "Global Trend of COVID-19", subtitle = "New Cases & New Deaths", x = "Date", y = '')
As we can see from the plot, the new cases and new deaths was increasing rapidly during the year 2020 to 2021, then slowly decreasing from 2022 to 2023.
Question 2: What are the top 5 countries with the HIGHEST COVID-19 fatality rate?
country_total %>%
slice_max(deaths_per_million, n = 5)
## # A tibble: 5 × 6
## Country_Region cases deaths cases_per_million deaths_per_m…¹ Popul…²
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Peru 4487553 219539 136103. 6658. 3.30e7
## 2 Bulgaria 1297523 38228 186736. 5502. 6.95e6
## 3 Hungary 2196804 48762 227404. 5048. 9.66e6
## 4 Bosnia and Herzegovina 401729 16280 122448. 4962. 3.28e6
## 5 North Macedonia 346852 9662 166485. 4638. 2.08e6
## # … with abbreviated variable names ¹deaths_per_million, ²Population
Pero, Bulgaria, Hungary, Bosnia and Herzegovina, North Macedonia are the top 5 countries with the highest COVID-19 fatality rate.
Question 3: What are the top 5 countries with the LOWEST COVID-19 fatality rate?
country_total %>%
slice_min(deaths_per_million, n = 5)
## # A tibble: 5 × 6
## Country_Region cases deaths cases_per_million deaths_per_million Population
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Holy See 29 0 35847. 0 809
## 2 Tuvalu 2828 0 239824. 0 11792
## 3 Korea, North 1 6 0.0388 0.233 25778815
## 4 Burundi 53631 38 4510. 3.20 11890781
## 5 Chad 7679 194 467. 11.8 16425859
The five countries with the lowest COVID-19 fatality rate are Holy See, Tuvalu, North Korea, Burundi, Chad.
Question4: How many total deaths caused by COVID-19 in every country?
ggplot(country_map) +
geom_map(dat = world_map, map = world_map, aes(map_id = region),
fill = "white", color = "lightgrey", linewidth = 0.25) +
geom_map(map = world_map, aes(map_id = Country_Region, fill = deaths), linewidth = 0.25) +
scale_fill_gradient(low = "lightyellow", high = "darkorange", labels = comma) +
expand_limits(x = world_map$long, y = world_map$lat) +
labs(title = "Global COVID-19 Deaths", subtitle = "2020-01-22 ~ 2023-03-09") +
theme(legend.title = element_blank())
There are a few countries showing 0 deaths in the map. For example, Turkmenistan. I checked on WHO website, it says “In Turkmenistan, from 3 January 2020 to 8:14am CEST, 26 July 2023, there have been 0 confirmed cases of COVID-19 with 0 deaths, reported to WHO”.
We’ll build a simple linear model to see the correlation between death_per_million and cases_per_million
mod <- lm(deaths_per_million ~ cases_per_million, data = country_total)
summary(mod)
##
## Call:
## lm(formula = deaths_per_million ~ cases_per_million, data = country_total)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2411.4 -606.1 -390.0 465.4 5544.4
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.239e+02 1.104e+02 5.649 5.68e-08 ***
## cases_per_million 3.601e-03 4.364e-04 8.250 2.37e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1146 on 194 degrees of freedom
## Multiple R-squared: 0.2597, Adjusted R-squared: 0.2559
## F-statistic: 68.07 on 1 and 194 DF, p-value: 2.369e-14
The p-value shows the statistical significance. But the adjusted R-squared is only 0.2559, meaning the model can only explained 25.59% of the variance. Apparently, there are many more other factors affecting the fatality rate, such as vaccination rate, country policy, medical care, and so on.
From the data, we have three findings as following:
The COVID-19 pandemic has had different degrees of impact on all of us. Therefore, when analyzing the data of global cases and deaths, it is inevitable that we may bring some personal emotions, especially when I saw the huge number of deaths. The awareness reminded me to put away the bias, just look into the data, let the data tell the story.
We only take into account of the population factor in this analysis, which is far from enough. The factors influencing COVID-19 fatality rate are complex and diverse. This can be considered bias at the data level.
sessionInfo()
## R version 4.2.1 (2022-06-23)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur ... 10.16
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.2/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.2/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] scales_1.2.1 maps_3.4.1 lubridate_1.8.0 forcats_0.5.2
## [5] stringr_1.5.0 dplyr_1.0.10 purrr_1.0.1 readr_2.1.3
## [9] tidyr_1.3.0 tibble_3.1.8 ggplot2_3.4.1 tidyverse_1.3.2
##
## loaded via a namespace (and not attached):
## [1] assertthat_0.2.1 digest_0.6.29 utf8_1.2.3
## [4] R6_2.5.1 cellranger_1.1.0 backports_1.4.1
## [7] reprex_2.0.2 evaluate_0.17 highr_0.9
## [10] httr_1.4.4 pillar_1.8.1 rlang_1.0.6
## [13] curl_4.3.3 googlesheets4_1.0.1 readxl_1.4.1
## [16] rstudioapi_0.14 jquerylib_0.1.4 rmarkdown_2.20
## [19] labeling_0.4.2 googledrive_2.0.0 bit_4.0.5
## [22] munsell_0.5.0 broom_1.0.1 compiler_4.2.1
## [25] modelr_0.1.9 xfun_0.37 pkgconfig_2.0.3
## [28] htmltools_0.5.4 tidyselect_1.2.0 fansi_1.0.4
## [31] crayon_1.5.2 tzdb_0.3.0 dbplyr_2.2.1
## [34] withr_2.5.0 grid_4.2.1 jsonlite_1.8.4
## [37] gtable_0.3.1 lifecycle_1.0.3 DBI_1.1.3
## [40] magrittr_2.0.3 cli_3.6.0 stringi_1.7.8
## [43] vroom_1.6.1 cachem_1.0.6 farver_2.1.1
## [46] fs_1.5.2 xml2_1.3.3 bslib_0.4.0
## [49] ellipsis_0.3.2 generics_0.1.3 vctrs_0.5.2
## [52] tools_4.2.1 bit64_4.0.5 glue_1.6.2
## [55] hms_1.1.2 parallel_4.2.1 fastmap_1.1.0
## [58] colorspace_2.0-3 gargle_1.2.1 rvest_1.0.3
## [61] knitr_1.40 haven_2.5.1 sass_0.4.2