Тест Кохрана-Армитиджа для линейных трендов в Python

Тест Кохрана-Армитиджа является мощным инструментом в статистике Python для выявления линейных трендов в категориальных данных, особенно при анализе дозо-зависимых реакций. С помощью функции test_ordinal_association из библиотеки statsmodels вы можете эффективно проверять гипотезу о наличии линейной ассоциации между упорядоченными категориями, такими как разные уровни дозировки и соответствующие ответы.

Содержание

• Введение в тест Кохрана-Армитиджа
• Подготовка данных в Python
• Использование функции test_ordinal_association
• Интерпретация результатов
• Практические примеры анализа дозо-зависимых реакций
• Источники
• Заключение

Введение в тест Кохрана-Армитиджа и его применение для выявления линейных трендов

Тест Кохрана-Армитиджа (Cochran-Armitage test) — это статистический метод, предназначенный для выявления линейных трендов в таблицах сопряженности, где одна переменная имеет упорядоченные категории. В контексте анализа дозо-зависимых реакций этот тест позволяет определить, существует ли линейная связь между уровнем дозировки (упорядоченной переменной) и частотой наблюдаемого отклика.

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

В статистике Python этот тест реализован в библиотеке statsmodels через функцию test_ordinal_association, которая предоставляет удобный интерфейс для проведения анализа и получения статистически обоснованных выводов о наличии линейного тренда в ваших данных.

Математическая основа теста

Математически тест Кохрана-Армитиджа основан на оценке линейного контраста в таблице сопряженности. Статистика теста рассчитывается как:

$$Z = \frac{\sum_{j=1}^{c} w_j (p_j - \bar{p})}{\sqrt{\bar{p}(1-\bar{p}) \sum_{j=1}^{c} \frac{n_j w_j^2}{n_{j\cdot}}}}$$

где:

• $w_j$ — веса для категорий (обычно используются номера категорий)
• $p_j$ — пропорция положительных откликов в j-й категории
• $\bar{p}$ — общая пропорция положительных откликов
• $n_j$ — объем выборки в j-й категории
• $n_{j\cdot}$ — общее количество наблюдений в j-й категории

Подготовка данных в Python для анализа линейных трендов

Перед применением теста Кохрана-Армитиджа важно правильно подготовить ваши данные в Python. В статистике Python существует несколько подходов к организации данных для анализа линейных трендов, но для теста Кохрана-Армитиджа данные должны быть представлены в виде таблицы сопряженности.

Формат данных

Ваши данные должны быть организованы в виде двумерной таблицы, где:

• Строки представляют различные уровни дозировки (упорядоченные категории)
• Столбцы представляют возможные отклики (обычно бинарные: “позитивный” и “негативный”)
• Ячейки содержат количество наблюдений для каждой комбинации дозировки и отклика

Пример такой таблицы для анализа дозо-зависимых реакций:

Создание таблицы сопряженности в Python

В Python вы можете создать таблицу сопряженности с помощью библиотеки pandas или numpy. Вот пример кода:

import pandas as pd
import numpy as np

# Создание таблицы сопряженности
data = np.array([[10, 40], # Низкая дозировка
 [25, 25], # Средняя дозировка
 [35, 15]]) # Высокая дозировка

# Создание DataFrame для удобства
df = pd.DataFrame(data, 
 columns=['Позитивный отклик', 'Негативный отклик'],
 index=['Низкая', 'Средняя', 'Высокая'])

print(df)

Проверка предпосылок теста

Перед применением теста Кохрана-Армитиджа убедитесь, что выполнены следующие предпосылки:

1. Данные представляют собой независимые наблюдения
2. Таблица имеет форму r×2, где r — количество уровней дозировки
3. Уровни дозировки имеют естественный порядок
4. Ожидаемые частоты не слишком малы (обычно не менее 5 в большинстве ячеек)

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

Использование функции test_ordinal_association из библиотеки statsmodels

