Nieklasyczne metody statystyki

Regresja kwantylowa

Zadanie

Stwórz model objaśniający średnie wybagrodzenia.

Wykonaj:

1. porównanie KMNK oraz regresji kwantylowej (różno-poziomowej) dla zmiennej “earnings” - wynagrodzenia.
2. Dobierz i przetestuj predyktory, kwantyle dla modeli.
3. Wykonaj testy różnic współczynnikow dla finalnych modeli.

Dane z CPSSW9298:

• earnings - średnie wynagrodzenie w $/godzine - (average hourly earnings (sum of annual pretax wages, salaries, tips, and bonuses, divided by the number of hours worked annually).)
• degree - stopień edukacji - szkoła średnia lub licencjat (factor indicating highest educational degree (“bachelor” or”highschool”).)
• age - wiek w latach

Wynagrodzenia zostaną zamodelowane oddzielnie dla roku 1992 (ankieta1) oraz 1998 (ankieta2).

ankieta1 <- subset(CPSSW9298, year=="1992")
ankieta2 <- subset(CPSSW9298, year=="1998")

MODELE - PORÓWNANIE

Model KNMK

Na podstawie modelu KMNK stwierdzamy, że w 1992 i 1998 wszystkie zmienne objaśniające - stopień edukacji, płeć oraz wiek istotnie wpływały na poziom wynagrodzeń.

• Wraz z zwiększeniem poziomu swojej edukacji i zdobyciem licencjatu, przeciętne wynagrodzenia zwiększały się średnio w 1992 roku o 37% , a w 1998 o 38%.

• W 1992 roku Kobiety względem mężczyzn w przepadanej próbie, zarabiały przeciętnie mniej o 17%. Natomiast w 1998 ta różnica w zarobkach się zwiększyła i średnie wynagrodzenie kobiet było przeciętnie niższe o 18% od wynagrodzenia męzczyzn.

• W 1992 roku Wraz z wiekiem, przeciętne wynagrodzenia zwiększają się o 26%. W 1998 r. wraz ze starzeniem się tzn. z każdnym kolejnym rokiem, ankietowani zarabiali śrendio o 21% więcej.

dla 1992 roku

lm1 <- lm(log(earnings) ~ degree + gender + age, data = ankieta1)  # const nie jest istotna czyli rozkład jest skośny, więcej ludzi zarabia więcej niż mniej, DLATEGO  zmienna earnings jest w logarytmie !!
summary(lm1)

## 
## Call:
## lm(formula = log(earnings) ~ degree + gender + age, data = ankieta1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.80573 -0.26371  0.02737  0.29516  1.53710 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     1.486349   0.053186   27.95   <2e-16 ***
## degreebachelor  0.375044   0.010162   36.91   <2e-16 ***
## genderfemale   -0.167263   0.010009  -16.71   <2e-16 ***
## age             0.026398   0.001764   14.97   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4305 on 7586 degrees of freedom
## Multiple R-squared:  0.1872, Adjusted R-squared:  0.1869 
## F-statistic: 582.3 on 3 and 7586 DF,  p-value: < 2.2e-16

dla 1998 roku

## 
## Call:
## lm(formula = log(earnings) ~ degree + gender + age, data = ankieta2)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.96583 -0.27644  0.02536  0.30209  1.50215 
## 
## Coefficients:
##                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     1.78845    0.06230   28.71   <2e-16 ***
## degreebachelor  0.38277    0.01173   32.64   <2e-16 ***
## genderfemale   -0.18003    0.01182  -15.23   <2e-16 ***
## age             0.02142    0.00207   10.35   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4457 on 5907 degrees of freedom
## Multiple R-squared:  0.1828, Adjusted R-squared:  0.1824 
## F-statistic: 440.5 on 3 and 5907 DF,  p-value: < 2.2e-16

Model regresji kwantylowej

