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

Sytuacja

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

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 $k_1$ 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.

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ąć)?

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 na$x=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.7057 -1.1003 \times 9.6 = 11.1428$. Ponieważ druga linia przechodzi przez punkt$(9.6, 11.428)$ ze spadkiem $0.2757$, jego przechwycenie to $11.1428 - 0.2757 \times 9.6 = 8.496$. Oczywiście, możesz połączyć te kroki razem, co upraszcza aż do przecięcia dla drugiego segmentu =$\beta_0 - \beta_2 k_1 = 21.7057 - 1.3760 \times 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^\top\hat\beta$, gdzie $a$ składa się z 1 i 0, wariancja wynosi $a^\top\text{Var}(\hat\beta)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.

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.