Введение в ANOVA (в сочетании со случаями COVID-19)

машинное обучение

Автор|ГОСТ БЛОГ Компилировать|ВКонтакте Источник | Аналитика Видья

вводить

"Факт — это простое утверждение, в которое все верят. То есть факт не является неверным, пока не будет доказано, что он ошибочен. Предположим, что существует предположение, в которое никто не хочет верить, и оно не может стать фактом, пока не будет установлено, что оно неверно. быть действительным».- Эдвард Тейлор

Мы имеем дело с пандемией беспрецедентного масштаба. Исследователи по всему миру лихорадочно пытаются разработать вакцину или лечение от COVID-19, а врачи пытаются остановить распространение пандемии по всему миру.

Недавно у меня появилась идея применить свои статистические знания к этому огромному количеству данных о COVID.

Рассмотрим сценарий, в котором у врача есть четыре медицинских подхода к лечению пациента. Как только мы получим результаты анализов, лучшим подходом будет лечение, которое займет наименьшее количество времени, чтобы вылечить пациента.

Но что, если некоторые из этих пациентов уже частично вылечены или их уже лечат другие препараты?

Чтобы принять уверенное и надежное решение, нам нужны доказательства, подтверждающие нашу практику. Здесь в игру вступает концепция ANOVA.

В этой статье я познакомлю вас с тестами ANOVA и различными типами, которые они используют для принятия лучших решений. Я продемонстрирую каждый тип теста ANOVA (дисперсионный анализ) в Python, чтобы визуализировать их и обработать данные о COVID-19.

Примечание. Для понимания этой темы необходимо знать основы статистики. Лучше всего знать о t-тестах и ​​проверке гипотез.

Что такое тест дисперсионного анализа (ANOVA)

ANOVA или ANOVA можно рассматривать как расширение t-критерия для более чем двух групп. Независимые t-тесты использовались для сравнения условных средних между двумя группами. ANOVA используется, когда мы хотим сравнить среднее значение состояния более чем двух групп пациентов.

ANOVA проверяет, есть ли разница в средних значениях где-то в модели (проверяет, есть ли общий эффект), но не может сказать нам, где разница (если она есть). Чтобы найти разницу между двумя группами, мы должны провести апостериорный тест.

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

  • Нулевая гипотеза – нет существенной разницы между группами
  • Альтернативная гипотеза — существенные различия между группами

По сути, дисперсионный анализ выполняется путем сравнения двух типов вариации: вариации между выборочными средними и вариации внутри каждой выборки. Следующие формулы представляют статистику одностороннего теста Anova.

Результат формулы ANOVA, статистика F (также известная как отношение F), позволяет анализировать несколько наборов данных для определения изменчивости между образцами и внутри них.

Формула однофакторного дисперсионного анализа может быть записана как:

Когда мы рисуем таблицу ANOVA, все вышеперечисленные компоненты могут выглядеть так:

В общем, если значение p, связанное с F, меньше 0,05, нулевая гипотеза будет отвергнута, а альтернативная гипотеза будет поддержана. Если нулевая гипотеза отвергнута, мы можем заключить, что средние значения всех групп не равны.

Примечание. Если между тестируемыми группами нет реальной разницы, так называемой нулевой гипотезы, то статистика F-отношения для ANOVA будет близка к 1.

Гипотезы для тестов ANOVA

Прежде чем делать ANOVA, нам нужно сделать некоторые предположения:

  1. Получайте наблюдения независимо и случайным образом из популяции, определяемой уровнями факторов
  2. Данные для каждого уровня факторов нормально распределены.
  3. Независимость от случая: примеры случаев должны быть независимы друг от друга
  4. Однородность дисперсии: Однородность означает, что дисперсии между группами должны быть примерно равными.

Предположение об однородности дисперсий можно проверить с помощью теста Левена или теста Брауна-Форсайта. Нормальность дробных распределений можно проверить с помощью гистограмм, значений асимметрии и эксцесса, а также графиков Шапиро-Уилка, Колмогорова-Смирнова или Q-Q. Предположение о независимости может быть определено на основе плана исследования.

