Вступление

Python - невероятно универсальный язык, полезный для множества задач в самых разных дисциплинах. Одной из таких дисциплин является статистический анализ наборов данных, и наряду с SPSS Python является одним из наиболее распространенных инструментов для статистики.

Дружественный и интуитивно понятный характер Python упрощает выполнение статистических тестов и внедрение аналитических методов, особенно с помощью библиотеки statsmodels .

Представляем библиотеку statsmodels в Python

Библиотека statsmodels - это модуль для Python, который обеспечивает легкий доступ к различным статистическим инструментам для выполнения статистических тестов и изучения данных. Библиотека предоставляет доступ к ряду статистических тестов и функций, включая обычные регрессии методом наименьших квадратов (OLS), обобщенные линейные модели, логит-модели, модели анализа главных компонентов (PCA) и модели авторегрессионного интегрированного скользящего среднего (ARIMA).

Результаты моделей постоянно проверяются по сравнению с другими статистическими пакетами, чтобы гарантировать точность моделей. В сочетании с SciPy и Pandas легко визуализировать данные, запускать статистические тесты и проверять взаимосвязи на предмет значимости.

Выбор набора данных

Прежде чем мы сможем практиковать статистику с Python, нам нужно выбрать набор данных. Мы будем использовать набор данных, собранный Gapminder Foundation.

Набор данных Gapminder® отслеживает множество переменных, используемых для оценки общего состояния здоровья и благополучия населения в странах по всему миру. Мы будем использовать набор данных, потому что он очень хорошо документирован, стандартизирован и полон. Нам не нужно много делать в плане предварительной обработки, чтобы использовать ее.

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

Отправной точкой нашего статистического анализа набора данных Gapminder является исследовательский анализ данных. Мы будем использовать некоторые функции построения графиков и построения графиков из Matplotlib и Seaborn, чтобы визуализировать некоторые интересные отношения и получить представление о том, какие отношения переменных мы можем захотеть изучить.

Исследовательский анализ данных и предварительная обработка

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

Мы импортируем эти две и любые другие библиотеки, которые будем здесь использовать:

 import statsmodels.formula.api as smf 
 import statsmodels.stats.multicomp as multi 
 import scipy 
 from scipy.stats import pearsonr 
 import pandas as pd 
 from seaborn import regplot 
 import matplotlib.pyplot as plt 
 import numpy as np 
 import seaborn as sns 

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

 # Check for missing data 
 def check_missing_values(df, cols): 
 
 for col in cols: 
 print("Column {} is missing:".format(col)) 
 print((df[col].values == ' ').sum()) 
 print() 
 
 # Convert to numeric 
 def convert_numeric(dataframe, cols): 
 
 for col in cols: 
 dataframe[col] = pd.to_numeric(dataframe[col], errors='coerce') 
 
 df = pd.read_csv("gapminder.csv") 
 
 print("Null values:") 
 print(df.isnull().values.any()) 
 
 cols = ['lifeexpectancy', 'breastcancerper100th', 'suicideper100th'] 
 norm_cols = ['internetuserate', 'employrate'] 
 
 df2 = df.copy() 
 
 check_missing_values(df, cols) 
 check_missing_values(df, norm_cols) 
 
 convert_numeric(df2, cols) 
 convert_numeric(df2, norm_cols) 

Вот результаты:

 Null values: 
 Column lifeexpectancy is missing: 
 22 
 
 Column breastcancerper100th is missing: 
 40 
 
 Column suicideper100th is missing: 
 22 
 
 Column internetuserate is missing: 
 21 
 
 Column employrate is missing: 
 35 

Есть несколько пропущенных значений, но наше числовое преобразование должно превратить их в NaN , позволяя проводить исследовательский анализ данных в наборе данных.