Библиотека statsmodels предоставляет удобную функцию test_ordinal_association для проведения теста Кохрана-Армитиджа в Python. Эта функция реализует алгоритм выявления линейных трендов в категориальных данных и возвращает статистику теста и p-value.

Установка и импорт необходимых библиотек

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

pip install statsmodels numpy pandas

Импортируйте необходимые модули:

import numpy as np
import pandas as pd
from statsmodels.stats.contingency_tables import test_ordinal_association

Пример применения функции

Вот как использовать функцию test_ordinal_association для анализа ваших данных:

# Таблица сопряженности (доза vs ответ)
table = np.array([[10, 40], # Низкая дозировка
 [25, 25], # Средняя дозировка
 [35, 15]]) # Высокая дозировка

# Проведение теста Кохрана-Армитиджа
result = test_ordinal_association(table)

# Вывод результатов
print("Статистика теста:", result.statistic)
print("p-value:", result.pvalue)
print("Z-оценка:", result.zscore)

Параметры функции

Функция test_ordinal_association принимает следующие параметры:

• table (array-like): Таблица сопряженности 2D-массив или DataFrame
• weights (array-like, optional): Веса для категорий (по умолчанию используются номера категорий 0, 1, 2,…)
• method (str, optional): Метод расчета (‘score’ по умолчанию или ‘wald’)
• alternative (str, optional): Альтернативная гипотеза (‘two-sided’ по умолчанию, ‘increasing’ или ‘decreasing’)

Пример с указанием весов

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

# Собственные веса для дозировок (например, 1, 2, 4 для низкой, средней, высокой)
custom_weights = np.array([1, 2, 4])

# Проведение теста с пользовательскими весами
result = test_ordinal_association(table, weights=custom_weights)

print("Результаты с пользовательскими весами:")
print("Статистика теста:", result.statistic)
print("p-value:", result.pvalue)

Обработка результатов

Функция test_ordinal_association возвращает объект с атрибутами:

• statistic: Значение статистики теста
• pvalue: p-value для проверки гипотезы
• zscore: Z-оценка теста
• weights: Использованные веса

Сравните полученное p-value с выбранным уровнем значимости (обычно 0.05). Если p-value < α, отвергайте нулевую гипотезу о отсутствии линейного тренда.

Интерпретация результатов теста Кохрана-Армитиджа в Python

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

Статистическая значимость

Основной показатель, на который следует обратить внимание — это p-value. В статистике Python стандартные уровни значимости:

• p-value < 0.05: Существует статистически значимый линейный тренд
• 0.05 ≤ p-value < 0.1: Слабые доказательства наличия тренда
• p-value ≥ 0.1: Недостаточно доказательств для заключения о наличии тренда

Пример интерпретации:

if result.pvalue < 0.05:
 print("Статистически значимый линейный тренд обнаружен (p =", round(result.pvalue, 4), ")")
else:
 print("Статистически значимый линейный тренд не обнаружен (p =", round(result.pvalue, 4), ")")

Направление тренда

Знак Z-оценки указывает на направление линейного тренда:

• Z > 0: Положительный линейный тренд (с увеличением дозы увеличивается частота позитивных откликов)
• Z < 0: Отрицательный линейный тренд (с увеличением дозы уменьшается частота позитивных откликов)

Пример интерпретации направления:

if result.zscore > 0:
 print("Положительный линейный тренд: с увеличением дозы увеличивается частота позитивных откликов")
elif result.zscore < 0:
 print("Отрицательный линейный тренд: с увеличением дозы уменьшается частота позитивных откликов")
else:
 print("Линейный тренд не обнаружен")

Размер эффекта

Для оценки силы обнаруженного линейного тренда можно рассчитать меру ассоциации, такую как коэффициент контингенции или коэффициент корреляции. В Python это можно сделать следующим образом:

