MATH2349 Data Wrangling

Assignment 2

Required packages

These are the packages required for this analysis.

library(readr)
library(tidyr)
library(dplyr)
library(editrules)
library(outliers)
library(forecast)

Executive Summary

Data pre-processing was done on three dataset to allow the exploration of the effect of gender on life expectancy across countries and years.

• Data: Data were firstly obtained from the World Bank and imported to RStudio for the data pre-processing. As the data was originally untidy, the data were tidied to transform the data from wide to long format as the original data has variable values as the column names. The data was cleaned prior to merging them together.

• Understand: After merging, the structure and metadata of the dataset are then checked to allow me to better understand the dataset. Data type conversion were also done to convert variables into the appropriate data type.

• Tidy and Manipulate: In addition to the tidying that were done prior to merging the data, a new variable called ‘life expectancy gap’ was also created to allow the investigation of gender differences on life expectancy across different countries and age.

• Scan: The data are then scanned for missing values, special values, or errors. Several observations with missing values were dropped because they were deemed not important for the analysis. No special values or errors were found in the data. The data was also scanned for outliers. Several outliers were found in the data but were not removed because they can provide valuable information and insight for the analysis.

• Transform: To address the potential impact or biases that may happen because the outliers were not removed, transformation was done on the data. After the transformation using cube and square root transformation, the previously skewed distributions are now normally distributed.

Data

Dataset 1: Female Life Expectancy at Birth