В частности, мы могли бы попытаться проанализировать взаимосвязь между уровнем использования Интернета и ожидаемой продолжительностью жизни или между уровнем использования Интернета и уровнем занятости. Давайте попробуем построить отдельные графики некоторых из этих отношений, используя Seaborn и Matplotlib:

 sns.lmplot(x="internetuserate", y="breastcancerper100th", data=df2, fit_reg=False) 
 plt.title("Internet Use Rate and Breast Cancer Per 100k") 
 plt.show() 
 
 sns.lmplot(x="internetuserate", y="lifeexpectancy", data=df2, fit_reg=False) 
 plt.title("Internet Use Rate and Life Expectancy") 
 plt.show() 
 
 sns.lmplot(x="internetuserate", y="employrate", data=df2, fit_reg=False) 
 plt.title("Internet Use Rate and Employment Rate") 
 plt.show() 

Вот результаты графиков:

{.ezlazyload}

{.ezlazyload}

{.ezlazyload}

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

Также кажется, что существует довольно сильная, хотя и менее линейная зависимость между продолжительностью жизни и уровнем использования Интернета.

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

Выбор подходящей гипотезы

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

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

• ANOVA
• Хи-квадрат тесты
• Регрессии.

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

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

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

Дисперсионный анализ (ANOVA)

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

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

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

Мы можем преобразовать данные из непрерывных в количественные, выбрав категорию и разделив соответствующую переменную на процентили. Независимая переменная будет преобразована в категориальную переменную, а зависимая переменная останется непрерывной. Мы можем использовать qcut() в Pandas, чтобы разделить фрейм данных на бункеры:

 def bin(dataframe, cols): 
 # Create new columns that store the binned data 
 for col in cols: 
 new_col_name = "{}_bins".format(col) 
 dataframe[new_col_name] = pd.qcut(dataframe[col], 10, labels=["1=10%", "2=20%", "3=30%", "4=40%", "5=50%", "6=60%", "7=70%", "8=80", "9=90%", "10=100%"]) 
 
 df3 = df2.copy() 
 
 # This creates new columns filled with the binned column data 
 bin(df3, cols) 
 bin(df3, norm_cols) 

После того, как переменные были преобразованы и готовы к анализу, мы можем использовать statsmodel для выполнения дисперсионного анализа выбранных функций. Мы распечатаем результаты ANOVA и проверим, является ли связь между двумя переменными статистически значимой:

 anova_df = df3[['lifeexpectancy', 'internetuserate_bins', 'employrate_bins']].dropna() 
 
 relate_df = df3[['lifeexpectancy', 'internetuserate_bins']] 
 
 anova = smf.ols(formula='lifeexpectancy ~ C(internetuserate_bins)', data=anova_df).fit() 
 
 print(anova.summary()) 
 
 # We may also want to check the mean and standard deviation for the groups 
 mean = relate_df.groupby("internetuserate_bins").mean() 
 sd = relate_df.groupby("internetuserate_bins").std() 
 print(mean) 
 print(sd) 

Вот результат работы модели:

 OLS Regression Results 
 ============================================================================== 
 Dep. Variable: lifeexpectancy R-squared: 0.689 
 Model: OLS Adj. R-squared: 0.671 
 Method: Least Squares F-statistic: 38.65 
 Date: Mon, 11 May 2020 Prob (F-statistic): 1.71e-35 
 Time: 17:49:24 Log-Likelihood: -521.54 
 No. Observations: 167 AIC: 1063. 
 Df Residuals: 157 BIC: 1094. 
 Df Model: 9 
 Covariance Type: nonrobust 
 ====================================================================================================== 
 coef std err t P>|t| [0.025 0.975] 
 ------------------------------------------------------------------------------------------------------ 
 Intercept 56.6603 1.268 44.700 0.000 54.157 59.164 
 C(internetuserate_bins)[T.2=20%] 1.6785 1.870 0.898 0.371 -2.015 5.372 
 C(internetuserate_bins)[T.3=30%] 5.5273 1.901 2.907 0.004 1.772 9.283 
 C(internetuserate_bins)[T.4=40%] 11.5693 1.842 6.282 0.000 7.932 15.207 
 C(internetuserate_bins)[T.5=50%] 14.6991 1.870 7.860 0.000 11.005 18.393 
 C(internetuserate_bins)[T.6=60%] 16.7287 1.870 8.946 0.000 13.035 20.422 
 C(internetuserate_bins)[T.7=70%] 17.8802 1.975 9.052 0.000 13.978 21.782 
 C(internetuserate_bins)[T.8=80] 19.8302 1.901 10.430 0.000 16.075 23.586 
 C(internetuserate_bins)[T.9=90%] 23.0723 1.901 12.135 0.000 19.317 26.828 
 C(internetuserate_bins)[T.10=100%] 23.3042 1.901 12.257 0.000 19.549 27.060 
 ============================================================================== 
 Omnibus: 10.625 Durbin-Watson: 1.920 
 Prob(Omnibus): 0.005 Jarque-Bera (JB): 11.911 
 Skew: -0.484 Prob(JB): 0.00259 
 Kurtosis: 3.879 Cond. No. 10.0 
 ============================================================================== 
 
 Warnings: 
 [1] Standard Errors assume that the covariance matrix of the errors is correctly specified. 

