Vprašanje2
Tibor Kovačič
2024-01-08
Za raziskovalno vprašanje sem si izbral sledeče: Kako Ginijev koeficient, brezposelnost celotne delovne sile, indeks zivljenskih stroskov in clanstvo v EU vplivajo na bruto domači proizvod na prebivalca v različnih državah v Evropi? Naredil bom dva modela, da vidim, ali ima članstvo v EU vpliv.
Vse podatke razen indeksa življenjskih stroškov sem dobil na spletni strani World Bank, podatke o indeksu življenjskih stroškov pa na spletni strani Numbeo.
V prvem koraku sem prebral podatke s funkcijo read.table
VPRASANJE2 <- read.table("~/Faks/Magistrski študij/1. letnik, 1. semester/Metode in tehnike raziskovalnega dela/Domača naloga R/PodatkiVprašanje2.csv", header=TRUE, sep=";", dec=",")BDP na prebivalca: Je moja odvisna spremenljivka, merjena je v dolarjih.
Ginijev koeficient: Merjeno v odstotkih. 0 pomeni popolno enakost, 100 pomeni popolno neenakost.
Stopnja brezposelnosti: Merjeno v odstotkih in sicer odstotek celotne delovne sile.
Indeks zivljenskih stroskov: Relativne cene potrošniških dobrin glede na cene v New York City-ju. Država z indeksom 120 bi bila za 20% dražja od New Yorka.
Članica EU: 0-> ni članica, 1-> je članica.
head(VPRASANJE2)## Drzava BDPnaPrebivalcav. GinijevKoeficient.
## 1 Albania 5343.04 29.4
## 2 Austria 48789.50 29.8
## 3 Belarus 6542.86 24.4
## 4 Belgium 45609.00 26.0
## 5 Bulgaria 10148.34 40.5
## 6 Croatia 14269.91 29.5
## Brezposelnost..CelotneDelovneSile. IndeksZivljenjskihStroskov ClanicaEU
## 1 13.07 36.4 0
## 2 5.20 70.4 1
## 3 4.05 34.7 0
## 4 5.56 71.8 1
## 5 5.13 36.7 1
## 6 7.51 49.7 1Kategorialni spremenljivki ClanicaEU sem pripisal vsebinsko vrednost, nato pa sem s summary predstavil opisno statistiko spremenljivk.
VPRASANJE2$ClanicaEUF <- factor(VPRASANJE2$ClanicaEU,
levels = c(0,1),
labels = c("Ne", "Da"))summary(VPRASANJE2[c(-1,-6)])## BDPnaPrebivalcav. GinijevKoeficient. Brezposelnost..CelotneDelovneSile.
## Min. : 3752 Min. :24.00 Min. : 1.210
## 1st Qu.: 13047 1st Qu.:27.10 1st Qu.: 4.990
## Median : 23595 Median :30.70 Median : 5.640
## Mean : 31326 Mean :30.73 Mean : 7.095
## 3rd Qu.: 46750 3rd Qu.:34.60 3rd Qu.: 8.290
## Max. :116905 Max. :40.50 Max. :16.600
## IndeksZivljenjskihStroskov ClanicaEUF
## Min. : 31.60 Ne: 9
## 1st Qu.: 39.20 Da:24
## Median : 53.40
## Mean : 55.98
## 3rd Qu.: 70.30
## Max. :100.50Lahko vidimo, da sta spremenljivki BDP na prebivalca in Brezposelnost rahlo asimetrični v desno.
Nato prevermi povezavo med našimi spremenljivkami s korelacijsko matriko.
library(car)## Loading required package: carDatascatterplotMatrix(VPRASANJE2[, c(-1, -6, -7)],
smooth= FALSE)
Kot vidimo, nebomo imeli nobenih težav s linearnostjo.
library(Hmisc)##
## Attaching package: 'Hmisc'## The following objects are masked from 'package:base':
##
## format.pval, unitsrcorr(as.matrix(VPRASANJE2[, c(-1, -6, -7)]))## BDPnaPrebivalcav. GinijevKoeficient.
## BDPnaPrebivalcav. 1.00 -0.17
## GinijevKoeficient. -0.17 1.00
## Brezposelnost..CelotneDelovneSile. -0.20 0.35
## IndeksZivljenjskihStroskov 0.85 -0.27
## Brezposelnost..CelotneDelovneSile.
## BDPnaPrebivalcav. -0.20
## GinijevKoeficient. 0.35
## Brezposelnost..CelotneDelovneSile. 1.00
## IndeksZivljenjskihStroskov -0.18
## IndeksZivljenjskihStroskov
## BDPnaPrebivalcav. 0.85
## GinijevKoeficient. -0.27
## Brezposelnost..CelotneDelovneSile. -0.18
## IndeksZivljenjskihStroskov 1.00
##
## n= 33
##
##
## P
## BDPnaPrebivalcav. GinijevKoeficient.
## BDPnaPrebivalcav. 0.3570
## GinijevKoeficient. 0.3570
## Brezposelnost..CelotneDelovneSile. 0.2618 0.0463
## IndeksZivljenjskihStroskov 0.0000 0.1309
## Brezposelnost..CelotneDelovneSile.
## BDPnaPrebivalcav. 0.2618
## GinijevKoeficient. 0.0463
## Brezposelnost..CelotneDelovneSile.
## IndeksZivljenjskihStroskov 0.3228
## IndeksZivljenjskihStroskov
## BDPnaPrebivalcav. 0.0000
## GinijevKoeficient. 0.1309
## Brezposelnost..CelotneDelovneSile. 0.3228
## IndeksZivljenjskihStroskovSpremenljivke so tako pozitivno kot tudi negativno povezane, kot vidimo iz zgornje matrike. Ko gledamo p vrednosti, vidimo, da sta najpomembnejši povezava med BDP na prebivalca in indeksom življenskih stroškov ter povezava med Brezposelnostjo in Ginijevim koeficientom.
Spodaj naredimo linearni model. Preverimo še VIF, če je pri kateri spremenljivki višji od 5, ga moramo dati strani.
fit1 <- lm(BDPnaPrebivalcav. ~ GinijevKoeficient. + Brezposelnost..CelotneDelovneSile. + IndeksZivljenjskihStroskov, data = VPRASANJE2)
vif(fit1)## GinijevKoeficient. Brezposelnost..CelotneDelovneSile.
## 1.199115 1.148939
## IndeksZivljenjskihStroskov
## 1.087038mean(vif(fit1))## [1] 1.145031Nimamo problema, vsi so dosti nižji.
Spodaj v tabelo dodamo dva stolpca, s katerima prevberimo standardne ostanke in cookovo razdalja. S tem se bomo poizkušali znebiti morebitnih osamelcev.
VPRASANJE2$StdOstanki <- round(rstandard(fit1), 3)
VPRASANJE2$CooksD <- round(cooks.distance(fit1), 3)
hist(VPRASANJE2$StdOstanki,
xlab="Standardizirani ostanki",
ylab="Frekvenca",
main="Histogram standardiziranih ostankov")
Imamo enega osamelca.
shapiro.test(VPRASANJE2$StdOstanki)##
## Shapiro-Wilk normality test
##
## data: VPRASANJE2$StdOstanki
## W = 0.77753, p-value = 1.238e-05Zavrnemo H0, ker je p<0,001 in ugotavljamo, da ni normalne porazdelitve.
head(VPRASANJE2[order(VPRASANJE2$StdOstanki),], 4)## Drzava BDPnaPrebivalcav. GinijevKoeficient.
## 33 Iceland 58848.40 24.8
## 19 Malta 29597.64 31.4
## 15 Italy 31922.92 35.2
## 11 France 39179.74 30.7
## Brezposelnost..CelotneDelovneSile. IndeksZivljenjskihStroskov ClanicaEU
## 33 5.50 100.5 0
## 19 4.35 67.5 1
## 15 9.16 67.3 1
## 11 8.01 74.1 1
## ClanicaEUF StdOstanki CooksD
## 33 Ne -2.009 0.324
## 19 Da -1.327 0.031
## 15 Da -1.112 0.033
## 11 Da -1.026 0.019Imamo enega osamelca, in sicer je to Islandija.
head(VPRASANJE2[order(-VPRASANJE2$StdOstanki),], 4)## Drzava BDPnaPrebivalcav. GinijevKoeficient.
## 18 Luxembourg 116905.37 33.4
## 14 Ireland 85973.09 29.2
## 31 North Macedonia 5965.50 33.5
## 28 Sweden 52837.90 28.9
## Brezposelnost..CelotneDelovneSile. IndeksZivljenjskihStroskov ClanicaEU
## 18 6.77 81.9 1
## 14 5.62 75.9 1
## 31 16.60 31.6 0
## 28 8.29 69.8 1
## ClanicaEUF StdOstanki CooksD
## 18 Da 4.095 0.636
## 14 Da 2.304 0.100
## 31 Ne 0.680 0.044
## 28 Da 0.485 0.004Tukaj pa imamo dva osamelca, in sicer Luksemburg in Irsko.
hist(VPRASANJE2$CooksD,
xlab="Cookove razdalje",
ylab="Frekvenca",
main="Histogram Cookovih razdalj")
head(VPRASANJE2[order(-VPRASANJE2$CooksD),], 6)## Drzava BDPnaPrebivalcav. GinijevKoeficient.
## 18 Luxembourg 116905.37 33.4
## 33 Iceland 58848.40 24.8
## 14 Ireland 85973.09 29.2
## 12 Greece 17617.29 33.6
## 31 North Macedonia 5965.50 33.5
## 15 Italy 31922.92 35.2
## Brezposelnost..CelotneDelovneSile. IndeksZivljenjskihStroskov ClanicaEU
## 18 6.77 81.9 1
## 33 5.50 100.5 0
## 14 5.62 75.9 1
## 12 15.90 55.7 1
## 31 16.60 31.6 0
## 15 9.16 67.3 1
## ClanicaEUF StdOstanki CooksD
## 18 Da 4.095 0.636
## 33 Ne -2.009 0.324
## 14 Da 2.304 0.100
## 12 Da -0.828 0.047
## 31 Ne 0.680 0.044
## 15 Da -1.112 0.033Tudi tukaj imamo tri osamelce, enake kot prej, zato lahko vse te tri odstranimo.
library(dplyr)##
## Attaching package: 'dplyr'## The following objects are masked from 'package:Hmisc':
##
## src, summarize## The following object is masked from 'package:car':
##
## recode## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionVPRASANJE2 <- VPRASANJE2 %>%
dplyr::filter(!Drzava %in% c("Luxembourg", "Iceland", "Ireland"))Ocenimo model
fit1 <- lm(BDPnaPrebivalcav. ~ GinijevKoeficient. + Brezposelnost..CelotneDelovneSile. + IndeksZivljenjskihStroskov, data = VPRASANJE2)
vif(fit1)## GinijevKoeficient. Brezposelnost..CelotneDelovneSile.
## 1.168919 1.140285
## IndeksZivljenjskihStroskov
## 1.055027VPRASANJE2$StdOstanki <- round(rstandard(fit1), 3)
VPRASANJE2$CooksD <- round(cooks.distance(fit1), 3)Preverimo homoskedastičnost in linearnost.
VPRASANJE2$StdOcenaVred <- scale(fit1$fitted.values)
library(car)
scatterplot(y=VPRASANJE2$StdOstanki, x=VPRASANJE2$StdOcenaVred,
ylab="Standardizirani ostanki",
xlab="Standardizirane ocenjene vrednosti",
regLine=FALSE,
smooth=FALSE)
S Breusch Paganovim testom preverimo homoskedastičnost. Glede linearnosti ni problema.
library(olsrr)##
## Attaching package: 'olsrr'## The following object is masked from 'package:datasets':
##
## riversols_test_breusch_pagan(fit1)##
## Breusch Pagan Test for Heteroskedasticity
## -----------------------------------------
## Ho: the variance is constant
## Ha: the variance is not constant
##
## Data
## ---------------------------------------------
## Response : BDPnaPrebivalcav.
## Variables: fitted values of BDPnaPrebivalcav.
##
## Test Summary
## -----------------------------
## DF = 1
## Chi2 = 3.859761
## Prob > Chi2 = 0.04945736Vidimo, da smo kršili predpostavko homoskedastičnosti pri p<0,05. Preden ocenimo naslednji model, moramo narediti popravek s funkcijo lm_robust. HC1 samo pomeni, je bodo te napake delo avtorja Huberja-Whita.
library(estimatr)fit1 <- lm_robust(BDPnaPrebivalcav. ~ GinijevKoeficient. + Brezposelnost..CelotneDelovneSile. + IndeksZivljenjskihStroskov, se_type="HC1", data = VPRASANJE2)summary(fit1)##
## Call:
## lm_robust(formula = BDPnaPrebivalcav. ~ GinijevKoeficient. +
## Brezposelnost..CelotneDelovneSile. + IndeksZivljenjskihStroskov,
## data = VPRASANJE2, se_type = "HC1")
##
## Standard error type: HC1
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -23816.04 6510.73 -3.6580 1.133e-03
## GinijevKoeficient. -79.88 167.39 -0.4772 6.372e-01
## Brezposelnost..CelotneDelovneSile. -353.30 288.41 -1.2250 2.316e-01
## IndeksZivljenjskihStroskov 1030.04 66.12 15.5773 1.061e-14
## CI Lower CI Upper DF
## (Intercept) -37199.0 -10433.0 26
## GinijevKoeficient. -423.9 264.2 26
## Brezposelnost..CelotneDelovneSile. -946.1 239.5 26
## IndeksZivljenjskihStroskov 894.1 1166.0 26
##
## Multiple R-squared: 0.9117 , Adjusted R-squared: 0.9015
## F-statistic: 95.3 on 3 and 26 DF, p-value: 3.773e-14Statistično značilna je pojasnjevalna spremenljivka Indeks zivljenskih stroskov. Primer razlage: Če se indeks življenjskih stroškov poveča za 1, se BDP na prebivalca poveča za 1030,04 dolarja pri p<0,0001. R na kvadrat je multipli koeficient determinacije in nam pove, da lahko 91,17% vse variabilnosti BDP-ja na prebivalca pojasnimo z linearnim vplivom vseh pojasnjevalnih spremenljivk.
V spodnjem modelu dodajamo še podatek o članstvu v EU.
fit2 <- lm_robust(BDPnaPrebivalcav. ~ GinijevKoeficient. + Brezposelnost..CelotneDelovneSile. + IndeksZivljenjskihStroskov + ClanicaEUF, se_type="HC1", data = VPRASANJE2)summary(fit2)##
## Call:
## lm_robust(formula = BDPnaPrebivalcav. ~ GinijevKoeficient. +
## Brezposelnost..CelotneDelovneSile. + IndeksZivljenjskihStroskov +
## ClanicaEUF, data = VPRASANJE2, se_type = "HC1")
##
## Standard error type: HC1
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -22664.2 6599.60 -3.4342 2.082e-03
## GinijevKoeficient. -114.3 171.53 -0.6662 5.114e-01
## Brezposelnost..CelotneDelovneSile. -336.0 301.28 -1.1151 2.754e-01
## IndeksZivljenjskihStroskov 1008.5 72.22 13.9643 2.610e-13
## ClanicaEUFDa 1265.5 1616.17 0.7830 4.410e-01
## CI Lower CI Upper DF
## (Intercept) -36256.4 -9072.1 25
## GinijevKoeficient. -467.5 239.0 25
## Brezposelnost..CelotneDelovneSile. -956.5 284.5 25
## IndeksZivljenjskihStroskov 859.7 1157.2 25
## ClanicaEUFDa -2063.1 4594.0 25
##
## Multiple R-squared: 0.9124 , Adjusted R-squared: 0.8984
## F-statistic: 70.31 on 4 and 25 DF, p-value: 3.122e-13Še vedno je statistično značilen samo indeks življenjskih stroškov. Tukaj R na kvadrat predpostavlja da lahko 91,24% vse variabilnosti BDP-ja na prebivalca pojasnimo z linearnim vplivom vseh pojasnjevalnih spremenljivk.
Ker smo modela ocenili z robustnimi standardnimi napakami, uporabimo waldtest
library(lmtest)## Loading required package: zoo##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numericwaldtest(fit1, fit2, test="F")## Wald test
##
## Model 1: BDPnaPrebivalcav. ~ GinijevKoeficient. + Brezposelnost..CelotneDelovneSile. +
## IndeksZivljenjskihStroskov
## Model 2: BDPnaPrebivalcav. ~ GinijevKoeficient. + Brezposelnost..CelotneDelovneSile. +
## IndeksZivljenjskihStroskov + ClanicaEUF
## Res.Df Df F Pr(>F)
## 1 26
## 2 25 1 0.6131 0.441Dobimo vrednost večjo od 0,05 in lahko rečemo, da drugi model, ki vključuje tudi članstvo v EU, bolje pojasnjuje.
Edina statistično značilna spremenljivka v prvem in drugem modelu je Indeks življenskih stroškov, lahko pa rečemo, da je drugi model boljši kot prvi.