Вычисление p-value в R

В задаче 4 из листочка к семинару 15 мы столкнулись с ситуацией, когда точно рассчитать p-value без помощи R у нас не получилось. По табличке мы смогли лишь оценить, что p-value (при двусторонней альтернативе) точно больше 0.2. Как посчитать вероятность вида \(P(T<t)\) для случайной величины \(T\), имеющей распределение Стьюдента с некоторым числом степеней свободы \(\text{df}\)? Для этого есть функция pt(). Например, посчитаем \(P(T<0)\) для распределения Стьюдента с \(\text{df}=10\):

pt(0, df = 10)
## [1] 0.5

Довольно логично, потому что плотность распределения Стьюдента с любым числом степеней свободы симметрична относительно 0. А теперь посчитаем \(P(T<2)\) для того же распределения:

pt(2, df = 10)
## [1] 0.963306

При желании всегда сможем вычесть результат из 1 и посчитать, например, \(P(T > 2)\):

1 - pt(2, df = 10)
## [1] 0.03669402

Теперь посчитаем \(P(T<1.11)\) для нашей задачи:

pt(1.11, df = 10)
## [1] 0.8535115

И, наконец, точное \(\text{p-value} = P(T>1.11) + P(T<-1.11)= 2 \cdot P(T>1.11)\):

(1 - pt(1.11, df = 10)) * 2
## [1] 0.2929771

Мы не ошиблись, когда на качественном уровне прикинули, что p-value точно больше 0.2.

Проверка гипотезы о равенстве доли числу

Теперь перейдем к проверке гипотез в R. Загрузим файл с данными, сохраним таблицу в переменную chipdale и удалим строки с пропущенными значениями:

chipdale <- read.csv("https://allatambov.github.io/psms/chip-n-dale.csv")
chipdale <- na.omit(chipdale)

Посмотрим на описательные статистики:

summary(chipdale)
##        X             group        female       spec        fav        assoc   
##  Min.   : 1.00   191    :16   Min.   :0.0000   PA:24   Chip  :15   Chip  :22  
##  1st Qu.:17.25   192    :17   1st Qu.:0.0000   PG:42   Dale  :10   Dale  :12  
##  Median :34.50   193    :10   Median :1.0000           Gadget:29   Gadget:16  
##  Mean   :34.23   194    :13   Mean   :0.5606           Monty : 6   Monty : 7  
##  3rd Qu.:50.75   195    : 8   3rd Qu.:1.0000           Other : 3   Other : 5  
##  Max.   :67.00   teacher: 2   Max.   :1.0000           Zipper: 3   Zipper: 4  
##      cheer            grump            cheeze            tech       
##  Min.   :  0.00   Min.   :  5.00   Min.   :  0.00   Min.   :  0.00  
##  1st Qu.: 42.50   1st Qu.: 50.00   1st Qu.: 50.00   1st Qu.: 25.00  
##  Median : 70.00   Median : 70.00   Median : 77.00   Median : 60.00  
##  Mean   : 61.67   Mean   : 68.47   Mean   : 70.59   Mean   : 51.21  
##  3rd Qu.: 80.00   3rd Qu.: 90.00   3rd Qu.:100.00   3rd Qu.: 72.25  
##  Max.   :100.00   Max.   :101.00   Max.   :100.00   Max.   :100.00  
##     reaction     
##  Min.   :  0.00  
##  1st Qu.: 50.00  
##  Median : 70.00  
##  Mean   : 60.54  
##  3rd Qu.: 78.53  
##  Max.   :100.00

Допустим, нас интересует доля студентов, которые хотели бы поступить на профиль «политический анализ», и мы хотим на 5%-ном уровне значимости проверить следующую нулевую гипотезу, выдвинув двустороннюю альтернативу:

\[ H_0: p = 0.25 \] \[ H_1: p \ne 0.25 \] Для проверки гипотезы о доле нам нужна функция prop.test(). На первом месте в этой функции указывается число успехов (24 человека хотят на «политический анализ»), на втором месте — общее число опрошенных (66 человек), а затем формулируется нулевая гипотеза.

prop.test(24, 66, p = 0.25)
## 
##  1-sample proportions test with continuity correction
## 
## data:  24 out of 66, null probability 0.25
## X-squared = 3.9596, df = 1, p-value = 0.0466
## alternative hypothesis: true p is not equal to 0.25
## 95 percent confidence interval:
##  0.2513999 0.4918234
## sample estimates:
##         p 
## 0.3636364