Мы видим, что модель дает очень маленькое P-значение ( вероятность F-статистика ) 1,71e-35 . Это намного меньше обычного порога значимости 0,05 , поэтому мы заключаем, что существует значимая взаимосвязь между продолжительностью жизни и уровнем использования Интернета.

Поскольку значение P корреляции действительно кажется значительным, и поскольку у нас есть 10 различных категорий, мы захотим запустить апостериорный тест, чтобы проверить, что разница между средними значениями все еще значительна даже после того, как мы проверим тип-1. ошибки. Мы можем проводить апостериорные тесты с помощью multicomp , используя тест Tukey Honestly Significant Difference (Tukey HSD):

 multi_comparison = multi.MultiComparison(anova_df["lifeexpectancy"], anova_df["internetuserate_bins"]) 
 results = multi_comparison.tukeyhsd() 
 print(results) 

Вот результаты теста:

 Multiple Comparison of Means - Tukey HSD, FWER=0.05 
 ======================================================= 
 group1 group2 meandiff p-adj lower upper reject 
 ------------------------------------------------------- 
 10=100% 1=10% -23.3042 0.001 -29.4069 -17.2015 True 
 10=100% 2=20% -21.6257 0.001 -27.9633 -15.2882 True 
 10=100% 3=30% -17.7769 0.001 -24.2097 -11.344 True 
 10=100% 4=40% -11.7349 0.001 -17.9865 -5.4833 True 
 10=100% 5=50% -8.6051 0.001 -14.9426 -2.2676 True 
 10=100% 6=60% -6.5755 0.0352 -12.913 -0.238 True 
 10=100% 7=70% -5.4241 0.2199 -12.0827 1.2346 False 
 10=100% 8=80 -3.4741 0.7474 -9.9069 2.9588 False 
 10=100% 9=90% -0.2319 0.9 -6.6647 6.201 False 
 1=10% 2=20% 1.6785 0.9 -4.3237 7.6807 False 
 1=10% 3=30% 5.5273 0.1127 -0.5754 11.6301 False 
 1=10% 4=40% 11.5693 0.001 5.6579 17.4807 True 
 1=10% 5=50% 14.6991 0.001 8.6969 20.7013 True 
 1=10% 6=60% 16.7287 0.001 10.7265 22.7309 True 
 1=10% 7=70% 17.8801 0.001 11.5399 24.2204 True 
 1=10% 8=80 19.8301 0.001 13.7274 25.9329 True 
 1=10% 9=90% 23.0723 0.001 16.9696 29.1751 True 
 2=20% 3=30% 3.8489 0.6171 -2.4887 10.1864 False 
 2=20% 4=40% 9.8908 0.001 3.7374 16.0443 True 
 2=20% 5=50% 13.0206 0.001 6.7799 19.2614 True 
 2=20% 6=60% 15.0502 0.001 8.8095 21.291 True 
 2=20% 7=70% 16.2017 0.001 9.6351 22.7683 True 
 2=20% 8=80 18.1517 0.001 11.8141 24.4892 True 
 2=20% 9=90% 21.3939 0.001 15.0563 27.7314 True 
 3=30% 4=40% 6.042 0.0678 -0.2096 12.2936 False 
 3=30% 5=50% 9.1718 0.001 2.8342 15.5093 True 
 3=30% 6=60% 11.2014 0.001 4.8638 17.5389 True 
 3=30% 7=70% 12.3528 0.001 5.6942 19.0114 True 
 3=30% 8=80 14.3028 0.001 7.87 20.7357 True 
 3=30% 9=90% 17.545 0.001 11.1122 23.9778 True 
 4=40% 5=50% 3.1298 0.8083 -3.0237 9.2833 False 
 4=40% 6=60% 5.1594 0.1862 -0.9941 11.3129 False 
 4=40% 7=70% 6.3108 0.0638 -0.1729 12.7945 False 
 4=40% 8=80 8.2608 0.0015 2.0092 14.5124 True 
 4=40% 9=90% 11.503 0.001 5.2514 17.7546 True 
 5=50% 6=60% 2.0296 0.9 -4.2112 8.2704 False 
 5=50% 7=70% 3.181 0.8552 -3.3856 9.7476 False 
 5=50% 8=80 5.131 0.2273 -1.2065 11.4686 False 
 5=50% 9=90% 8.3732 0.0015 2.0357 14.7108 True 
 6=60% 7=70% 1.1514 0.9 -5.4152 7.718 False 
 6=60% 8=80 3.1014 0.8456 -3.2361 9.439 False 
 6=60% 9=90% 6.3436 0.0496 0.0061 12.6812 True 
 7=70% 8=80 1.95 0.9 -4.7086 8.6086 False 
 7=70% 9=90% 5.1922 0.2754 -1.4664 11.8508 False 
 8=80 9=90% 3.2422 0.8173 -3.1907 9.675 False 
 ------------------------------------------------------- 

