Uzyskanie wzoru na granice predykcji w modelu liniowym (tj .: przedziały prognozowania)

18

Weźmy następujący przykład:

set.seed(342)
x1 <- runif(100)
x2 <- runif(100)
y <- x1+x2 + 2*x1*x2 + rnorm(100)
fit <- lm(y~x1*x2)

Tworzy to model y oparty na x1 i x2 przy użyciu regresji OLS. Jeśli chcemy przewidzieć y dla danego x_vec, moglibyśmy po prostu użyć wzoru, który otrzymujemy z summary(fit).

Co jednak, jeśli chcemy przewidzieć dolną i górną prognozę y? (dla danego poziomu ufności).

Jak więc zbudowalibyśmy formułę?

Tal Galili
źródło
Przedział ufności na nowe obserwacje części tej strony mogą pomóc.
GaBorgulya,
@ Tal Przepraszam, ale tak naprawdę nie jest dla mnie jasne, co tak naprawdę rozumiesz przez „przewidywanie dolnej i górnej prognozy y”. Czy ma to coś wspólnego z przedziałami predykcji lub tolerancji?
chl
@Tal - kilka zapytań. Kiedy mówisz „.. y na podstawie x1 i x2, używając regresji OLS”. , masz na myśli stworzenie modelu liniowego i oszacowanie parametrów za pomocą OLS . Czy mam rację? i pytanie @ chl - czy chcesz przewidzieć dolną i górną granicę przedziału prognoz?
suncoolsu,
@chl, przepraszam, że nie było bardziej jasne. Szukam dwóch formuł, które podadzą przedział, w którym „złapie” „rzeczywistą” wartość y 95% czasu. Mam na myśli to, jak używam definicji dla CI, gdy jest prawdopodobnie inny termin, którego powinienem użyć, przepraszam za to ...
Tal Galili
@suncoolsu - tak i tak.
Tal Galili

Odpowiedzi:

25

Będziesz potrzebował arytmetyki macierzowej. Nie jestem pewien, jak sobie z tym poradzi Excel. W każdym razie oto szczegóły.

Załóżmy, że regresja jest zapisana jako y=Xβ+e .

Niech będzie wektorem wiersza zawierającym wartości predyktorów dla prognoz (w tym samym formacie co X ). Wtedy prognoza jest przez r = X * β = X * ( X ' X ) - 1 X " Y w skojarzony wariancji σ 2 [ 1 + X * ( X ' X ) - 1 ( X * ) " ] .XX

y^=Xβ^=X(XX)1XY
σ2[1+X(XX)1(X)].
Następnie 95% przedział predykcji może być obliczone (zakładając, że rozkład normalny błąd), a Y ± 1,96 Ď Uwzględnia to niepewność związaną ze składnikiem błędueoraz niepewność w oszacowaniach współczynników. Jednak ignoruje wszelkie błędy w X. Jeśli więc przyszłe wartości predyktorów nie są pewne, wówczas przedział przewidywania obliczony przy użyciu tego wyrażenia będzie zbyt wąski.
y^±1.96σ^1+X(XX)1(X).
eX
Rob Hyndman
źródło
1
+1, doskonała odpowiedź. Powinienem jednak zauważyć, że model regresji zawsze szacuje oczekiwanie warunkowe, więc jest tak dobry, jak jego regresory. Ostatni komentarz, choć jest bardzo dobry, nie jest absolutnie konieczny, ponieważ jeśli zbudujesz model regresji, musisz zaufać regresorom.
mpiktas
y^=Xβ+X(XX)1Xevary^=varX(XX)1Xe=σ2X(XX)1(X) ?
mpiktas,
y^
N×N
X
7

Czy jesteś przypadkiem po różnych typach przedziałów prognoz? Strona predict.lmpodręcznika ma

 ## S3 method for class 'lm'
 predict(object, newdata, se.fit = FALSE, scale = NULL, df = Inf, 
         interval = c("none", "confidence", "prediction"),
         level = 0.95, type = c("response", "terms"),
         terms = NULL, na.action = na.pass,
         pred.var = res.var/weights, weights = 1, ...)

i

Ustawienie „przedziałów” określa obliczanie przedziałów ufności lub przewidywania (tolerancji) na określonym „poziomie”, czasem określanym jako przedziały wąskie vs. szerokie.

Czy o to ci chodziło?

Dirk Eddelbuettel
źródło
Cześć Dirk, to jest rzeczywiście to, co chciałbym znaleźć, ale chcę, aby górne i dolne wiązania miały postać formuły (aby później wdrożyć w jakiejś niskiej formie oprogramowania statystycznego, na przykład excel ...)
Tal Galili,
ps: Widzę teraz, że w tytule mojego pytania pojawiła się zmiana, która mogła sprawić, że pomyślałem, że pytam o parametr przedziału predykcji.lm (którego nie jestem) :)
Tal Galili
8
Nadużywacie tutaj terminologii. Excel nie jest oprogramowaniem statystycznym.
Dirk Eddelbuettel
1
Masz rację, moja oferta, a może „aplikacja do arkuszy kalkulacyjnych”?
Tal Galili,
3
Moge z tym zyc; woła diabła po imieniu ;-)
Dirk Eddelbuettel
6

@Tal: Mogę zasugerować Kutnera i in. Jako wspaniałe źródło modeli liniowych.

Istnieje różnica między 1) prognozą Y z indywidualnej nowej obserwacji X_vec, 2) oczekiwaną wartością Y uwarunkowaną na X_vec, mi(Y|Xvmido) oraz 3) Y z kilku wystąpień x_vec - wszystkie szczegółowo omówione w tekście.

Myślę, że szukasz wzoru na przedział ufności mi(Y|Xvmido) and that is Y^ ± t(1-α /2)s{Y^} where t has n-2 d.f. and s{Y^} is the standard error of Y^, σ2n+(XvecX¯)2σ2(XiX¯)2

B_Miner
źródło
1
(+1) for making the distinction. However, I believe the OP is asking for (1), not (2) (and I have edited the question's title accordingly). Also note that your formula appears to assume the regression depends only on one variable.
whuber