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.
who_df <- read.csv('https://raw.githubusercontent.com/catfoodlover/DATA605/main/who.csv')
glimpse(who_df)
## Rows: 190
## Columns: 10
## $ Country <chr> "Afghanistan", "Albania", "Algeria", "Andorra", "Angola…
## $ LifeExp <int> 42, 71, 71, 82, 41, 73, 75, 69, 82, 80, 64, 74, 75, 63,…
## $ InfantSurvival <dbl> 0.835, 0.985, 0.967, 0.997, 0.846, 0.990, 0.986, 0.979,…
## $ Under5Survival <dbl> 0.743, 0.983, 0.962, 0.996, 0.740, 0.989, 0.983, 0.976,…
## $ TBFree <dbl> 0.99769, 0.99974, 0.99944, 0.99983, 0.99656, 0.99991, 0…
## $ PropMD <dbl> 0.000228841, 0.001143127, 0.001060478, 0.003297297, 0.0…
## $ PropRN <dbl> 0.000572294, 0.004614439, 0.002091362, 0.003500000, 0.0…
## $ PersExp <int> 20, 169, 108, 2589, 36, 503, 484, 88, 3181, 3788, 62, 1…
## $ GovtExp <int> 92, 3128, 5184, 169725, 1620, 12543, 19170, 1856, 18761…
## $ TotExp <int> 112, 3297, 5292, 172314, 1656, 13046, 19654, 1944, 1907…
plot(who_df$TotExp, who_df$LifeExp)
le.lm <- lm(LifeExp ~ TotExp, data = who_df)
summary(le.lm)
##
## Call:
## lm(formula = LifeExp ~ TotExp, data = who_df)
##
## 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 and p-value tells us whether our model explains the response variable better than a model with just the intercept. The p-value is 7.714e-14 which is smaller than 0.05 so we can reject the null hypothesis that our model is no better than the intercept.
Article on Interpreting the F-statistic
The R-squared is 0.2577 which means that ~26% of the variation in Life Expectancy is explained by the TotExp predictor.
The Standard Error is the standard error of the estimate. It tells you on average how far an observation is from the regression line. You would like the Standard Error to be 5-10 times smaller than the coefficient. The Standard Error for the TotExp is 8 times smaller than the coefficient.
plot(le.lm)
I’m not sure that the assumptions of linear regression were met. The plot doesn’t show a linear relationship. As the residuals are not normal. The magnitude of min and max are not the same and the median deviates from zero. The residuals are not evenly distributed and the Q-Q looks skewed.
who_df$LifeExp_4 <- who_df$LifeExp^4.6
who_df$TotExp_log <- who_df$TotExp^0.06
plot(who_df$TotExp_log, who_df$LifeExp_4)
le_trans.lm <- lm(LifeExp_4 ~ TotExp_log, data = who_df)
summary(le_trans.lm)
##
## Call:
## lm(formula = LifeExp_4 ~ TotExp_log, data = who_df)
##
## 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_log 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 p-value for the F-statistic is less than 0.05 so we can reject the null hypothesis that this model is no better than the model just containing the intercept.
The R-squared is 0.7298 which means that ~ 73% of the variation in response variable can be explained by the predictor variable.
The Standard Error is 22 times smaller than the Beta Estimate which is good thing. We were hoping for 5-10 and we got 22.
The scatter plot looks linear but the residuals are hard to interpret. There is some skew to the residuals but the residuals look evenly distributed and the Q-Q looks good.
plot(le_trans.lm)
3. Using the results from 2, forecast life expectancy when TotExp^.06 =1.5. Then forecast life expectancy when TotExp^.06=2.5.
new_obs = data.frame(TotExp_log = 1.5)
life_4 <- predict(le_trans.lm, newdata = new_obs)
life_4^(1/4.6)
## 1
## 63.31153
new_obs = data.frame(TotExp_log = 2.5)
life_4 <- predict(le_trans.lm, newdata = new_obs)
life_4^(1/4.6)
## 1
## 86.50645
new_model <- lm(LifeExp ~ PropMD + TotExp + PropMD*TotExp, data = who_df)
summary(new_model)
##
## Call:
## lm(formula = LifeExp ~ PropMD + TotExp + PropMD * TotExp, data = who_df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -27.320 -4.132 2.098 6.540 13.074
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.277e+01 7.956e-01 78.899 < 2e-16 ***
## PropMD 1.497e+03 2.788e+02 5.371 2.32e-07 ***
## TotExp 7.233e-05 8.982e-06 8.053 9.39e-14 ***
## PropMD:TotExp -6.026e-03 1.472e-03 -4.093 6.35e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 8.765 on 186 degrees of freedom
## Multiple R-squared: 0.3574, Adjusted R-squared: 0.3471
## F-statistic: 34.49 on 3 and 186 DF, p-value: < 2.2e-16
The p-value of the F-statistic is less that 0.05 so we can reject the null hypothesis that the added covariates perform just as well as a model limited to the intercept.
The R-squared is 0.3574, so ~36% of the variation in the response variable LifeExp can be explained by our model.
Our Standard Errors are between 4-8 times smaller than our beta estimates which are meeting the minimum requirements
Our p-values are all less than 0.05 rejecting the null hypothesis that the true betas are zero.
We have had models that explain more or the variation in our response variable. This model is okay but not great.
Let’s compare this model to our very first model. The AIC and BIC are better for the minimal model compared to this new model
library(performance)
## Warning: package 'performance' was built under R version 4.1.2
compare_performance(le.lm, new_model)
## # Comparison of Model Performance Indices
##
## Name | Model | AIC | AIC weights | BIC | BIC weights | R2 | R2 (adj.) | RMSE | Sigma
## -------------------------------------------------------------------------------------------------------
## le.lm | lm | 1393.483 | < 0.001 | 1403.224 | < 0.001 | 0.258 | 0.254 | 9.322 | 9.371
## new_model | lm | 1370.067 | 1.000 | 1386.302 | 1.000 | 0.357 | 0.347 | 8.673 | 8.765
new_obs = data.frame(PropMD = 0.3, TotExp = 14)
new_prod <- predict(new_model, newdata = new_obs)
new_prod
## 1
## 511.9966
summary(who_df$PropMD)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.0000196 0.0002444 0.0010474 0.0017954 0.0024584 0.0351290
summary(who_df$TotExp)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 13 584 5541 41696 26331 482750
We are using combinations of values that don’t ever occur in the data. We are using 10 times the max value for proportion of population that is are clinicians and the minimum for health care expenditure.