Стоит отметить, что дисперсионный анализ неустойчив к нарушениям предположения о независимости. То есть, даже если вы нарушаете предположение об однородности или нормальности, вы можете запустить тест и в основном доверять результатам.

Однако, если предположение о независимости нарушается, результаты дисперсионного анализа недействительны. В целом анализ считается достоверным, если имеются группы одинакового размера в случае нарушения однородности. Для нарушений нормальности, если размер выборки большой, обычно можно продолжить дисперсионный анализ.

Типы тестов ANOVA

  1. Односторонний дисперсионный анализ: Однофакторный дисперсионный анализ только с одной независимой переменной.

    • Например, различия в случаях коронавируса можно оценить по странам, и страна может сравнить 2, 20 или более различных категорий.
  2. Двусторонний дисперсионный анализ: Двухфакторный дисперсионный анализ (также называемый факторным дисперсионным анализом) относится к дисперсионному анализу с использованием двух независимых переменных.

    • Расширяя приведенный выше пример, двусторонний ANOVA может исследовать различия в случаях коронавируса (зависимая переменная) по возрастной группе (независимая переменная 1) и полу (независимая переменная 2). Двухфакторный дисперсионный анализ можно использовать для изучения взаимодействия между двумя независимыми переменными. Взаимодействие показывает, что различия между всеми категориями независимой переменной неравномерны.
    • Например, в группе пожилых людей в целом может быть больше коронарных случаев, чем в группе молодежи, но разница может быть больше (или меньше) в азиатских странах по сравнению с европейскими странами.
  3. N-факторный дисперсионный анализ: Исследователь также может использовать более двух независимых переменных, это N-факторный дисперсионный анализ (N — количество независимых переменных, которые у вас есть), также известный как тест MANOVA.

    • Например, потенциальные различия в случаях заражения коронавирусом могут быть изучены одновременно по стране, полу, возрастной группе, этнической принадлежности и т. д.
    • ANOVA даст вам одномерное f-значение, тогда как ANOVA даст вам многомерное f-значение

С копированием или без

Вы можете часто слышать об ANOVA с повторением и без него. Давайте разберемся, что это такое:

  1. Двусторонний дисперсионный анализ с повторением: Две группы и члены этих групп делают более одной вещи

    • Например, если предположить, что вакцина против COVID-19 еще не разработана, врачи пробуют два разных метода лечения, чтобы вылечить две группы пациентов с COVID-19.
  2. Двусторонний ANOVA (без повторения):Используется, когда имеется только одна группа и для одной и той же группы выполняется двойной тест.

    • Например, предположим, что вакцина против COVID-19 была разработана, и исследователи тестируют группу добровольцев до и после вакцинации, чтобы увидеть, работает ли она.

посмертный

Когда мы проводим дисперсионный анализ, мы пытаемся определить, есть ли статистически значимая разница между группами. Если мы находим разницу, нам нужно исследовать расположение групповой разницы.

По сути, апостериорные тесты сообщают исследователю, какие группы отличаются друг от друга.

На этом этапе вы можете запустить вскрытие, т.е.т-тест,Используется для проверки различий в средних значениях между группами. Для контроля частоты ошибок типа I можно провести несколько сравнительных тестов, включая тесты Бонферрони, Шеффе, Даннета и Тьюки.

Теперь давайте разберемся с каждым типом теста ANOVA с некоторыми реальными данными, используя Python.

Односторонний тест ANOVA в Python

Я загрузил эти данные из текущего соревнования Kaggle:Вместе.

В этом тесте мы попытаемся проанализировать взаимосвязь между плотностью региона или штата и количеством корональных случаев. Поэтому мы составим карту каждого штата в зависимости от плотности его населения.

Сначала импортируйте все необходимые библиотеки и данные:

import pandas as pd
import numpy as np

import scipy.stats as stats
import os
import random

import statsmodels.api as sm
import statsmodels.stats.multicomp

from statsmodels.formula.api import ols
from statsmodels.stats.anova import anova_lm
import matplotlib.pyplot as plt
from scipy import stats
import seaborn as sns

