We choose WHO Life Expectancy Dataset as our main dataset in our study. The data can be accessed here. The code is read as following:
suppressMessages(library(tidyverse))
suppressMessages(library(dbplyr))
options(tibble.print_min = 5)
#Get life expectancy data
LifeExpectancy = read_csv("data/who_life_exp.csv",col_types = cols(.default = col_guess())) %>%
janitor::clean_names()Life expectancy datasets contains a total of 183
countries with year ranging from 2000 to 2016.
There are a total of 32 variables. Out of which, we
consider the variables and description is listed below:
*country: country name
*country_code: Three letter of the country id
*region: The region of the country
*year: year of all values stored
*life_expect: Life expectancy at birth measured in year.
It is continuous variable.
*life_exp60: Life expectancy at age 60 measured in year.
It is continuous variable.
*adult_mortality: Adult Mortality Rates of both sexes
(probability of dying between 15 and 60 years per 1000 population)
*infant_mort: Death rate up to age 1
*age1_4mort: Death rate between ages 1 and 4
*alcohol: Alcohol, recorded per capita (15+) consumption
(in litres of pure alcohol)
*bmi: Mean BMI (kg/\(m^2\)) (18+) (age-standardized
estimate)
*age5_19thinness: Prevalence of thinness among children
and adolescents. This is measured in a crude estimate percentage for
children with BMI < (median - 2 s.d.)
*age5_19obesity: Prevalence of obesity among children
and adolescents. This is measured in a crude estimate percentage for
children with BMI < (median - 2 s.d.)
*hepatitis: Hepatitis B (HepB) immunization coverage
among 1-year-olds (%)
*measles: Measles-containing-vaccine first-dose (MCV1)
immunization coverage among 1-year-olds (%)
*polio: Polio (Pol3) immunization coverage among
1-year-olds (%)
*diphtheria: Diphtheria tetanus toxoid and pertussis
(DTP3) immunization coverage among 1-year-olds (%)
*basic_water: Population using at least basic
drinking-water services
*doctors: Number of medical doctors (per 10,000)
*hospitals: Total density of hospitals per 100 000
population
*gni_capita: Gross national income per capita in
dollars. This is measured from GHO server
*gghe_d: Domestic general government health expenditure
(GGHE-D) as percentage of gross domestic product (GDP). This is measured
from GHO server
*che_gdp: Current health expenditure (CHE) as percentage
of gross domestic product (GDP) (%)
*une_pop: Population (thousands)
*une_infant: Mortality rate, infant (per 1,000 live
births). This is measured from GHO server
*une_life: (Our response variable) Life expectancy at
birth, total (years). This is measured from GHO server. It is contains
less missing value than life expectancy.
*une_hiv: Prevalence of HIV, total (% of population ages
15-49). This is measured from GHO server
*une_gni: GNI per capita, PPP (current international $).
This is measured from UNESCO server
*une_poverty: Government expenditure on education as a
percentage of GDP (%). This is measured from UNESCO server
*une_edu_spend: Adult literacy rate, population 15+
years, both sexes (%). This is measured from UNESCO server
*une_literacy: Mean years of schooling (ISCED 1 or
higher), population 25+ years, both sexes. This is measured from UNESCO
server
Country code dataset contains the status of countries, such as the information of the independence status and its capital city. The data is extracted from here. We include the variable specifies whether the country is a developed or developing country. The data is read in the following code and is merged by the code.
country_code =
read_csv("data/country-codes.csv", show_col_types = FALSE) %>%
select(`ISO3166-1-Alpha-3`, `Developed / Developing Countries`)
merged_data_PM = merge(LifeExpectancy, country_code, by.x = "country_code", by.y = "ISO3166-1-Alpha-3" )
PM2.5 dataset is a dataset that specifies the percentage of total population exposed to levels exceeding WHO guideline value. Such value exists over year 1960 to 2021 but most of the value before 2010 and after 2017 are missing. The dataset can be accessed here. The following code read and clean the dataset.
merge the data with dataset contain pm 2.5 information:
#read pm2.5 data
PM_dataset =
read_csv("data/pm2.5.csv", show_col_types = FALSE, skip = 4)
PM_dataset_clean =
PM_dataset %>%
janitor::clean_names() %>%
select(-(x1960:x1999),-(x2017:x2021)) %>%
pivot_longer(
x2000:x2016,
names_to = "year",
names_prefix = "x",
values_to = "PM_value"
) %>%
select(country_code, year, PM_value) %>%
mutate(year = as.numeric(year))The following code is used to merge with all dataset.
#merged data with developed/ developing countries
merged_data_d =
merged_data_PM %>%
left_join(PM_dataset_clean, by = c("country_code","year"))
Income level data is a dataset specifies the income status accessed here. It is a factor variable with four levels, namely high income, low income, lower middle income and upper middle income. The dataset is cleaned and merged in the following code.
#read income level data
Income_dataset = read_csv("data/income_level.csv", show_col_types = FALSE) %>%
janitor::clean_names() %>%
select(country_code,income_group)
#merged data with income level
merged_data_i = merged_data_d %>%
left_join( Income_dataset, by = c("country_code")) %>%
select(country_code:une_school,PM_value,income_group,`Developed / Developing Countries`)
Health expenditure dataset measures the health expenditure from government with respect to different countries. The dataset is retrieved from here. The dataset is measured in percentage of GDP or dollar value per capita. We choose percentage of GDP as our unit because the dataset measured in percentage of GDP has more complete values. Here is the coding for read, clean and merge the Health expenditure dataset to original dataset.
health_exp = read_csv("data/government_compulsory_health expenditure_from_1970_2020.csv") %>%
janitor::clean_names() %>%
filter(measure=="PC_GDP") %>%
filter(subject=="TOT") %>%
mutate("country_code"=location,
"year" = time) %>%
select(country_code,year, value)
merged_data_h =
merged_data_i %>%
left_join(health_exp, by = c("country_code","year")) %>%
mutate("health_exp" = value) %>%
select(country_code:une_school,PM_value,income_group,health_exp,`Developed / Developing Countries`)
Parliament information data is a dataset that measured the percentage of women in national parliament around 2015~2018 for different countries. Since different countries have different election dates, and there is rare data before 2015, we choose the dataset that is around 2015~2018. Then, we adjust the date to year 2015, and add the percentage of women in lower house of national parliament to the last column of the original dataset. The data origins from here and the following code cleans and merge the data:
women_in_parliament <-
read_csv("data/around_2019_women_percent_in_national_parliaments.csv", skip = 5) %>%
janitor::clean_names() %>%
select(x2, elections_3, percent_w_6) %>%
mutate("country" = x2,
"percent_w_in_lower_house" = percent_w_6,
"year" = as.double(2015)) %>%
select(country, year, percent_w_in_lower_house)
merged_data =
merged_data_h %>%
left_join(women_in_parliament, by = c("country","year")) %>%
select(country_code:une_school,PM_value,health_exp,income_group,`Developed / Developing Countries`, percent_w_in_lower_house)Store the data
write.csv(merged_data , file = "data/Merged_expectation.csv")
For the simplicity and accuracy purpose, imputation using weighted k nearest neighbour algorithm (w-kNN) is considered in my study. It is built on simple k nearest neighbour(kNN) but unlike kNN, it considers the weight of each dimension so that it circumvents the problem of neglecting the correlation between missing variables and other variables according to Ling et. al.. Additionally, such a non-parametric method does not make any assumption on the distribution of the input vector hence it is suitable for correlation analysis by mukid’s paper. In my method, we choose Euclidean distance as my measurement for the similarity between authorities as it performs better than the other metrics, namely Manhattan distance. Thus, any missing value \(x_{i,j}\) at time t can be imputed using an inverse of Euclidean distance as a weighting factor, where the subscript i corresponds to the \(i_{th}\) country and subscript j represents the \(j_{th}\) sub-indicator (or the \(j_{th}\) variable in the dataset). \
To obtain the imputation for missing value \(x_{i,j}\), euclidean distances between all underlying pairs of country \((a_i,a_k)\) at year t without missing value need to be calculated. The pair with the largest euclidean distance $d_t(a_i,a_k) $, however, will not be considered, and instead it is treated as the denominator for standardization of other available euclidean distance values. For a pair of country \((a_i,a_k)\) at year t, the euclidean distance $d_t(a_i,a_k) $ can be obtained by the following formula: \[ \begin{equation} \mathbf{d}_t(a_i,a_k)=||a_i,a_k||=\sqrt{\sum_{j=1}^J w(a_{ij}-a_{kj})^2} \end{equation} \] where J is the total number of variables or non missing data and w is the weighting function by Ling et. al.. As my available data may have different length due to different number of missing values, I choose \(w=\frac{1}{J}\) as my weight in the Euclidean distance. Then, I can normalize the Euclidean distance by dividing the largest euclidean distance in the available data set for country i. Without loss of generality, let \(\mathbf{d}_t(a_i,a_l)\) be the maximum euclidean distance for country i. Then for \(k\neq l\): \[ \begin{equation} \mathbf{D}_t(a_i,a_k)=\frac{\mathbf{d}_t(a_i,a_k)}{\mathbf{d}_t(a_i,a_{l})} \end{equation} \]
To balance the size of the imputed value for missing value \(x_{i,j}\), I consider the weight as the inverse of Euclidean distance. Such an approach is validated in Laumann et. al’s work where they believe the imputed value will be similar to countries with sub-indicator that has small Euclidean distance. With all other things being equal, countries are replaced by authorities and sub-indicators are variables in my data. Hence, having computed euclidean distances between all possible pairs and years, the missing value is measured as:
\[\begin{equation} x^t_{i,j}=\frac{1}{K}\sum_k \frac{1}{|\mathbf{D}_t(a_i,a_k)|}x^t_{k,j} \end{equation}\] where K is the total number of countries that have values at year t for variable j.
This is nicely packed in filling package in R. In our case, we choose k to be 4. This is according to the income level of countries. More detail is justified in clustering part in our analysis.
raw_mat =
merged_data %>%
select(-country, -country_code, -region,-income_group, -`Developed / Developing Countries`,percent_w_in_lower_house) %>%
as.matrix()
imputed_mat = filling::fill.KNNimpute(raw_mat, k = 2)#colnames_exp = colnames(LifeExpectancy)
Imputed_expectancy = merged_data
for (i in 1:(length(merged_data)-6)) {
Imputed_expectancy[i+3][[1]] = imputed_mat[[1]][,i]
}
head(Imputed_expectancy %>% select(-percent_w_in_lower_house))## country_code country region year life_expect life_exp60
## 1 AFG Afghanistan Eastern Mediterranean 2000 55.89618 15.14620
## 2 AFG Afghanistan Eastern Mediterranean 2001 56.52523 15.20886
## 3 AFG Afghanistan Eastern Mediterranean 2002 57.43997 15.24703
## 4 AFG Afghanistan Eastern Mediterranean 2003 57.97100 15.35566
## 5 AFG Afghanistan Eastern Mediterranean 2004 58.42978 15.44257
## 6 AFG Afghanistan Eastern Mediterranean 2005 58.90924 15.52284
## adult_mortality infant_mort age1_4mort alcohol bmi age5_19thinness
## 1 316.0496 0.098245 0.011050 0.4562055 21.7 20.6
## 2 307.2416 0.095925 0.010625 3.5850006 21.8 20.4
## 3 292.3430 0.093330 0.010130 2.8226608 21.9 20.2
## 4 286.4569 0.090470 0.009655 1.4059016 22.0 20.0
## 5 281.8943 0.087595 0.009210 1.2471070 22.1 19.8
## 6 277.1813 0.084630 0.008785 0.0162300 22.2 19.6
## age5_19obesity hepatitis measles polio diphtheria basic_water doctors
## 1 0.6 76.13455 27 24 24 27.77190 0.8403400
## 2 0.7 71.38447 37 35 33 27.79726 1.8990000
## 3 0.8 59.58910 35 36 36 29.90076 0.8168644
## 4 0.9 45.72630 39 41 41 32.00507 0.7011373
## 5 1.0 72.57173 48 50 50 34.12623 0.4960309
## 6 1.1 62.91274 50 58 58 36.26526 0.8882858
## hospitals gni_capita gghe_d che_gdp une_pop une_infant une_life une_hiv
## 1 0.4254508 802.1334 0.7777427 4.704074 20779.95 90.6 55.841 0.1
## 2 0.4229862 1060.0000 0.3869379 7.475627 21606.99 88.0 56.308 0.1
## 3 0.4191697 870.0000 0.0841800 9.443390 22600.77 85.4 56.784 0.1
## 4 2.0438096 920.0000 0.6509600 8.941260 23680.87 82.8 57.271 0.1
## 5 0.4081089 920.0000 0.5429300 9.808470 24726.68 80.1 57.772 0.1
## 6 0.4073613 1020.0000 0.5291800 9.948290 25654.28 77.4 58.290 0.1
## une_gni une_poverty une_edu_spend une_literacy une_school PM_value
## 1 931.0142 75.35382 3.831208 23.86542 1.486406 100
## 2 4217.6150 50.62891 3.657184 47.14927 1.458467 100
## 3 1084.1678 50.52512 2.576543 30.36625 2.665628 100
## 4 1431.1906 56.79214 2.701828 38.24146 2.370850 100
## 5 1524.5696 57.40101 3.567654 38.10863 2.240961 100
## 6 1157.3764 14.05812 1.623665 77.07397 2.163486 100
## health_exp income_group Developed / Developing Countries
## 1 4.262374 Low income Developing
## 2 4.170745 Low income Developing
## 3 4.134256 Low income Developing
## 4 4.057298 Low income Developing
## 5 4.060042 Low income Developing
## 6 4.066611 Low income Developingwrite.csv(Imputed_expectancy , file = "data/Imputed_expectation.csv")