Príklad s údajmi

# install.packages("quantreg")  
library(quantreg)
library(ggplot2)
library(dplyr)
df <- read_csv("smart_logistics_dataset.csv")
# Nastavenie moho datasetu

Scatterplot údajov

ggplot(df, aes(x = Demand_Forecast, y = Waiting_Time)) + # Vytvorenie základného bodového grafu, kde sa definuje závislá a nezávislá premenná.
  geom_point(alpha = 0.6) + # Zobrazenie jednotlivých pozorovaní
  labs(
    title = "Vzťah medzi dopytom a čakacím časom",
    x = "Demand Forecast",
    y = "Waiting Time"
  ) + # Nastavenie grafu
  theme_minimal()

Odhad OLS a kvantilovej regresie

We first estimate the standard OLS model, and then quantile regressions for:

m_ols <- lm(Waiting_Time ~ Demand_Forecast, data = df)
m_q10 <- rq(Waiting_Time ~ Demand_Forecast, tau = 0.10, data = df)
m_q50 <- rq(Waiting_Time ~ Demand_Forecast, tau = 0.50, data = df)
m_q90 <- rq(Waiting_Time ~ Demand_Forecast, tau = 0.90, data = df)
# Odhad kvantilovej regresie pre 10, 50 a 90 kvantil, aby sme videli, ako sa vplyv mení v rôznych častiach rozdelenia.

Výsledky OLS

summary(m_ols)
## 
## Call:
## lm(formula = Waiting_Time ~ Demand_Forecast, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -25.5649 -12.5455  -0.3536  13.4133  25.6110 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     36.508501   1.590785  22.950   <2e-16 ***
## Demand_Forecast -0.007258   0.007645  -0.949    0.343    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 14.48 on 998 degrees of freedom
## Multiple R-squared:  0.0009025,  Adjusted R-squared:  -9.86e-05 
## F-statistic: 0.9015 on 1 and 998 DF,  p-value: 0.3426

Výsledky kvantilovej regresie

summary(m_q10)
## 
## Call: rq(formula = Waiting_Time ~ Demand_Forecast, tau = 0.1, data = df)
## 
## tau: [1] 0.1
## 
## Coefficients:
##                 coefficients lower bd upper bd
## (Intercept)     17.57471     14.02385 20.40110
## Demand_Forecast -0.01149     -0.02458  0.01392
summary(m_q50)
## 
## Call: rq(formula = Waiting_Time ~ Demand_Forecast, tau = 0.5, data = df)
## 
## tau: [1] 0.5
## 
## Coefficients:
##                 coefficients lower bd upper bd
## (Intercept)     35.00000     32.57501 40.70094
## Demand_Forecast  0.00000     -0.03103  0.00868
summary(m_q90)
## 
## Call: rq(formula = Waiting_Time ~ Demand_Forecast, tau = 0.9, data = df)
## 
## tau: [1] 0.9
## 
## Coefficients:
##                 coefficients lower bd upper bd
## (Intercept)     52.65116     50.87986 56.88108
## Demand_Forecast  0.01163     -0.00738  0.02191

Porovnanie odhadnutých koeficientov

coef_table <- rbind(
  OLS = coef(m_ols),
  Q10 = coef(m_q10),
  Q50 = coef(m_q50),
  Q90 = coef(m_q90)
)
# Spojenie všetkých odhadnutých koeficientov do jednej tabuľky pre porovnanie rozdielov medzi OLS a kvantilmi.
coef_table
##     (Intercept) Demand_Forecast
## OLS    36.50850   -7.258491e-03
## Q10    17.57471   -1.149425e-02
## Q50    35.00000   -4.449020e-17
## Q90    52.65116    1.162791e-02
# Zobrazenie tabuľky

Interpretácia

Ak sa koeficienty sklonu líšia medzi kvantilmi, znamená to, že vplyv vysvetľujúcej premennej nie je rovnaký v celom podmienenom rozdelení závislej premennej.

V mojem prípade vidíme, že:

v nižšom kvantile (Q10) je vzťah medzi Demand_Forecast a Waiting_Time negatívny, v mediáne (Q50) je efekt prakticky nulový, v hornom kvantile (Q90) sa efekt mení na pozitívny.

To znamená, že vplyv Demand_Forecast na Waiting_Time je heterogénny.

Grafické porovnanie regresných priamok

