350-anova.qmd

# Дисперсионный анализ (ANOVA) {#sec-anova}

Дисперсионный анализ или ANOVA[^350-anova-1] - один из самых распространенных методов статистического анализа в психологии, биологии, медицине и многих других дисциплинах. Дисперсионный анализ очень хорошо подходит для анализа данных, полученных в эксперименте - методе организации исследования, при котором исследователь напрямую управляет уровнями независимой переменной. Терминологическая связь между дисперсионным анализом и планированием эксперимента настолько тесная, что многие термины пересекаются, поэтому нужно быть осторожными. Как и в случае с линейной регрессией, если мы что-то называем "независимой переменной" (или "фактором"), это не порождает никакой каузальной связи.

[^350-anova-1]: ANOVA от ANalysis Of VAriance, по-русски часто читается как "АНОВА".

> Еще одна важная вещь, которую нужно понимать про дисперсионный анализ, это то, что у этого метода очень запутывающее название: из названия кажется, что этот статистический метод для сравнения дисперсий. Нет, это не так (хотя такие статистические тесты тоже есть, и они нам сегодня пригодятся - см. \@ref(aova)). Нет, это просто сравнение средних в случае, если есть больше, чем 2 группы для сравнения.

У дисперсионного анализа очень много разновидностей, для которых придумали множество названий. "Обычная" ANOVA называется One-Way ANOVA, она же межгрупповая ANOVA, это аналог независимого т-теста для нескольких групп.