Теперь мы лучше понимаем, какие группы в нашем сравнении имеют статистически значимые различия.

Если reject имеет метку False , мы знаем, что рекомендуется отклонить нулевую гипотезу и предположить, что существует значительная разница между двумя сравниваемыми группами.

Проверка независимости хи-квадрат

ANOVA подходит для случаев, когда одна переменная является непрерывной, а другая категориальной. Теперь мы посмотрим, как провести тест независимости по критерию хи-квадрат .

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

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

Подобно тому, как мы преобразовали нашу независимую переменную в категориальную переменную (путем ее объединения), для теста ANOVA нам нужно сделать обе переменные категориальными, чтобы выполнить тест хи-квадрат. Наша гипотеза для этой проблемы такая же, как и гипотеза в предыдущей задаче, что существует значительная взаимосвязь между продолжительностью жизни и уровнем использования Интернета.

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

Мы будем проводить апостериорное сравнение, чтобы защититься от ошибок типа 1 (ложных срабатываний), используя подход, называемый корректировкой Бонферрони . Для этого вы можете выполнить сравнения для различных возможных пар вашей переменной ответа, а затем проверить их скорректированную значимость.

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

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

 def half_bin(dataframe, cols): 
 
 for col in cols: 
 new_col_name = "{}_bins_2".format(col) 
 dataframe[new_col_name] = pd.qcut(dataframe[col], 2, labels=["1=50%", "2=100%"]) 
 
 half_bin(df3, ['internetuserate']) 
 
 # Recoding scheme 
 recode_2 = {"3=30%": "3=30%", "7=70%": "7=70%"} 
 recode_3 = {"2=20%": "2=20%", "8=80": "8=80"} 
 recode_4 = {"6=60%": "6=60%", "9=90%": "9=90%"} 
 recode_5 = {"4=40%": "4=40%", "7=70%": "7=70%"} 
 
 # Create the new features 
 df3['Comp_3v7'] = df3['lifeexpectancy_bins'].map(recode_2) 
 df3['Comp_2v8'] = df3['lifeexpectancy_bins'].map(recode_3) 
 df3['Comp_6v9'] = df3['lifeexpectancy_bins'].map(recode_4) 
 df3['Comp_4v7'] = df3['lifeexpectancy_bins'].map(recode_5) 

