library(zoo)
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(tseries)
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
library(lmtest)
library(sandwich)
library(car)
## Loading required package: carData
rm(list=ls())

Využívame databázu hotel_bookings.cv. Zameriame sa na údaje, kde je počet zrušených a počet nezrušených rezervácií aspoň 1, aby sa nám jednoduchšie robila regresia.

udaje <- read.csv("Ekonometria/hotel_bookings.csv", header = TRUE, sep = ",", dec = ".")
udaje1 <- udaje[udaje$previous_cancellations > 0,]
udaje.new <- udaje1[udaje1$previous_bookings_not_canceled > 0, c("lead_time", "arrival_date_week_number", "arrival_date_day_of_month", "stays_in_weekend_nights", "stays_in_week_nights", "is_repeated_guest", "previous_cancellations", "previous_bookings_not_canceled")]
head(udaje.new)

Pozeráme sa na vzťah medzi lead_time, teda čas urobenia rezervácie pred začiatkom pobytu, a ostatných premenných. Model lineárnej regresie teda vyzerá nasledovne:

attach(udaje.new)
model <- lm(lead_time ~ +1 + arrival_date_week_number + arrival_date_day_of_month + stays_in_weekend_nights + stays_in_week_nights + is_repeated_guest + previous_cancellations + previous_bookings_not_canceled, data = udaje.new)
summary(model)
## 
## Call:
## lm(formula = lead_time ~ +1 + arrival_date_week_number + arrival_date_day_of_month + 
##     stays_in_weekend_nights + stays_in_week_nights + is_repeated_guest + 
##     previous_cancellations + previous_bookings_not_canceled, 
##     data = udaje.new)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -100.373  -21.442   -8.445    4.564  309.289 
## 
## Coefficients:
##                                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                      0.52933    7.51078   0.070 0.943837    
## arrival_date_week_number         0.45626    0.13330   3.423 0.000659 ***
## arrival_date_day_of_month        0.04418    0.22907   0.193 0.847137    
## stays_in_weekend_nights         13.21855    2.78788   4.741 2.62e-06 ***
## stays_in_week_nights             9.57520    1.56698   6.111 1.73e-09 ***
## is_repeated_guest              -14.67823    5.39897  -2.719 0.006732 ** 
## previous_cancellations           7.21382    0.77572   9.299  < 2e-16 ***
## previous_bookings_not_canceled  -0.94936    0.15481  -6.132 1.52e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 49.03 on 639 degrees of freedom
## Multiple R-squared:  0.2465, Adjusted R-squared:  0.2383 
## F-statistic: 29.87 on 7 and 639 DF,  p-value: < 2.2e-16

Vidíme, že hoco väčšina veličín je pre náš model štatisticky významná, model nie je optimálny. Podľa našich výsledkov napríklad deň v mesiaci nemá významný vplyv na to, ako skoro človek robí rezevráciu.

Na zlepšenie nášho modelu využijeme test RESET.

1. Test RESET (test chyby špecifikácie Ramseyho regresnej rovnice - Ramsey Reset Test)

Ide o najštandardnejší formálny test nesprávnej špecifikácie funkčnej formy ale dá sa použiť aj pre prípady, ak sme nešpecifikovali všetky vysvetľujúce premenné.

Myšlienka: Nech pôvodný model má tvar \[y_t = \beta_0 + \beta_1 x_{t10} + \dots +\beta_k x_{tk} + u_t\] Ak je náš model správne špecifikovaný, potom pridaním mocnín vyrovnaných hodnôt (napr. \(\hat y_t^2\), \(\hat{y}_t^3\)) by sa pôvodný model nemal podstatne zlepšiť, teda budeme testovať pôvodný model uvedený vyššie a model

\[y_t = \beta_0 + \beta_1 x_{t10} + \dots +\beta_k x_{tk} + \gamma_2\hat y_t^2 + \gamma_3\hat{y}_t^3 + u_t\] Budeme testovať hypotézu

\(H_0:\) model je správne špecifikovaný (\(\gamma_2 = \gamma_3 = 0\))

oproti

