About me

Background

  • Economics (cand.oecon)
  • PhD (epidemiology and biostatistics)

Interests

  • Data visualization, Statistics and epidemiology
  • Health economics
  • R/SQL/SAS

Research

  • Observational research (registry data)
  • Cancer (hematology/oncology)
  • Collaboration with pharma industry

Agenda today

  • Reproducible research
  • OECD API in R
  • Life expectancy in OECD countries
  • Healthcare funding in OECD countries
  • Combination

Reproducible research


Non-reproducible workflow

  • Download data as file from a website
  • Unstructured data handeling: copy-paste, delete outliers, comments in data
  • Unstructured analysis: Results from mix of sources
  • Hard to repeat analysis when data are updated, and to check for errors

Examples:

Rainhart & Rogoff - Growth in a Time of Debt (2010)

The replication crisis

Reproducible workflow

  • Data are fetched directly from the source (e.g. API)
    • Raw data is not edited
  • Every step is documented in code
  • Update analysis by re-running the script
  • Easy to share, review and reuse

Unfortunately incentives doesn’t always align

Important initiatives:

  • Open source tools (R, Python, Julia, LaTex, etc.)
  • Non-commercial repositories/directories (Zenodo, arXiv, ORCID)

What is an API?

API

An Application Programming Interface lets code ask a server for data directly - no clicking in a web page.

Many institutions issue public facing API

  • OECD
  • Statistics Denmark
  • DMI/CVR/DAWA
  • EuroStat
  • IMF/WB
  • Stock exchanges

Increase in use of API keys in public api - to manage access

  • AWS/Azure/Google cloud

OECD data explorer

OECD Data Explorer See if you can find life expectancy in Denmark in 2015


What do you think the lift expectancy for DK, USA and Mexico was in 2015 (and why?)


Health > Health status > Life expectancy > DK: 80.7 years, US: 78.7, MEX: 75.1

OECD API in R

We will be using the unofficial R package OECD Github

#install.packages("pak")
#install.packages("tidyverse")

library(tidyverse)
library(pak)
pak("expersso/OECD")
library(OECD)

The server returns the data, and the R OECD package turns it into a data.frame. The OECD uses the SDMX standard, so all datasets are organised in the same way: dimensions (e.g. country, age, sex) and observations (the value).

Your turn - try the API!

life <- get_dataset(dataset = "OECD.ELS.HD,DSD_HEALTH_STAT@DF_LE,1.1")


Tasks:

  1. Find out what was the expected life expectancy of a 65 year old woman in Denmark (DNK) in 2015
  2. Plot the life expectancy af a 0, 40, 60, 65 and 80 year old in Denmark (DNK) in 2015
  3. Stratify by country

Task 1.1

Use the data explorer in R or:

life |> 
  filter(REF_AREA    == "DNK", 
         MEASURE     == "LFEXP", 
         SEX         == "F", 
         AGE         == "Y65", 
         TIME_PERIOD == 2015) |> 
  select(ObsValue)
# A tibble: 1 × 1
  ObsValue
  <chr>   
1 20.6

Task 1.2

life |> 
  filter(REF_AREA    == "DNK", 
         MEASURE     == "LFEXP", 
         SEX         == "F", 
         TIME_PERIOD == 2015) |> 
  ggplot(aes(x = AGE, y = as.numeric(ObsValue))) + geom_point()

Task 1.3

life |> 
  filter(MEASURE     == "LFEXP", 
         SEX         == "F", 
         TIME_PERIOD == 2015) |> 
  ggplot(aes(x = AGE, y = as.numeric(ObsValue))) + geom_point() + facet_wrap(~REF_AREA, nrow=5)

Life expectancy in OECD countries

Example - stratification by sex

Code 
life2 <- life |> 
  select(REF_AREA, MEASURE, AGE, TIME_PERIOD, SEX, ObsValue) |> 
  filter(TIME_PERIOD == 2015,
         MEASURE == "LFEXP",
         REF_AREA %in% c("MEX","SWE","NOR","DNK","DEU","USA")) |> 
  mutate(AGE2 = parse_number(AGE),
         SEX2 = case_when(SEX=="M"~"Male",
                          SEX=="F"~"Female",
                          SEX=="_T"~"Total"),
         ObsValue = as.numeric(ObsValue))