dla 1992 roku

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## 
## Call: rq(formula = log(earnings) ~ degree + gender + age, tau = kwantyle1, 
##     data = ankieta1)
## 
## tau: [1] 0.25
## 
## Coefficients:
##                Value     Std. Error t value   Pr(>|t|) 
## (Intercept)      1.32454   0.07633   17.35368   0.00000
## degreebachelor   0.41185   0.01357   30.35236   0.00000
## genderfemale    -0.14206   0.01416  -10.03079   0.00000
## age              0.02206   0.00249    8.87098   0.00000
## 
## Call: rq(formula = log(earnings) ~ degree + gender + age, tau = kwantyle1, 
##     data = ankieta1)
## 
## tau: [1] 0.5
## 
## Coefficients:
##                Value     Std. Error t value   Pr(>|t|) 
## (Intercept)      1.41968   0.06021   23.57894   0.00000
## degreebachelor   0.39768   0.01159   34.31471   0.00000
## genderfemale    -0.18271   0.01200  -15.23177   0.00000
## age              0.02955   0.00205   14.38757   0.00000
## 
## Call: rq(formula = log(earnings) ~ degree + gender + age, tau = kwantyle1, 
##     data = ankieta1)
## 
## tau: [1] 0.75
## 
## Coefficients:
##                Value     Std. Error t value   Pr(>|t|) 
## (Intercept)      1.69431   0.06590   25.70880   0.00000
## degreebachelor   0.35528   0.01217   29.18287   0.00000
## genderfemale    -0.18447   0.01360  -13.56140   0.00000
## age              0.02978   0.00210   14.17697   0.00000

dla 1998 roku

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## 
## Call: rq(formula = log(earnings) ~ degree + gender + age, tau = kwantyle2, 
##     data = ankieta2)
## 
## tau: [1] 0.25
## 
## Coefficients:
##                Value     Std. Error t value   Pr(>|t|) 
## (Intercept)      1.53975   0.09458   16.27957   0.00000
## degreebachelor   0.39850   0.01674   23.80190   0.00000
## genderfemale    -0.18133   0.01728  -10.49616   0.00000
## age              0.02026   0.00312    6.50131   0.00000
## 
## Call: rq(formula = log(earnings) ~ degree + gender + age, tau = kwantyle2, 
##     data = ankieta2)
## 
## tau: [1] 0.5
## 
## Coefficients:
##                Value     Std. Error t value   Pr(>|t|) 
## (Intercept)      1.71413   0.07286   23.52507   0.00000
## degreebachelor   0.39563   0.01418   27.91005   0.00000
## genderfemale    -0.19526   0.01392  -14.02939   0.00000
## age              0.02479   0.00242   10.25211   0.00000
## 
## Call: rq(formula = log(earnings) ~ degree + gender + age, tau = kwantyle2, 
##     data = ankieta2)
## 
## tau: [1] 0.75
## 
## Coefficients:
##                Value     Std. Error t value   Pr(>|t|) 
## (Intercept)      1.94914   0.07408   26.30978   0.00000
## degreebachelor   0.36898   0.01402   26.32117   0.00000
## genderfemale    -0.21198   0.01445  -14.66584   0.00000
## age              0.02666   0.00250   10.65240   0.00000

Interpretacja oszacowań (0.25,0.5,0.75)

## 
## Wyniki regresji kwantylowych w 1992r.
## ============================================
##                     Dependent variable:     
##                -----------------------------
##                        log(earnings)        
##                   (1)       (2)       (3)   
## --------------------------------------------
## degreebachelor 0.412***  0.398***  0.355*** 
##                 (0.015)   (0.012)   (0.011) 
##                                             
## genderfemale   -0.142*** -0.183*** -0.184***
##                 (0.015)   (0.012)   (0.011) 
##                                             
## age            0.022***  0.030***  0.030*** 
##                 (0.003)   (0.002)   (0.002) 
##                                             
## Constant       1.325***  1.420***  1.694*** 
##                 (0.080)   (0.062)   (0.058) 
##                                             
## --------------------------------------------
## Observations     7,590     7,590     7,590  
## ============================================
## Note:            *p<0.1; **p<0.05; ***p<0.01

Wszystkie zmienne istotnie wpływają na przeciętne wynagrodzenia w 1992r.

Regresja kwantylowa pozwola dostrzec różnice między zależnościami dla zarabiających najmniej i najwięcej.

• W przypadku wpływu płci na zarobki, widzimy, że dla zarabiających w pierwszej ćwiartce, kobiety zarabiały średnio mniej o 14,2% w stostunku do mężczyz. Natomiast rozpartując grupe zarabiających najwięcej, tzn. ankietowanych z ostatniej ćwiatki, kobiety zarabiały przeciętnie mniej o 18,4% niż mężczyźni. Zatem im wyższe zarobki tym różnice w wynagrodzeniach kobiet i mężczyzn są większe. Inaczej można powiedzieć, że przy niższych stawkach dyskrymiacja ze względu na płeć jest mniejsza.