Загрузить данные из каталога:

StatewiseTestingDetails=pd.read_csv('./StatewiseTestingDetails.csv')
population_india_census2011=pd.read_csv('./population_india_census2011.csv')

StatewiseTestingDetails содержит информацию об общем количестве положительных и отрицательных случаев в каждом штате за день. Принимая во внимание, что human_india_census2011 содержит информацию о плотности каждого штата и другую соответствующую информацию о населении.

population_india_census2011.head() 
StatewiseTestingDetails.head() #了解数据
StatewiseTestingDetails['Positive'].sort_values().head() #排序

Из приведенного выше фрагмента кода видно, что в нескольких штатах за день зарегистрировано 0 или совсем нет случаев заражения коронавирусом. Итак, давайте посмотрим на такие состояния:

StatewiseTestingDetails['State'][StatewiseTestingDetails['Positive']==1].unique()

Мы видели, что в Нагаленде и Сиккиме также не было случаев короны за один день. С другой стороны, в Аруначале и Мизораме происходит только один случай короны в день.

Вменить пропущенные значения: Замечаем, что в столбце «Положительные» много пропущенных значений. Итак, давайте вменим этим пропущенным значениям положительную медиану для каждого состояния:

stateMedianData=StatewiseTestingDetails.groupby('State')[['Positive']].median().reset_index().rename(columns={'Positive':'Median'})
stateMedianData.head()

for index,row in StatewiseTestingDetails.iterrows():

    if pd.isnull(row['Positive']):

        StatewiseTestingDetails['Positive'][index]=int(stateMedianData['Median'][stateMedianData['State']==row['State']])

#合并StatewiseTestingDetails & population_india_census2011
data=pd.merge(StatewiseTestingDetails,population_india_census2011,on='State')

Теперь мы можем написать функцию, которая создает корзину группы плотности на основе плотности каждого состояния, где Dense1

def densityCheck(data):
    data['density_Group']=0
    for index,row in data.iterrows():
        status=None
        i=row['Density'].split('/')[0]
        try:
            if (',' in i):
                i=int(i.split(',')[0]+i.split(',')[1])
            elif ('.' in i):
                i=round(float(i))
            else:
                i=int(i)
        except ValueError as err:
            pass
        try:
            if (0<i<=300):
                status='Dense1'
            elif (300<i<=600):
                status='Dense2'
            elif (600<i<=900):
                status='Dense3'
            else:
                status='Dense4'
        except ValueError as err:
            pass
        data['density_Group'].iloc[index]=status
    return data

Теперь сопоставьте каждое состояние с группой плотности. Мы можем экспортировать эти данные для последующего использования в двухстороннем тесте ANOVA:

data=densityCheck(data)
#导出
stateDensity=data[['State','density_Group']].drop_duplicates().sort_values(by='State')

Давайте изменим набор данных, который можно использовать для тестирования ANOVA:

df=pd.DataFrame({'Dense1':data[data['density_Group']=='Dense1']['Positive'],
                 'Dense2':data[data['density_Group']=='Dense2']['Positive'],
                 'Dense3':data[data['density_Group']=='Dense3']['Positive'],
                 'Dense4':data[data['density_Group']=='Dense4']['Positive']})

Одно из предположений нашего теста ANOVA состоит в том, что выборка должна быть выбрана случайным образом и что выборка должна быть близка к распределению Гаусса. Итак, давайте выберем 10 случайных выборок из каждого фактора или уровня:

np.random.seed(1234)
dataNew=pd.DataFrame({'Dense1':random.sample(list(data['Positive'][data['density_Group']=='Dense1']), 10),
'Dense2':random.sample(list(data['Positive'][data['density_Group']=='Dense1']), 10),
'Dense3':random.sample(list(data['Positive'][data['density_Group']=='Dense1']), 10),
'Dense4':random.sample(list(data['Positive'][data['density_Group']=='Dense1']), 10)})

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

Мы можем ясно видеть, что данные не следуют распределению Гаусса.