life2 |> filter(SEX!="_T") |> 
  ggplot(aes(x=AGE2, y=ObsValue, color=REF_AREA)) + 
  geom_line(linewidth = 1) + 
  scale_x_continuous(breaks=c(0,40,60,65,80)) + 
  scale_color_brewer(palette="Paired")+
  facet_wrap(~SEX) +
  theme(panel.grid = element_blank()) +
  labs(title="Life expectancy 2015", y="years",x="Age")

Example - stratification by country

Code 
life2 |> filter(SEX!="_T") |> 
  ggplot(aes(x=AGE2, y=ObsValue, color=SEX)) + 
  geom_line(linewidth = 1) + 
  scale_x_continuous(breaks=c(0,40,60,65,80)) + 
  scale_color_brewer(palette="Paired") +
  facet_wrap(~REF_AREA) +
  theme(panel.grid = element_blank()) +
  labs(title="Expected lifespwn", y="years",x="Age")

Example - ratio male-female

Code 
life2 |> filter(SEX!="_T") |> 
  select(-SEX) |> 
  pivot_wider(names_from = SEX2, values_from = ObsValue) |> 
  mutate(rate = Male/(Female)) |> 
  ggplot(aes(x=AGE2, y=rate)) + 
  geom_line(linewidth = 1) + 
  scale_x_continuous(breaks=c(0,40,60,65,80)) + 
  scale_color_brewer(palette="Paired")+
  facet_wrap(~REF_AREA) +
  theme(panel.grid = element_blank()) +
  labs(title="Expected lifespwn", y="years",x="Age")

Example - ratio male-female through time

Code 
life3 <- life |> 
  select(REF_AREA, MEASURE, AGE, TIME_PERIOD, SEX, ObsValue) |> 
  filter(TIME_PERIOD %in% c(1990,2000,2010,2020),
         MEASURE == "LFEXP",
         REF_AREA %in% c("MEX","SWE","NOR","DNK","DEU","USA")) |> 
  mutate(AGE2 = parse_number(AGE),
         SEX2 = case_when(SEX=="M"~"Male",
                          SEX=="F"~"Female",
                          SEX=="_T"~"Total"),
         ObsValue = as.numeric(ObsValue))

life3 |> filter(SEX!="_T") |> 
  select(-SEX) |> 
  pivot_wider(names_from = SEX2, values_from = ObsValue) |> 
  mutate(rate = Male/(Female)) |> 
  ggplot(aes(x=AGE2, y=rate, linetype=TIME_PERIOD)) + 
  geom_line(linewidth = 0.5) + 
  scale_x_continuous(breaks=c(0,40,60,65,80)) + 
  scale_color_brewer(palette="Paired")+
  facet_wrap(~REF_AREA) +
  theme(panel.grid = element_blank()) +
  labs(title="Expected lifespwn", y="years",x="Age")

Example - ratio male-female through time

Code 
life3 <- life |> 
  select(REF_AREA, MEASURE, AGE, TIME_PERIOD, SEX, ObsValue) |> 
  filter(MEASURE == "LFEXP",
         REF_AREA %in% c("MEX","SWE","NOR","DNK","DEU","USA")) |> 
  mutate(AGE2 = parse_number(AGE),
         SEX2 = case_when(SEX=="M"~"Male",
                          SEX=="F"~"Female",
                          SEX=="_T"~"Total"),
         ObsValue = as.numeric(ObsValue))

life3 |> filter(SEX!="_T") |> 
  select(-SEX) |> 
  pivot_wider(names_from = SEX2, values_from = ObsValue) |> 
  mutate(rate = Male/(Female)) |> 
  ggplot(aes(x=AGE2, y=rate, color=TIME_PERIOD)) + 
  geom_line(linewidth = 0.5) + 
  scale_x_continuous(breaks=c(0,40,60,65,80)) + 
  facet_wrap(~REF_AREA) +
  theme(panel.grid = element_blank()) +
  labs(title="Expected lifespwn", y="years",x="Age")

Example - expected lifetime through time