Давайте начнем сразу с проведения теста. Мы будем использовать [данные с курса по статистике Университета Шеффилда про эффективность диет](https://www.sheffield.ac.uk/polopoly_fs/1.570199!/file/stcp-Rdataset-Diet.csv), на которых мы разбирались с t-тестом (\@ref(dep_ttest)).

```{r}
library(tidyverse)
diet <- read_csv("data/stcp-Rdataset-Diet.csv")
```

Сделаем небольшой препроцессинг данных. Создадим дополнительные "факторные" переменные, создадим переменную, в которой будет разница массы "до" и "после", удалим `NA`.

```{r}
diet <- diet %>%
  mutate(weight.loss = weight6weeks - pre.weight,
         Dietf = factor(Diet, labels = LETTERS[1:3]),
         Person = factor(Person)) %>%
  drop_na()
```

### Функция aov() {#sec-aov}

Попробуем сразу провести дисперсионных анализ с помощью функции `aov()`:

```{r}
aov_model <- aov(weight.loss ~ Dietf, diet)
aov_model
summary(aov_model)
```

Мы получили что-то похожее на результат применения функции `lm()`. Правда, лаконичнее, но с новыми столбцами `Sum Sq`, `Mean Sq` и новой статистикой *F* вместо *t*. Что будет, если с теми же данными с той же формулой запустить `lm()` вместо `aov()`?

```{r}
summary(lm(weight.loss ~ Dietf, diet))
```

`lm()` превратил `Dietf` в две переменные, но *F* и p-value у двух моделей одинаковые! Кроме того, функция `aov()` является, по сути, просто "оберткой" над `lm()`:

> This provides a wrapper to lm for fitting linear models to balanced or unbalanced experimental designs.

## Тестирование значимости нулевой гипотезы в ANOVA. {#sec-anova_nhst}

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

1.  **Формулирование нулевой и альтернативной гипотезы.** Нулевая гипотеза говорит, что между средними в генеральной совокупности нет различий:

$$H_0:\mu_1 = \mu_2 = ... = \mu_n$$ Можно было бы предположить, что ненулевая гипотеза звучит как "все средние не равны", но вообще-то это не так. Альтернативная гипотеза в дисперсионном анализе звучит так:

$$H_1: \text{Не все средние равны}$$

2.  **Подсчет статистики.** Как мы уже видели раньше, в дисперсионном анализе используется новая для нас статистика *F*. Впрочем, мы ее видели, когда смотрели на аутпут функции `lm()`, когда делали линейную регрессию. Чтобы считать *F* (если вдруг мы хотим сделать это вручную), нужно построить талбицу ANOVA (ANOVA table).

| Таблица ANOVA   | Степени свободы |      Суммы квадратов      |        Средние квадраты         |       F-статистика        |
|:----------------|:---------------:|:-------------------------:|:-------------------------------:|:-------------------------:|
| Межгрупповые    |    $df_{b}$     |         $SS_{b}$          | $MS_{b} =\frac{SS_{b}}{df_{b}}$ | $F=\frac{MS_{b}}{MS_{w}}$ |
| Внутригрупповые |    $df_{w}$     |         $SS_{w}$          | $MS_{w} =\frac{SS_{w}}{df_{w}}$ |                           |
| Общие           |    $df_{t}$     | $SS_{t}= SS_{b} + SS_{w}$ |                                 |                           |

Именно эту таблицу мы видели, когда использовали функцию `aov()`:

```{r}
summary(aov_model)
```

Вот как это все считается:

| Таблица ANOVA | Степени свободы |                                    Суммы квадратов                                    |        Средние квадраты         |       F-статистика        |
|:--------------|:---------------:|:-------------------------------------------------------------------------------------:|:-------------------------------:|:-------------------------:|
| Между         |  $df_{b}=J-1$   | $SS_{b}= \sum\limits_{j=1}^J \sum\limits_{i=1}^{n_j} (\overline{x_j}-\overline{x})^2$ | $MS_{b} =\frac{SS_{b}}{df_{b}}$ | $F=\frac{MS_{b}}{MS_{w}}$ |
| Внутри        |  $df_{w}=N-J$   |    $SS_{w}= \sum\limits_{j=1}^J \sum\limits_{i=1}^{n_j} (x_{ij}-\overline{x_j})^2$    | $MS_{w} =\frac{SS_{w}}{df_{w}}$ |                           |
| Общие         |  $df_{t}=N-1$   |     $SS_{t}= \sum\limits_{j=1}^J \sum\limits_{i=1}^{n_j} (x_{ij}-\overline{x})^2$     |                                 |                           |

$J$ означает количество групп, $N$ - общее количество наблюдений во всех группах, $n_j$ означает количество наблюдений в группе j, а $x_{ij}$ - наблюдение под номером $i$ в группе $j$.

Вариабельность обозначается $SS$ и означает "сумму квадратов" (sum of squares) - это то же, что и дисперсия, только мы не делим вме в конце на количество наблюдений (или количество наблюдений минус один): $$SS = \sum\limits_{i=1}^{n_j} (x_{i}-\overline{x})^2$$

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

3.  **Подсчет p-value.** В t-тесте мы смотрели, как статистика распределена при условии верности нулевой гипотезы. То есть что будет, если нулевая гипотеза верна, мы будем повторять эксперимент с точно таким же дизайном (и размером выборок) бесконечное количество раз и считать *F*.

```{r, fig.cap="\\label{fig:fig1}F-распределение при верности нулевой гипотезы (см. детали в тексте)", out.width = '60%'}

betweendf <- 2
withindf <- 73
f <- summary(aov_model)[[1]]$F[1]

v <- seq(0.1,10, 0.01)
fdist <- data.frame(fvalues = v, pdf = df(v, betweendf, withindf))

library(ggplot2)

label <- paste0("F(", betweendf, ", ", withindf, ") = ", round(f, 3))

ggplot(fdist, aes(x = fvalues, y = pdf))+
  geom_line()+
  geom_vline(xintercept = f)+
  annotate("text", x = f+1, y = 0.2, label = label)+
  scale_y_continuous(expand=c(0,0)) + 
  theme_minimal()+
theme(axis.line.y = element_blank(),
      axis.ticks.y = element_blank(),
      axis.text.y = element_blank(),
      axis.title.y = element_blank())  
```

Заметьте, распределение *F* несимметричное[^350-anova-2]. Это значит, что мы всегда считаем считаем площадь от *F* до плюс бесконечности (без умножения на 2, как мы это делали в t-тесте):

[^350-anova-2]: Форма *F*-распределения будет сильно зависеть от числа степеней свободы. Но оно всегда определено от 0 до плюс бесконечности: в числителе и знаменателе всегда неотрицательные числа.

```{r}
1 - pf(f, betweendf, withindf)
```

Это и есть наш p-value!

4.  **Сравнение p-value с уровнем** $\alpha$. Самый простой этап: если наш p-value меньше, чем $\alpha$ (который обычно равен 0.05), то мы отвергаем нулевую гипотезу. Если нет - не отвергаем.

В нашем случае это `r 1 - pf(f, betweendf, withindf)`, что, очевидно, меньше, чем 0.05. Отвергаем нулевую гипотезу (о том, что нет различий), принимаем ненулевую (о том, что различия есть). Все!

## Post-hoc тесты {#sec-anova_posthoc}

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

Post-hoc переводится с латыни как "после этого". Post-hoc тесты или просто "пост-хоки" проводятся, если в результате ANOVA была отвергнута нулевая гипотеза. Собственно, пост-хоки никак не связаны с дисперсионным анализом на уровне расчетов - это абсолютно независимые тесты, но исторически так сложилось, что они известны именно как дополнительный этап ANOVA.

Самый простой вариант пост-хок теста - это попарные т-тесты [^350-anova-3] с поправками на множественные сравнения:

[^350-anova-3]: По умолчанию функция `pairwise.t.test()` использует объединенное стандартное отклонение для всех условий при расчетах. Это дает немного другие результаты, чем просто попарные t-тесты, хотя разница обычно незначительная. Чтобы отключить такое поведение функции `pairwise.t.test()`, можно поставить `FALSE` для параметра `pool.sd()`

```{r}
pairwise.t.test(diet$weight.loss, diet$Dietf)
```

Второй подход связан с использованием специализированных тестов, таких как тест Тьюки (Tukey Honest Significant Differences = Tukey HSD). Для этого в R есть функция `TukeyHSD()`, которую нужно применять на объект `aov`:

```{r}
TukeyHSD(aov_model)
```

## ANOVA и т-тест как частные случаи линейной регрессии {#sec-aov_as_lm}

Как мы уже видели, если применить `lm()` или `aov()` на одних и тех же данных с одной и той же формулой, то результат будет очень похожим. Но есть одно но: `lm()` создает из одного фактора две переменных-предиктора:

```{r}
summary(lm(weight.loss ~ Dietf, diet))
```

Дело в том, что мы не можем просто так загнать номинативную переменную в качестве предиктора в линейную регрессию. Мы можем это легко сделать, если у нас всего два уровня в номинативном предикторе. Тогда один из уровней можно обозначить за 0, другой - за 1. Такие переменные иногда называются "бинарными". Тогда это легко использовать в линейной регрессии:

```{r}
summary(lm(weight.loss ~ gender, diet))
```

Можно ли так делать? Вполне! Допущения линейной регрессии касаются остатков, а не переменных самих по себе. Разве что это немного избыточно: линейная регрессия с бинарным предиктором - это фактически независимый t-тест:

```{r}
t.test(weight.loss ~ gender, diet, var.equal = TRUE)
```

Как видите, p-value совпадают! А *t* статистика в квадрате - это *F* (при двух группах):

```{r}
t.test(weight.loss ~ gender, diet, var.equal = TRUE)$statistic^2
```

Более того, те же самые результаты можно получить и с помощью коэффициента корреляции Пирсона:

```{r}
cor.test(diet$gender, diet$weight.loss)
```

Теперь должно быть понятно, почему все эти функции делают вроде бы разные статистические тесты, но выдают такой похожий результат - это фактически один и тот же метод! Все эти методы (и некоторые из тех, что будем рассматривать далее) можно рассматривать как разновидности множественной линейной регрессии. [^350-anova-4]

[^350-anova-4]: Обобщением множественной линейной регрессии (вернее, одним из) можно считать общую линейную модель (general linear model). Общая линейная модель может предсказывать не одну, а сразу несколько объясняемых переменных в отличие от множественной линейной регрессии. Следующим этапом обобщения служит обобщенная линейная модель (generalized linear model). Фишка последней в том, что можно использовать не только модели с нормально распределенными остатками, но и, например, логистическую и пуассоновскую регрессию.

## Dummy coding {#sec-dummy}

Тем не менее, вопрос остается открытым: как превратить номинативную переменную в количественную и загнать ее в регрессию? Для этого можно использовать **"фиктивное кодирование" (dummy coding):**

```{r}
diet <- diet %>%
  mutate(isA = as.numeric(Dietf == "A"),
         isB = as.numeric(Dietf == "B"),
         isC = as.numeric(Dietf == "C"))

diet %>%
  group_by(Dietf) %>%
  slice(1:2) %>%
  select(Dietf, isA:isC)
```

Заметьте, что такое кодирование избыточно. Если мы знаем, что диет 3, а данная диета - это не диета В и не диета С, то это диета А. Значит, одна из созданных нами колонок - "лишняя":

```{r}
diet$isA <- NULL
```

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

```{r}
summary(lm(weight.loss ~ isB + isC, diet))
summary(lm(weight.loss ~ Dietf, diet))
```

То же самое!

## Допущения ANOVA {#sec-aova}

1.  **Нормальность распределения ошибок:**

```{r}
hist(residuals(aov_model))
```

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

2.  **Гомогенность дисперсий.**

То есть их равенство. Можно посмотреть на распределение остатков:

```{r}
diet$residuals <- residuals(aov_model)
ggplot(diet, aes(x = Dietf, y = residuals))+ geom_jitter(width = 0.1, alpha = 0.5)
```

Все выглядит неплохо: нет какой-то одной группы, у которой разброс сильно больше или меньше. Есть и более формальные способы проверить равенство дисперсий. Например, с помощью теста Ливиня (Levene's test). Для того, чтобы его провести, мы воспользуемся новым пакетом `ez` (читать как "easy"). Этот пакет сильно упрощает проведение дисперсионного анализа, особенно для более сложных дизайнов.

```{r, eval = FALSE}
install.packages("ez")
```

Синтаксис довольно простой: нужно указать, данные, зависимую переменную, переменную с ID, факторы. Необходимо прописать фактор в `between =` или `within =`. В данном случае - в `between =`.

```{r}
library(ez)
ez_model <- ezANOVA(data = diet,
        dv= weight.loss,
        wid = Person, 
        between = Dietf,
        detailed = T, 
        return_aov = T)
ez_model
```

Если при проведении теста Ливиня мы получаем *p \< .05*, то мы отбрасываем нулевую гипотезу о равенстве дисперсий. В данном случае мы не можем ее отбросить и поэтому принимаем [^350-anova-5]

[^350-anova-5]: Вообще-то эта логика не совсем корректна. Тест Ливиня - это такой же статистический тест, как и остальные. Поэтому считать, что допущения соблюдаются на основании того, что p-value больше допустимого уровня $\alpha$, - это неправильно. Но для проверки допущений такая не очень корректная практика считается допустимой.

Полученный объект (если поставить `return_aov = T`) содержит еще и объект `aov()` - на случай, если у Вас есть функции, которые работают с этим классом:

```{r}
TukeyHSD(ez_model$aov)
```

3.  **Примерно одинаковое количество испытуемых в разных группах.** Здесь у нас все в порядке:

```{r}
diet %>%
  count(Dietf)
```

Небольшие различия в размерах групп - это ОК, тем более, что на практике такое очень часто случается: кого-то пришлось выкинуть из анализа, для какой-то строчки были потеряны данные и т.д. Однако больших различий в размерах групп стоит избегать. Самое плохое, когда группы различаются значительно по размеру (более чем в 2 раза) и вариабельность внутри групп отличается значительно (более чем в 2 раза).

## Многофакторный дисперсионный анализ (Factorial ANOVA) {#sec-fact_aov}

На практике можно встретить One-Way ANOVA (однофакторную ANOVA) довольно редко. Обычно в исследованиях встречается многофакторный дисперсионный анализ, в котором проверяется влияние сразу нескольких факторов. В научных статьях это обозначается примерно так: "3х2 ANOVA". Это означает, что был проведен двухфакторный дисперсионный анализ, причем в одном факторе было три уровня, во втором - два. В нашем случае это будут факторы "Диета" и "Пол". Это означает, что у нас две гипотезы: о влиянии диеты на потерю веса и о влиянии пола на потерю веса. Кроме того, появляется гипотеза о взаимодействии факторов - то есть о том, что разные диеты по разному влияют на потерю веса для разных полов.

Взаимодействие двух факторов хорошо видно на графике с линиями: если две линии параллельны, то взаимодействия нет. Если они не параллельны (пересекаются, сходятся, расходятся), то взаимодействие есть.

```{r}
diet <- diet %>%
  mutate(genderf = factor(gender, labels = c("ж", "м"))) #превращаем в бинарную переменную в фактор

sem <- function(x) sd(x)/sqrt(length(x)) #пишем функцию для стандартной ошибки
pd = position_dodge(0.05) #немного раздвигаем положение точек на будущем графике

diet %>%
  group_by(Dietf, genderf) %>%
  summarise(meanloss = mean(weight.loss),
            se = sem(weight.loss)) %>%
  ggplot(aes(x = Dietf, 
             y = meanloss, 
             colour = genderf)) +
  geom_line(aes(group = genderf), position = pd) +
  geom_pointrange(aes(ymin = meanloss - se, 
                      ymax = meanloss + se), position = pd) +
  theme_minimal()
```

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

```{r}
ezANOVA(data = diet,
        dv= weight.loss,
        wid = Person, 
        between = .(Dietf, gender),
        detailed = T, 
        return_aov = T)
```

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

## Дисперсионный анализ с повторными измерениями (Repeated-measures ANOVA) {#sec-rm_aov}

Если обычный дисперсионный анализ - это аналог независимого т-теста для нескольких групп, то дисперсионный анализ с повторными измерениями - это аналог зависимого т-теста. В функции `ezANOVA()` для проведения дисперсионного анализа с повторными измерениями нужно просто поставить нужным параметром внутригрупповую переменную. Это означает, что в данном случае мы должны иметь данные в длинном формате, для чего мы воспользуемся функцией `pivot_longer()`:

```{r}
dietlong <- diet %>%
  pivot_longer(cols = c(pre.weight, weight6weeks),
               names_to = "time",
               values_to = "weight")
dietlongC <- dietlong %>%
          filter(Dietf == "C") %>%
          droplevels()
```

```{r}
ezANOVA(dietlongC,
        dv = weight, 
        wid = Person,
        within = time)
```

## Смешанный дисперсионный анализ (Mixed between-within subjects ANOVA) {#sec-mixed_aov}

Нам никто не мешает совмещать и внутригруппоdые, и межгрупповые факторы вместе. В `ezANOVA()` это делается просто с помощью прописывания разных факторов в нужных переменных: `between =` и `within =`.

```{r}
ezANOVA(dietlong,
        dv = weight, 
        wid = Person,
        within = time,
        between = Dietf)
```

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

## Непараметрические аналоги ANOVA {#sec-nonparam_aov}

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

### Тест Краскела-Уоллеса {#sec-krask_aov}

Это тест Краскела-Уоллеса - обобщение теста Манна-Уитни на несколько выборок (т.е. аналог межгруппового ANOVA):

```{r}
kruskal.test(weight.loss ~ Dietf, diet)
```

### Тест Фридмана {#sec-friedman_aov}

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

```{r}
friedman.test(weight ~ time | Person, dietlongC)
```

## Заключение {#sec-aov_final}

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

-   многомерный ANOVA (Multivariate ANalysis Of Variance; MANOVA) - расширение ANOVA для ситуации нескольких зависимых переменных. Это довольно редкая разновидность ANOVA, который берет во внимание ковариации между зависимыми переменными.

-   тестирование сферичности для дисперсионного анализа с повторноми измерениями с помощью теста сферичности Моучли (Mauchly's sphericity test). Этот тест проверяет использует матрицу ковариаций разниц каждого условия с каждым: дисперсии разниц между условиями должны быть примерно одинаковыми. Если нулевая гипотеза о сферичности может быть отброшена (p-value \< $\alpha$), то нужно проводить специальные поправки, обычно это поправки Гринхауса-Гейсера (Greenhouse-Geisser corrections). Мы этого не делали, потому что в ситуации RM-ANOVA с всего двумя условиями эта сферичность никогда не нарушается: у нас всего одна дисперсия разниц между условиями, которую просто-напросто не с чем сравнивать. Тест сферичности Моучли вместе с поправками Гринхауса-Гейсера проводятся автоматически для RM-ANOVA с тремя или более группами при использовании функции `ezANOVA()`, так что особо париться над этим не стоит. Правда, нужно помнить, что как и все подобные статистические тесты допущений, они обладают проблемами, связанными с тем, что это статистические тесты: на маленьких выборках они не заметят даже серьезных отклонений от сферичности, на больших - даже маленькие отклонения, а $p-value < 0.05$, по сути, не может интерпретироваться как верность нулевой гипотезы. Тем не менее, это довольно стандартная процедура.

-   Как правильно репортить результаты дисперсионного анализа. Здесь все, конечно, зависит от стиля, используемого конкретным журналом. В психологии и близких к ней дисциплинам фактическим lingua franca является стиль Американской Психологической Ассоциации (APA). И тут у меня есть для Вас хорошие новости: есть куча пакетов в R, которые позволяют репортить результаты статистических тестов в APA-стиле! Спасибо дотошным авторам руководства APA по офромлению статей, что этот стиль настолько точно прописывает, как нужно описывать результаты исследований, что это можно запрограммировать. Я лично пользуюсь пакетом `apa`, он весьма удобен:

```{r, eval = FALSE}
install.packages("apa")
```

```{r, message = FALSE}
library(apa)
anova_apa(ez_model)
```

В тексте это будет выглядеть это будет вот так:

Dietf: *F*(2, 73) = 5.38, *p* = .007, $\eta^2_p$ = .13

Еще есть пакеты `apaStyle` и `papaja`, которые могут даже сразу делать весь документ в APA-формате! Если же Вы описываете результаты самостоятельно вручную, то нужно помнить: ни в коем случае не описывайте только p-value. Обязательно прописывайте значение $F$ и степени свободы, желательно с размером эффекта. Для post-hoc теста часто репортятся только p-value (зачастую только для статистически значимых сравнений), но обязательно нужно прописывать какие именно post-hoc тесты проводились, какой показатель размера эффекта использовался (если использовался), применялись ли тест сферичности Моучли вместе с поправками Гринхауса-Гейсера для дисперсионного анализа с повторными измерениями.

-   Модели со смешанными эффектами (mixed-effects models) / иерархическая регрессия (hierarchical regression) / многоуровневое моделирование (multilevel modelling). очень популярный нынче метод, которому повезло иметь много названий - в зависимости от области, в которой он используется. В экспериментальной психологии обычно он называется "модели со смешанными эффектами" и позволяет включать в линейную регрессию не только фиксированные эффекты (fixed effects), но и случайные эффекты (random effects). Для экспериментальной психологии это интересно тем, что в таких моделях можно не усреднять показатели по испытуемым, а учитывать влияние группирующей переменной "испытуемый" как случайный эффект. Подобные модели используются в самых разных областях. Для их использования в R есть два известных пакета: `nlme` и `lme4`.