Существуют различные методы преобразования данных, позволяющие приблизить данные к распределению Гаусса. Выполняем преобразование Бокса-Кокса:

dataNew['Dense1'],fitted_lambda = stats.boxcox(dataNew['Dense1'])
dataNew['Dense2'],fitted_lambda = stats.boxcox(dataNew['Dense2'])
dataNew['Dense3'],fitted_lambda = stats.boxcox(dataNew['Dense3'])
dataNew['Dense4'],fitted_lambda = stats.boxcox(dataNew['Dense4'])

Теперь давайте снова построим их распределение, чтобы проверить:

Метод 1: Односторонний ANOVA с использованием модуля statsmodels

В Python есть два способа выполнить тест ANOVA. Одним из них является метод stats.f_oneway():

F, p = stats.f_oneway(dataNew['Dense1'],dataNew['Dense2'],dataNew['Dense3'],dataNew['Dense4'])
## 看看整体模型是否重要
print('F-Statistic=%.3f, p=%.3f' % (F, p))

Мы нашли p-значение

Метод 2: Односторонний ANOVA с моделью OLS

Как мы знаем из регрессии, мы можем регрессировать каждую входную переменную и исследовать ее влияние на целевую переменную. Итак, мы собираемся следовать тому же подходу, что и в линейной регрессии.

model = ols('Count ~ C(density_Group)', newDf).fit()
model.summary()

## 看看整体模型是否重要
print(f"Overall model F({model.df_model: .0f},{model.df_resid: .0f}) = {model.fvalue: .3f}, p = {model.f_pvalue: .4f}")

#创建方差分析表
res = sm.stats.anova_lm(model, typ= 2)
res

Как видно из приведенного выше вывода, значение p меньше 0,05. Следовательно, мы можем отвергнуть нулевую гипотезу об отсутствии различий между различными группами плотности.

F-статистика = 5,817, p-значение = 0,002, что указывает на то, чтоdensity_Group оказывает общее значительное влияние на случаи с положительным результатом на коронавирус. Однако мы пока не знаем, в чем разница между desnity_groups. Следовательно, исходя из p-значения, мы можем отклонить H0, нет существенных различий с точки зрения плотности площади и количества корональных случаев.

апостериорный сравнительный тест

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

Если обнаружены различия, необходимо изучить расположение групповых различий. Поэтому мы будем использоватьТест Тьюки HSDопределить разницу:

mc = statsmodels.stats.multicomp.MultiComparison(newDf['Count'],newDf['density_Group'])
mc_results = mc.tukeyhsd()
print(mc_results)

Тест Tuckey HSD ясно показывает, что существуют значительные различия между группами Group1 – Group3, Group1 – Group4, Group2 – Group3 и Group3 – Group4.

Это показывает, что все попарные сравнения числа коронарных случаев, за исключением двух вышеупомянутых групп, отвергают нулевую гипотезу без статистически значимой разницы.

Проверка гипотез/диагностика модели

Проверка гипотезы нормального распределения

При использовании моделей линейной регрессии и ANOVA допущения относятся к остаткам, а не к самим переменным.

Метод 1: тест Шапиро-Уилка:
### 正态性假设检查
w, pvalue = stats.shapiro(model.resid)
print(w, pvalue)

Из фрагмента кода выше мы видим, что все группы плотности имеют p-значения больше 0,05. Следовательно, мы можем сделать вывод, что они следуют распределению Гаусса.

Метод 2: тест графика Q-Q:

Мы можем использовать график Q-Q, чтобы проверить эту гипотезу:

res = model.resid
fig = sm.qqplot(res, line='s')
plt.show()

На приведенном выше графике мы видим, что все точки данных близки к линии 45 градусов, поэтому мы можем сделать вывод, что он следует нормальному распределению.

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

Однородность предположений о дисперсии следует проверять для каждого уровня категориальных переменных. Мы можем использовать тест Левена для проверки равной дисперсии между группами.

w, pvalue = stats.bartlett(newDf['Count'][newDf['density_Group']=='Dense1'], newDf['Count'][newDf['density_Group']=='Dense2']
                           , newDf['Count'][newDf['density_Group']=='Dense3'], newDf['Count'][newDf['density_Group']=='Dense4'])