# Расчет коэффициента контингенции
n = np.sum(table)
phi = result.statistic / np.sqrt(n * (n - 1)) * np.sqrt(np.sum(table[0, :]))

print("Коэффициент контингенции (φ):", round(phi, 4))

Визуализация результатов

Для лучшего понимания результатов рекомендуется визуализировать данные. В Python это можно сделать с помощью matplotlib или seaborn:

import matplotlib.pyplot as plt
import seaborn as sns

# Расчет пропорций позитивных откликов по дозировкам
proportions = table[:, 0] / np.sum(table, axis=1)

# Создание графика
plt.figure(figsize=(10, 6))
sns.lineplot(x=['Низкая', 'Средняя', 'Высокая'], y=proportions, marker='o')
plt.title('Зависимость частоты позитивных откликов от уровня дозировки')
plt.xlabel('Уровень дозировки')
plt.ylabel('Пропорция позитивных откликов')
plt.ylim(0, 1)
plt.grid(True, alpha=0.3)
plt.show()

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

Практические примеры анализа дозо-зависимых реакций с помощью теста

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

Пример 1: Анализ эффективности дозы лекарства

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

import numpy as np
from statsmodels.stats.contingency_tables import test_ordinal_association

# Данные: [улучшилось, не улучшилось] для каждой дозы
dose_data = np.array([[8, 32], # Плацебо
 [15, 25], # Низкая доза (10 мг)
 [22, 18], # Средняя доза (20 мг)
 [28, 12]]) # Высокая доза (40 мг)

# Проведение теста
result = test_ordinal_association(dose_data)

print("Пример 1: Анализ эффективности дозы лекарства")
print("Статистика теста:", round(result.statistic, 4))
print("p-value:", round(result.pvalue, 6))
print("Z-оценка:", round(result.zscore, 4))

# Интерпретация
if result.pvalue < 0.05:
 print("\nВывод: Обнаружен статистически значимый линейный тренд")
 if result.zscore > 0:
 print("С увеличением дозы эффективность лекарства повышается")
 else:
 print("С увеличением дозы эффективность лекарства снижается")
else:
 print("\nВывод: Статистически значимый линейный тренд не обнаружен")

Пример 2: Токсикологический анализ

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

# Данные: [реакция отсутствует, реакция присутствует] для каждой дозы
toxic_data = np.array([[45, 5], # Контроль (0 мг/кг)
 [38, 12], # Низкая доза (5 мг/кг)
 [25, 25], # Средняя доза (10 мг/кг)
 [10, 40]]) # Высокая доза (20 мг/кг)

# Проведение теста
result = test_ordinal_association(toxic_data)

print("\nПример 2: Токсикологический анализ")
print("Статистика теста:", round(result.statistic, 4))
print("p-value:", round(result.pvalue, 6))
print("Z-оценка:", round(result.zscore, 4))

# Визуализация
import matplotlib.pyplot as plt

proportions = toxic_data[:, 1] / np.sum(toxic_data, axis=1)
plt.figure(figsize=(10, 6))
plt.plot(['Контроль', '5 мг/кг', '10 мг/кг', '20 мг/кг'], proportions, 'o-', linewidth=2, markersize=8)
plt.title('Зависимость частоты реакций от дозы токсического вещества')
plt.xlabel('Доза (мг/кг)')
plt.ylabel('Пропорция реакций')
plt.ylim(0, 1)
plt.grid(True, alpha=0.3)
plt.show()

Пример 3: Анализ дозо-зависимого эффекта в экологических исследованиях

В экологических исследованиях часто изучают влияние концентрации загрязняющего вещества на выживаемость организмов:

# Данные: [выжило, не выжило] для каждой концентрации
ecotox_data = np.array([[50, 0], # Контроль (0 мг/л)
 [45, 5], # Низкая концентрация (1 мг/л)
 [35, 15], # Средняя концентрация (5 мг/л)
 [20, 30], # Высокая концентрация (10 мг/л)
 [5, 45]]) # Очень высокая концентрация (20 мг/л)

