--- title: "Занятие 2" author: "Тушавин В. А." date: '15 октября 2016 г ' output: html_document --- ## Цели, задачи и процесс анализа данных За последние несколько десятилетий арсенал прикладной статистики значительно пополнился за счет развития новых методов, краткий список которых включает в себя: линейную регрессию, обобщённые линейные и смешанные модели (generalized linear models), обобщённые аддитивные модели (generalized additive models), модели со смешанными эффектами (mixed-effects models), регрессионные и классификационные деревья (regression and classification trees), нейронные сети (neural networks), анализ выживаемости (survival analysis), многомерный анализ (multivariate analysis) с такими его методами, как: анализ главных компонент (principal component analysis), канонический анализ соответствий (canonical correspondence analysis), различные временные ряды (time series) и т.д. Часть из них появилась достаточно давно, но появление быстрых современных компьютеров и свободного программного обеспечения, такого как GNU R, сделало все эти требующие вычислительных ресурсов методы доступными практически для каждого исследователя. Однако такая доступность ещё больше обострила хорошо известную проблему всех статистических методов, которую на английском языке часто описывают как "rubbish in, rubbish out"", т.е. "мусор на входе - мусор на выходе". Речь здесь идет о следующем: чудес не бывает, и если мы не будем уделять должного внимания тому, как тот или иной метод работает и какие требования предъявляет к анализируемым данным, то получаемые с его помощью результаты нельзя будет воспринимать всерьез. Поэтому каждый раз исследователю следует начинать свою работу с тщательного ознакомления со свойствами полученных данных и проверки необходимых условий применимости соответствующих статистических методов. Этот начальный этап анализа называют разведочным (Exploratory Data Analysis). ## Выявление значений-промахов При обработке реальных экономических, финансовых, технических и иных данных, полученных в процессе наблюдений, измерений, расчётов, иногда один или несколько результатов наблюдений резко выделяются, т.е. далеко отстоят от основной массы данных. В рассматриваемых случаях возникает естественная мысль о том, что подобные наблюдения не относятся к изучаемой совокупности, поскольку содержат грубую погрешность, а получены они в результате ошибки, промаха. В справочнике по метрологии об этом явлении говорится так: "Грубая погрешность (промах) - это случайная погрешность результата отдельного наблюдения, входящего в ряд измерений, которая для данных условий резко отличается от остальных результатов этого ряда. Они, как правило, возникают из-за ошибок или неправильных действий оператора (его психофизиологического состояния, неверного отсчета, ошибок в записях или вычислениях, неправильного включения приборов или сбоев в их работе и др.). Возможной причиной возникновения промахов также могут быть кратковременные резкие изменения условий проведения измерений. Если промахи обнаруживаются в процессе измерений, то результаты, их содержащие, отбрасывают. Однако чаще всего промахи выявляют только при окончательной обработке результатов измерений с помощью специальных критериев" [1, C. 88]. Чувствительность разных статистических методов к наличию аномальных наблюдений в данных неодинакова. Так, при использовании обобщённой линейной модели для анализа зависимой переменной, распределённой по закону Пуассона, наличие промахов может вызвать избыточную дисперсию, что сделает модель неприменимой [2]. В то же время при использовании непараметрического многомерного шкалирования, основанного на индексе Жаккара, все исходные данные переводятся в номинальную шкалу с двумя значениями (1/0), и наличие промахов никак не сказывается на результат анализа [3]. Исследователь должен чётко понимать эти различия между разными методами и при необходимости выполнять проверку на наличие аномальных наблюдений в данных. Обычно для выявления промахов используют диаграммы размахов [4], на которых отображаются медианы и разбросы значений анализируемых совокупностей. Диаграммы размахов, или "ящики с усами"" (box-whisker plots), получили свое название за характерный вид: точку или линию, соответствующую медиане или средней арифметической, окружает прямоугольник ("ящик""), длина которого соответствует одному из показателей разброса или точности оценки генерального параметра. Дополнительно от этого прямоугольника отходят "усы"", также соответствующие по длине одному из показателей разброса или точности. Графики этого типа очень популярны, поскольку позволяют дать очень полную статистическую характеристику анализируемой совокупности. Кроме того, диаграммы размаха можно использовать для визуальной экспресс-оценки разницы между двумя и более группами (например, между датами отбора проб, экспериментальными группами, участками пространства, и т.п.).Обычно медиана в этом графике представлена в виде горизонтальной линии. Эта линия окаймляется прямоугольником, длина которого соответствует интерквартильному размаху (ИКР). От прямоугольника отходят линии, ограничивающие интервал длиной 1,5ИКР. Наблюдения, находящиеся за пределами интервала 1,5ИКР, рассматриваются как потенциальные промахи. Однако всегда следует внимательно относиться к такого рода нестандартным наблюдениям. поскольку поскольку эти "промахи"" вполне могут оказаться в порядке вещей для исследуемой совокупности. Так в нормальном распределении шансы случайной величины превысить два стандартных отклонения составляют 1 к 21, а шансы выйти за три стандартных отклонения составляют 1 к 369. В R для построения диаграмм размаха используется функция `boxplot`. Подключим стандартные данные, имеющиеся в R: ```{r} data("cars") head(cars) ``` В данном случае мы имеем дело с данными зависимости тормозного пути автомобиля от скорости. Подробно сведени можно получить с помощью команды `?cars` ``` cars {datasets} R Documentation Speed and Stopping Distances of Cars Description The data give the speed of cars and the distances taken to stop. Note that the data were recorded in the 1920s. Usage cars Format A data frame with 50 observations on 2 variables. [,1] speed numeric Speed (mph) [,2] dist numeric Stopping distance (ft) Source Ezekiel, M. (1930) Methods of Correlation Analysis. Wiley. References McNeil, D. R. (1977) Interactive Data Analysis. Wiley. ``` Построим диаграммы для обоих колонок значений ```{r} boxplot(cars$speed) boxplot(cars$dist) ``` Обратите внимание, что на второй диаграмме явно виден потенциальный выброс. Другим очень полезным, но, к сожалению, недостаточно используемым графическим средством выявления промахов является точечная диаграмма Кливленда. На таком графике по оси ординат откладывают порядковые номера отдельных наблюдений, а по оси абсцисс - значения этих наблюдений. Наблюдения, "значительно"" выделяющиеся из основного облака точек, потенциально могут быть промахами. Для построения точечных диаграмм Кливленда в R служит функция `dotchart`. ```{r} dotchart(cars$dist) points(cars[49,2],49,pch=19,col="red") ``` Обратите внимание на точку в правом верхнем углу, специально покрашенную второй командой в красный цвет. Возможно, что это выброс. Другим не менее полезным инструментом визуализации яляется диаграмма типа «скрипка» (violin plot), совмещенная с классической диаграммой размахов (box-and-whisker plot). Эта комбинированная графическая форма представления дает больше информации о характере распределения, чем диаграмма размахов или гистограмма, т.к. помимо данных о медиане и квартилях, отражает еще и показатели ядерной плотности распределения [5]. К сожалению, выбросы она не показывает. ```{r} library(vioplot) vioplot(cars$dist) ``` Поэтому можно использовать другой график, `beanplot` из пакета `beanplot` который следует предварительно установить из меню RStudio Tools. ```{r} library(beanplot) beanplot(cars$dist,cars$speed) ``` Подробнее об этом графике смотри https://cran.r-project.org/web/packages/beanplot/vignettes/beanplot.pdf До этого момента мы называли "промахом"" или "выбросом" наблюдение, которое существенно отличается от большинства других наблюдений в исследуемой совокупности. Однако более строгий подход к определению промахов состоит в оценке того, какое влияние эти необычные наблюдения оказывают на результаты анализа. При этом следует делать различие между необычными наблюдениями для зависимых и независимых переменных (предикторов). Решение отбросить какие-либо данные в конечном счете всегда субъективно. Наиболее строгий подход к проблеме анализа подобных аномальных данных - повторить измерения много раз. Если подобная аномалия повторится снова, вероятно, можно будет выяснить её причину. Если же аномальное показание не повторится, то даже оставляя его, оно мало скажется на дальнейших вычислениях из-за большого объёма измерений. Однако при небольшом объёме опытных данных и отсутствии явных причин появления грубых ошибок целесообразно использовать статистические методы выявления промахов [6, c.245]. Рассмотрим несколько статистических тестов на наличие аномальных наблюдений, имеющихся в R. ### Q тест Диксона (Dixon's Q test) применим для малых выборок (меньше 30 наблюдений). Он применим только для данных, распределённых нормально. Если тест выявляет выбросы, то после их устранения, повторно этот тест использовать нельзя [7]. Для его проведения используют функцию `dixon.test` из пакета `outliers` [8]. Нулевая гипотеза для этого и последующего тестов на наличие аномальных значений заключается в том, что наименьшее или наибольшее отклонение от среднего в выборке не является промахом. Поскольку у нас больше 30 значений, возьмём для примера другую случайную выборку и проведем этот тест. ```{r} library(outliers) sample<-c(28.28,29.10,30.95,31.49,32.49,32.79,32.79,31.99,31.99) dixon.test(sample) ``` Как мы видим, наименьшее значение 28.28 не является выбросом (p-value > 0.05). ### Тест хи-квадрат для промахов Другой тест, тест хи-квадрат для промахов, требует знаний о дисперсии генеральной совокупности. В противном случае его использовать нельзя. Тест применим для тестирования остатков моделей, входит в пакет `outliers` и вызывается командой `chisq.out.test`. ```{r} chisq.out.test(cars$dist) ``` ### Тест Граббса Для проверки наличия промахов в нормально распределённых данных также используется тест Граббса (Grubbs test). Тест также входит в пакет `outliers` и вызывается командой `grubbs.test`. Вызов команды имеет вид `grubbs.test( sample, type = 10, opposite = FALSE, two.sided = TRUE)`. Параметр `sample` должен содержать исследуемые данные, параметр `type` принимает одно из трех значений: 10 - для одного промаха, 11 - для двух промахов с противоположных концов распределения, 20 - для двух промахов с одной стороны распределения. ```{r} grubbs.test(cars$dist) ``` Если удаление необычных значений предиктора по тем или иным причинам не представляется возможным, помочь может определенное преобразование этого предиктора (например, логарифмирование). С необычными значениями зависимой переменной бороться сложнее, особенно при построении регрессионных моделей. Преобразование путём, например, логарифмирования, может помочь, но поскольку зависимая переменная представляет особый интерес при построении регрессионных моделей, лучше попытаться подобрать метод анализа, который основан на распределении вероятностей, допускающем больший разброс значений для больших средних величин (например, гамма-распределение для непрерывных переменных или распределение Пуассона для дискретных количественных переменных). Такой подход позволит работать с исходными значениями зависимой переменной. В конечном счёте решение об удалении необычных значений из анализа принимает сам исследователь. При этом он должен помнить о том, что причины для возникновения таких наблюдений могут быть разными. Так, удаление промахов, возникших из-за неудачного планирования эксперимента, может быть вполне оправданным. Оправданным будет также удаление промахов, явно возникших из-за ошибок при выполнении измерений. В то же время необычные наблюдения среди значений зависимой переменной могут потребовать более тонкого подхода, особенно если они отражают естественную вариабельность этой переменной. В этой связи важно вести подробное документирование условий, при которых происходит экспериментальная часть исследования - это может помочь интерпретировать промахи в ходе анализа данных. Независимо от причин возникновения необычных наблюдений, в итоговом научном отчёте (например, в статье) важно сообщить читателю как о самом факте выявления таких наблюдений, так и о принятых в их отношении мерах. ## Экспресс-обзор регрессионный анализ Поскольку windows XP продолжло радовать мы перешли к регрессии, немного забежав вперед. Выделим следующие этапы построения эконометрической модели для срезов данных: 1. Спецификация модели: * определение цели и назначения модели (зависимой переменной Y) * определение потенциальных объясняющих (независимых) переменных * предварительное сокращение потенциальных переменных * выбор аналитической формы модели * формулирование гипотезы моделирования 2. Оценивание структурных параметров модели. 3. Верификация эконометрической модели * оценивание существенности влияния конкретных объясняющих переменных на зависимую переменную, t𝑢-тест Стьюдента — исключение a posteriori, F-тест Снедекера; * оценивание степени соответствия модели эмпирическим данным; * оценивание нормальности распределения остатков; * оценивание однородности дисперсии остатков — тест на гетероскедатичность; * оценивание линейности аналитической формы модели. 4. Формулирование выводов и их интерпретация. В R линейная регрессионная модель классическим методом наименьшиъ квадратов строится ледующим образом: ```{r} lm1<-lm(dist~speed,data=cars) summary(lm1) ``` Интерпретация результатов, полученных выше с помощью команды summary, осуществляется по следующему алгоритму. Рассмотрим полученные результаты по порядку сверху вниз, за исключением наиболее важной статистики, F-статистики, которая находится в самом низу. ``` F-statistic: 89.57 on 1 and 48 DF, p-value: 1.49e-12 ``` F-статистика вычисляется по формуле \(F=\frac{\frac{Q_R}{k-1}}{\frac{Q_e}{n-k}}\), далее результат сравнивается с квантилью распределения Фишера с данными (1 и 48) степенями свободы и вычисляется вероятность того, что F превышает соответсвующий квантиль распределения Фишера. При необходимость, тот же результат можно получить используя команду `pf` c параметром `lower.tail=F`. Команда `pf` возвращает значение функции распределения Фишера, т.е. вероятность того, что значения функции не превосходят заданную величину. Для решения противоположной задачи используется параметр `lower.tail=T`. ```{r} pf(89.57,1,48,lower.tail=F) ``` Полученное p-значение *p-value* показывает вероятность того, что модель является незначимой. Иными словами, если оно получилось меньше чем 0,05, то это означает, что наиболее вероятен факт, что модель значима. Если же, напротив, она превышает 0,05, то модель, скорее всего, назначаемая. В приведённом выше примере `p-value: < 1.49e-12`, откуда следует, что модель, скорее всего, является значимой. Следует отметить, что большинство людей в первую очередь смотрят на \(R^2\), считая его самой важной статистикой, однако, если модель является незначимой, то все остальные её показатели, какими бы хорошими они не были, не имеют никакого значения. Самые первые строки, которые выдаёт команда `summary`,содержат в себе описание модели, которая была построена с помощью команды `lm`. ``` Call: lm(formula = dist ~ speed, data = cars) ``` Следующие за этим строки содержат данные об остатках модели (Residuals). Как известно, классический метод наименьших квадратов (КМНК) предполагает, что математическое ожидание остатков модели равно нулю. В идеальной ситуации, остатки модели должны быть распределены нормально. Как мы видим из следующих строк ``` Residuals: Min 1Q Median 3Q Max -29.069 -9.525 -2.272 9.215 43.201 ``` это условие не соблюдается. Медиана (Median) смещена влево относительно нуля, первый квартиль (1Q) не совсем симметричен относительно третьего (3Q), с незначительным смещением влево, в то время как для нормального распределения характерна симметрия. Максимумы (Max) и минимумы (Min) же позволяют достаточно быстро оценить наличие потенциальных промахов в исходных данных, поскольку промахи в исходных данных порождают промахи в остатках модели. Следующий блок содержит информацию о коэффициентах модели и их значимости. ``` Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) -17.5791 6.7584 -2.601 0.0123 * speed 3.9324 0.4155 9.464 1.49e-12 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 ``` Колонка **Estimate** содержит подобранные коэффициенты регрессии \(\tilde{\beta}_j\), колонка **Std. Error** содержит их среднеквадратичные отклонения \(s_j\). Колонка **t value** рассчитывается как их отношение (без учёта модуля, который будет учтен далее). В колонке **Pr(>|t|)** показана вероятность того, что данный коэффициент является незначимым. Для удобства интерпретации значимость коэффициентов также отображается символами справа в соответствии с легендой под чертой. Дальше идут знакомые вам Эр-квардрат и скорретированный Эр-квадрат. Поговорим о статистике R-квадрат или, как её иногда называют, коэффициента детерминации. Она показывает, насколько условная дисперсия модели отличается от дисперсии реальных значений Y. Более точно — это единица минус доля необъяснённой дисперсии (дисперсии случайной ошибки модели, или условной по факторам дисперсии зависимой переменной) в дисперсии зависимой переменной. Его рассматривают как универсальную меру зависимости одной случайной величины от множества других. В частном случае линейной зависимости \(R^2\) является квадратом так называемого множественного коэффициента корреляции между зависимой переменной и объясняющими переменными. В частности, для модели парной линейной регрессии коэффициент детерминации равен квадрату обычного коэффициента корреляции между y и x. Если этот коэффициент близок к 1, то условная дисперсия модели достаточно мала и весьма вероятно, что модель неплохо описывает данные.Для приемлемых моделей предполагается, что коэффициент детерминации должен быть хотя бы не меньше 50% (в этом случае коэффициент множественной корреляции превышает по модулю 70%). Модели с коэффициентом детерминации выше 80% можно признать достаточно хорошими (коэффициент корреляции превышает 90%). Значение коэффициента детерминации 1 означает функциональную зависимость между переменными. Однако, у статистики R-квадрат есть один серьезный недостаток: при увеличении числа предикторов эта статистика может только возрастать. Поэтому, может показаться, что модель с большим количеством предикторов лучше, чем модель с меньшим, даже если все новые предикторы никак не влияют на зависимую переменную. Тут можно вспомнить про принцип бритвы Оккама. Следуя ему, по возможности, стоит избавляться от лишних предикторов в модели, поскольку она становится более простой и понятной. Для этих целей была придумана статистика скорректированный R-квадрат (скорректированный коэффициент детерминации). Она представляет собой обычный R-квадрат, но со штрафом за большое количество предикторов. Основная идея: если новые независимые переменные дают большой вклад в качество модели, значение этой статистики растет, если нет — то наоборот уменьшается. Для построения диагностических графиков остатков используется команда `plot`. ```{r} plot(lm1) ``` Для построения доверительного интервала коэффициентов команда `confint` ```{r} confint(lm1) ``` В качестве иллюстрации команд построим модель без свободного члена и выведем параметры модели, причем сделаем это одной командой ```{r} summary(lm2<-lm(dist~speed-1,data=cars)) ``` И регрессионную зависимость от квадрата скорости ```{r} summary(lm3<-lm(dist~I(speed^2)-1,data=cars)) plot(lm3) ``` Устраним 49 измерение,содержащее выброс и построим регрессию заново. ```{r} cars2<-cars[-49,] summary(lm4<-lm(dist~I(speed^2)-1,data=cars2)) plot(lm4) ``` ## Список источников 1. Сергеев А. Г., Крохин В. В. Метрология: Учеб. пособие для вузов.— М.: Логос, 2001. — 408 с. — ISBN: 5-94010-039-2. 2. Hilbe J. M. Negative Binomial Regression. – 1st. edition. – Cambridge Universit Press, 2007. – 251 p. – ISBN: 978-0521857727. 3.Legendre P., Legendre L. F. J. Numerical Ecology (Developments in EnvironmentalbModelling) [Kindle Edition]. – Elsevier Science, 1998. – 870 p. – ASIN : B0092L8KFA. 4. Graphical Methods for Data Analysis (Wadsworth & Brooks/Cole Statistics/Probability Series) / John M. Chambers, William S. Cleveland, Paul A. Tukey, Beat Kleiner. – Duxbury Press, 1983. – 395 p. – ISBN: 978-0534980528. 5. Hintze J.L., Nelson R.D. Violin Plots: A Box Plot-Density Trace Synergism / Jerry L. Hintze, Ray D. Nelson // The American Statistician. – 1998. – Vol. 52. – P. 181-184 6. Куликов Е. И. Прикладной статитический анализ. Учебное пособие для вузов. — 2-е. изд. — М. : Горячая линия–Телеком, 2008. — 464 с. — ISBN: 978-5-99120021-9. 7.Lewis N. D. 100 Statistical Tests in R [Kindle Edition]. – Heather Hills Press, 2013. – 496 p. – ASIN : B00C143F14. 8.Komsta L. – outliers: Tests for outliers, 2011. – R package version 0.14. URL: http://CRAN.R-project.org/package=outliers