Code 
life3 |> filter(SEX=="M") |> 
  ggplot(aes(x=AGE2, y=ObsValue, color=TIME_PERIOD)) + 
  geom_line(linewidth = 0.5) + 
  scale_x_continuous(breaks=c(0,40,60,65,80)) + 
  facet_wrap(~REF_AREA) +
  theme(panel.grid = element_blank()) +
  labs(title="Expected lifespwn", y="years",x="Age")

Healthcare funcing in OECD countries

New API request - healthcare funding

funding <- get_dataset(dataset = "OECD.ELS.HD,DSD_SHA@DF_SHA", 
                  filter = "SWE+NLD+GBR+USA+DEU+DNK+FRA.A.EXP_HEALTH.PT_B1GQ.HF122+HF121+HF3+HF2+HF11.._T._T._T..._Z")

Copy the code and do the following tasks.

Tasks:

  1. Find out what programs (FINANCING_SCHEMES) financed danish healthcare in 2015
  2. Plot the healthcare funcing for Denmark through time
  3. Stratify by country
code Financing Scheme
HF11 Government schemes
HF121 Social health insurance schemes
HF122 Compulsory private insurance schemes
HF2 Voluntary healthcare payment schemes
HF3 Household out-of-pocket payment

Task 2.1

funding |> 
  filter(REF_AREA    == "DNK", 
         TIME_PERIOD == 2015) |> 
  select(FINANCING_SCHEME, ObsValue)
# A tibble: 3 × 2
  FINANCING_SCHEME ObsValue
  <chr>            <chr>   
1 HF11             8.726   
2 HF2              0.237   
3 HF3              1.396

Task 2.2

funding |> 
  filter(REF_AREA == "DNK") |> 
  ggplot(aes(x=as.numeric(TIME_PERIOD), y=as.numeric(ObsValue), fill=FINANCING_SCHEME)) + 
  geom_area()

Task 2.3

funding |> 
  ggplot(aes(x=as.numeric(TIME_PERIOD), y=as.numeric(ObsValue), fill=FINANCING_SCHEME)) + 
  geom_area() + facet_wrap(~REF_AREA)

Nice plot healthcare expenditures

Code 
plot <- funding |> 
  mutate(obsTime = as.numeric(TIME_PERIOD)) |> 
  mutate(ObsValue = as.numeric(ObsValue)) |> 
  filter(between(obsTime, 2001, 2025)) |> 
  mutate(date = as.Date(paste0(obsTime,"-01-01"))) |> 
  select(date, ObsValue, FINANCING_SCHEME, REF_AREA) |> 
  add_row(date = as.Date("2013-01-01"), ObsValue = 0, FINANCING_SCHEME = "HF122", REF_AREA = "USA") |> 
  add_row(date = as.Date("2005-01-01"), ObsValue = 0, FINANCING_SCHEME = "HF122", REF_AREA = "NLD") |> 
  add_row(date = as.Date("2015-01-01"), ObsValue = 0, FINANCING_SCHEME = "HF122", REF_AREA = "FRA") |> 
  add_row(date = as.Date("2005-01-01"), ObsValue = 0, FINANCING_SCHEME = "HF121", REF_AREA = "FRA") |> 
  add_row(date = as.Date("2008-01-01"), ObsValue = 0, FINANCING_SCHEME = "HF122", REF_AREA = "DEU") |> 
  mutate(REF_AREA = fct_relevel(REF_AREA, "DNK","SWE","GBR")) |> 
  ggplot(aes(date, ObsValue, fill = FINANCING_SCHEME)) + 
  geom_area() + 
  facet_wrap(~REF_AREA, nrow = 1) + 
  scale_y_continuous(breaks = c(seq(0,20, by = 4)),
                     label = unit_format(unit = "%", sep = "", accuracy = 1)) +
  scale_fill_manual("",values = palblue(5), 
                    labels = c("Government\nschemes", 
                               "Social health\ninsurance schemes",
                               "Compulsory private\ninsurance schemes", 
                               "Voluntary healthcare\npayment schemes",
                               "Household out-of-pocket\npayment")) + 
  scale_x_date(breaks = seq(ymd("2002-01-01"), ymd("2023-01-01"), by = '3 years'),
               date_labels = "%y", 
               limits = as.Date(c("2001-01-01","2024-01-01"))) +
  labs(title = NULL, 
       caption = "Source: OECD - SHA (Health expenditure and financing)", 
       x = NULL, y = "GDP") + 
  theme(legend.key.height = unit(2,"line"), 
        legend.key = element_blank(), 
        legend.position = "top",
        panel.grid.major = element_blank(), 
        panel.grid.minor = element_blank())