ggplot(df, aes(x = Demand_Forecast, y = Waiting_Time)) + #Zobrazenie pôvodných dát
  geom_point(alpha = 0.5) +
  geom_abline(
    intercept = coef(m_ols)[1],
    slope = coef(m_ols)[2],
    linewidth = 1
  ) + # Pridanie regresnej priamky OLS
  geom_abline(
    intercept = coef(m_q10)[1],
    slope = coef(m_q10)[2],
    linetype = "dashed",
    color = "red",
    linewidth = 1
  ) + # Pridanie kvantilovej regresie pre 10 % kvantil
  geom_abline(
    intercept = coef(m_q50)[1],
    slope = coef(m_q50)[2],
    linetype = "dotdash",
    color = "blue",
    linewidth = 1
  ) + # Pridanie mediánovej regresie 50 % kvantil
  geom_abline(
    intercept = coef(m_q90)[1],
    slope = coef(m_q90)[2],
    linetype = "longdash",
    color = "green",
    linewidth = 1
  ) + # Pridanie 90 % kvantilovej regresie
  labs(
    title = "Porovnanie OLS a kvantilových regresií",
    subtitle = "čierna = OLS, červená = Q10, modrá= Q50, zelená  = Q90",
    x = "Demand Forecast",
    y = "Waiting Time"
  ) + #Názva grafa
  theme_minimal() #Vizuálny štýl

Odhadovanie viacerých kvantilov naraz

Odhadnem celú sadu kvantilov. Napríklad od 10. do 90. percentilu.

taus <- seq(0.1, 0.9, by = 0.1)
# Vytvorenie postupnosti kvantilov od 10 % do 90 %.
m_many <- rq(Waiting_Time ~ Demand_Forecast,tau = taus, data = df)
summary(m_many)
## 
## Call: rq(formula = Waiting_Time ~ Demand_Forecast, tau = taus, data = df)
## 
## tau: [1] 0.1
## 
## Coefficients:
##                 coefficients lower bd upper bd
## (Intercept)     17.57471     14.02385 20.40110
## Demand_Forecast -0.01149     -0.02458  0.01392
## 
## Call: rq(formula = Waiting_Time ~ Demand_Forecast, tau = taus, data = df)
## 
## tau: [1] 0.2
## 
## Coefficients:
##                 coefficients lower bd upper bd
## (Intercept)     25.59483     22.44877 27.57589
## Demand_Forecast -0.02586     -0.03637 -0.01301
## 
## Call: rq(formula = Waiting_Time ~ Demand_Forecast, tau = taus, data = df)
## 
## tau: [1] 0.3
## 
## Coefficients:
##                 coefficients lower bd upper bd
## (Intercept)     27.91228     23.15614 33.57068
## Demand_Forecast -0.01754     -0.04968  0.00628
## 
## Call: rq(formula = Waiting_Time ~ Demand_Forecast, tau = taus, data = df)
## 
## tau: [1] 0.4
## 
## Coefficients:
##                 coefficients lower bd upper bd
## (Intercept)     30.00000     27.59142 34.69634
## Demand_Forecast  0.00000     -0.02478  0.01538
## 
## Call: rq(formula = Waiting_Time ~ Demand_Forecast, tau = taus, data = df)
## 
## tau: [1] 0.5
## 
## Coefficients:
##                 coefficients lower bd upper bd
## (Intercept)     35.00000     32.57501 40.70094
## Demand_Forecast  0.00000     -0.03103  0.00868
## 
## Call: rq(formula = Waiting_Time ~ Demand_Forecast, tau = taus, data = df)
## 
## tau: [1] 0.6
## 
## Coefficients:
##                 coefficients lower bd upper bd
## (Intercept)     42.40260     34.77767 46.72730
## Demand_Forecast -0.01299     -0.04286  0.01984
## 
## Call: rq(formula = Waiting_Time ~ Demand_Forecast, tau = taus, data = df)
## 
## tau: [1] 0.7
## 
## Coefficients:
##                 coefficients lower bd upper bd
## (Intercept)     46.74000     41.32850 49.63927
## Demand_Forecast -0.00667     -0.02832  0.04255
## 
## Call: rq(formula = Waiting_Time ~ Demand_Forecast, tau = taus, data = df)
## 
## tau: [1] 0.8
## 
## Coefficients:
##                 coefficients lower bd upper bd
## (Intercept)     49.46597     47.79816 57.67740
## Demand_Forecast  0.00524     -0.00956  0.01467
## 
## Call: rq(formula = Waiting_Time ~ Demand_Forecast, tau = taus, data = df)
## 
## tau: [1] 0.9
## 
## Coefficients:
##                 coefficients lower bd upper bd
## (Intercept)     52.65116     50.87986 56.88108
## Demand_Forecast  0.01163     -0.00738  0.02191
# Odhad kvantilovej regresie

