Oto mały przykład:
MyDf<-data.frame(x=c(1,2,3,4), y=c(1.2, .7, -.5, -3))
Teraz z base::lm
:
> lm(y~x, data=MyDf) %>% summary
Call:
lm(formula = y ~ x, data = MyDf)
Residuals:
1 2 3 4
-0.47 0.41 0.59 -0.53
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 3.0500 0.8738 3.491 0.0732 .
x -1.3800 0.3191 -4.325 0.0495 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.7134 on 2 degrees of freedom
Multiple R-squared: 0.9034, Adjusted R-squared: 0.8551
F-statistic: 18.71 on 1 and 2 DF, p-value: 0.04952
Teraz, spróbuj to samo z biglm
z biglm
pakietu:
XX<-biglm(y~x, data=MyDf)
print(summary(XX), digits=5)
Large data regression model: biglm(y ~ x, data = MyDf)
Sample size = 4
Coef (95% CI) SE p
(Intercept) 3.05 1.30243 4.79757 0.87378 0.00048
x -1.38 -2.01812 -0.74188 0.31906 0.00002
Zauważ, że potrzebujemy print
i, digits
aby zobaczyć wartość p. Współczynniki i błędy standardowe są takie same, ale wartości p są bardzo różne. Dlaczego tak jest?
r
regression
p-value
linear-model
Jan Paweł
źródło
źródło
pt(-3.491, 2)*2
dopnorm(-3.491)*2
, na przykład.Odpowiedzi:
Aby zobaczyć, które wartości p są poprawne (jeśli jedno z nich), powtórzmy obliczenia dla danych symulowanych, w których hipoteza zerowa jest prawdziwa. W obecnym ustawieniu obliczenia są dopasowane do danych (x, y) metodą najmniejszych kwadratów, a hipotezą zerową jest to, że nachylenie wynosi zero. W pytaniu są cztery wartości x 1,2,3,4, a szacowany błąd wynosi około 0,7, więc uwzględnijmy to w symulacji.
Oto konfiguracja napisana, aby była zrozumiała dla wszystkich, nawet tych, którzy się jej nie znają
R
.Symulacja generuje niezależne błędy, dodaje je
y.expected
, wywołujelm
w celu dopasowania isummary
obliczenia wartości p. Chociaż jest to nieefektywne, sprawdza rzeczywisty użyty kod. Nadal możemy wykonać tysiące iteracji w ciągu sekundy:Prawidłowo obliczone wartości p będą działać jak jednolite liczby losowe od do10 1 gdy hipoteza zerowa jest prawdziwa. Histogram tych wartości p pozwoli nam to sprawdzić wizualnie - czy wygląda mniej więcej poziomo - a test jednolitości chi-kwadrat pozwoli na bardziej formalną ocenę. Oto histogram:
a dla tych, którzy mogą sobie wyobrazić, że to nie jest wystarczająco jednolite, oto test chi-kwadrat:
Duża wartość p w tym teście pokazuje, że wyniki te są zgodne z oczekiwaną jednorodnością. Innymi słowy,
lm
jest poprawny.Skąd zatem biorą się różnice w wartościach p? Sprawdźmy prawdopodobne formuły, które można wywołać w celu obliczenia wartości p. W każdym razie statystyki testowe będą
równa rozbieżności między oszacowanym współczynnikiem a hipotetyczną (i poprawną wartością) , wyrażoną jako wielokrotność błędu standardowego oszacowania współczynnika. W pytaniu są to wartości beta=0β^ β=0
dla oszacowania przechwytywania i
do oszacowania nachylenia. Zwykle byłyby one porównywane z rozkładem Studenta którego parametr stopni swobody wynosi (ilość danych) minus (liczba oszacowanych współczynników). Obliczmy to dla przechwytywania:4 2t 4 2
(To obliczenie zwielokrotnia prawdopodobieństwo -Studenta po lewej stronie przez ponieważ jest to test stosunku do dwustronnej alternatywy ). Zgadza się z wynikiem.2t 2 H A : β ≠ 0H0:β=0 HA:β≠0
lm
Alternatywne obliczenia wykorzystują standardowy rozkład normalny do przybliżenia rozkładu Studenta . Zobaczmy, co produkuje:t
Rzeczywiście:t
biglm
zakłada, że rozkład zerowy statystyki jest standardowy Normalny. Ile to jest błędu? Ponowne uruchomienie poprzedniej symulacji zamiast zamiast tego daje histogram wartości p:biglm
lm
Prawie 18% tych wartości p jest mniejszych niż , co stanowi standardowy próg „istotności”. To ogromny błąd.0.05
Oto niektóre lekcje, których możemy się nauczyć z tego małego dochodzenia:
Nie używaj przybliżeń pochodzących z analiz asymptotycznych (takich jak standardowy rozkład normalny) z małymi zestawami danych.
Poznaj swoje oprogramowanie.
źródło