Выполнение теста хи-квадрат и апостериорного сравнения включает в себя сначала построение сравнительной таблицы с перекрестными таблицами. Таблица сравнения перекрестных таблиц показывает процент встречаемости переменной ответа для разных уровней объясняющей переменной.

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

 # Get table of observed counts 
 count_table = pd.crosstab(df3['internetuserate_bins_2'], df3['lifeexpectancy_bins']) 
 print(count_table) 

 lifeexpectancy_bins 1=10% 2=20% 3=30% 4=40% ... 7=70% 8=80 9=90% 10=100% 
 internetuserate_bins_2 ... 
 1=50% 18 19 16 14 ... 4 4 1 0 
 2=100% 0 0 1 4 ... 15 11 16 19 

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

Теперь нам нужно вычислить кросс-таблицы для различных пар, которые мы создали выше, так как это то, что мы запускаем в тесте хи-квадрат:

 count_table_3 = pd.crosstab(df3['internetuserate_bins_2'], df3['Comp_3v7']) 
 count_table_4 = pd.crosstab(df3['internetuserate_bins_2'], df3['Comp_2v8']) 
 count_table_5 = pd.crosstab(df3['internetuserate_bins_2'], df3['Comp_6v9']) 
 count_table_6 = pd.crosstab(df3['internetuserate_bins_2'], df3['Comp_4v7']) 

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

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

 def chi_sq_test(table): 
 
 print("Results for:") 
 print(str(table)) 
 
 # Get column percentages 
 col_sum = table.sum(axis=0) 
 col_percents = table/col_sum 
 print(col_percents) 
 
 chi_square = scipy.stats.chi2_contingency(table) 
 print("Chi-square value, p-value, expected_counts") 
 print(chi_square) 
 
 print() 
 
 print("Initial Chi-square:") 
 chi_sq_test(count_table) 
 print(" ") 
 
 chi_sq_test(count_table_3) 
 chi_sq_test(count_table_4) 
 chi_sq_test(count_table_5) 
 chi_sq_test(count_table_6) 

