Протокол разведочного анализа данных: выявление выбросов

За последние несколько десятилетий арсенал прикладной статистики значительно пополнился за счет развития новых методов, таких как, например, обобщенные линейные модели (Generalized Linear Models), обобщенные аддитивные модели (Generalized Additive Models), модели со смешанными эффектами (Mixed-Effects Models), деревья принятия решений (Regression and Classification Trees), анализ выживаемости (Survival Analysis), многомерное шкалирование (Multidimensional Scaling), и др. Более того, появление быстрых современных компьютеров и свободного программного обеспечения (вроде R) сделало все эти требующие вычислительных ресурсов методы доступными практически для каждого исследователя. Однако такая доступность еще больше обостряет хорошо известную проблему всех статистических методов, которую на английском языке часто описывают как "rubbish in, rubbish out", т.е. "мусор на входе - мусор на выходе". Речь здесь идет о следующем: чудес не бывает, и если мы не будем уделять должного внимания тому, как тот или иной метод работает и какие требования предъявляет к анализируемым данным, то получаемые с его помощью результаты нельзя будет воспринимать всерьез. Поэтому каждый раз исследователю следует начинать свою работу с тщательного ознакомления со свойствами полученных данных и проверки необходимых условий применимости соответствующих статистических методов. Этот начальный этап анализа называют разведочным (Exploratory Data Analysis).

В литературе по статистике можно найти немало рекомендаций по выполнению разведочного анализа данных (РДА). Два года назад в журнале Methods in Ecology and Evolution была опубликована отличная статья, в которой эти рекомендации сведены в единый протокол по выполнению РДА: Zuur A. F., Ieno E. N., Elphick C. S. (2010) A protocol for data exploration to avoid common statistical problems. Methods in Ecology and Evolution 1(1): 3-14. Несмотря на то, что статья написана для биологов (в частности, для экологов), изложенные в ней принципы, безусловно, верны и в отношении других научных дисциплин. В этом и последующих сообщениях блога я приведу выдержки из работы Zuur et al. (2010) и опишу предложенный авторами РДА-протокол. Подобно тому, как это сделано в оригинальной статье, описание отдельных шагов протокола будет сопровождаться краткими рекомендациями по использованию соответствующих функций и пакетов системы R.

Предлагаемый протокол включает следующие основные элементы:

1. Формулировка исследовательской гипотезы. Выполнение экспериментов/наблюдений для сбора данных.
2. Разведочный анализ данных:
   
   • Выявление точек-выборосов
   • Проверка однородности дисперсий
   • Проверка нормальности распределения данных
   • Выявление избыточного количества нулевых значений
   • Выявление коллинеарных переменных
   • Выявление характера связи между анализируемыми переменными
   • Выявление взаимодействий между переменными-предикторами
   • Выявление пространственно-временных корреляций среди значений зависимой переменной
3. Применение соответствующего ситуации статистического метода (модели).

Zuur et al. (2010) отмечают, что РДА наиболее эффективен при использовании разнообразных графических средств, поскольку графики часто позволяют лучше понять структуру и свойства анализируемых данных, чем формальные статистические тесты.

Рассмотрение приведенного РДА-протокола начнем с выявления точек-выбросов. Чувствительность разных статистических методов к наличию выбросов в данных неодинакова. Так, при использовании обобщенной линейной модели для анализа зависимой переменной, распределенной по закону Пуассона (например, количество случаев какого-либо заболевания в разных городах), наличие выбросов может вызвать избыточную дисперсию, что сделает модель неприменимой. В то же время при использовании непараметрического многомерного шкалирования, основанного на индексе Жаккара, все исходные данные переводятся в номинальную шкалу с двумя значениями (1/0), и наличие выбросов никак не сказывается на результат анализа. Исследователь должен четко понимать эти различия между разными методами и при необходимости выполнять проверку на наличие выборосов в данных. Дадим рабочее определение: под "выбросом" мы будем понимать наблюдение, которое "слишком" велико или "слишком" мало по сравнению с большинством других имеющихся наблюдений.

Обычно для выявления выбросов используют диаграммы размахов. В R при построении диаграмм размахов используются устойчивые (робастные) оценки центральной тенденции (медиана) и разброса (интерквартильный размах, ИКР). Верхний "ус" простирается от верхней границы "ящика" до наибольшего выборочного значения, находящегося в пределах расстояния 1.5 х ИКР от этой границы. Аналогично, нижний "ус" простирается от нижней границы "ящика" до наименьшего выборочного значения, находящегося в пределах расстояния 1.5 х ИКР от этой границы. Наблюдения, находящиеся за пределами "усов", рассматриваются как потенциальные выбросы (Рисунок 1).

Рисунок 1. Строение диаграммы размахов.

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

• Базовая функция boxplot() (подробнее см. здесь).
• Пакет ggplot2:  геометрический объект ("geom") boxplot. Например:
  

  p <- ggplot(mtcars, aes(factor(cyl), mpg))

  p + geom_boxplot()

  # или:

  qplot(factor(cyl), mpg, data = mtcars, geom = "boxplot") 

Другим очень полезным, но, к сожалению, недостаточно используемым графическим средством выявления выборосов является точечная диаграмма Кливленда. На таком графике по оси ординат откладывают порядковые номера отдельных наблюдений, а по оси абсцисс - значения этих наблюдений. Наблюдения, "значительно" выделяющиеся из основного облака точек, потенциально могут быть выбросами (Рисунок 2).

Рисунок 2. Точечная диаграмма Кливленда, изображающая данные о длине крыла у 1295 воробьев (Zuur et al. 2010). В этом примере данные предварительно были упорядочены в соответствии с весом птиц, и поэтому облако точек имеет примерно S-образную форму.

На Рисунке 2 хорошо выделяется точка, соответствующая длине крыла 68 мм. Однако это значение длины крыла не следует рассматривать в качестве выброса, поскольку оно лишь незначительно отличается от других значений длины. Эта точка выделяется на общем фоне лишь потому, что исходные значения длины крыла были упорядочены по весу птиц. Соответственно, выброс скорее стоит искать среди значений веса (т.е. очень высокое значение длины крыла (68 мм) было отмечено у необычно мало весящего для этого воробья).

Для построения точечных диаграмм Кливленда в R служит функция dotchart(). Я подробно описывал работу с ней в одном из предыдущих сообщений.

До этого момента мы называли "выбросом" наблюдение, которое "значительно" отличается от большинства других наблюдений в исследуемой совокупности. Однако более строгий подход к определению выбросов состоит в оценке того, какое влияние эти необычные наблюдения оказывают на результаты анализа. При этом следует делать различие между необычными наблюдениями для зависимых и независимых переменных (предикторов). Например, при изучении зависимости численности какого-либо биологического вида от температуры большинство значений температуры может лежать в пределах от 15 до 20 °С, и лишь одно значение может оказаться равным 25 °С. Такой план эксперимента, мягко говоря, неидеален, поскольку диапазон температур от 20 до 25 °С будет исследован неравномерно. Однако при проведении реальных полевых исследований возможность выполнить измерения для высокой температуры может представиться только однажды. Что же тогда делать с этим необычным измерением, выполненным при 25 °С? При большом объеме наблюдений подобные редкие наблюдения можно исключить из анализа. Однако при относительно небольшом объеме данных еще большее его уменьшение может быть нежелательным с точки зрения статистической значимости получаемых результатов. Если удаление необычных значений предиктора по тем или иным причинам не представляется возможным, помочь может определенное преобразование этого предиктора (например, логарифмирование).

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

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