Przedziały ufności w R: t.test(), confint() i prop.test()

Co naprawdę znaczy "95% ufności"

Przedział ufności zamienia estymację punktową ("średnia z próby wynosi 5.61") w uczciwy zakres ("prawdziwa średnia prawdopodobnie leży między 5.29 a 5.93"). Ale wyrażenie 95% ufności to jedno z najczęściej źle rozumianych w statystyce, więc wyjaśnijmy je na początku.

95% opisuje procedurę, a nie przedział. Wyobraź sobie, że powtarzasz badanie w kółko: za każdym razem nowa próba i nowy przedział. Przedziały skakałyby w różne strony, a około 95% z nich obejmowałoby prawdziwą wartość; 5% by chybiło. Twój jeden rzeczywisty przedział to pojedyncze losowanie z tego procesu. Czego 95% nie oznacza: "jest 95% prawdopodobieństwa, że prawdziwa średnia leży między tymi liczbami". Prawdziwa średnia to stała (nieznana) liczba, która nie wędruje do przedziałów i z nich; zmieniają się przedziały.

W praktyce przyjazne odczytanie jest w porządku: przedział to zakres wartości zgodnych z twoimi danymi. Po prostu wiedz, jakie twierdzenie wolno ci postawić, gdy ktoś zacznie dopytywać.

Przedział ufności dla średniej, łatwy sposób

Każdy t.test() zawiera przedział ufności; możesz uruchomić test tylko po to, żeby go zebrać:

times <- c(5.1, 4.9, 6.2, 5.8, 5.5, 5.9, 6.1, 5.2, 5.4, 6.0)

result <- t.test(times)
print(result$conf.int)
print(result$estimate)

Średnia 5.61, 95% przedział ufności mniej więcej od 5.29 do 5.93. Przedział robi to, czego nie zrobi sama średnia: pokazuje, ile precyzji naprawdę dają dziesięć obserwacji.

Ten sam przedział ręcznie

Wzór stojący za tym przedziałem warto raz zbudować samemu, bo odczarowuje każdy przedział ufności, jaki kiedykolwiek przeczytasz: estymacja ± wartość krytyczna × błąd standardowy.

times <- c(5.1, 4.9, 6.2, 5.8, 5.5, 5.9, 6.1, 5.2, 5.4, 6.0)

n <- length(times)
m <- mean(times)
se <- sd(times) / sqrt(n)          # standard error of the mean
t_crit <- qt(0.975, df = n - 1)    # t critical value, 2.5% in each tail

margin <- t_crit * se
print(t_crit)
print(margin)
print(c(m - margin, m + margin))

Porównaj z wynikiem t.test(): ten sam przedział co do ostatniej cyfry. Trzy ruchome części:

• Błąd standardowy sd(x)/sqrt(n) mierzy, jak bardzo chwieje się średnia z próby (różnicę między sd a SE opisują statystyki opisowe).
• qt(0.975, df) to wartość krytyczna t: dla pokrycia 95% zostawiasz 2.5% w każdym ogonie, stąd 0.975. Przy df = 9 wynosi około 2.26, więcej niż 1.96 dla rozkładu normalnego. To podatek za małą próbę, bo sd szacujesz z tych samych danych.
• Margines błędu to ich iloczyn, czyli liczba "±" cytowana w nagłówkach.

Zmiana poziomu: 90%, 99% i kompromis

conf.level kontroluje pokrycie, a razem z nim szerokość:

times <- c(5.1, 4.9, 6.2, 5.8, 5.5, 5.9, 6.1, 5.2, 5.4, 6.0)

print(t.test(times, conf.level = 0.90)$conf.int)
print(t.test(times, conf.level = 0.95)$conf.int)
print(t.test(times, conf.level = 0.99)$conf.int)

Więcej ufności kosztuje więcej szerokości: 99% przedział musi być na tyle szeroki, żeby trafiać 99 razy na 100, więc się rozciąga; 90% przedział jest węższy, ale chybia dwa razy częściej niż 95%. Nie ma nic za darmo, jest tylko pokrętło. 95% to czysta konwencja: rozsądna wartość domyślna, jeśli nie masz powodu jej zmieniać, a nie prawo natury.

Przedział ufności dla proporcji: prop.test()

Załóżmy, że 47 ze 120 ankietowanych użytkowników zaczęło korzystać z nowej funkcji. Jaki jest wiarygodny zakres prawdziwego odsetka?

result <- prop.test(47, 120)
print(result$estimate)
print(result$conf.int)

