Regresja
Gabriela Malec
Nazwa uczelni: Akademia Górniczo-Hutnicza im. Stanisława Staszica w Krakowie
Kierunek: Geoinformacja II rok
Przedmiot: Przetwarzanie danych środowiskowych
1. Dane
Ze strony https://www.kaggle.com/unsdsn/world-happiness pobrano wyniki badan dotyczacych wplywu roznych czynnikow na szczescie ludzi. Dane zawieraja takie zmienne jak wspolczynnik szczescia (wynik), GDP, pomoc socjalna, dlugosc zycia w zdrowiu, wolnosc podejmowania decyzji, hojnosc, poziom korupcji.
dat <- read_csv("C:/Users/admin/Desktop/STUDIA/IV semestr/PDŚ/cwiczenia/z2/regresja/2019.csv")
View(dat)
colnames(dat) <- c("Pozycja", "Kraj", "Wynik", "GDP", "Pomoc_socjalna", "Zdrowie", "Wolnosc", "Hojnosc", "Korupcja")
dat <- as.data.frame(dat)
dat = dat[, -c(1,2)]Podstawowe informacje:
summary(dat)## Wynik GDP Pomoc_socjalna Zdrowie
## Min. :2.853 Min. :0.0000 Min. :0.000 Min. :0.0000
## 1st Qu.:4.545 1st Qu.:0.6028 1st Qu.:1.056 1st Qu.:0.5477
## Median :5.380 Median :0.9600 Median :1.272 Median :0.7890
## Mean :5.407 Mean :0.9051 Mean :1.209 Mean :0.7252
## 3rd Qu.:6.184 3rd Qu.:1.2325 3rd Qu.:1.452 3rd Qu.:0.8818
## Max. :7.769 Max. :1.6840 Max. :1.624 Max. :1.1410
## Wolnosc Hojnosc Korupcja
## Min. :0.0000 Min. :0.0000 Min. :0.0000
## 1st Qu.:0.3080 1st Qu.:0.1087 1st Qu.:0.0470
## Median :0.4170 Median :0.1775 Median :0.0855
## Mean :0.3926 Mean :0.1848 Mean :0.1106
## 3rd Qu.:0.5072 3rd Qu.:0.2482 3rd Qu.:0.1412
## Max. :0.6310 Max. :0.5660 Max. :0.4530hist(dat$Wynik, main="Histogram",col="pink")
shapiro.test(dat$Wynik)##
## Shapiro-Wilk normality test
##
## data: dat$Wynik
## W = 0.9872, p-value = 0.1633P-value jest większe od 0.05, więc mozna przyjac hipotezę zerową, czyli wykres jest podobny do rozkładu normalnego.
Dane zostaly podzielone na dwie grupy - treningowa i testowa.
sample_split <- sample.split(Y = dat$Wynik, SplitRatio = 0.5)
train_set <- subset(x = dat, sample_split == TRUE)
test_set <- subset(x = dat, sample_split == FALSE)
datatable(train_set, caption = "Grupa treningowa", options = list(scrollX = TRUE))datatable(test_set, caption = "Grupa testowa", options = list(scrollX = TRUE))2. Korelacja
M=cor(dat)
corrplot(M, method="circle")
pairs.panels(dat,
method = "pearson",
hist.col = "#376485",
density = FALSE,
ellipses = FALSE)
3. Modele
Model - wspolczynnik szczescia w zaleznosci od wszystkich zmiennych
m_all=lm(Wynik~.,data=train_set, na.action = na.exclude)
summary(m_all)##
## Call:
## lm(formula = Wynik ~ ., data = train_set, na.action = na.exclude)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.37238 -0.34375 0.05006 0.37800 1.01900
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.1232 0.2970 7.148 6.28e-10 ***
## GDP 1.1070 0.3276 3.380 0.00118 **
## Pomoc_socjalna 1.0664 0.3653 2.919 0.00470 **
## Zdrowie 0.3510 0.4952 0.709 0.48083
## Wolnosc 1.3732 0.5919 2.320 0.02321 *
## Hojnosc 1.0567 0.7283 1.451 0.15119
## Korupcja 0.4417 0.7214 0.612 0.54233
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5305 on 71 degrees of freedom
## Multiple R-squared: 0.7714, Adjusted R-squared: 0.7521
## F-statistic: 39.94 on 6 and 71 DF, p-value: < 2.2e-16Model - wspolczynnik szczescia w zaleznosci od wszystkich zmiennych z ich interakcjami
m_all_with=lm(Wynik~.*.,data=train_set, na.action = na.exclude)
summary(m_all_with)##
## Call:
## lm(formula = Wynik ~ . * ., data = train_set, na.action = na.exclude)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.23565 -0.21870 -0.01137 0.19863 0.85799
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.9188 1.0224 4.811 1.17e-05 ***
## GDP 1.7082 1.8013 0.948 0.347025
## Pomoc_socjalna -4.6884 1.2988 -3.610 0.000655 ***
## Zdrowie -2.6965 2.6660 -1.011 0.316141
## Wolnosc 9.4119 2.2808 4.127 0.000123 ***
## Hojnosc -0.4355 2.9506 -0.148 0.883186
## Korupcja -2.8038 4.4631 -0.628 0.532420
## GDP:Pomoc_socjalna -1.0338 1.4523 -0.712 0.479529
## GDP:Zdrowie -2.5134 0.9782 -2.570 0.012874 *
## GDP:Wolnosc 3.0437 3.0814 0.988 0.327508
## GDP:Hojnosc 4.0684 4.4488 0.914 0.364388
## GDP:Korupcja 2.7799 4.8428 0.574 0.568241
## Pomoc_socjalna:Zdrowie 9.2339 2.3861 3.870 0.000287 ***
## Pomoc_socjalna:Wolnosc -3.8917 2.2530 -1.727 0.089624 .
## Pomoc_socjalna:Hojnosc 6.3105 3.7813 1.669 0.100724
## Pomoc_socjalna:Korupcja 9.3486 3.9407 2.372 0.021136 *
## Zdrowie:Wolnosc -7.3952 5.3915 -1.372 0.175641
## Zdrowie:Hojnosc -8.3856 6.4316 -1.304 0.197632
## Zdrowie:Korupcja -13.2792 7.5175 -1.766 0.082773 .
## Wolnosc:Hojnosc -8.2513 5.8583 -1.408 0.164518
## Wolnosc:Korupcja 1.6720 6.2205 0.269 0.789083
## Hojnosc:Korupcja -2.8442 7.5837 -0.375 0.709046
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4024 on 56 degrees of freedom
## Multiple R-squared: 0.8962, Adjusted R-squared: 0.8573
## F-statistic: 23.03 on 21 and 56 DF, p-value: < 2.2e-16Regresja krokowa wsteczna
mAIC=stepAIC(m_all, direction = "backward")## Start: AIC=-92.24
## Wynik ~ GDP + Pomoc_socjalna + Zdrowie + Wolnosc + Hojnosc +
## Korupcja
##
## Df Sum of Sq RSS AIC
## - Korupcja 1 0.1055 20.084 -93.829
## - Zdrowie 1 0.1413 20.120 -93.690
## <none> 19.979 -92.240
## - Hojnosc 1 0.5924 20.571 -91.961
## - Wolnosc 1 1.5148 21.493 -88.540
## - Pomoc_socjalna 1 2.3983 22.377 -85.397
## - GDP 1 3.2138 23.192 -82.605
##
## Step: AIC=-93.83
## Wynik ~ GDP + Pomoc_socjalna + Zdrowie + Wolnosc + Hojnosc
##
## Df Sum of Sq RSS AIC
## - Zdrowie 1 0.1694 20.253 -95.174
## <none> 20.084 -93.829
## - Hojnosc 1 0.7957 20.880 -92.799
## - Pomoc_socjalna 1 2.3206 22.405 -87.301
## - Wolnosc 1 2.3702 22.454 -87.128
## - GDP 1 3.4995 23.584 -83.301
##
## Step: AIC=-95.17
## Wynik ~ GDP + Pomoc_socjalna + Wolnosc + Hojnosc
##
## Df Sum of Sq RSS AIC
## <none> 20.253 -95.174
## - Hojnosc 1 0.8447 21.098 -93.987
## - Wolnosc 1 2.4525 22.706 -88.258
## - Pomoc_socjalna 1 2.6427 22.896 -87.608
## - GDP 1 9.4666 29.720 -67.261summary(mAIC)##
## Call:
## lm(formula = Wynik ~ GDP + Pomoc_socjalna + Wolnosc + Hojnosc,
## data = train_set, na.action = na.exclude)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.39565 -0.32552 0.06633 0.38389 1.06908
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.1594 0.2898 7.452 1.50e-10 ***
## GDP 1.3180 0.2256 5.841 1.34e-07 ***
## Pomoc_socjalna 1.0457 0.3388 3.086 0.00286 **
## Wolnosc 1.5594 0.5245 2.973 0.00399 **
## Hojnosc 1.2117 0.6945 1.745 0.08522 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5267 on 73 degrees of freedom
## Multiple R-squared: 0.7683, Adjusted R-squared: 0.7556
## F-statistic: 60.51 on 4 and 73 DF, p-value: < 2.2e-16Regresja LEAPS
mLEAPS=regsubsets(Wynik~., data=train_set, nbest=5)
plot(mLEAPS, scale="r2")
mLEAPS2=lm(Wynik~GDP+Zdrowie, data=train_set, na.action = na.exclude)Algorytm EARTH
mEarth=earth(Wynik~.,train_set)
plot(mEarth)
ev=evimp(mEarth)
plot(ev)
mEarth2=lm(Wynik~GDP+Pomoc_socjalna+Zdrowie+Wolnosc, data=train_set, na.action = na.exclude)Algorytm Boruta
mBoruta=Boruta(Wynik~.,train_set)
boruta_signif <- names(mBoruta$finalDecision[mBoruta$finalDecision %in% c("Confirmed", "Tentative")])
plot(mBoruta, cex.axis=.9, las=2, xlab="")
mBoruta2=lm(Wynik~Zdrowie+GDP+Pomoc_socjalna, data=train_set, na.action = na.exclude)Wnioski
Z powyższych wykresów wynika, że przy budowaniu modeli najwieksze znaczenie mialy parametry, takie jak zdrowie i pomoc socjalna, a w dalszej kolejnosci równiez wskaznik GDP.
4. Walidacja krzyżowa
Control <- trainControl(method = "cv", number = 10, verboseIter = F)
mCV <- caret::train(Wynik~.,
data=train_set,
method = "lm",
metric = 'MAE',
trControl = Control,
na.action = na.omit,
tuneLength = 15
)
mCV## Linear Regression
##
## 78 samples
## 6 predictor
##
## No pre-processing
## Resampling: Cross-Validated (10 fold)
## Summary of sample sizes: 70, 70, 70, 71, 70, 71, ...
## Resampling results:
##
## RMSE Rsquared MAE
## 0.5526091 0.739669 0.450548
##
## Tuning parameter 'intercept' was held constant at a value of TRUEsummary(mCV)##
## Call:
## lm(formula = .outcome ~ ., data = dat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.37238 -0.34375 0.05006 0.37800 1.01900
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.1232 0.2970 7.148 6.28e-10 ***
## GDP 1.1070 0.3276 3.380 0.00118 **
## Pomoc_socjalna 1.0664 0.3653 2.919 0.00470 **
## Zdrowie 0.3510 0.4952 0.709 0.48083
## Wolnosc 1.3732 0.5919 2.320 0.02321 *
## Hojnosc 1.0567 0.7283 1.451 0.15119
## Korupcja 0.4417 0.7214 0.612 0.54233
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5305 on 71 degrees of freedom
## Multiple R-squared: 0.7714, Adjusted R-squared: 0.7521
## F-statistic: 39.94 on 6 and 71 DF, p-value: < 2.2e-165. Porównanie modeli
Model - wspolczynnik szczescia w zaleznosci od wszystkich zmiennych
aa=predict(m_all)
plot(train_set$Wynik,aa)
ggof(aa, train_set$Wynik)
dane2=data.frame(train_set$Wynik,aa)
ggplot(dane2, aes(x=train_set$Wynik,y=aa))+
geom_point()+
stat_smooth(method = "lm", col="red")
Model - wspolczynnik szczescia w zaleznosci od wszystkich zmiennych z ich interakcjami
bb=predict(m_all_with)
plot(train_set$Wynik,bb)
ggof(bb, train_set$Wynik)
dane3=data.frame(train_set$Wynik,bb)
ggplot(dane3, aes(x=train_set$Wynik,y=bb))+
geom_point()+
stat_smooth(method = "lm", col="red")
Regresja krokowa wsteczna
cc=predict(mAIC)
plot(train_set$Wynik,cc)
ggof(cc, train_set$Wynik)
dane4=data.frame(train_set$Wynik,cc)
ggplot(dane4, aes(x=train_set$Wynik,y=cc))+
geom_point()+
stat_smooth(method = "lm", col="red")
Regresja LEAPS
dd=predict(mLEAPS2)
plot(train_set$Wynik,dd)
ggof(dd, train_set$Wynik)
dane5=data.frame(train_set$Wynik,dd)
ggplot(dane5, aes(x=train_set$Wynik,y=dd))+
geom_point()+
stat_smooth(method = "lm", col="red")
Algorytm EARTH
ee=predict(mEarth2)
plot(train_set$Wynik,ee)
ggof(ee, train_set$Wynik)
dane6=data.frame(train_set$Wynik,ee)
ggplot(dane6, aes(x=train_set$Wynik,y=ee))+
geom_point()+
stat_smooth(method = "lm", col="red")
Algorytm Boruta
ff=predict(mBoruta2)
plot(train_set$Wynik,ff)
ggof(ff, train_set$Wynik)
dane7=data.frame(train_set$Wynik,ff)
ggplot(dane7, aes(x=train_set$Wynik,y=ff))+
geom_point()+
stat_smooth(method = "lm", col="red")
Ogólne porównanie
Kryterium AIC
datatable(AIC(m_all,m_all_with,mAIC,mLEAPS2,mEarth2,mBoruta2))a <- cor(train_set$Wynik,aa)
b <- cor(train_set$Wynik,bb)
c <- cor(train_set$Wynik,cc)
d <- cor(train_set$Wynik,dd)
e <- cor(train_set$Wynik,ee)
f <- cor(train_set$Wynik,ff)Korelacja dane - predykcja
x <- data.frame("Model" = c("m_all","m_all_with","mAIC","mLEAPS2","mEarth2","mBoruta2"), "Cor" = c(a,b,c,d,e,f))
x## Model Cor
## 1 m_all 0.8783081
## 2 m_all_with 0.9466963
## 3 mAIC 0.8765158
## 4 mLEAPS2 0.8180104
## 5 mEarth2 0.8724186
## 6 mBoruta2 0.8558757Zaleznosc liniowa
Model - wspolczynnik szczescia w zaleznosci od wszystkich zmiennych
ggplot(m_all, aes(.fitted, .resid)) +
modelr::geom_ref_line(h = 0) +
geom_point() +
geom_smooth(se = FALSE) +
ggtitle("Residuals vs Fitted")
***
Model - wspolczynnik szczescia w zaleznosci od wszystkich zmiennych z ich interakcjami
ggplot(m_all_with, aes(.fitted, .resid)) +
modelr::geom_ref_line(h = 0) +
geom_point() +
geom_smooth(se = FALSE) +
ggtitle("Residuals vs Fitted")
***
Regresja krokowa wsteczna
ggplot(mAIC, aes(.fitted, .resid)) +
modelr::geom_ref_line(h = 0) +
geom_point() +
geom_smooth(se = FALSE) +
ggtitle("Residuals vs Fitted")
***
MSE
Model - wspolczynnik szczescia w zaleznosci od wszystkich zmiennych
# test MSE
t_s <- test_set %>%
add_predictions(m_all) %>%
summarise(MSE = mean((Wynik - pred)^2))
# training MSE
t_s_2 <- train_set %>%
add_predictions(m_all) %>%
summarise(MSE = mean((Wynik - pred)^2))
y <- data.frame("1" = c("Testowa",t_s,"Treningowa",t_s_2 ))
y## X1..Testowa. X1.MSE X1..Treningowa. X1.MSE.1
## 1 Testowa 0.324418 Treningowa 0.2561347Model - wspolczynnik szczescia w zaleznosci od wszystkich zmiennych z ich interakcjami
# test MSE
t_s1 <- test_set %>%
add_predictions(m_all_with) %>%
summarise(MSE = mean((Wynik - pred)^2))
# training MSE
t_s_21 <- train_set %>%
add_predictions(m_all_with) %>%
summarise(MSE = mean((Wynik - pred)^2))
y <- data.frame("2" = c("Testowa",t_s1,"Treningowa",t_s_21 ))
y## X2..Testowa. X2.MSE X2..Treningowa. X2.MSE.1
## 1 Testowa 0.4340285 Treningowa 0.1162775Regresja krokowa wsteczna
# test MSE
t_s3 <- test_set %>%
add_predictions(mAIC) %>%
summarise(MSE = mean((Wynik - pred)^2))
# training MSE
t_s_23 <- train_set %>%
add_predictions(mAIC) %>%
summarise(MSE = mean((Wynik - pred)^2))
y <- data.frame("3" = c("Testowa",t_s3,"Treningowa",t_s_23 ))
y## X3..Testowa. X3.MSE X3..Treningowa. X3.MSE.1
## 1 Testowa 0.3516567 Treningowa 0.2596596. Wnioski
Z przeprowadzonej analizy można wywnioskowac, że modele ze wszystkimi danymi maja przewage nad pozostalymi. Najskuteczniejszym z nich okazal m_all_with, czyli model ze wszystkimi danymi wraz z ich interakcjami. Biorac pod uwage pozostale modele, dobrze sprawuje sie polaczenie Zdrowie+GDP+Pomoc_socjalna.