Graf koeficientov sklonu naprieč kvantilmi

coef_many <- data.frame(
  tau = taus,
  intercept = coef(m_many)[1, ],
  slope = coef(m_many)[2, ]
)

coef_many
##          tau intercept         slope
## tau= 0.1 0.1  17.57471 -1.149425e-02
## tau= 0.2 0.2  25.59483 -2.586207e-02
## tau= 0.3 0.3  27.91228 -1.754386e-02
## tau= 0.4 0.4  30.00000  0.000000e+00
## tau= 0.5 0.5  35.00000 -4.449020e-17
## tau= 0.6 0.6  42.40260 -1.298701e-02
## tau= 0.7 0.7  46.74000 -6.666667e-03
## tau= 0.8 0.8  49.46597  5.235602e-03
## tau= 0.9 0.9  52.65116  1.162791e-02
# Extrahovanie interceptu a sklonu pre každý kvantil do tabuľky
ggplot(coef_many, aes(x = tau, y = slope)) +
  geom_line(linewidth = 1) +
  geom_point(size = 2) +
  geom_hline(yintercept = coef(m_ols)[2], linetype = "dashed") +
  labs(
    title = "Sklon koeficientu naprieč kvantilmi",
    subtitle = "čiarkovaná horizontálna čiara = OLS sklon",
    x = "Quantile (tau)",
    y = "Odhadnuté sklony regresných koeficientov"
  ) +
  theme_minimal()

# Vytvorenie grafu, ktorý ukazuje, ako sa mení sklon (vplyv x) podľa kvantilu (10,50,90) a pridanie OLS sklonu ako referencie pre porovnanie.

Interpretácia koeficientov

Z výsledkov kvantilovej regresie vidíme, že tento vplyv nie je konštantný:

v nižších kvantiloch (napr. τ = 0.1 až 0.3) je odhadovaný sklon záporný, v stredných kvantiloch (τ ≈ 0.4–0.6) je efekt blízky nule, v horných kvantiloch (τ ≈ 0.8–0.9) sa efekt mení na mierne pozitívny.

To znamená, že vplyv Demand_Forecast na Waiting_Time sa mení naprieč celým p odmieneným rozdelením závislej premennej. Pre τ = 0.90 je odhadovaný sklon približne 0.0116 Zvýšenie Demand_Forecast o jednu jednotku je spojené so zvýšením približne o 0.0116 jednotky v podmienenom 90. percentile Waiting_Time.

Voliteľný empirický príklad

m_ols <- lm(Waiting_Time ~ Demand_Forecast, data = df)
m_q25 <- rq(Waiting_Time ~ Demand_Forecast, tau = 0.25, data = df)
m_q50 <- rq(Waiting_Time ~ Demand_Forecast, tau = 0.50, data = df)
m_q75 <- rq(Waiting_Time ~ Demand_Forecast, tau = 0.75, data = df)
# Odhad kvantilovej regresie pre vybrané kvantily 25, 50, 75, ktorá umožňuje analyzovať vplyv premennej v rôznych častiach rozdelenia závislej premennej.
rbind(
  OLS = coef(m_ols),
  Q25 = coef(m_q25),
  Q50 = coef(m_q50),
  Q75 = coef(m_q75)
)
##     (Intercept) Demand_Forecast
## OLS    36.50850   -7.258491e-03
## Q25    27.00000   -2.439024e-02
## Q50    35.00000   -4.449020e-17
## Q75    47.45833    5.208333e-03
# Spojenie odhadnutých koeficientov do jednej tabuľky

Záver

V mojej analýze som zistil, že vplyv premennej Demand_Forecast na Waiting_Time nie je konštantný. V nižších kvantiloch je efekt negatívny, v stredných kvantiloch je takmer nulový a v horných kvantiloch sa mení na pozitívny. To znamená, že vzťah medzi premennými je heterogénny a závisí od úrovne čakacej doby.

Zároveň graf koeficientov naprieč kvantilmi ukazuje, že tento efekt sa postupne mení v celom rozdelení, čo by OLS regresia nedokázala zachytiť, pretože poskytuje iba priemerný efekt.

Celkovo môžem konštatovať, že kvantilová regresia je vhodným nástrojom na analýzu tohto typu dát, pretože umožňuje lepšie pochopiť rozdielne správanie systému pri nízkych, stredných a vysokých hodnotách závislej premennej.