Proporcja w próbie 0.39, 95% przedział ufności mniej więcej od 0.30 do 0.49, więc "około 40% użytkowników" jest uczciwe tylko z dołączoną aureolą ±9 punktów. prop.test() używa lepszego przybliżenia niż podręcznikowy wzór p ± 1.96 × sqrt(p(1−p)/n) (to przedział typu Wilsona z poprawką na ciągłość), co ma największe znaczenie blisko 0 lub 1: przedział podręcznikowy potrafi wyjść poza [0, 1], ten nie. Przy bardzo małych liczebnościach binom.test(47, 120)$conf.int daje wersję dokładną.

Przedziały ufności dla współczynników modelu: confint()

Dopasowane modele dostają przedziały przez jedną ogólną funkcję:

fit <- lm(mpg ~ wt, data = mtcars)
print(confint(fit))
print(confint(fit, level = 0.99))

Każdy wiersz obejmuje jeden współczynnik regresji liniowej: 95% przedział dla nachylenia wt biegnie od około −6.5 do −4.2 mpg na 1000 funtów. To więcej informacji niż p-value z podsumowania: efekt jest nie tylko niezerowy, ale wynosi co najmniej około 4 mpg, a może nawet około 6.5.

Dla logistycznego glm() funkcja confint(fit) działa tak samo, ale zwraca granice w log-szansach; zastosuj funkcję wykładniczą, żeby dostać przedziały dla ilorazów szans (exp(confint(fit))), i pamiętaj, że wartością odniesienia "brak efektu" staje się 1 zamiast 0.

Liczebność próby: prawo pierwiastka

Szerokość maleje z sqrt(n), co ma łatwą do zapamiętania konsekwencję: cztery razy więcej danych, dwa razy węższy przedział. Zobacz to na symulowanych danych: ten sam rozkład, jedna próba 4× większa (ziarno zapewnia powtarzalność; zobacz liczby losowe):

set.seed(42)
small <- rnorm(50, mean = 100, sd = 15)
large <- rnorm(200, mean = 100, sd = 15)   # 4 times the n

ci_small <- t.test(small)$conf.int
ci_large <- t.test(large)$conf.int

print(diff(ci_small))   # width with n = 50
print(diff(ci_large))   # width with n = 200 - about half

Szerokości wychodzą blisko proporcji 2:1 (szum losowania sprawia, że nie dokładnie). Przez pierwiastek precyzja drożeje: pierwsze 100 obserwacji zawęża przedział bardziej niż kolejne 300 razem, a zmniejszenie o połowę już wąskiego przedziału zawsze kosztuje 4× tyle, ile kosztowało wszystko do tej pory.

Przedział ufności i p-value: dwa widoki jednego testu

Przedział ufności i test hipotezy to ta sama informacja w innym ubraniu. 95% przedział ufności zawiera dokładnie te wartości parametru, których dwustronny test na poziomie 0.05 nie odrzuciłby. Zatem:

• Przedział dla różnicy średnich nie obejmuje 0 ⇔ test t daje p < 0.05.
• Przedział dla ilorazu szans nie obejmuje 1 ⇔ p współczynnika < 0.05.

Gdy masz przedział, zwykle masz lepsze podsumowanie: daje ten sam werdykt o istotności oraz wielkość efektu w prawdziwych jednostkach. "p = 0.03" mówi, że różnica istnieje; "95% CI: 0.2 do 7.6" mówi, że istnieje i może być błaha albo ogromna, a często to właśnie jest najważniejszy wniosek.

Co warto zapamiętać

• 95% opisuje długoterminową skuteczność procedury, a nie prawdopodobieństwo, że ten jeden przedział złapał prawdę.
• Średnia: t.test(x)$conf.int albo ręcznie średnia ± qt(0.975, n−1) × SE.
• Proporcja: prop.test(x, n)$conf.int; współczynniki modelu: confint(fit) (dla ilorazów szans z glm zastosuj exp).
• Wyższa ufność = szerszy przedział; 95% to konwencja, nie prawo.
• Szerokość maleje z sqrt(n): 4× więcej danych to dwa razy węższy przedział.
• Przedział zawiera każdą wartość, której test na poziomie 0.05 by nie odrzucił, i w przeciwieństwie do p-value pokazuje wielkość efektu.

Dalej: zestaw narzędzi do symulacji stojący za tym wszystkim, czyli liczby losowe z rnorm, runif, sample i set.seed.