# Проведение теста
result = test_ordinal_association(ecotox_data)

print("\nПример 3: Анализ дозо-зависимого эффекта в экологических исследованиях")
print("Статистика теста:", round(result.statistic, 4))
print("p-value:", round(result.pvalue, 6))
print("Z-оценка:", round(result.zscore, 4))

# Расчет LC50 (концентрация, вызывающая гибель 50% организмов)
cumulative = np.cumsum(ecotox_data[:, 1])
total_organisms = np.sum(ecotox_data)
lc50_dose = None

for i, cum in enumerate(cumulative):
 if cum >= total_organisms * 0.5:
 lc50_dose = [1, 5, 10, 20][i] # Дозы в примере
 break

print(f"\nОценочная LC50: {lc50_dose} мг/л")

Пример 4: Сравнение линейного и нелинейного трендов

Интересно сравнить результаты теста Кохрана-Армитиджа с другими методами выявления трендов, такими как линейная регрессия Python:

from scipy import stats

# Данные с нелинейным трендом (U-образная зависимость)
nonlinear_data = np.array([[30, 10], # Низкая доза
 [20, 20], # Средне-низкая доза
 [15, 25], # Средняя доза
 [20, 20], # Средне-высокая доза
 [30, 10]]) # Высокая доза

# Тест Кохрана-Армитиджа
result_cochran = test_ordinal_association(nonlinear_data)
print("\nПример 4: Сравнение методов")
print("Тест Кохрана-Армитиджа:")
print("p-value:", round(result_cochran.pvalue, 6))

# Линейная регрессия Python
x = np.arange(len(nonlinear_data)) # Дозы как числовые значения
y = nonlinear_data[:, 0] / np.sum(nonlinear_data, axis=1) # Пропорции

slope, intercept, r_value, p_value, std_err = stats.linregress(x, y)
print("\nЛинейная регрессия Python:")
print("p-value для наклона:", round(p_value, 6))
print("R-квадрат:", round(r_value**2, 4))

# Хи-квадрат тест на общую зависимость
from scipy.stats import chi2_contingency
chi2, p_chi2, dof, expected = chi2_contingency(nonlinear_data)
print("\nХи-квадрат тест на общую зависимость:")
print("p-value:", round(p_chi2, 6))

Этот пример показывает, что тест Кохрана-Армитиджа может не обнаружить тренд, если он не является линейным, в то время как общий хи-квадрат тест может показать наличие зависимости, а линейная регрессия Python может не выявить значимый наклон.

Источники

1. Statsmodels Documentation — Официальная документация по функции test_ordinal_association: https://www.statsmodels.org/stable/generated/statsmodels.stats.contingency_tables.test_ordinal_association.html
2. Cochran-Armitage Trend Test — Математическое обоснование теста Кохрана-Армитиджа: https://en.wikipedia.org/wiki/Cochran–Armitage_test_for_trend
3. Ordinal Association Testing — Методы проверки порядковой ассоциации в статистике: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC2935870/
4. Python for Statistical Analysis — Использование Python для статистического анализа: https://realpython.com/python-statistics/
5. Dose-Response Analysis — Методы анализа дозо-зависимых реакций: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC4086060/

Заключение

Тест Кохрана-Армитиджа является мощным инструментом в статистике Python для выявления линейных трендов в категориальных данных, особенно при анализе дозо-зависимых реакций. Функция test_ordinal_association из библиотеки statsmodels предоставляет удобный и надежный способ проверки гипотез о наличии линейной ассоциации между упорядоченными категориями.

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

1. Правильная подготовка данных в виде таблицы сопряженности
2. Учет предпосылок теста, особенно ожидаемых частот
3. Интерпретация как статистической значимости (p-value), так и направления тренда (Z-оценка)
4. Визуализация результатов для лучшего понимания обнаруженных закономерностей

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

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