library(readxl)
district <- read_excel("district.xls")
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.1 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.1
## ✔ purrr 1.0.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(pastecs)
##
## Attaching package: 'pastecs'
##
## The following objects are masked from 'package:dplyr':
##
## first, last
##
## The following object is masked from 'package:tidyr':
##
## extract
library(lmtest)
## Loading required package: zoo
##
## Attaching package: 'zoo'
##
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
library(MASS)
##
## Attaching package: 'MASS'
##
## The following object is masked from 'package:dplyr':
##
## select
#Linear model "lm()" from the variables.
data_multiple <- lm(DDH00A001S22R~DPETBILP+DPSTBIFP+DPFPABILP,data = district)
#Assumptions: Linearity(plot and raintest)
plot(data_multiple,which = 1)

raintest(data_multiple)
##
## Rainbow test
##
## data: data_multiple
## Rain = 0.87303, df1 = 599, df2 = 594, p-value = 0.9513
#Independence of errors (durbin-watson)
library(car)
## Loading required package: carData
##
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
##
## recode
## The following object is masked from 'package:purrr':
##
## some
durbinWatsonTest(data_multiple)
## lag Autocorrelation D-W Statistic p-value
## 1 0.01915531 1.961555 0.5
## Alternative hypothesis: rho != 0
#Homoscedasticity(plot,bptest)
plot(data_multiple,which = 3)

bptest(data_multiple)
##
## studentized Breusch-Pagan test
##
## data: data_multiple
## BP = 4.4104, df = 3, p-value = 0.2204
#Normality of residuals (QQ plot, shapiro test)
plot(data_multiple,which = 2)

shapiro.test(data_multiple$residuals)
##
## Shapiro-Wilk normality test
##
## data: data_multiple$residuals
## W = 0.89531, p-value < 2.2e-16
#No multicollinearity (VIF, cor)
vif(data_multiple)
## DPETBILP DPSTBIFP DPFPABILP
## 1.989514 1.459168 1.598077
data_multiple_vars <- district %>% dplyr::select(DPETBILP,DPSTBIFP,DPFPABILP)
cor(data_multiple_vars)
## DPETBILP DPSTBIFP DPFPABILP
## DPETBILP 1 NA NA
## DPSTBIFP NA 1 NA
## DPFPABILP NA NA 1
#Log transformation - Normality of Residuals
district_drop <- district %>% drop_na(DDH00A001S22R,DPETBILP,DPSTBIFP,DPFPABILP) %>% filter(DDH00A001S22R>0 & DPETBILP>0)
data_multiple_log <- lm(log(DDH00A001S22R)~log(DPETBILP)+DPSTBIFP+DPFPABILP,data = district_drop)
bptest(data_multiple_log)
##
## studentized Breusch-Pagan test
##
## data: data_multiple_log
## BP = 0.79281, df = 3, p-value = 0.8512
shapiro.test(data_multiple_log$residuals)
##
## Shapiro-Wilk normality test
##
## data: data_multiple_log$residuals
## W = 0.86179, p-value < 2.2e-16
summary(data_multiple)
##
## Call:
## lm(formula = DDH00A001S22R ~ DPETBILP + DPSTBIFP + DPFPABILP,
## data = district)
##
## Residuals:
## Min 1Q Median 3Q Max
## -74.914 -5.757 0.802 7.769 27.088
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 73.91439 0.47738 154.833 < 2e-16 ***
## DPETBILP -0.22250 0.03451 -6.448 1.65e-10 ***
## DPSTBIFP 0.09047 0.07542 1.200 0.231
## DPFPABILP 0.44518 0.34458 1.292 0.197
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 12.5 on 1193 degrees of freedom
## (10 observations deleted due to missingness)
## Multiple R-squared: 0.04516, Adjusted R-squared: 0.04276
## F-statistic: 18.81 on 3 and 1193 DF, p-value: 6.377e-12
#Assumptions met-
#The Homoscedasticity test shows that the model has a mostly straight line.
#The model also meets the vif command. All the variables show a weak correlation of less than 5.
#The durbinWatsonTest shows a p-value greater than 0.05 and the D-W Statistic is very close to 2.
#The bptest shows that the model's p-value is more that 0.05 which is not significant, meaning the model is homoscedastic, not heteroscedastic.
#The rainbow test shows that the model is very linear with a p-value of 0.95.
#Assumptions violated-
#The Shapiro-Wilk test shows a significant p-value, meaning the residuals are not normal. A log transformation was performed for the normality of residuals.