\(H_1:\) model je nesprávne špecifikovaný (\(\gamma_2 \ne 0 \quad \text{alebo} \quad \gamma_3 \ne 0\))

# Suppose your model is:
model <- lm(lead_time ~ +1 + arrival_date_week_number + arrival_date_day_of_month + stays_in_weekend_nights + stays_in_week_nights + is_repeated_guest + previous_cancellations + previous_bookings_not_canceled, data = udaje.new)

# RESET test from 'lmtest' package:
library(lmtest)
resettest(model)
## 
##  RESET test
## 
## data:  model
## RESET = 11.036, df1 = 2, df2 = 637, p-value = 1.943e-05

p-value nám v modeli vychádza < 0.05, teda zamietame hypotézu \(H_0:\). Náš model je pravdepodobne zle špecifikovaný, chýbajú mu niektoré dodatočné vysvetľujúce premenné, alebo je potrebné urobiť nelineárnu transformáciu používaných premenných.

2. Grafická analýza

Graf Residuals vs. Fitted

Grafická analýza vzťahu medzi vyrovnanými hodnotami náhodnej premennej a rezíduami:

plot(model, which = 1)

Model nám vykazuje nekonštantné reziduá, krivka najprv lineárne klesá a pri vyšších hodnotách stúpa. Z toho nám vyplýva, že model by mohol benefitovať z nelineárnych členov.

Grafy C+R **

Táto analýza nám môže pomôcť pri hľadaní odpovede na otázku, ktorú premennú by sme mali transfomovať pomocou nejakej známej funkcie. Vychádzajme z pôvodnej rovnice

\[y_t = \beta_0 + \beta_1 x_{t1} + \dots +\beta_k x_{tk} + u_t\] Túto rovnicu najprv odhadneme a potom vykresľujeme grafy, kde výraz component+residual (C+R) plot vykresľuje na zvislej osi \(\hat{\beta}_ix_{ti}+e_t\) a na vodorovnej osi vykresľuje hodnoty \(x_{ti}\)

Tieto grafy pomáhajú identifikovať nelineárne vzťahy pre každý regresor:

car::crPlots(model)

Vidíme, že najviac sa nám od priamky odkláňajú veličiny “stays_in_weekend_nights”, “stays_in_week_nights” a “previous_cancellations”. Tieto premenné skúsime nejako nelineárne v modeli upraviť.

3. Nelineárna špecifikácia

Častokrát môžeme aj zložitejšie nelineárne vzťahy modelovať s pomocou ich aproximácie s polynómom, teda v prípade kvadratických členov

\[y_t = \beta_0 + \beta_1 x_{t10} + \dots +\beta_k x_{tk} + \dots + \gamma_i\hat x_{ik}^2 + \dots + \gamma_j\hat x_{jk}^2 + \dots + u_t\]

3. Porovnanie základného a modifikovaného modelu

Predpokladajme, že sa pri nelineárnych úpravách pôvodnej rovnice dostaneme k zavedeniu kvadrátu premennej stays_in_weekend_nights, stays_in_week_nights a previous_cancellations, pri ktorých sa krivka najviac líši od priamky.

model <- lm(lead_time ~ +1 + arrival_date_week_number + arrival_date_day_of_month + stays_in_weekend_nights + stays_in_week_nights + is_repeated_guest + previous_cancellations + previous_bookings_not_canceled, data = udaje.new)
model_schooling_kvadr <- lm(lead_time ~ +1 + arrival_date_week_number + arrival_date_day_of_month + stays_in_weekend_nights + stays_in_week_nights + is_repeated_guest + previous_cancellations + previous_bookings_not_canceled + I(stays_in_weekend_nights^2) + I(stays_in_week_nights^2) + I(previous_cancellations^2), data = udaje.new)

