Standardowy błąd nachyleń w częściowej regresji liniowej ze znanymi punktami przerwania

9

Sytuacja

Mam zestaw danych z jednym zależnym y i jedna zmienna niezależna x. Chcę dopasować do ciągłej częściowej regresji liniowejk znane / ustalone punkty przerwania występujące w (a1,a2,,ak). Breakpoins są znane bez wątpliwości, więc nie chcę ich szacować. Następnie dopasowuję regresję (OLS) formularza

yi=β0+β1xi+β2max(xia1,0)+β3max(xia2,0)++βk+1max(xiak,0)+ϵi
Oto przykład w R
set.seed(123)
x <- c(1:10, 13:22)
y <- numeric(20)
y[1:10] <- 20:11 + rnorm(10, 0, 1.5)
y[11:20] <- seq(11, 15, len=10) + rnorm(10, 0, 2)

Załóżmy, że punkt przerwania k1 występuje o 9.6:

mod <- lm(y~x+I(pmax(x-9.6, 0)))
summary(mod)

Coefficients:
                    Estimate Std. Error t value Pr(>|t|)    
(Intercept)          21.7057     1.1726  18.511 1.06e-12 ***
x                    -1.1003     0.1788  -6.155 1.06e-05 ***
I(pmax(x - 9.6, 0))   1.3760     0.2688   5.120 8.54e-05 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1   1

Punkt przecięcia i nachylenie dwóch segmentów to: 21.7 i 1.1 po raz pierwszy i 8.5 i 0.27 odpowiednio dla drugiego.

Punkt przerwania


pytania

  1. Jak łatwo obliczyć przecięcie i nachylenie każdego segmentu? Czy model można ponownie sparametryzować, aby zrobić to w jednym obliczeniu?
  2. Jak obliczyć błąd standardowy każdego nachylenia każdego segmentu?
  3. Jak sprawdzić, czy dwa sąsiednie zbocza mają takie same zbocza (tj. Czy punkt przerwania można pominąć)?
COOLSerdash
źródło

Odpowiedzi:

7
  1. Jak łatwo obliczyć przecięcie i nachylenie każdego segmentu?

Nachylenie każdego segmentu oblicza się, dodając wszystkie współczynniki do aktualnej pozycji. Oszacuj nachylenie nax=15 jest 1.1003+1.3760=0.2757.

Punkt przecięcia jest nieco trudniejszy, ale jest to liniowa kombinacja współczynników (obejmująca węzły).

W twoim przykładzie druga linia spełnia pierwszą o x=9.6, więc czerwony punkt znajduje się w pierwszej linii na 21.70571.1003×9.6=11.1428. Ponieważ druga linia przechodzi przez punkt(9.6,11.428) ze spadkiem 0.2757, jego przechwycenie to 11.14280.2757×9.6=8.496. Oczywiście, możesz połączyć te kroki razem, co upraszcza aż do przecięcia dla drugiego segmentu =β0β2k1=21.70571.3760×9.6.

Czy model można sparametryzować, aby zrobić to w jednym obliczeniu?

Cóż, tak, ale prawdopodobnie łatwiej jest po prostu obliczyć to na podstawie modelu.

2. Jak obliczyć błąd standardowy każdego nachylenia każdego segmentu?

Ponieważ oszacowanie jest liniową kombinacją współczynników regresji aβ^, gdzie a składa się z 1 i 0, wariancja wynosi aVar(β^)a. Błąd standardowy jest pierwiastkiem kwadratowym z tej sumy wariancji i warunków kowariancji.

np. w twoim przykładzie standardowy błąd nachylenia drugiego segmentu to:

Sb <- vcov(mod)[2:3,2:3]
sqrt(sum(Sb))

alternatywnie w formie matrycy:

Sb <- vcov(mod)
a <- matrix(c(0,1,1),nr=3)
sqrt(t(a) %*% Sb %*% a)

3. Jak sprawdzić, czy dwa sąsiednie zbocza mają takie same zbocza (tj. Czy punkt przerwania można pominąć)?

Sprawdza się to, patrząc na współczynnik w tabeli tego segmentu. Zobacz tę linię:

I(pmax(x - 9.6, 0))   1.3760     0.2688   5.120 8.54e-05 ***

To zmiana nachylenia przy 9,6. Jeśli ta zmiana jest różna od 0, dwa zbocza nie są takie same. Zatem wartość p dla testu, że drugi segment ma takie samo nachylenie jak pierwszy, znajduje się na końcu tej linii.

Glen_b - Przywróć Monikę
źródło
(+1) Dzięki Glen za odpowiedź. Tylko jedno małe pytanie dotyczące nr 2: W moim przykładzie potrzebowałbym macierzy wariancji-kowariancji xi I(pmax(x-9.6,0))czy to prawda?
COOLSerdash
Nie. Zredagowałem wyraźny przykład na podstawie twojego przykładu. Jeśli chcesz uzyskać więcej informacji, zapytaj.
Glen_b
Bardzo dziękuję za edycję, która trochę mnie wyjaśnia. Czy więc rozumiem to poprawnie: błąd standardowy jest taki sam dla każdego zbocza?
COOLSerdash
1
Nie. Procedura jest taka sama, ale wartość nie jest. Standardowy błąd nachylenia pierwszego segmentu znajduje się w tabeli regresji (0,1788). Standardowy błąd nachylenia drugiego segmentu wynosi 0,1160. Gdybyśmy mieli trzeci segment, wymagałoby to większej liczby warunków wariancji-kowariancji w jego sumie (przed obliczeniem pierwiastka kwadratowego).
Glen_b
6

Moje naiwne podejście, które odpowiada na pytanie 1:

mod2 <- lm(y~I((x<9.6)*x)+as.numeric((x<9.6))+
             I((x>=9.6)*x)+as.numeric((x>=9.6))-1)
summary(mod2)

#                        Estimate Std. Error t value Pr(>|t|)    
# I((x < 9.6) * x)        -1.1040     0.2328  -4.743 0.000221 ***
# as.numeric((x < 9.6))   21.7188     1.3099  16.580 1.69e-11 ***
# I((x >= 9.6) * x)        0.2731     0.1560   1.751 0.099144 .  
# as.numeric((x >= 9.6))   8.5442     2.6790   3.189 0.005704 ** 

Ale nie jestem pewien, czy statystyki (w szczególności stopnie swobody) są wykonywane poprawnie, jeśli zrobisz to w ten sposób.

Roland
źródło
(+1) Dziękuję bardzo za odpowiedź. Zapewnia bardzo wygodny sposób bezpośredniego obliczania przechwyceń i nachyleń, dzięki!
COOLSerdash