• Rozważając wpływ poziomu edukacji, W przypadku osób zarabiających najmniej (pierwsza ćwiartka), jeżeli ankietowany posiadał stopień licencjata śrendie zarobki były przeciętnie wyższe o 41% w porówaniu do osób o średnim wykształceniu. Dla osób zarabiających najwięcej ta różnica w zależności od posiomu wykształcenia była mniejsza, gdyż stopnień licenccjata dawał przeciętnie wyższe zarobki o 35%. Zatem można wnioskować, że gdy się zarabia więcej, stopień wykształcenia w mniejszym stopniu przyczynia się do wzrostu zaraobków niż na poziomie niskich wynagrodzeń.

• Rozpatrując wpływ wieku na zarobki widzimy, że dla osób młodych każdny kolejny rok kiedy stają się starsi, zarobki były większe o 2,2%. Natomiast dla grupy zarabiajacych w przedziale środkowym (30lat i więcej), wraz z wiekiem wynagrodzenia były wyższe o 3%. Zatem można dostrzec, że dla grupy osób najmłodszych wśród badanych, wzrost dochodów był mniejszy wraz wiekiem, niż dla osób które miały więcej niż 30 lat.

## [1] 30

## 
## Wyniki regresji kwantylowych w 1998r.
## ============================================
##                     Dependent variable:     
##                -----------------------------
##                        log(earnings)        
##                   (1)       (2)       (3)   
## --------------------------------------------
## degreebachelor 0.398***  0.396***  0.369*** 
##                 (0.016)   (0.014)   (0.015) 
##                                             
## genderfemale   -0.181*** -0.195*** -0.212***
##                 (0.016)   (0.014)   (0.015) 
##                                             
## age            0.020***  0.025***  0.027*** 
##                 (0.003)   (0.002)   (0.003) 
##                                             
## Constant       1.540***  1.714***  1.949*** 
##                 (0.085)   (0.075)   (0.079) 
##                                             
## --------------------------------------------
## Observations     5,911     5,911     5,911  
## ============================================
## Note:            *p<0.1; **p<0.05; ***p<0.01

Tabela. Porównanie współczynników

Porównanie współczynników w modelu KMNK i regresji kwantylowej (0.1, 0.25, 0.5, 0.75, 0.90, 0.95))

Oszacowania z KMNK oraz quantreg dla roku 1992

	OLS	\(\tau_{0.10}\)	\(\tau_{0.25}\)	\(\tau_{0.50}\)	\(\tau_{0.75}\)	\(\tau_{0.90}\)	\(\tau_{0.95}\)
(Intercept)	1.4863490	1.1233290	1.3245402	1.4196786	1.6943127	1.9324262	2.1211185
degreebachelor	0.3750441	0.3897862	0.4118489	0.3976830	0.3552754	0.3369236	0.3226939
genderfemale	-0.1672629	-0.1140067	-0.1420630	-0.1827112	-0.1844718	-0.1827728	-0.1965213
age	0.0263975	0.0186965	0.0220623	0.0295486	0.0297762	0.0296363	0.0280745

Oszacowania z KMNK oraz quantreg dla 1998 roku

	OLS	\(\tau_{0.10}\)	\(\tau_{0.25}\)	\(\tau_{0.50}\)	\(\tau_{0.75}\)	\(\tau_{0.90}\)	\(\tau_{0.95}\)
(Intercept)	1.7884485	1.1233290	1.3245402	1.4196786	1.6943127	1.9324262	2.1211185
degreebachelor	0.3827709	0.3897862	0.4118489	0.3976830	0.3552754	0.3369236	0.3226939
genderfemale	-0.1800279	-0.1140067	-0.1420630	-0.1827112	-0.1844718	-0.1827728	-0.1965213
age	0.0214183	0.0186965	0.0220623	0.0295486	0.0297762	0.0296363	0.0280745

TESTY RÓŻNIC WSPÓŁCZYNNIKÓW

Badanie statystycznej różnicy między 25. i 50. kwantylem warunkowym - Test Walda.

testowanie wspólne (joint=TRUE)