summary(model_schooling_kvadr)
## 
## Call:
## lm(formula = lead_time ~ +1 + arrival_date_week_number + arrival_date_day_of_month + 
##     stays_in_weekend_nights + stays_in_week_nights + is_repeated_guest + 
##     previous_cancellations + previous_bookings_not_canceled + 
##     I(stays_in_weekend_nights^2) + I(stays_in_week_nights^2) + 
##     I(previous_cancellations^2), data = udaje.new)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -103.598  -21.541   -5.584    6.815  300.568 
## 
## Coefficients:
##                                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                      8.581430   8.370285   1.025 0.305647    
## arrival_date_week_number         0.539925   0.131477   4.107 4.54e-05 ***
## arrival_date_day_of_month        0.011556   0.224529   0.051 0.958970    
## stays_in_weekend_nights         27.234390   6.264363   4.348 1.60e-05 ***
## stays_in_week_nights            12.574622   3.456371   3.638 0.000297 ***
## is_repeated_guest              -14.007796   5.296977  -2.644 0.008383 ** 
## previous_cancellations          -8.654769   3.467859  -2.496 0.012823 *  
## previous_bookings_not_canceled  -0.458599   0.189878  -2.415 0.016006 *  
## I(stays_in_weekend_nights^2)    -7.144605   3.043809  -2.347 0.019218 *  
## I(stays_in_week_nights^2)        0.009575   0.550682   0.017 0.986133    
## I(previous_cancellations^2)      1.276130   0.275876   4.626 4.52e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 48.01 on 636 degrees of freedom
## Multiple R-squared:  0.2809, Adjusted R-squared:  0.2696 
## F-statistic: 24.84 on 10 and 636 DF,  p-value: < 2.2e-16

Výsledok kvadratického modelu nám ukazuje, že arrival_date_day_of_month stále nie je štatisticky významný pre model. Urobíme teda úpravu, kde túto veličinu z modelu vyhodíme, a taktiež vyhodíme druhú mocninu stays_in_week_nights, ktorá na model tiež nemá významný vplyv.

model_schooling_kvadr_2 <- lm(lead_time ~ +1 + arrival_date_week_number + stays_in_weekend_nights + stays_in_week_nights + is_repeated_guest + previous_cancellations + previous_bookings_not_canceled + I(stays_in_weekend_nights^2) + I(previous_cancellations^2), data = udaje.new)

summary(model_schooling_kvadr_2)
## 
## Call:
## lm(formula = lead_time ~ +1 + arrival_date_week_number + stays_in_weekend_nights + 
##     stays_in_week_nights + is_repeated_guest + previous_cancellations + 
##     previous_bookings_not_canceled + I(stays_in_weekend_nights^2) + 
##     I(previous_cancellations^2), data = udaje.new)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -103.428  -21.587   -5.511    6.797  300.723 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                      8.7255     7.3361   1.189  0.23473    
## arrival_date_week_number         0.5408     0.1300   4.161 3.60e-05 ***
## stays_in_weekend_nights         27.2067     6.0532   4.495 8.28e-06 ***
## stays_in_week_nights            12.6230     1.8670   6.761 3.09e-11 ***
## is_repeated_guest              -14.0229     5.2771  -2.657  0.00807 ** 
## previous_cancellations          -8.6645     3.4577  -2.506  0.01246 *  
## previous_bookings_not_canceled  -0.4583     0.1888  -2.428  0.01545 *  
## I(stays_in_weekend_nights^2)    -7.1175     2.6032  -2.734  0.00643 ** 
## I(previous_cancellations^2)      1.2766     0.2753   4.637 4.30e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 47.94 on 638 degrees of freedom
## Multiple R-squared:  0.2809, Adjusted R-squared:  0.2719 
## F-statistic: 31.15 on 8 and 638 DF,  p-value: < 2.2e-16

V tomto modeli sa nám všetky veličiny ukazujú viac alebo menej štatisticky významné. Porovnáme si tento nový kvadratický model s pôvodným lineárnym.

anova(model, model_schooling_kvadr_2)
resettest(model_schooling_kvadr_2)
## 
##  RESET test
## 
## data:  model_schooling_kvadr_2
## RESET = 12.353, df1 = 2, df2 = 636, p-value = 5.456e-06

Z výsledkov anova testu vidíme, že nový kvadratický model lepšie popisuje správanie hodnôt, ako pôvodný lineárny. Podľa RESET testu však ani nový kvadratický model nie je špecifikovaný správne a zišli by sa mu nejaké úpravy.