Вот результаты:

 Initial Chi-square: 
 Results for: 
 lifeexpectancy_bins 1=10% 2=20% 3=30% 4=40% ... 7=70% 8=80 9=90% 10=100% 
 internetuserate_bins_2 ... 
 1=50% 18 19 16 14 ... 4 4 1 0 
 2=100% 0 0 1 4 ... 15 11 16 19 
 
 [2 rows x 10 columns] 
 lifeexpectancy_bins 1=10% 2=20% 3=30% ... 8=80 9=90% 10=100% 
 internetuserate_bins_2 ... 
 1=50% 1.0 1.0 0.941176 ... 0.266667 0.058824 0.0 
 2=100% 0.0 0.0 0.058824 ... 0.733333 0.941176 1.0 
 
 [2 rows x 10 columns] 
 Chi-square value, p-value, expected_counts 
 (102.04563740451277, 6.064860600653971e-18, 9, array([[9.45251397, 9.97765363, 8.9273743 , 9.45251397, 9.45251397, 
 9.97765363, 9.97765363, 7.87709497, 8.9273743 , 9.97765363], 
 [8.54748603, 9.02234637, 8.0726257 , 8.54748603, 8.54748603, 
 9.02234637, 9.02234637, 7.12290503, 8.0726257 , 9.02234637]])) 
 ----- 
 
 Results for: 
 Comp_3v7 3=30% 7=70% 
 internetuserate_bins_2 
 1=50% 16 4 
 2=100% 1 15 
 Comp_3v7 3=30% 7=70% 
 internetuserate_bins_2 
 1=50% 0.941176 0.210526 
 2=100% 0.058824 0.789474 
 Chi-square value, p-value, expected_counts 
 (16.55247678018576, 4.7322137795376575e-05, 1, array([[ 9.44444444, 10.55555556], 
 [ 7.55555556, 8.44444444]])) 
 ----- 
 Results for: 
 Comp_2v8 2=20% 8=80 
 internetuserate_bins_2 
 1=50% 19 4 
 2=100% 0 11 
 Comp_2v8 2=20% 8=80 
 internetuserate_bins_2 
 1=50% 1.0 0.266667 
 2=100% 0.0 0.733333 
 Chi-square value, p-value, expected_counts 
 (17.382650301643437, 3.0560286589975315e-05, 1, array([[12.85294118, 10.14705882], 
 [ 6.14705882, 4.85294118]])) 
 ----- 
 Results for: 
 Comp_6v9 6=60% 9=90% 
 internetuserate_bins_2 
 1=50% 6 1 
 2=100% 13 16 
 Comp_6v9 6=60% 9=90% 
 internetuserate_bins_2 
 1=50% 0.315789 0.058824 
 2=100% 0.684211 0.941176 
 Chi-square value, p-value, expected_counts 
 (2.319693757720874, 0.12774517376836148, 1, array([[ 3.69444444, 3.30555556], 
 [15.30555556, 13.69444444]])) 
 ----- 
 Results for: 
 Comp_4v7 4=40% 7=70% 
 internetuserate_bins_2 
 1=50% 14 4 
 2=100% 4 15 
 Comp_4v7 4=40% 7=70% 
 internetuserate_bins_2 
 1=50% 0.777778 0.210526 
 2=100% 0.222222 0.789474 
 Chi-square value, p-value, expected_counts 
 (9.743247922437677, 0.0017998260000241526, 1, array([[8.75675676, 9.24324324], 
 [9.24324324, 9.75675676]])) 
 ----- 

Если мы смотрим только на результаты для полной таблицы подсчета, похоже, что существует P-значение 6.064860600653971e-18 .

Однако, чтобы выяснить, как разные группы отличаются друг от друга, нам необходимо выполнить тест хи-квадрат для разных пар в нашем фреймворке данных. Мы проверим, есть ли статистически значимая разница для каждой из выбранных пар. Обратите внимание, что P-значение, которое указывает на значительный результат, изменяется в зависимости от того, сколько сравнений вы проводите, и, хотя мы не будем рассматривать это в этом руководстве, вам необходимо помнить об этом.

Сравнение 6 и 9 дает нам P-значение 0.127 , что выше 0.05 , что указывает на то, что разница для этой категории может быть незначительной. Видение различий в сравнениях помогает нам понять, почему нам нужно сравнивать разные уровни друг с другом.

Корреляции Пирсона

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

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

Тест корреляции Пирсона используется для анализа силы взаимосвязи между двумя предоставленными переменными, причем обе являются количественными по своей природе. Значение или сила корреляции Пирсона будет между +1 и -1 .

Корреляция 1 указывает на идеальную связь между переменными, а корреляция может быть положительной или отрицательной. Коэффициенты корреляции около 0 указывают на очень слабые, почти несуществующие корреляции. Хотя есть и другие способы измерения корреляции между двумя переменными, такие как корреляция Спирмена или ранговая корреляция Кендалла, корреляция Пирсона, вероятно, является наиболее часто используемым корреляционным тестом.

Поскольку набор данных Gapminder имеет свои функции, представленные с помощью количественных переменных, нам не нужно выполнять какое-либо категориальное преобразование данных перед запуском для них корреляции Пирсона. Обратите внимание, что предполагается, что обе переменные имеют нормальное распределение и в наборе данных не так много значительных выбросов. Нам понадобится доступ к SciPy, чтобы провести корреляцию Пирсона.

