The attached who.csv dataset contains real-world data from 2008. The variables included follow.
Country: name of the country LifeExp: average life expectancy for the country in years
* InfantSurvival: proportion of those surviving to one year or more
* Under5Survival: proportion of those surviving to five years or more
* TBFree: proportion of the population without TB.
* PropMD: proportion of the population who are MDs
* PropRN: proportion of the population who are RNs
* PersExp: mean personal expenditures on healthcare in US dollars at average exchange rate
* GovtExp: mean government expenditures per capita on healthcare, US dollars at average exchange rate
* TotExp: sum of personal and government expenditures.
* Load data from CSV file
getwd() # Path where the file will be downloaded
## [1] "C:/MSDS_Course/Spring_2022/DATA_605/Week_12"
who_data <- read.csv("who.csv")
head(who_data)
## Country LifeExp InfantSurvival Under5Survival TBFree PropMD
## 1 Afghanistan 42 0.835 0.743 0.99769 0.000228841
## 2 Albania 71 0.985 0.983 0.99974 0.001143127
## 3 Algeria 71 0.967 0.962 0.99944 0.001060478
## 4 Andorra 82 0.997 0.996 0.99983 0.003297297
## 5 Angola 41 0.846 0.740 0.99656 0.000070400
## 6 Antigua and Barbuda 73 0.990 0.989 0.99991 0.000142857
## PropRN PersExp GovtExp TotExp
## 1 0.000572294 20 92 112
## 2 0.004614439 169 3128 3297
## 3 0.002091362 108 5184 5292
## 4 0.003500000 2589 169725 172314
## 5 0.001146162 36 1620 1656
## 6 0.002773810 503 12543 13046
Scatterplot of Life Expectancy and Total Expenditure
plot(who_data$TotExp, who_data$LifeExp, xlab = "Total Expenditures", ylab = "Life Expectancy")
Linear regression of Life Expectancy with Total Expenditure
LifeExp_lm <- lm(LifeExp ~ TotExp, data = who_data)
summary(LifeExp_lm)
##
## Call:
## lm(formula = LifeExp ~ TotExp, data = who_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -24.764 -4.778 3.154 7.116 13.292
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.475e+01 7.535e-01 85.933 < 2e-16 ***
## TotExp 6.297e-05 7.795e-06 8.079 7.71e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 9.371 on 188 degrees of freedom
## Multiple R-squared: 0.2577, Adjusted R-squared: 0.2537
## F-statistic: 65.26 on 1 and 188 DF, p-value: 7.714e-14
The F-statistic is 65.26 with a p-value of 7.714e-14. The p-value is very small, which suggests that we can reject the null hypothesis (regression model is same as intercept-only model) and accept the alternative that this simple linear regression model provides a better fit than the intercept-only model. The R2=0.2577 value is not strong, which tells us that the model only accounts for 25.77% of the variation in the data.
Create new columns with \(LifeExp^{4.6}\) and \(TotalExp^{4.6}\) values
who_data$LifeExp_4.6 <- (who_data$LifeExp)^4.6
who_data$TotExp_0.06 <- (who_data$TotExp)^0.06
Plot total expenditures vs life expectancy for who data set
plot(who_data$TotExp_0.06, who_data$LifeExp_4.6, xlab = "Total Expenditures ^ 0.06", ylab = "Life Expectancy ^ 4.6")
Run simple linear regression on transformed LifeExp and TotExp
data_trans_model <- lm(LifeExp_4.6 ~ TotExp_0.06, data = who_data)
Run summary function and get the F-statistics, R-squared, standard error, and P-values
summary(data_trans_model)
##
## Call:
## lm(formula = LifeExp_4.6 ~ TotExp_0.06, data = who_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -308616089 -53978977 13697187 59139231 211951764
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -736527910 46817945 -15.73 <2e-16 ***
## TotExp_0.06 620060216 27518940 22.53 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 90490000 on 188 degrees of freedom
## Multiple R-squared: 0.7298, Adjusted R-squared: 0.7283
## F-statistic: 507.7 on 1 and 188 DF, p-value: < 2.2e-16
The F-statistic for the transformed model is 507.7 (same degrees of freedom as model from 1) is much better vs. the the previous model. The P-value is very low which is a good thing. In the transformed model the \(R^2\)=0.7298 (72.98%), which is a lot better than the previous model (25.77%).
Transformed model
forecast_LifeExp <- function(TotExp_0.06_value){
return (736527910 + 620060216 * TotExp_0.06_value)^(1/4.6)
}
Forecast life expectancy when \(TotExp^{.06}\) =1.5
forecast_LifeExp(1.5)
## [1] 1666618234
Then forecast life expectancy when \(TotExp^{.06}\)=2.5
forecast_LifeExp(2.5)
## [1] 2286678450
Build multiple regression model using transformed variables
who_data$PropMD_TotExp_0.06 <- who_data$PropMD * who_data$TotExp_0.06
LifeExp_multiple_model <- lm(who_data$LifeExp_4.6 ~ who_data$PropMD + who_data$TotExp_0.06 + who_data$PropMD_TotExp_0.06)
Run summary function to provide the F-statistics, R-squared, Standard Error and P-values.
summary(LifeExp_multiple_model)
##
## Call:
## lm(formula = who_data$LifeExp_4.6 ~ who_data$PropMD + who_data$TotExp_0.06 +
## who_data$PropMD_TotExp_0.06)
##
## Residuals:
## Min 1Q Median 3Q Max
## -296470018 -47729263 12183210 60285515 212311883
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -7.244e+08 5.083e+07 -14.253 <2e-16 ***
## who_data$PropMD 4.727e+10 2.258e+10 2.094 0.0376 *
## who_data$TotExp_0.06 6.048e+08 3.023e+07 20.005 <2e-16 ***
## who_data$PropMD_TotExp_0.06 -2.121e+10 1.131e+10 -1.876 0.0622 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 88520000 on 186 degrees of freedom
## Multiple R-squared: 0.7441, Adjusted R-squared: 0.74
## F-statistic: 180.3 on 3 and 186 DF, p-value: < 2.2e-16
The F-statistic from our multiple regression model is 180.3. The P-value is strong for the model, and strong for all of the factors except PropMD x TotExp0.06 (0.0622). Finally, the \(R^2\)=0.7441, which means that the model accounts for 74.41% of variability in the data.
Co-efficients of the model
LifeExp_multiple_model$coefficients
## (Intercept) who_data$PropMD
## -724418697 47273338389
## who_data$TotExp_0.06 who_data$PropMD_TotExp_0.06
## 604795792 -21214671638
forecast_LifeExp_2 <- function(PropMD, TotExp_0.06){
LifeExp_4.6 <- 724418697 + 47273338389 * PropMD + 604795792 * TotExp_0.06 + -21214671638 * PropMD * TotExp_0.06
return(LifeExp_4.6^(1/4.6))
}
Forecast LifeExp when PropMD=.03 and TotExp = 14
forecast_LifeExp_2(.03, 14^0.06)
## [1] 106.3698
The forecast for proportion of population that are medical doctors is 3% and total expenditure of $14 is about 82.57 years. Looking at the data below for propMD that falls around 3%, I do see that 80 to 82 years of life expectancy is present; however, the total expenditure is a lot higher than 14. So, I do not think that this is a realistic forecast.
Life expentancy of propMD around 3%
LifeExp_propMD_3 <- subset(who_data, who_data$PropMD >.02 & who_data$PropMD < .04, select = c("LifeExp", "PropMD", "TotExp"))
head(LifeExp_propMD_3)
## LifeExp PropMD TotExp
## 45 80 0.03322813 40749
## 146 82 0.03512903 281653
min(who_data$TotExp)
## [1] 13
To get an idea of the life expectancy of countries with total expenditure close to $14, below is a list of entries with total expenditures under $100.
There is only one entry with total expenditure that is closest to $14.
The life expectancy of this country is 49 years.
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
## LifeExp PropMD TotExp
## 1 49 0.000024500 13
## 2 56 0.000023900 70
## 3 47 0.000096100 71
## 4 62 0.000194783 80
## 5 63 0.000274894 87
## 6 53 0.000107505 87
## 7 63 0.000045800 88
## 8 42 0.000021500 94
Based on actual data, the life expectancy of 82 years for a country with propMD of 3% and total expenditure of 14 is not realistic. Countries with propMD of about 3% tend to have total expenditure that is way above $14.