• Dla roku 1992 (dane=ankieta1) Na podstawie poniższego testu na istotność różnicy między 25 i 50 kwantylem warunkowym, odrzucamy hipotezę zerową, na rzecz H1, ponieważ p-value< alfa=0.05. Zatem różnica pomiędzy współczynnikami dla 25 i 50 kwantylu jest statystycznie istotna.
• Dla roku 1998 (dane=ankieta2) Na podstwie testu Walda stwierdzamy brak podstaw do odrzucenia hipotezy zerowej, ponieważ p-value=0.210 > alfa=0.05. Co jest równoznanczne z tym, że współczynniki dla oszacowań w modelu dla 25 i 50 kwantyla spełniają pewne liniowe ograniczenia. Brak isotności statystycznej.

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Quantile Regression Analysis of Deviance Table
## 
## Model: log(earnings) ~ degree + gender + age
## Joint Test of Equality of Slopes: tau in {  0.25 0.5  }
## 
##   Df Resid Df F value    Pr(>F)    
## 1  3    15177  8.1014 2.178e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Quantile Regression Analysis of Deviance Table
## 
## Model: log(earnings) ~ degree + gender + age
## Joint Test of Equality of Slopes: tau in {  0.25 0.5  }
## 
##   Df Resid Df F value Pr(>F)
## 1  3    11819  1.5073 0.2104

Badanie statystycznej różnicy między 25, 50 i 75 kwantylem warunkowym - Test Walda

testowanie wspólne (joint=TRUE)

• dla 1992r

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Quantile Regression Analysis of Deviance Table
## 
## Model: log(earnings) ~ degree + gender + age
## Joint Test of Equality of Slopes: tau in {  0.25 0.5 0.75  }
## 
##   Df Resid Df F value    Pr(>F)    
## 1  6    22764  7.1177 1.356e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

• dla 1998r

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Quantile Regression Analysis of Deviance Table
## 
## Model: log(earnings) ~ degree + gender + age
## Joint Test of Equality of Slopes: tau in {  0.25 0.5 0.75  }
## 
##   Df Resid Df F value  Pr(>F)  
## 1  6    17727  2.0933 0.05064 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Uzyskano wyniki podobne co dla testów dla dwoch kwantyli.

Dla roku 1992 odrzucamy Ho, na rzecz hipotezy alternatywnej (p-value<0.05), zatem istnieją istotne różnice pomiędzy współczynnikami dla osób zarabiających w pierwszym kwantylu 0.25, medianie i tych zarabiających najwięcej. Natomiast takie różnice w roku 1998 “osłabiły się” i występlują przy założeniu istotności alfa=10% (p-value=0.05064).

Test WALDA - Testowanie oddzielne (joint=FALSE)

Test dla różnic między dwoma kwantylami 25 i 50

## Quantile Regression Analysis of Deviance Table
## 
## Model: log(earnings) ~ degree + gender + age
## Tests of Equality of Distinct Slopes: tau in {  0.25 0.5  }
## 
##                Df Resid Df F value    Pr(>F)    
## degreebachelor  1    15179  1.2313 0.2671772    
## genderfemale    1    15179 10.4838 0.0012068 ** 
## age             1    15179 11.1764 0.0008304 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

## Quantile Regression Analysis of Deviance Table
## 
## Model: log(earnings) ~ degree + gender + age
## Tests of Equality of Distinct Slopes: tau in {  0.25 0.5  }
## 
##                Df Resid Df F value  Pr(>F)  
## degreebachelor  1    11821  0.0416 0.83830  
## genderfemale    1    11821  0.9879 0.32028  
## age             1    11821  3.3290 0.06809 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Test dla różnic między trzema kwantylami 25,50 i 75

## Quantile Regression Analysis of Deviance Table
## 
## Model: log(earnings) ~ degree + gender + age
## Tests of Equality of Distinct Slopes: tau in {  0.25 0.5 0.75  }
## 
##                Df Resid Df F value    Pr(>F)    
## degreebachelor  2    22768  9.1631 0.0001052 ***
## genderfemale    2    22768  5.4138 0.0044606 ** 
## age             2    22768  5.7401 0.0032191 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

## Quantile Regression Analysis of Deviance Table
## 
## Model: log(earnings) ~ degree + gender + age
## Tests of Equality of Distinct Slopes: tau in {  0.25 0.5 0.75  }
## 
##                Df Resid Df F value Pr(>F)
## degreebachelor  2    17731  2.1488 0.1167
## genderfemale    2    17731  1.4965 0.2239
## age             2    17731  2.2269 0.1079