Мы построим график зависимости между ожидаемой продолжительностью жизни и уровнем использования Интернета, а также уровнем использования Интернета и уровнем занятости, просто чтобы увидеть, как может выглядеть другой график корреляции. После создания графической функции мы будем использовать personr() из SciPy, чтобы провести корреляцию и проверить результаты:

 df_clean = df2.dropna() 
 df_clean['incomeperperson'] = df_clean['incomeperperson'].replace('', np.nan) 
 
 def plt_regression(x, y, data, label_1, label_2): 
 
 reg_plot = regplot(x=x, y=y, fit_reg=True, data=data) 
 plt.xlabel(label_1) 
 plt.ylabel(label_2) 
 plt.show() 
 
 plt_regression('lifeexpectancy', 'internetuserate', df_clean, 'Life Expectancy', 'Internet Use Rate') 
 plt_regression('employrate', 'internetuserate', df_clean, 'Employment Rate', 'Internet Use Rate') 
 
 print('Assoc. - life expectancy and internet use rate') 
 print(pearsonr(df_clean['lifeexpectancy'], df_clean['internetuserate'])) 
 
 print('Assoc. - between employment rate and internet use rate') 
 print(pearsonr(df_clean['employrate'], df_clean['internetuserate'])) 

Вот результаты:

{.ezlazyload}

{.ezlazyload}

 Assoc. - life expectancy and internet use rate 
 (0.77081050888289, 5.983388253650836e-33) 
 Assoc. - between employment rate and internet use rate 
 (-0.1950109538173115, 0.013175901971555317) 

Первое значение - это направление и сила корреляции, а второе - P-значение. Цифры указывают на довольно сильную корреляцию между продолжительностью жизни и уровнем использования Интернета, которая не является случайной. Между тем, существует более слабая, но все же значимая корреляция между уровнем занятости и уровнем использования Интернета.

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

Модераторы и статистическое взаимодействие

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

Модерация - это когда третья (или более) переменная влияет на силу связи между независимой переменной и зависимой переменной.

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

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

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

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

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

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

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

 def income_groups(row): 
 if row['incomeperperson'] <= 744.23: 
 return 1 
 elif row['incomeperperson'] <= 942.32: 
 return 2 
 else: 
 return 3 
 
 
 # Apply function and set the new features in the dataframe 
 df_clean['income_group'] = df_clean.apply(lambda row: income_groups(row), axis=1) 
 
 # Create a few subframes to try test for moderation 
 subframe_1 = df_clean[(df_clean['income_group'] == 1)] 
 subframe_2 = df_clean[(df_clean['income_group'] == 2)] 
 subframe_3 = df_clean[(df_clean['income_group'] == 3)] 
 
 print('Assoc. - life expectancy and internet use rate for low income countries') 
 
 print(pearsonr(subframe_1['lifeexpectancy'], subframe_1['internetuserate'])) 
 
 print('Assoc. - life expectancy and internet use rate for medium income countries') 
 
 print(pearsonr(subframe_2['lifeexpectancy'], subframe_2['internetuserate'])) 
 
 print('Assoc. - life expectancy and internet use rate for high income countries') 
 
 print(pearsonr(subframe_3['lifeexpectancy'], subframe_3['internetuserate'])) 

Вот результаты:

 Assoc. - life expectancy and internet use rate for low income countries 
 (0.38386370068495235, 0.010101223355274047) 
 
 Assoc. - life expectancy and internet use rate for medium income countries 
 (0.9966009508278395, 0.05250454954743393) 
 
 Assoc. - life expectancy and internet use rate for high income countries 
 (0.7019997488251704, 6.526819886007788e-18) 

Еще раз, первое значение - это направление и сила корреляции, а второе - это P-значение.

Заключение

statsmodels - чрезвычайно полезная библиотека, которая позволяет пользователям Python анализировать данные и запускать статистические тесты для наборов данных. Вы можете выполнять ANOVA, тесты хи-квадрат, корреляции Пирсона и тестировать на модерацию.

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