This dataset was downloaded as a csv file from The World Bank (URL: https://data.worldbank.org/indicator/SP.DYN.LE00.FE.IN?view=chart). The dataset consisted of the female life expectancy at birth across 264 countries from the year 1960 to 2018. The dataset contains 63 variables, but only 61 variables are relevant for this analysis. The variables are:

• Country Name: the name of the country
• Country Code: the ISO code used to represent the country
• 59 variables for the female life expectancy for each year from 1960 to 2018. For example, the variable ‘1960’ for the female life expectancy in 1960

As the dataset is in csv file format, the data was imported into R as ‘female_life_exp’ using the read_csv function from the readr library. The variables that were not relevant for this analysis were removed and then the head function was used to print first 5 rows of the dataset.

# Import the female life expectancy data from the current directory, remove irrelevant variables, and print the first 5 rows
female_life_exp <- read_csv("API_SP.DYN.LE00.FE.IN_DS2_en_csv_v2_1347350.csv",
                            skip = 4)

female_life_exp <- female_life_exp %>% 
  select(-(3:4))

head(female_life_exp, 5)

Dataset 2: Male Life Expectancy at Birth

This dataset was downloaded as a csv file from The World Bank (URL: https://data.worldbank.org/indicator/SP.DYN.LE00.MA.IN?view=chart). The dataset consisted of the male life expectancy at birth across 264 countries from the year 1960 to 2018. The dataset contains 63 variables, but only 61 variables are relevant for this analysis. The variables are:

• Country Name: the name of the country
• Country Code: the ISO code used to represent the country
• 59 variables for the male life expectancy for each year from 1960 to 2018. For example, the variable ‘1960’ for the female life expectancy in 1960

As the dataset is in csv file format, the data was imported into R as ‘male_life_exp’ using the read_csv function from the readr library. The variables that were not relevant for this analysis were removed and then the head function was used to print first 5 rows of the dataset.

# Import the male life expectancy data from the current directory, remove irrelevant variables, and print the first 5 rows
male_life_exp <- read_csv("API_SP.DYN.LE00.MA.IN_DS2_en_csv_v2_1346754.csv", 
                          skip = 4)

male_life_exp <- male_life_exp %>%
  select(-(3:4))

head(male_life_exp, 5)

Dataset 3: Countries Metadata

This dataset was obtained as part of the csv file from the zip file of the first and second dataset (female and male life expectancy dataset) downloaded from The World Bank. The dataset consisted metadata information of 264 countries. The dataset contains 5 variables, but only 3 variables are relevant for this analysis. The variables are:

• Country Code: the ISO code used to represent the country
• Region: the region of the country
• Income Group: the income group of the country from ‘Low income’, ‘Lower middle income’, ‘Upper middle income’, to ‘High income’

As the dataset is in csv file format, the data was imported into R as ‘metadata’ using the read_csv function from the readr library. The variables that were not relevant for this analysis were removed and then the head function was used to print first 5 rows of the dataset.

# Import the countries metadata from the current directory, remove irrelevant variables, and print the first 5 rows
metadata <- read_csv("Metadata_Country_API_SP.DYN.LE00.FE.IN_DS2_en_csv_v2_1347350.csv")

metadata <- metadata %>%
  select(1:3)

head(metadata, 5)

As dataset 1 and 2 both have the years as the column names and are untidy, merging them together now would make it more difficult to clean the data. Hence, both dataset 1 and 2 would first be tidied and manipulated (Tidy & Manipulate Data I) before merging all three dataset together.

Tidy & Manipulate Data I

Tidying the Female and Male Life Expectancy Data

The data for both the female and male life expectancy are untidy because the column names 1960 to 2018 are not the names of the variables but it represent the year of each life expectancy values, and should be the values of a variable called ‘year’ instead. Hence, the data is tidied by gathering and transforming these column names into a variable called ‘year’. The gather() function was used to transform the data from wide to long format. The head() function is also used to show the data after the transformation.

# Transforming the female and male life expectancy data from wide to long format
female_life_exp_long <- female_life_exp %>% 
  gather(key="Year", value="Female_Life_Expectancy", 3:61)

head(female_life_exp_long, 5)

male_life_exp_long <- male_life_exp %>% 
  gather(key="Year", value="Male_Life_Expectancy", 3:61)

head(male_life_exp_long, 5)

Merging the data

After tidying the data, the three dataset would then be combined together using left_join() function. As the female and male life expectancy data shares the same variable other than the life expectancy (‘Female_Life_Expectancy’ and ‘Male_Life_Expectancy’), the data would be joined using ‘Country Name’, ‘Country Code’, and ‘Year’ as the key variables. The new combined data would then be joined again with the metadata using ‘Country Code’ as the key variable. The head() function is then used to show the data after the merge

# Merging the three dataset together
life_exp_comb <- left_join(female_life_exp_long, male_life_exp_long, by = c("Country Name", "Country Code", "Year"))

life_exp_comb <- left_join(life_exp_comb, metadata, by='Country Code')

# Show the data after the merge
head(life_exp_comb, 5)

Renaming the Columns

Several columns in the data were renamed to make a consistent format amongst the column names and to replace the spaces between words with the underscore symbol (_).

# Renaming the columns
life_exp_comb <- rename(life_exp_comb,
                        'Country_Name'='Country Name',
                        'Country_Code'='Country Code',
                        'Income_Group'='IncomeGroup')

#Show the data after renaming
head(life_exp_comb, 5)

Understand

The str() function is used to check the structure of the data and the type of variables in the data. The str() function showed that there are 15,576 observations and 7 variables in the dataframe.

# Check the structure of the data
str(life_exp_comb)

tibble [15,576 x 7] (S3: tbl_df/tbl/data.frame)
 $ Country_Name          : chr [1:15576] "Aruba" "Afghanistan" "Angola" "Albania" ...
 $ Country_Code          : chr [1:15576] "ABW" "AFG" "AGO" "ALB" ...
 $ Year                  : chr [1:15576] "1960" "1960" "1960" "1960" ...
 $ Female_Life_Expectancy: num [1:15576] 67.1 33.3 38.8 63.2 NA ...
 $ Male_Life_Expectancy  : num [1:15576] 64.1 31.7 36.3 61.3 NA ...
 $ Region                : chr [1:15576] "Latin America & Caribbean" "South Asia" "Sub-Saharan Africa" "Europe & Central Asia" ...
 $ Income_Group          : chr [1:15576] "High income" "Low income" "Lower middle income" "Upper middle income" ...

The str() function also showed the data type of each variables. The variables in the dataset were in the character and numeric data type. It was found that the variable ‘Year’ is in character variable type which is incorrect as it should be integer data type instead. Additionally, ‘Region’ and ‘Income_Group’ are also in character variable type even though it would be more appropriate for these variables to be in factor data type. Hence, the following data type conversion were done:

• ‘Year’ variable is converted into integer data type using as.integer function.
• ‘Region’ variable is converted into labeled factor data type using factor function with the levels set as each of the regions.
• ‘Income_Group’ variable is converted into ordered factor data type using factor function with the levels ordered from low income, lower middle income, upper middle income, to high income.

After doing data type conversion for the three variables, the str() function was used again to re-check the structure of the variables.

# Converting 'Year' variable into integer
life_exp_comb$Year <- as.integer(life_exp_comb$Year) 

# Converting 'Region' variable into labeled factor
life_exp_comb$Region <- factor(life_exp_comb$Region,
                               levels=c('Latin America & Caribbean',
                                        'South Asia',
                                        'Sub-Saharan Africa',
                                        'Europe & Central Asia',
                                        'Middle East & North Africa',
                                        'East Asia & Pacific',
                                        'North America'))

# Converting 'Income_Group' variable into ordered factor
life_exp_comb$Income_Group<- factor(life_exp_comb$Income_Group,
                                    levels=c('Low income',
                                             'Lower middle income',
                                             'Upper middle income',
                                             'High income'),
                                    labels = c('Low',
                                               'Lower middle',
                                               'Upper middle',
                                               'High'),
                                    ordered=TRUE)

# Checking the structure of the data again
str(life_exp_comb)

tibble [15,576 x 7] (S3: tbl_df/tbl/data.frame)
 $ Country_Name          : chr [1:15576] "Aruba" "Afghanistan" "Angola" "Albania" ...
 $ Country_Code          : chr [1:15576] "ABW" "AFG" "AGO" "ALB" ...
 $ Year                  : int [1:15576] 1960 1960 1960 1960 1960 1960 1960 1960 1960 1960 ...
 $ Female_Life_Expectancy: num [1:15576] 67.1 33.3 38.8 63.2 NA ...
 $ Male_Life_Expectancy  : num [1:15576] 64.1 31.7 36.3 61.3 NA ...
 $ Region                : Factor w/ 7 levels "Latin America & Caribbean",..: 1 2 3 4 4 NA 5 1 4 6 ...
 $ Income_Group          : Ord.factor w/ 4 levels "Low"<"Lower middle"<..: 4 1 2 3 4 NA 4 3 3 3 ...

All columns are now in the proper data type, and the dataset now have characters, numeric, and factors variables. Next, the metadata of the dataset is checked using the attributes() function. The metadata are the column names, row names, and class of the dataframe object.

# Checking the metadata of the data
attributes(life_exp_comb)

$row.names
   [1]    1    2    3    4    5    6    7    8    9   10   11   12
  [13]   13   14   15   16   17   18   19   20   21   22   23   24
  [25]   25   26   27   28   29   30   31   32   33   34   35   36
  [37]   37   38   39   40   41   42   43   44   45   46   47   48
  [49]   49   50   51   52   53   54   55   56   57   58   59   60
  [61]   61   62   63   64   65   66   67   68   69   70   71   72
  [73]   73   74   75   76   77   78   79   80   81   82   83   84
  [85]   85   86   87   88   89   90   91   92   93   94   95   96
  [97]   97   98   99  100  101  102  103  104  105  106  107  108
 [109]  109  110  111  112  113  114  115  116  117  118  119  120
 [121]  121  122  123  124  125  126  127  128  129  130  131  132
 [133]  133  134  135  136  137  138  139  140  141  142  143  144
 [145]  145  146  147  148  149  150  151  152  153  154  155  156
 [157]  157  158  159  160  161  162  163  164  165  166  167  168
 [169]  169  170  171  172  173  174  175  176  177  178  179  180
 [181]  181  182  183  184  185  186  187  188  189  190  191  192
 [193]  193  194  195  196  197  198  199  200  201  202  203  204
 [205]  205  206  207  208  209  210  211  212  213  214  215  216
 [217]  217  218  219  220  221  222  223  224  225  226  227  228
 [229]  229  230  231  232  233  234  235  236  237  238  239  240
 [241]  241  242  243  244  245  246  247  248  249  250  251  252
 [253]  253  254  255  256  257  258  259  260  261  262  263  264
 [265]  265  266  267  268  269  270  271  272  273  274  275  276
 [277]  277  278  279  280  281  282  283  284  285  286  287  288
 [289]  289  290  291  292  293  294  295  296  297  298  299  300
 [301]  301  302  303  304  305  306  307  308  309  310  311  312
 [313]  313  314  315  316  317  318  319  320  321  322  323  324
 [325]  325  326  327  328  329  330  331  332  333  334  335  336
 [337]  337  338  339  340  341  342  343  344  345  346  347  348
 [349]  349  350  351  352  353  354  355  356  357  358  359  360
 [361]  361  362  363  364  365  366  367  368  369  370  371  372
 [373]  373  374  375  376  377  378  379  380  381  382  383  384
 [385]  385  386  387  388  389  390  391  392  393  394  395  396
 [397]  397  398  399  400  401  402  403  404  405  406  407  408
 [409]  409  410  411  412  413  414  415  416  417  418  419  420
 [421]  421  422  423  424  425  426  427  428  429  430  431  432
 [433]  433  434  435  436  437  438  439  440  441  442  443  444
 [445]  445  446  447  448  449  450  451  452  453  454  455  456
 [457]  457  458  459  460  461  462  463  464  465  466  467  468
 [469]  469  470  471  472  473  474  475  476  477  478  479  480
 [481]  481  482  483  484  485  486  487  488  489  490  491  492
 [493]  493  494  495  496  497  498  499  500  501  502  503  504
 [505]  505  506  507  508  509  510  511  512  513  514  515  516
 [517]  517  518  519  520  521  522  523  524  525  526  527  528
 [529]  529  530  531  532  533  534  535  536  537  538  539  540
 [541]  541  542  543  544  545  546  547  548  549  550  551  552
 [553]  553  554  555  556  557  558  559  560  561  562  563  564
 [565]  565  566  567  568  569  570  571  572  573  574  575  576
 [577]  577  578  579  580  581  582  583  584  585  586  587  588
 [589]  589  590  591  592  593  594  595  596  597  598  599  600
 [601]  601  602  603  604  605  606  607  608  609  610  611  612
 [613]  613  614  615  616  617  618  619  620  621  622  623  624
 [625]  625  626  627  628  629  630  631  632  633  634  635  636
 [637]  637  638  639  640  641  642  643  644  645  646  647  648
 [649]  649  650  651  652  653  654  655  656  657  658  659  660
 [661]  661  662  663  664  665  666  667  668  669  670  671  672
 [673]  673  674  675  676  677  678  679  680  681  682  683  684
 [685]  685  686  687  688  689  690  691  692  693  694  695  696
 [697]  697  698  699  700  701  702  703  704  705  706  707  708
 [709]  709  710  711  712  713  714  715  716  717  718  719  720
 [721]  721  722  723  724  725  726  727  728  729  730  731  732
 [733]  733  734  735  736  737  738  739  740  741  742  743  744
 [745]  745  746  747  748  749  750  751  752  753  754  755  756
 [757]  757  758  759  760  761  762  763  764  765  766  767  768
 [769]  769  770  771  772  773  774  775  776  777  778  779  780
 [781]  781  782  783  784  785  786  787  788  789  790  791  792
 [793]  793  794  795  796  797  798  799  800  801  802  803  804
 [805]  805  806  807  808  809  810  811  812  813  814  815  816
 [817]  817  818  819  820  821  822  823  824  825  826  827  828
 [829]  829  830  831  832  833  834  835  836  837  838  839  840
 [841]  841  842  843  844  845  846  847  848  849  850  851  852
 [853]  853  854  855  856  857  858  859  860  861  862  863  864
 [865]  865  866  867  868  869  870  871  872  873  874  875  876
 [877]  877  878  879  880  881  882  883  884  885  886  887  888
 [889]  889  890  891  892  893  894  895  896  897  898  899  900
 [901]  901  902  903  904  905  906  907  908  909  910  911  912
 [913]  913  914  915  916  917  918  919  920  921  922  923  924
 [925]  925  926  927  928  929  930  931  932  933  934  935  936
 [937]  937  938  939  940  941  942  943  944  945  946  947  948
 [949]  949  950  951  952  953  954  955  956  957  958  959  960
 [961]  961  962  963  964  965  966  967  968  969  970  971  972
 [973]  973  974  975  976  977  978  979  980  981  982  983  984
 [985]  985  986  987  988  989  990  991  992  993  994  995  996
 [997]  997  998  999 1000
 [ reached getOption("max.print") -- omitted 14576 entries ]

$names
[1] "Country_Name"           "Country_Code"          
[3] "Year"                   "Female_Life_Expectancy"
[5] "Male_Life_Expectancy"   "Region"                
[7] "Income_Group"          

$class
[1] "tbl_df"     "tbl"        "data.frame"

Tidy & Manipulate Data II

To allow the investigation of gender difference on life expectancy across different countries and years, the mutate() function will be used to create a new variable called ‘Life_Expectancy_Gap’ by subtracting the ‘Female_Life_Expectancy’ with the ‘Male_Life_Expectancy’. The head() function will be used to show the data after creating the new variable.

# Create a new variable for the life expectancy gap between female and male
life_exp_df <- mutate(life_exp_comb, 
                      Life_Expectancy_Gap = Female_Life_Expectancy - Male_Life_Expectancy)

# Show the data with the newly created variable
head(life_exp_df, 5)

Scan I

Checking and Handling Missing Values

The data was firstly scanned to check for missing values. The function colSums() and is.na() were used to check the total number of missing values in each columns.

# Check the total number of missing values for each columns
colSums(is.na(life_exp_df))

          Country_Name           Country_Code                   Year 
                     0                      0                      0 
Female_Life_Expectancy   Male_Life_Expectancy                 Region 
                  1331                   1331                   2773 
          Income_Group    Life_Expectancy_Gap 
                  2773                   1331 

Several columns were found to have missing values. Inspection of the rows with missing values in the ‘Female_Life_Expectancy’ column would be done first. Subsetting of the data for rows with missing values in the ‘Female_Life_Expectancy’ column were done using the subset and is.na functions. The head function is then used to display the subsetted rows with missing values.

# Subset the rows with missing value in the female life expectancy column
subset_female_NA <- life_exp_df %>% subset(is.na(Female_Life_Expectancy))

# Show the subsetted data
head(subset_female_NA, 5)

Inspection of the subsetted data above showed that rows with missing values in the ‘Female_Life_Expectancy’ column are also the rows with missing values in the ‘Male_Life_Expectancy’ and ‘Life_Expectancy_Gap’ columns. Further inspection were done to check the total number of rows with missing life expectancy data for each country. This inspection was done using the group_by function to categorise the subsetted data based on the country name, and then counting the total number of rows with missing life expectancy value (the ‘Female_Life_Expectancy’ column was used in the code but previous inspection showed that rows with missing values for female, male , and life expectancy gap are all the same). The data from the inspection was then printed.

# Group the data by country and calculate the total missing life expectancy values 
subset_female_NA_count <- subset_female_NA %>% 
  group_by(Country_Name) %>% 
  summarise(Count_NA = sum(is.na(Female_Life_Expectancy)))

# Print the calculated data
subset_female_NA_count

The inspection showed that there are 29 countries with missing life expectancy values. The ‘Count_NA’ represents the number of rows with missing life expectancy values for a particular country. As the rows represents the year from 1960 to 2018, the ‘Count_NA’ showed the number of years that the countries have missing data for their citizens life expectancy values. For countries with ‘Count_NA’ of 59, this suggested that the World Bank may not have life expectancy data for these countries at all. While for countries with ‘count_NA’ smaller than 59, this suggested that the World Bank may not have the life expectancy data for these countries only for certain number of years.

As using the life expectancy value of the countries from another year to replace the missing values on other rows may cause the data to be inaccurate, biased, and unable to depict the real changes of life expectancy of countries across different gender and years, it was decided that these rows with missing values would be dropped instead. The complete.cases function was then used to drop the rows with missing life expectancy values. The colSums() function was then used again to check the number of missing values on the life expectancy columns after dropping the rows.

# Drop the rows with missing life expectancy values
life_exp_df <- life_exp_df[complete.cases(life_exp_df$Female_Life_Expectancy),]

# Check for missing values again after dropping the rows
colSums(is.na(life_exp_df))

          Country_Name           Country_Code                   Year 
                     0                      0                      0 
Female_Life_Expectancy   Male_Life_Expectancy                 Region 
                     0                      0                   2714 
          Income_Group    Life_Expectancy_Gap 
                  2714                      0 

Next, inspection of the rows with missing values in the ‘Region’ column would be done. Subsetting of the data for rows with missing values in the ‘Region’ column were done using the subset and is.na functions. The head function is then used to display the subsetted rows with missing values.

# Subset the rows with missing value in the region column
subset_region_NA <- life_exp_df %>% subset(is.na(Region))

# Show the subsetted data
head(subset_region_NA, 5)

Inspection of the subsetted data above showed that rows with missing values in the ‘Region’ column are also the rows with missing values in the ‘Income_Group’ column. Further inspection were done to check the total number of rows with missing region and income group data for each country. This inspection was done using the group_by function to categorise the subsetted data based on the country name, and then counting the total number of rows with missing region and income group data (the ‘Region’ column was used in the code but previous inspection showed that rows with missing values for region and income group are the same). The data from the inspection was then printed.

# Group the data by country and calculate the total missing region and income group
subset_region_NA_count <- subset_region_NA %>% 
  group_by(Country_Name) %>% 
  summarise(Count_NA = sum(is.na(Region)))

# Print the calculated data
subset_region_NA_count

The inspection showed that there are 46 countries with missing region and income group. The ‘Count_NA’ represents the number of rows with missing region and income group for a particular country. Inspection of the countries in the ‘Country_Name’ column suggested that these countries are not actual countries and that they are groups of several countries instead. For example, ‘Arab World’ is a group of 22 Arab Countries (e.g., Algeria, Bahrain, Comoros, etc.) that are members of the Arab league.

As we are only interested in examining changes of life expectancy between gender across different countries and years, it was decided that these rows with missing values would be dropped instead. The complete.cases function was then used to drop the rows with missing region and income group. The colSums() function was then used again to check the number of missing values on the life expectancy columns after dropping the rows.

# Drop the rows with missing life expectancy values
life_exp_df <- life_exp_df[complete.cases(life_exp_df$Region),]

# Check for missing values again after dropping the rows
colSums(is.na(life_exp_df))

          Country_Name           Country_Code                   Year 
                     0                      0                      0 
Female_Life_Expectancy   Male_Life_Expectancy                 Region 
                     0                      0                      0 
          Income_Group    Life_Expectancy_Gap 
                     0                      0 

Checking and Handling Special Values

The data was then checked for special values: infinite or Not a Number (NaN) values. A function called is.special was created to check every cell numerical cell whether they have infinite or NaN values. As the is.special function only accept vectorial input, sapply function was used alongside the is.special function to allow checking of the pressence of special values in the dataframe. Additionally, another function was written inside the sapply to calculate the total special values for each column in the dataframe. The inspection showed that there are no special values in the data.

# Create the function to check for special (infinite or NaN) values
is.special <- function(x){
  if (is.numeric(x)) (is.infinite(x) | is.nan(x))
}

# Show the total special values of each columns
sapply(life_exp_df, function(x) sum(is.special(x)))

          Country_Name           Country_Code                   Year 
                     0                      0                      0 
Female_Life_Expectancy   Male_Life_Expectancy                 Region 
                     0                      0                      0 
          Income_Group    Life_Expectancy_Gap 
                     0                      0 

Checking and Handling Errors

The data was then checked for errors, particularly negative value in the female and male life expectancy columns. The life expectancy gap column is allowed to be negative because negative values would indicate that the male life expectancy is higher than the female life expectancy by x amount of years. While positive values would indicate that the female life expectancy is higher than the male life expectancy by x amount of years.

The editset function from the editrules package was used to define the rule for the female and male life expectancy columns. The data was then checked against the rule previously set using the violatedEdits function. The summary function was then used to display the results of the inspection. The results showed that there was no error in the data (NULL), showing that there was no rows that violated the rules and that life expectancy for female and male are all positive values.

# Define the rule that female and male life expectancy cannot be negative
(Rule <- editset(c("Female_Life_Expectancy >= 0", "Male_Life_Expectancy >= 0")))

Edit set:
num1 : 0 <= Female_Life_Expectancy
num2 : 0 <= Male_Life_Expectancy 

# Check the data against the rule
Error_Check <- violatedEdits(Rule, life_exp_df)

# Show the results of the check
summary(Error_Check)

NULL

Scan II

The data was also checked for outliers. The three life expectancy variables (female and male life expectancy, and life expectancy gap) were the only variables that should be checked for outliers as they are the only variables that may have outliers in the data. Firstly, boxplot function was used to create box plot visualisations to detect outliers.

The boxplots indicated the presence of several outliers in the data for all three variables. Inspection of the female and male life expectancy boxplots indicated that the outliers were all below the lower outlier fence. Additionally, the boxplots for female and male life expectancy seem to have a longer lower whisker as compared to the upper whisker, suggesting that the female and male life expectancy are not normally distributed and are left skewed. Further inspection of the data normality would be done on the next section (Transform).

# Create boxplots for female and male life expectancy
par(mfrow=c(1,2))
life_exp_df$Female_Life_Expectancy %>%  boxplot(main="Female Life Expectancy", ylab="Life Expectancy (year)", col = "grey")
life_exp_df$Male_Life_Expectancy %>%  boxplot(main="Male Life Expectancy", ylab="Life Expectancy (year)", col = "grey")

For the life expectancy gap, the boxplots indicated that there are presence of outliers below the lower outlier fence and above the upper outlier fence. More outliers were found in the upper outlier fence as compared to the lower outlier fence. Inspection of the boxplot seem to suggest that the life expectancy gap is approximately normally distributed.

# Create boxplot for life expectancy gap
par(mfrow=c(1,1))
life_exp_df$Life_Expectancy_Gap %>%  boxplot(main="Life Expectancy Gap", ylab="Life Expectancy (year)", col = "grey")

Further inspection of the outliers were done using the z-score method. The scores function from the outlier package was used to calculate the z-score of each variables. The summary function was then used to display the z-score summary statistics.

# Calculate and show the z-score summary statistics for female life expectancy
z_female <- life_exp_df$Female_Life_Expectancy %>%  scores(type = "z")
z_female %>% summary()

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
-3.6653 -0.7089  0.2766  0.0000  0.7634  1.7564 

# Calculate and show the z-score summary statistics for male life expectancy
z_male <- life_exp_df$Male_Life_Expectancy %>%  scores(type = "z")
z_male %>% summary()

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
-4.1489 -0.6818  0.2456  0.0000  0.7433  2.0251 

# Calculate and show the z-score summary statistics for life expectancy gap
z_gap <- life_exp_df$Life_Expectancy_Gap %>%  scores(type = "z")
z_gap %>% summary()

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
-2.9937 -0.7396 -0.1389  0.0000  0.6888  4.5238 

When using z-score to detect outlier, an observation is considered as an outlier if the absolute value of the z-score is greater than 3. From the summary output, it can be seen that the minimum z score for the female and male life expectancy have an absolute value that are greater than 3. While for the life expectancy gap, it was only only found that the maximum z score value have an absolute value that are greater than 3, while the minimum value still has an absolute value of smaller than 3. This was in contrast with the boxplot that indicated there are values that were below the lower outlier fence.

The length and which functions are then used to get the total number of outliers for each variables according to the z-score method. It was found that there are 13 outliers in the female life expectancy column, 18 outliers in the male life expectancy column, and 76 outliers in the life expectancy gap column. Further inspection of the outliers would be done for each variables.

# Get the total number of outliers for female life expectancy
length(which(abs(z_female) > 3))

[1] 13

# Get the total number of outliers for male life expectancy
length(which(abs(z_male) > 3))

[1] 18

# Get the total number of outliers for life expectancy gap
length(which(abs(z_gap) > 3))

[1] 76

Investigation of the outliers in the ‘Female_Life_Expectancy’ column will be done first. The outliers in the female life expectancy column were selected using which function to get the rows with z-score that are higher than 3. These rows are then printed for investigation. It was found that the countries with the outliers were mostly low income or lower middle income countries.

# Select the rows that are outliers based on the z-score value
female_outliers <- life_exp_df[which(abs(z_female) > 3),]

# Print the outliers
female_outliers

Investigation of the outliers in the ‘Male_Life_Expectancy’ column will be done next. The outliers in the male life expectancy column were selected using which function to get the rows with z-score that are higher than 3. These rows are then printed for investigation. Similar to the female life expectancy column, it was found that the countries with the outliers were mostly low income or lower middle income countries.

# Select the rows that are outliers based on the z-score value
male_outliers <- life_exp_df[which(abs(z_male) > 3),]

# Print the outliers
male_outliers

Lastly, investigation of the outliers in the ‘Life_Expectancy_Gap’ column will be done. The outliers in the life expectancy gap column were selected using which function to get the rows with z-score that are higher than 3. These rows are then printed for investigation. Unlike the male and female life expectancy columns, it was found that the countries with the outliers in the life expectancy gap are were mostly high or upper middle income countries, with the exception of Syria being the only outliers in the low income group.

# Select the rows that are outliers based on the z-score value
gap_outliers <- life_exp_df[which(abs(z_gap) > 3),]

# Print the outliers
gap_outliers

Several different countries were found to be outliers based on the z-score values. However, upon further inspection of each of the outliers, it was found that there may be a reasoning or story behind as to why they are outliers. For example, in the female and male life expectancy columns, it was found that most of the countries are lower and middle lower income country. This suggests that the economy of the country may be an underlying reason as to why the life expectancy of these countries are worse as compared to countries with better economy.

On the other hand, the life expectancy gap outliers were mostly in the higher and upper middle country. However, these outliers also suggest that the economy of the country may be reason as to why there are bigger gap between female and males, in which countries with higher economy may have a bigger life expectancy gap. Although these rows in the three life expectancy columns are considered as outliers, these outliers were considered to be valuable because they provide suggestion that there may be an effect of economy on life expectancy and a reason as to why these countries have lower life expectancy or higher life expectancy gap as compared to other countries.

Additionally, when looking at the country Cambodia in the female and male life expectancy, it was found that Cambodia started having life expectancy that are outliers from the year 1974 (for male) and 1975 (for female) till 1980. Investigation on the history of Cambodia showed that 1974 were when the Cambodian Civil War happened and males were recruited into the army for war. While 1975 was when the Cambodian genocide started; Cambodian genocide resulted in 1.5 to 2 million deaths from 1975 to 1979, which were approximated to be nearly a quarter of Cambodia’s population at that time. These outliers were actually telling a story that something happened in the 1974 that causes Cambodia life expectancy to suddenly change and become outliers. Hence, these outliers found in the data may be able to show an interesting and valuable story about the life expectancy of the countries.

Overall, although several rows in the data were considered as outliers according to the boxplot and z-score method. These rows were kept because they are valuable for analysis and able to show interesting story or insight about the data. Transformation (the next section) would be done on the three variables instead to reduce the potential impact of these outliers.

Transform

To address potential impact or bias that may occur because of the outliers, transformation will be done to the three life expectancy variables. Three histograms would be created first prior to the transformation to see how the initial distribution of the variables.

# Create the histogram of the female life expectancy
life_exp_df$Female_Life_Expectancy %>% hist(main = "Original", ylab = "Frequency", xlab = "Female Life Expectancy")

# Create the histogram of the male life expectancy
life_exp_df$Male_Life_Expectancy %>% hist(main = "Original", ylab = "Frequency", xlab = "Male Life Expectancy")

# Create the histogram of the life expectancy gap
life_exp_df$Life_Expectancy_Gap %>% hist(main = "Original", ylab = "Frequency", xlab = "Life Expectancy Gap")

We will firstly transform the female life expectancy variable. As the original histogram is left skewed, transformation using squares, cube, and fourth power would be done on the female life expectancy variable. The hist function would be used to create the original histogram of the female life expectancy alongside the histogram after using different transformation.

# Create the original and after transformation histograms
par(mfrow=c(2,2))

female_original <- life_exp_df$Female_Life_Expectancy
hist(female_original, main = "Original", ylab = "Frequency", xlab = "Female Life Expectancy")

female_square <- (life_exp_df$Female_Life_Expectancy)^2
hist(female_square, main = "Square Transformation", ylab = "Frequency", xlab = "(Female Life Expectancy)^2")

female_cube <- (life_exp_df$Female_Life_Expectancy)^3
hist(female_cube, main = "Cube Transformation", ylab = "Frequency", xlab = "(Female Life Expectancy)^3")

female_fourth <- (life_exp_df$Female_Life_Expectancy)^4
hist(female_fourth, main = "Fourth Power Transformation", ylab = "Frequency", xlab = "(Female Life Expectancy)^4")

From the transformation, it seems that cube transformation is the best in reducing the left skewness of the female life expectancy. The histogram and box plot before and after the transformation was created using the hist and boxplot functions to highlight the changes after the cube transformation.

# Create the histograms and boxplots before and after transformation
par(mfrow=c(2,2))

hist(female_original, main = "Original Histogram", ylab = "Frequency", xlab = "Female Life Expectancy")

boxplot(female_original, main="Original Boxplot", ylab="Female Life Expectancy", col = "grey")

hist(female_cube, main = "Histogram After Transformation", ylab = "Frequency", xlab = "Female Life Expectancy")

boxplot(female_cube, main="Boxplot After Transformation", ylab="Female Life Expectancy", col = "grey")

Next, we will transform the male life expectancy variable. Similar to the female life expectancy, the original histogram for male life expectancy is left skewed. Hence, squares, cube, and fourth power transformation would be done on the male life expectancy variable. The hist function would be used to create the original histogram of the male life expectancy alongside the histogram after using different transformation.

# Create the original and after transformation histograms
par(mfrow=c(2,2))

male_original <- life_exp_df$Male_Life_Expectancy
hist(male_original, main = "Original", ylab = "Frequency", xlab = "Male Life Expectancy")

male_square <- (life_exp_df$Male_Life_Expectancy)^2
hist(male_square, main = "Square Transformation", ylab = "Frequency", xlab = "(Male Life Expectancy)^2")

male_cube <- (life_exp_df$Male_Life_Expectancy)^3
hist(male_cube, main = "Cube Transformation", ylab = "Frequency", xlab = "(Male Life Expectancy)^3")

male_fourth <- (life_exp_df$Male_Life_Expectancy)^4
hist(male_fourth, main = "Fourth Power Transformation", ylab = "Frequency", xlab = "(Male Life Expectancy)^4")

From the transformation, it seems that cube transformation is the best in reducing the left skewness of the male life expectancy as well. The histogram and box plot before and after the transformation was created using the hist and boxplot functions to highlight the changes after the cube transformation.

# Create the histograms and boxplots before and after transformation
par(mfrow=c(2,2))

hist(male_original, main = "Original Histogram", ylab = "Frequency", xlab = "Male Life Expectancy")

boxplot(male_original, main="Original Boxplot", ylab="Male Life Expectancy", col = "grey")

hist(male_cube, main = "Histogram After Transformation", ylab = "Frequency", xlab = "Male Life Expectancy")

boxplot(male_cube, main="Boxplot After Transformation", ylab="Male Life Expectancy", col = "grey")

Lastly, we will transform the life expectancy gap variable. Unlike the female and male life expectancy, the original histogram for the life expectancy gap seem to be slightly skewed to the right. Hence, logarithm, natural logarithm, and square root transformation would be done on the life expectancy gap variable. Reciprocal transformation was not included because the histogram are only mildly skewed to the right, while reciprocal transformation is a very strong transformation with drastic effect on the distribution shape. The hist function would be used to create the original histogram of the male life expectancy alongside the histogram after using different transformation.

# Create the original and after transformation histograms
par(mfrow=c(2,2))

gap_original <- life_exp_df$Life_Expectancy_Gap
hist(gap_original, main = "Original", ylab = "Frequency", xlab = "Life Expectancy Gap")

gap_log <- log10(life_exp_df$Life_Expectancy_Gap)
hist(gap_log, main = "Log Transformation", ylab = "Frequency", xlab = "log10(Life Expectancy Gap)")

gap_nlog <- log(life_exp_df$Life_Expectancy_Gap)
hist(gap_nlog, main = "Natural Log Transformation", ylab = "Frequency", xlab = "log(Life Expectancy Gap)")

gap_sqrt <- sqrt(life_exp_df$Life_Expectancy_Gap)
hist(gap_sqrt, main = "Root Transformation", ylab = "Frequency", xlab = "sqrt(Life Expectancy Gap)")

From the transformation, it seems that the square root transformation is the best in reducing the right skewness of the life expectancy gap. The histogram and box plot before and after the transformation was created using the hist and boxplot functions to highlight the changes after the cube transformation.

# Create the histograms and boxplots before and after transformation
par(mfrow=c(2,2))

hist(gap_original, main = "Original Histogram", ylab = "Frequency", xlab = "Life Expectancy Gap")

boxplot(gap_original, main="Original Boxplot", ylab="Life Expectancy Gap", col = "grey")

hist(gap_sqrt, main = "Histogram After Transformation", ylab = "Frequency", xlab = "Life Expectancy Gap")

boxplot(gap_sqrt, main="Boxplot After Transformation", ylab="Life Expectancy Gap", col = "grey")

Although transformation was done on the life expectancy gap, the boxplot still show the pressence of several outliers. This was deemed to be alright because these outliers may be highlighting countries with very contrasting life expectancy between female and male, which may suggest that these countries should be further inspected by the analysts to see what causes such big life expectancy difference between female and males.

References

• The World Bank (2020). Life Expectancy at birth, female (years) | Data [Data Set]. The World Bank. https://data.worldbank.org/indicator/SP.DYN.LE00.FE.IN?view=chart

• The World Bank (2020). Life Expectancy at birth, male (years) | Data [Data Set]. The World Bank. https://data.worldbank.org/indicator/SP.DYN.LE00.MA.IN?view=chart