Критерий, который использует R, немного отличается от того, что мы использовали при расчетах руками, но это сейчас не главное. Что есть в выдаче? Во-первых, название теста. По первой строчке видно, что это критерий для доли и для одной выборки. Во-вторых, альтернативная гипотеза. Здесь true p is not equal to 0.25, и это неслучайно: R по умолчанию использует двустороннюю альтернативу. В-третьих, sample estimates, оценка доли \(\hat{p}\), полученная по выборке, равная 0.36. Наконец, само p-value, равное 0.0466.

Интерпретация выдачи. На 5% уровне значимости на имеющихся данных есть основания отвергнуть нулевую гипотезу, так как p-value меньше \(\alpha=5\)%. Доля студентов, которые планируют поступать на «политический анализ», не равна 0.25.

Если мы хотим сформулировать одностороннюю альтернативу, R позволит нам это сделать, в prop.test() можно указать аргумент alternative:

Представим теперь, что при проверке той же нулевой гипотезы мы выбрали правосторонню альтернативу:

prop.test(24, 66, p = 0.25, alternative = "greater")
## 
##  1-sample proportions test with continuity correction
## 
## data:  24 out of 66, null probability 0.25
## X-squared = 3.9596, df = 1, p-value = 0.0233
## alternative hypothesis: true p is greater than 0.25
## 95 percent confidence interval:
##  0.2665397 1.0000000
## sample estimates:
##         p 
## 0.3636364

Интерпретация у нас глобально не изменилась, нулевую гипотезу стоит отвергнуть.

Проверка гипотезы о равенстве среднего числу

Давайте теперь проверим гипотезу о равенстве среднего числу. Проверим гипотезу о равенстве средней степени любви к сыру 70:

\[ H_0: \mu = 70 \]

\[ H_1: \mu \ne 70 \]

Для проверки этой гипотезы нам потребуется критерий Стьюдента и функция t.test(), в которой мы укажем показатель cheeze, по которому мы считаем среднее, и саму гипотезу:

t.test(chipdale$cheeze, mu = 70)
## 
##  One Sample t-test
## 
## data:  chipdale$cheeze
## t = 0.16314, df = 65, p-value = 0.8709
## alternative hypothesis: true mean is not equal to 70
## 95 percent confidence interval:
##  63.35715 77.82467
## sample estimates:
## mean of x 
##  70.59091

Итак, видим, что это критерий Стьюдента для одной выборки. Наблюдаемое значение статистики критерия \(t_{набл} = 0.1631\). Это значение взято из распределения Стьюдента с числом степеней свободы \(\text{df} = n - 1 = 65\). Выборочное среднее \(\bar{x}\) равно 70.59. Значение p-value равно 0.8709.

Интерпретация выдачи. Так как p-value больше любого конвенционального уровня значимости (1%, 5%, 10%), на имеющихся данных у нас нет оснований отвергнуть нулевую гипотезу. Среднее значение любви к сыру можно считать равным 70.

Проверка гипотезы о равенстве двух средних

Теперь проверим гипотезу о равенстве средних в двух группах. Правда ли, что средние значения показателя любви к технике одинаковы для тех, кто хочет поступать на «политический анализ» и «политическое управление»?

\[ H_0: \mu_{PA} = \mu_{PG} \] \[ H_1: \mu_{PA} \ne \mu_{PG} \]

Воспользуемся критерием Стьюдента для двух выборок и в функции t.test() через ~ укажем переменную группировки — специализацию:

t.test(chipdale$tech ~ chipdale$spec)
## 
##  Welch Two Sample t-test
## 
## data:  chipdale$tech by chipdale$spec
## t = 1.2873, df = 47.498, p-value = 0.2042
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  -5.813436 26.492008
## sample estimates:
## mean in group PA mean in group PG 
##         57.79167         47.45238

Интерпретация выдачи. Несмотря на то, что разница между выборочными средними есть и примерно равна 10 (значения 57.79 и 47.45 в sample estimates), эта разница не является статистически значимой, то есть такой, чтобы нулевую гипотезу о равенстве средних можно было отвергнуть. На 5%-ном уровне значимости нет оснований отвергнуть нулевую гипотезу. Средний уровень любви к технике у студентов, собирающихся поступать на «политический анализ» и «политическое управление» можно считать одинаковым.