Assignment 12
The attached who.csv dataset contains real-world data from 2008.
1. Provide a scatterplot of LifeExp~TotExp, and run simple linear regression. Do not transform the variables. Provide and interpret the F statistics, R^2, standard error,and p-values only. Discuss whether the assumptions of simple linear regression met.
First, let’s load in the dataset
df <- read.csv("https://raw.githubusercontent.com/davidblumenstiel/CUNY-MSDS-Data-605/master/Week-12/who.csv")
head(df)## 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## Country LifeExp InfantSurvival Under5Survival
## Afghanistan : 1 Min. :40.00 Min. :0.8350 Min. :0.7310
## Albania : 1 1st Qu.:61.25 1st Qu.:0.9433 1st Qu.:0.9253
## Algeria : 1 Median :70.00 Median :0.9785 Median :0.9745
## Andorra : 1 Mean :67.38 Mean :0.9624 Mean :0.9459
## Angola : 1 3rd Qu.:75.00 3rd Qu.:0.9910 3rd Qu.:0.9900
## Antigua and Barbuda: 1 Max. :83.00 Max. :0.9980 Max. :0.9970
## (Other) :184
## TBFree PropMD PropRN PersExp
## Min. :0.9870 Min. :0.0000196 Min. :0.0000883 Min. : 3.00
## 1st Qu.:0.9969 1st Qu.:0.0002444 1st Qu.:0.0008455 1st Qu.: 36.25
## Median :0.9992 Median :0.0010474 Median :0.0027584 Median : 199.50
## Mean :0.9980 Mean :0.0017954 Mean :0.0041336 Mean : 742.00
## 3rd Qu.:0.9998 3rd Qu.:0.0024584 3rd Qu.:0.0057164 3rd Qu.: 515.25
## Max. :1.0000 Max. :0.0351290 Max. :0.0708387 Max. :6350.00
##
## GovtExp TotExp
## Min. : 10.0 Min. : 13
## 1st Qu.: 559.5 1st Qu.: 584
## Median : 5385.0 Median : 5541
## Mean : 40953.5 Mean : 41696
## 3rd Qu.: 25680.2 3rd Qu.: 26331
## Max. :476420.0 Max. :482750
## Below is a scatterplot of the data, along with a linear fit.
model <- lm(df$LifeExp~df$TotExp)
plot(df$LifeExp~df$TotExp)
abline(model$coefficients, col = "red")
This model obiously need’s work. Let’s examine some stats below.
##
## Call:
## lm(formula = df$LifeExp ~ df$TotExp)
##
## 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 ***
## df$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-14Let’s start with the F statistic: the ratio of the mean regression sum of squares over the mean error sum of squares. This, like the p-value derived from it, can be used to assess the overall significance of the model; unlike the t value, which only considers one variable at a time. If our F statistic here is greater than the critcal value of the F distribution at a certain signicance level, we can reject the null hypothesis (at that signifcance level) and say that the model fits the data better than an intercept only model. The critical value of the F distribution with 1 and 188 degrees of freedom (as shown above) and a significance of 0.05 (standard) is about 3.89. This is much lower than our F statistic for the model, meaning the model does fit the data better than an intercept only model.
The R-squared is the proportion of variance that can be explained by the model, with 0 being none of it, and 1 being a complete fit. The 0.26 we have isn’t that great; our model doesn’t fit the majority of the data very well, but it’s better than nothing. R-squared is calculated as the sum of squared residuals over the the total sum of squares.
The standard error (for the intercept and slope here) is an indication of how well the components of the model function. Small standard errors in comparison to the coefficients are good. Dividing a coefficient by it’s standard error gives us the t-value, which we can see above to the right of the coefficients. Generally speaking, having coefficients with an absolute value around 5-10 times more than the standard error is a good benchmark. The intercept standard error is around 86 times less than the coefficient, which bodes well for it, and the slope standard error is about 8 times less which is also good.
The p-values indicate whether or not the fits are significant. The intercept and slope fits have p-values much less than 0.05, which indicates that they are signicant. The overall p-value is also less than 0.05 (the same value as for the slope), which indicates that the model is significant.
As to whether or not the assumptions of linear regression are met, I can tell you right now just from looking at the scatter plot that this is not a linear relationship. Also, look above at the residual median and quartiles: those look super skewed. I’ll do Q-Q and residual plots below to further exemplify that
Those residuals are definitely not norally distributed, otherwise they would follow the Q-Q line above. This throws the residual normality and linear relationship assumptions out.
The residuals once again are definitely not normally distributed. Homoscedasticity may sill hold though, as although they are not normally distributed, they don’t appear to change
Overall, the assumptions for linear regression are not met.
2. Raise life expectancy to the 4.6 power (i.e., LifeExp^4.6). Raise total expenditures to the 0.06 power (nearly a log transform, TotExp^.06). Plot LifeExp^4.6 as a function of TotExp^.06, and re-run the simple regression model using the transformed variables. Provide and interpret the F statistics, R^2, standard error, and p-values. Which model is “better?”
dftrans <- df[,c(2,10)]
dftrans$LifeExp <- dftrans$LifeExp^4.6
dftrans$TotExp <- dftrans$TotExp^0.06
model2 <- lm(LifeExp~TotExp, data = dftrans)
plot(dftrans$LifeExp~dftrans$TotExp)
abline(model2$coefficients, col = "red")
Looks much better, right? Let’s look at the stats.
##
## Call:
## lm(formula = LifeExp ~ TotExp, data = dftrans)
##
## 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 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-16The F statistic here is much higher than the first model, which was alerady significant. The p-values are all at the lowest R will show, which is also better than the last, although those were also highly significant. The standard error is greater for the intercept in comparison to it’s coefficient, which is a step down although signifcant. The standard error for the slope is better though, now 22 times less than the coefficent compared to 8 times less for the origional model.
Most importantly however, the R-squared is much higher, suggesting the model is a much better fit and now accounts for about 73% of the variance. I’d definitely say this is the better model, although the residuals still aren’t totally normally distributed, as shown below. But, they are closer.
3. Using the results from 2, forecast life expectancy when TotExp^.06 =1.5. Then forecast life expectancy when TotExp^.06=2.5.
## 1 2
## 63.31153 86.50645It predicts life expectancies of around 63.3, and 86.5 years respectively.
4. Build the following multiple regression model and interpret the F Statistics, R^2, standard error, and p-values. How good is the model?
LifeExp = b0+b1 x PropMd + b2 x TotExp +b3 x PropMD x TotExp
To give you an idea of how this fits visually without using higher dimension figures, I’m going to plot predictions against actual values for the model. I find this useful for assessing multiple regression models.
model3 <- lm(LifeExp ~ PropMD + TotExp + (PropMD*TotExp), data = df)
plot(df$LifeExp ~ predict(model3, df))
abline(0,1, col = 'red')
So, with a perfect model we would expect all data to fall on the red line above, indicating that the preditions were the same as the actual values. This is obiously not the case, and this model clearly fails to capture trends as well as the transformation model did. Let’s examine some statistics below.
##
## Call:
## lm(formula = LifeExp ~ PropMD + TotExp + (PropMD * TotExp), data = 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-16Breifly, the F statistic (34.49) is higher than the critical F value (2.65), indicating significance, although not as much as the transformaton model. The overall P value is still very close to 0, also indicating overall signicance. Some of the p-values for the coefficients are greater or lower, but all are still very significant. The standard error for the intercept looks good, but the standard error for PropMD and PropMD*TotExp isn’t as good, being only 4-5 times lower than the coefficients. Most importantly here though, the R-squared value is 0.36, which indicates that this model fits the data better than the first model, but much worse than the transformation model.
5. Forecast LifeExp when PropMD=.03 and TotExp = 14. Does this forecast seem realistic? Why or why not?
## 1
## 107.696107.7 years doesn’t seem realistic. First, that’s a pretty long time. Second, because PropMD is much higher and TotExp is much lower than their average values: such values would be hard for this model to use because the data it was trained on really doesn’t have anything like that. Also, I suspect a high value for somthing like PropMd might not be able to make up for a low value like TotExp in reality; it wouldn’t matter how many doctors there were if you weren’t using them.