Mam przykładowy zestaw danych zawierający 31 wartości. Przeprowadziłem dwustronny test t, używając R, aby sprawdzić, czy prawdziwa średnia jest równa 10:
t.test(x=data, mu=10, conf.level=0.95) Wynik:
t = 11.244, df = 30, p-value = 2.786e-12 alternative hypothesis: true mean is not equal to 10 95 percent confidence interval: 19.18980 23.26907 sample estimates: mean of x 21.22944 Teraz próbuję zrobić to samo ręcznie:
t.value = (mean(data) - 10) / (sd(data) / sqrt(length(data))) p.value = dt(t.value, df=length(lengths-1)) Wartość t obliczona na podstawie tego jest taki sam jak wynik funkcji testu t. Jednak wartość p okazuje się wynosić 3,025803e-12.
Jakieś pomysły, co robię źle?
Dzięki!
EDYTUJ
Oto pełny kod R , w tym mój zbiór danych:
# Raw dataset -- 32 observations data = c(21.75, 18.0875, 18.75, 23.5, 14.125, 16.75, 11.125, 11.125, 14.875, 15.5, 20.875, 17.125, 19.075, 25.125, 27.75, 29.825, 17.825, 28.375, 22.625, 28.75, 27, 12.825, 26, 32.825, 25.375, 24.825, 25.825, 15.625, 26.825, 24.625, 26.625, 19.625) # Student t-Test t.test(x=data, mu=10, conf.level=0.95) # Manually calculate p-value t.value = (mean(data) - 10) / (sd(data) / sqrt(length(data))) p.value = dt(t.value, df=length(data) - 1) Odpowiedź
Użyj pt i uczyń go dwustronnym.
> 2*pt(11.244, 30, lower=FALSE) [1] 2.785806e-12 Komentarze
- I myślę, że brakuje jakiegoś szczegółu: Kiedy używać lower = F? – Zobacz moje pytanie tutaj: stats.stackexchange.com/questions/133091/…
- Wartość musi być dodatnia, więc jeśli używasz jej ze zmienną, która może być ujemna, zawiń
abs. - W przypadku testu dwustronnego ' ponownie poszukaj prawdopodobieństwa, że wartość jest mniejsza niż -11,244 lub większa niż +11,244. lower = F mówi R, aby obliczył prawdopodobieństwo, że wartość będzie większa niż pierwszy parametr. W przeciwnym razie daje prawdopodobieństwo, że wartość będzie mniejsza niż pierwszy parametr. W związku z tym możesz również zrobić 2 * pt (-11,244, 30). Osobiście zwykle robię 2 * pt (-abs (q), df = n-1), ponieważ R domyślnie jest niższe = T.
Odpowiedź
Opublikowałem to jako komentarz, ale kiedy chciałem dodać trochę więcej w edycji, stało się zbyt długie, więc przeniosłem to tutaj.
Edytuj : Twoje statystyki testowe i df są poprawne. Druga odpowiedź zwraca uwagę na problem z obliczeniem obszaru ogona w wywołaniu do pt() i podwojenie dla dwóch ogonów, które rozwiązuje twoją różnicę. Niemniej jednak zostawię moją wcześniejszą dyskusję / komentarz, ponieważ zawiera on bardziej ogólnie istotne punkty na temat wartości p w skrajnych ogonach:
Jest możliwe, że nie robisz nic złego i nadal uzyskasz różnicę, ale jeśli opublikujesz odtwarzalny przykład, może być możliwe dalsze zbadanie, czy masz jakiś błąd (powiedzmy w df).
Te rzeczy są obliczane na podstawie przybliżeń, które mogą nie być szczególnie dokładne w bardzo skrajnym ogonie .
Jeśli te dwie rzeczy nie używają identycznych przybliżeń, mogą nie zgadzać się ściśle, ale ten brak zgodności nie powinien mieć znaczenia (aby dokładny obszar ogonowy był tak dalece znaczącą liczbą, wymagane założenia trzeba trzymać się zdumiewających stopni dokładności). Czy naprawdę masz dokładną normalność, dokładną niezależność, dokładnie stałą wariancję?
Nie należy oczekiwać dużej dokładności, jeśli liczby i tak nic nie znaczą. W jakim stopniu ma znaczenie, czy obliczona przybliżona wartość p wynosi $ 2 \ times 10 ^ {- 12} $ czy $ 3 \ times 10 ^ {- 12} $? Żadna liczba nie jest miarą rzeczywistej wartości p twojej prawdziwej sytuacji. Nawet jeśli jedna z liczb reprezentuje prawdziwą wartość p twojej prawdziwej sytuacji, skoro spadła poniżej około 0,0001 $, dlaczego miałbyś się przejmować, jaka była ta wartość?
Odpowiedź
Najlepszym sposobem ręcznego obliczenia jest:
t.value = (mean(data) - 10) / (sd(data) / sqrt(length(data))) p.value = 2*pt(-abs(t.value), df=length(data)-1) Potrzebujesz abs (), ponieważ w przeciwnym razie ryzykujesz uzyskanie wartości p większych niż 1 $ (gdy średnia danych jest większa niż podana średnia)!
Odpowiedź
Bardzo podoba mi się odpowiedź udzielona przez @Aaron wraz z komentarzami abs. Uważam, że przydatnym potwierdzeniem jest uruchomienie
pt(1.96, 1000000, lower.tail = F) * 2
, co daje 0.04999607.
Tutaj „używamy dobrze znanej własności, że 95% obszaru pod rozkładem normalnym występuje przy ~ 1,96 odchylenia standardowego, więc wynik ~ 0,05 daje naszą wartość p. Użyłem 1000000 ponieważ kiedy N jest duże, rozkład t jest prawie taki sam jak rozkład normalny. Uruchomienie tego zapewniło mi komfort w rozwiązaniu @Aarona.