Na podstawie powyższych testów zauważamy, że dla oszacowań dla 1998 roku występuje brak istotnych różnic pomiędzy współczynnikami dla wybranych kwantyli, w przeciwieństwie do roku 1992, kiedy taki różnice były statystycznie istotne.

Dobroć dopasowania vs R-kwadrat

Porówanie R^2 i dobroci dopasowania w modelach KMNK i regresji kwantylowej dla danych z roku 1992

Dobroć dopasowania dla modelu kwantylowego (1992)

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## [1] 0.1598034

R-kwadrat (model KMNK)

## [1] 0.1871756

Lepiej wypada regresja liniowa niż regresja kwantylowa, ponieważ R-kwadrat jest wyższe.

Porówanie R^2 i dobroci dopasowania w modelach KMNK i regresji kwantylowej dla danych z roku 1998

Dobroć dopasowania dla modelu kwantylowego (1998)

## [1] 0.1410847

R-kwadrat (model KMNK)

## [1] 0.1828058

Dl adanych z roku 1998 również lepiej wypada regresja liniowa KMNK niż regresja kwantylowa, ponieważ R^2=18,28% jest większe od miernika dobrocji dopasowania.

Kryterium AIC

## Warning in AIC.default(model_lm_1992, model_kwantylowy_1992, model_lm_1998, :
## nie wszystkie modele są dopasowane do tej samej liczby obserwacji

##                       df      AIC
## model_lm_1992          5 8750.943
## model_kwantylowy_1992  4 9249.473
## model_lm_1998          5 7228.160
## model_kwantylowy_1998  4 7555.391

Kryterium AIC również potwierdza, że lepszym modelem w obu latach okazał się model regresji liniowej (KMNK).

DODATKOWE PORÓWNANIE MODELI

Model KMNK vs Model Medianowy

## Warning: pakiet 'sjPlot' został zbudowany w wersji R 4.3.3

## Learn more about sjPlot with 'browseVignettes("sjPlot")'.

Poniższe graficzne porównanie modeli pozawala na zobarazowanie różnic oszacowań parametrów. Dzięki temeu można zidentyfikować czy model liniowy niedoszacowywał lub przeszacowywał parametry względem mediany, która może być lepszą miarą w omawianym przykładzie.

śrendia, mediana, histogram wynagrodzeń

## [1] 11.65373

## [1] 10.57692

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

Graficzne prównanie jak objaśniane zjawisko różni się w zależności od decyla oraz jakie wnioski można wysnuć z porówanania wyników do śrendiej (KMNK).

Dla uproszczenia na osi X są zarobki w [$/h]. Oś OY odzwierciedlają współczynniki dla poszczególnych decyli. Czerwona ciągła linia to poziom średniej (OLS).

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

• Wykres poziomu edukacji (degree)

Na wykresie widzimy , że im więcej ankietowani zarabiali, tym wpłw stopnia edukacji był silniejszy dodatnio - Czarna linia jest rosnąca. Różnice pomiędzy zarobkami, a w wpływem posiadania licencjatu były większe przy wyższych poziomach płac.

Dodaktowo dla wynagrodzeń poniżej mediany estymacja modelem regresji liniowej jest przeszacowana, a dla zarbiających powyżej mediany model KMNK niedoszacowaywał paramtery. (czerwona prosta poniżej czarnej).

• Wykres wynagrodzenia~płeć

Regresja prosta, która wyznacza parametry modelu na podstawie wartości średnich (mean), przeszacowywały negatywny wpływ płci na zarobiki w grupie zarabaijących najmniej. Na wykresie można zaonserwowac, że w grupie 10% najmniej zarabiających, kobiety zarabiały jedynie średnio 0,5$/h mniej niż mężczyźni. Jednowcześnie KMNK niedoszacowywało negatywnego wpływu płci na zarobki, inaczej klasyczny model zaniżał dyskryminacje płciową.W grupie 10% osób najwiecej zarabiających, tam średnio kobiety zarabiały mniej o 3,4-3,6dolarów/h.

Nie ma takiego przedziału zarobków, w którym kobiety zarabiły by więcej niż mężczyźni.

• Wykres zależności wynagrodzeń i wieku

W oparciu o powyższy wykres dostrzegamy, że model KMNK zawyżał parametr dla zmiennej wiek dla grup ankietowanych zarabiających poniżej mediany, natomiast zaniżał ten parametr w przedziałach wynagrodzeń powyżej mediany.