Zmiana hipotezy zerowej w regresji liniowej

Mam pewne wysoce skorelowane dane. Jeśli przeprowadzę regresję liniową, otrzymam linię regresji o nachyleniu zbliżonym do jednego (= 0,93). Chciałbym przetestować, czy to nachylenie różni się znacznie od 1,0. Oczekuję, że tak nie jest. Innymi słowy, chciałbym zmienić hipotezę zerową regresji liniowej z nachylenia zera na nachylenie jednego. Czy to rozsądne podejście? Byłbym również bardzo wdzięczny za to, że w odpowiedzi mógłbyś podać kod R, aby móc wdrożyć tę metodę (lub lepszą, którą sugerujesz!). Dzięki.

set.seed(20); y = rnorm(20); x = y + rnorm(20, 0, 0.2) # generate correlated data
summary(lm(y ~ x))                  # original model
summary(lm(y ~ x, offset= 1.00*x))  # testing against slope=1
summary(lm(y-x ~ x))                # testing against slope=1

Wyjścia:

            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  0.01532    0.04728   0.324     0.75    
x            0.91424    0.04128  22.148 1.64e-14 ***

            Estimate Std. Error t value Pr(>|t|)  
(Intercept)  0.01532    0.04728   0.324   0.7497  
x           -0.08576    0.04128  -2.078   0.0523 .

            Estimate Std. Error t value Pr(>|t|)  
(Intercept)  0.01532    0.04728   0.324   0.7497  
x           -0.08576    0.04128  -2.078   0.0523 .

$R\beta=r$$\beta$$R$$r$

$\beta_1=0$$R=[0,1]$$r=1$

Do tego rodzaju hipotez można użyć linearHypothesisfunkcji z pakietu samochodu :

set.seed(20); y = rnorm(20); x = y + rnorm(20, 0, 0.2) # generate correlated data
mod <- lm(y ~ x))                  # original model


> linearHypothesis(mod,matrix(c(0,1),nrow=1),rhs=c(1))
Linear hypothesis test

Hypothesis:
x = 1

Model 1: restricted model
Model 2: y ~ x

  Res.Df     RSS Df Sum of Sq      F  Pr(>F)  
1     19 0.96022                              
2     18 0.77450  1   0.18572 4.3162 0.05234 .
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Wygląda na to, że nadal próbujesz odrzucić hipotezę zerową. Jest z tym mnóstwo problemów, z których nie mniej ważne jest to, że nie masz wystarczającej mocy, aby zobaczyć, że różnisz się od 1. Wygląda na to, że nie obchodzi cię, że nachylenie różni się od 0,07 od 1. Ale co, jeśli tak naprawdę nie możesz powiedzieć? Co jeśli faktycznie szacujesz nachylenie, które zmienia się dziko i może być całkiem dalekie od 1 z czymś w rodzaju przedziału ufności ± 0,4. Twoja najlepsza taktyka tutaj nie polega na zmianie hipotezy zerowej, ale właściwie na rozsądnym oszacowaniu przedziału. Jeśli zastosujesz polecenie confint () w swoim modelu, możesz uzyskać 95% przedział ufności wokół swojego nachylenia. Następnie możesz użyć tego do omówienia nachylenia. Jeśli 1 mieści się w przedziale ufności, możesz stwierdzić, że mieści się w zakresie wartości, które Twoim zdaniem mogą zawierać prawdziwą wartość. Ale co ważniejsze, możesz również określić, jaki jest ten zakres wartości.

Chodzi o to, że chcesz odrzucić hipotezę zerową, a nie ją potwierdzić. Fakt, że nie ma znaczącej różnicy, nie jest w żaden sposób dowodem braku znaczącej różnicy. W tym celu musisz określić, jaki rozmiar efektu uważasz za rozsądny, aby odrzucić wartość null.

$slope - 1$

set.seed(20); y = rnorm(20); x = y + rnorm(20, 0, 0.2)
model <- lm(y~x)

coefx <- coef(summary(model))[2,1]
seslope <- coef(summary(model))[2,2]
DF <- model$df.residual

# normal test
p <- (1 - pt(coefx/seslope,DF) )*2
# test whether different from 1
p2 <- (1 - pt(abs(coefx-1)/seslope,DF) )*2

Teraz powinieneś zdawać sobie sprawę z faktu, że wielkość efektu, dla którego różnica staje się znacząca, jest

> qt(0.975,DF)*seslope
[1] 0.08672358

pod warunkiem, że mamy przyzwoity estymator błędu standardowego na zboczu. Dlatego jeśli zdecydujesz, że istotną różnicę należy wykryć tylko od 0,1, możesz obliczyć niezbędny DF w następujący sposób:

optimize(
    function(x)abs(qt(0.975,x)*seslope - 0.1),
    interval=c(5,500)
) 
$minimum
[1] 6.2593

Pamiętaj, że to dość zależy od oszacowania seslope. Aby uzyskać lepsze oszacowanie na seslope, możesz dokonać ponownego próbkowania swoich danych. Naiwnym sposobem byłoby:

n <- length(y)
seslope2 <-
  mean(
    replicate(n,{
      id <- sample(seq.int(n),1)
      model <- lm(y[-id]~x[-id])
      coef(summary(model))[2,2]
    })
  )

umieszczenie seslope2 w funkcji optymalizacji zwraca:

$minimum
[1] 6.954609

Wszystko to powie Ci, że Twój zestaw danych zwróci znaczący wynik szybciej, niż uważasz za niezbędny, i że potrzebujesz tylko 7 stopni swobody (w tym przypadku 9 obserwacji), jeśli chcesz mieć pewność, że nieistotne oznacza to, czego chcesz. znaczy.