print(w, pvalue)

## Levene 方差测试
stats.levene(dataNew['Dense1'],dataNew['Dense2'],dataNew['Dense3'],dataNew['Dense4'])

Мы нашли значения p более 0,05 для всех групп плотности. Таким образом, можно сделать вывод, что группы имеют одинаковую дисперсию.

Двусторонний тест ANOVA в Python

Опять же, используя тот же набор данных, мы попытаемся понять, существует ли значительная связь между плотностью региона или штата, возрастом населения и количеством случаев заражения коронавирусом. Итак, мы составим карту каждого штата в зависимости от плотности населения, проживающего в нем.

Давайте импортируем данные и проверим двусмысленность данных:

individualDetails=pd.read_csv('./individualDetails.csv',parse_dates=[2])
stateDensity=pd.read_csv('./stateDensity.csv')

Из приведенного выше фрагмента кода видно, что записей о зараженных детях нет. Далее проверяем пропущенные значения в данных:

individualDetails.isna().sum()

print('Percentage of missing values in age & gender columns respectively :', \
      (individualDetails['age'].isna().sum()/individualDetails.shape[0])*100,'%',\
      (individualDetails['gender'].isna().sum()/individualDetails.shape[0])*100,'%')

Мы обнаружили, что более 91% и 80% записей отсутствовали в столбцах возраста и пола соответственно. Поэтому нам нужно придумать способ их оценки.

Здесь я оценю возраст по медианному показателю пола штата и медианному показателю пола штата. Поэтому вычислю медиану и моду:

ageMedianPerState=individualDetails[~individualDetails['age'].isna()]
ageMedianPerState['age']=ageMedianPerState['age'].astype(str).astype(int)
ageMedianPerState=ageMedianPerState.groupby('State')[['age']].median().reset_index()
ageMedianPerState['age']=ageMedianPerState['age'].apply(lambda x:math.ceil(x))

#通过COVID-19查找每个州的最常感染的性别
genderModePerState=individualDetails.groupby(['State'])['gender'].agg(pd.Series.mode).to_frame().reset_index()

#这没有获得有关总体性别的信息性别
genderModePerState=genderModePerState[genderModePerState['State']!='Arunachal Pradesh']
#现在在年龄和性别栏填充丢失的值
for index,row in individualDetails.iterrows():
    if row['State']=='Arunachal Pradesh':    
        individualDetails.drop([index],inplace=True)
        continue
    if pd.isnull(row['age']):
        individualDetails['age'][index]=list(ageMedianPerState['age'][ageMedianPerState['State']=='West Bengal'].values)[0]
    if pd.isnull(row['gender']):
        if len(genderModePerState['gender'][genderModePerState['State']==row['State']].values)>0:
            individualDetails['gender'][index]=(genderModePerState['gender'][genderModePerState['State']==row['State']].values[0])

Теперь давайте объединим кадры данных IndividualDetails и stateDensity, чтобы создать для нас общий набор данных:

data = pd.merge(individualDetails,stateDensity,on='State',how='left').reset_index(drop=True)

Теперь мы можем создать сегменты возрастных групп:

data.dropna(subset=['density_Group'],inplace=True)
data.reset_index(drop=True,inplace=True)

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

patient_Count=data.groupby(['diagnosed_date','density_Group'])[['diagnosed_date']].count().\

                        rename(columns={'diagnosed_date':'Count'}).reset_index()

data=pd.merge(data,patient_Count,on=['diagnosed_date','density_Group'],how='inner')

newData=data.groupby(['density_Group','age_Group'])['Count'].apply(list).reset_index()

newData.head()

np.random.seed(1234)

AnovaData=pd.DataFrame(columns=['density_Group','age_Group','Count'])

for index,row in newData.iterrows():
    count=17
    tempDf=pd.DataFrame(index=range(0,count),columns=['density_Group','age_Group','Count'])

    tempDf['age_Group']=newData['age_Group'][index]

    tempDf['density_Group']=newData['density_Group'][index]

    tempDf['Count']=random.sample(list(newData['Count'][index]),count)

    AnovaData=pd.concat([AnovaData,tempDf],axis=0)