plot
  • Healthcare is a superior good

Combination healthcare expenditures and expected lifetime

All ages

Code 
fund <- get_dataset(dataset = "OECD.ELS.HD,DSD_SHA@DF_SHA", 
                  filter = ".A.EXP_HEALTH.PT_B1GQ.HF1+HF2+HF3.._T._T._T..._Z") |> 
  select(FINANCING_SCHEME, REF_AREA, TIME_PERIOD, ObsValue) |> 
  group_by(TIME_PERIOD, REF_AREA) |> 
  summarize(healthcare = sum(as.numeric(ObsValue), na.rm=T))

combined <- life |> filter(MEASURE == "LFEXP", SEX == "F") |> 
  select(AGE, ObsValue, REF_AREA, TIME_PERIOD) |>
  left_join(fund) |> 
  mutate(expectedLife = as.numeric(ObsValue),
         year=as.numeric(TIME_PERIOD)) |> 
  filter(REF_AREA!="RUS", TIME_PERIOD<2025, !is.na(healthcare))

combined |> 
  ggplot(aes(x=healthcare, y=expectedLife, color=REF_AREA)) + 
  geom_point(size=0.5,show.legend = F) + facet_wrap(~AGE, scales="free_y") + 
  geom_path(show.legend = F) +
  geom_text(data=combined |> filter(TIME_PERIOD==2024), color="black",size=2,aes(label=REF_AREA))

Expected remaining lifetime 65 years

Code 
combined |> filter(AGE == "Y65") |> 
  ggplot(aes(x=healthcare, y=expectedLife, color=REF_AREA)) + 
  geom_point(size=0.5,show.legend = F) + facet_wrap(~AGE, scales="free_y") + 
  geom_path(show.legend = F) +
  geom_text(data=combined |> filter(TIME_PERIOD==2024,AGE == "Y65"), color="black",size=2,aes(label=REF_AREA))

Expected remaining lifetime 65 years (country)

Code 
combined |> filter(AGE == "Y65", year>1999) |> 
  ggplot(aes(x=healthcare, y=expectedLife, color=REF_AREA)) + 
  geom_point(size=0.5,show.legend = F) + facet_wrap(~REF_AREA, scales="free") + 
  geom_path(show.legend = F)+
  geom_text(data=combined |> filter(TIME_PERIOD==2024,AGE == "Y65"), color="black",size=2,aes(label=REF_AREA))

Code 
combined2 <- combined |> filter(AGE == "Y65") |> 
  filter(REF_AREA %in% c("DNK","USA","SWE","MEX","GER","DEU","FRA","JPN","POL","ITA")) |> 
  mutate(healthcare = if_else(REF_AREA=="ITA" & healthcare <4,NA,healthcare))

  combined2 |> 
  ggplot(aes(x=healthcare, y=expectedLife, color=REF_AREA)) + 
  geom_path(show.legend = F,size=0.7) +
  geom_point(size=1.2,show.legend = T, shape=21, fill="white",stroke = 1.3) + 
  scale_y_continuous(breaks=seq(0,26,2),labels = scales::label_number(suffix = " years"))+
  scale_x_continuous(breaks=seq(0,26,2),labels = scales::label_number(suffix = "%"))+
    scale_color_brewer(palette="Paired")+
  geom_text(data=combined2 |> filter(TIME_PERIOD==2024,AGE == "Y65"), color="black",size=3,aes(label=REF_AREA,y=expectedLife+0.2)) +
  labs(title = "Healthcare spending vs. expected remaining lifetime at 65 1970 - 2024", 
       caption = "Source: OECD - SHA (Health expenditure and financing) and LE (Life expectancy)", 
       x = "Healthcare spending as fraction of GDP", y = "Expected remaining lifetime at 65", 
       color="Country") + 
  theme(legend.key = element_blank(), 
        legend.position = "right",
        panel.grid.major = element_blank(), 
        panel.grid.minor = element_blank())

Thanks for your attention