Проверьте распределение столбца Count в данных и используйте метод boxplot для проверки выбросов в данных:

plt.hist(AnovaData['Count'])
plt.show() sns.kdeplot(AnovaData['Count'],cumulative=False,bw=2)

Мы обнаружили много выбросов в наших данных. Даже распределение переменной count не является гауссовым. Поэтому мы будем использовать метод преобразования Бокса-Кокса, чтобы справиться с этой ситуацией:

sns.boxplot(x='age_Group', y='Count',
                 data=AnovaData,
                 palette="colorblind")

AnovaData['Count']=AnovaData['Count'].astype(int)

## 数据转换
AnovaData['newCount'],fitted_lambda = stats.boxcox(AnovaData['Count'])

import matplotlib.pyplot as plt
sns.kdeplot(AnovaData['newCount'],cumulative=False,bw=2)

Теперь давайте воспользуемся моделью OLS для проверки нашей гипотезы:

## 拟合OlS模型-方法1
model2 = ols('newCount ~ C(age_Group)+ C(density_Group)', AnovaData).fit()

print(f"Overall model F({model2.df_model: },{model2.df_resid: }) = {model2.fvalue: }, p = {model2.f_pvalue: }")
model2.summary()

## 创建方差分析表
res2 = sm.stats.anova_lm(model2, typ=2)
res2

#检验残差的正态分布
res = model2.resid
fig = sm.qqplot(res, line='s')
plt.show()

Из приведенного выше графика Q-Q видно, что остатки распределены почти нормально (хотя в крайних конечных точках их можно не учитывать). Следовательно, мы можем сделать вывод, что он удовлетворяет предположению о нормальности теста ANOVA.

#方法2-检查组之间的交互
formula = 'newCount ~ C(age_Group) *C(density_Group)'
model = ols(formula, AnovaData).fit()
model.summary()

aov_table = anova_lm(model, typ=2)
print(aov_table.round(4))

Интерпретация результатов

Количество коронарных случаев, возрастная группа и группа плотности, а также значения P для взаимодействия, полученные с помощью анализа ANOVA, были статистически значимыми (P Делаем вывод, что типdensity_Group существенно влияет на результаты короны. ***

age_Group существенно влияет на результаты коронарных случаев, а взаимодействие age_Group иdensity_Group также существенно влияет на результаты коронарных случаев.

посмертный

Наконец, давайте определим, какие группы статистически различны. Мы будем использовать метод Tuckey HSD:

mc = statsmodels.stats.multicomp.MultiComparison(AnovaData['newCount'],AnovaData['density_Group'])
mc_results = mc.tukeyhsd()
print(mc_results)

mc = statsmodels.stats.multicomp.MultiComparison(AnovaData['newCount'],AnovaData['age_Group'])
mc_results = mc.tukeyhsd()
print(mc_results)

Из приведенных выше результатов теста Tuckey HSD мы можем ясно видеть, что существуют также значительные различия между группой 1 — группой 3, группой 1 — группой 4 в группе плотности и группой «молодой — взрослый и молодой — пожилой» в возрастной группе.

Таким образом, приведенные выше результаты Tukey HSD показывают, что все другие попарные сравнения числа коронарных случаев, за исключением вышеуказанных групп, отвергают нулевую гипотезу и не показывают статистически значимой разницы.

конец

Во время пандемии я попытался объяснить дисперсионный анализ с помощью коррелированного случая. Вы можете клонировать мой репозиторий Github, чтобы загрузить полный код и данные:GitHub.com/Praveen76/A…

Оригинальная ссылка:Woohoo.Со слов аналитиков vi.com/blog/2020/0…

Добро пожаловать на сайт блога Panchuang AI:panchuang.net/

sklearn машинное обучение китайские официальные документы:sklearn123.com/

Добро пожаловать на станцию ​​сводки ресурсов блога Panchuang:docs.panchuang.net/