Квадратичная подгонка для сгруппированных данных в R

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

|-----|-------|-------|-------|
| 0mg | 100mg | 200mg | 300mg |
|-----|-------|-------|-------|
| 25  |  16   |   6   |   8   |
| 19  |  15   |  14   |  18   |
| 22  |  19   |   9   |   9   |
| 15  |  11   |   5   |  10   |
| 16  |  14   |   9   |  12   |
| 20  |  23   |  11   |  13   |

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

Я пробовал следующее:

scores = c(25, 19, 22, 15, 16, 20,
           16, 15, 19, 11, 14, 23,
            6, 14,  9,  5,  9, 11,
            8, 18,  9, 10, 12, 13)

x_groups = rep(c(0,100, 200, 300), each = 6)
scores.quadratic = lm(scores ~ poly(x_groups, 2, raw = TRUE))

Затем я могу использовать функцию summary() для просмотра результатов. Я запутался в функции lm() и в том, как она должна соответствовать квадратичной функции. Насколько я понимаю, он возьмет каждый индекс в x_groups и возведет его в квадрат, а затем использует линейную подгонку с этим новым вектором, но мне это не кажется правильным.

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

Спасибо.


r model quadratic
person Cristopher Garduno    schedule 09.12.2017    source источник
comment
Квадратичная формула — это частный случай полиномиальной формулы, но она имеет порядок = 2. Полиномиальное соответствие с порядком = n переменной x будет соответствовать intercept + x + x^2 + x^3 + ... + x^n. Следовательно, квадратное число будет соответствовать intercept + x + x^2, и это именно те коэффициенты, которые вы получаете на выходе вашей модели. Похоже, вы ожидали, что это будет intercept + x^2.   -  person AntoniosK    schedule 09.12.2017

Ответы (1)


arrow_upward
2
arrow_downward

Давайте рассмотрим ваш образ мышления шаг за шагом. Во-первых, вы замечаете это падение по вашим числам для группы C. Лучший способ визуализировать это —

library(ggplot2)
library(dplyr)

scores = c(25, 19, 22, 15, 16, 20,
           16, 15, 19, 11, 14, 23,
           6, 14,  9,  5,  9, 11,
           8, 18,  9, 10, 12, 13)

x_groups = rep(c(0,100, 200, 300), each = 6)

# create dataset
d1 = data.frame(scores, x_groups) 

# calcuate average scores for each group
d2 = d1 %>% group_by(x_groups) %>% summarise(Avg = mean(scores))

# plot them
ggplot() + 
  geom_point(data = d1, aes(x_groups, scores)) +
  geom_line(data = d2, aes(x_groups, Avg), col="blue")

введите здесь описание изображения

Теперь вы действительно можете увидеть провал, и это тот паттерн, который вы хотите смоделировать.

Затем вы хотите подогнать свою квадратичную модель. Имейте в виду, что квадратичный - это частный случай полиномиальной формулы, но он имеет порядок = 2. Полиномиальное соответствие с порядком = n переменной x будет соответствовать intercept + x + x^2 + x^3 + ... + x^n. Следовательно, квадратное число будет соответствовать intercept + x + x^2, и это именно те коэффициенты, которые вы получаете на выходе вашей модели:

scores.quadratic = lm(scores ~ poly(x_groups, 2, raw = TRUE))
summary(scores.quadratic)

# Call:
#   lm(formula = scores ~ poly(x_groups, 2, raw = TRUE))
# 
# Residuals:
#   Min      1Q  Median      3Q     Max 
# -6.1250 -2.3333 -0.2083  1.8542  8.7917 
# 
# Coefficients:
#                                    Estimate Std. Error t value Pr(>|t|)    
#   (Intercept)                    20.2083333  1.5925328  12.689 2.58e-11 ***
#   poly(x_groups, 2, raw = TRUE)1 -0.0745833  0.0255747  -2.916  0.00825 ** 
#   poly(x_groups, 2, raw = TRUE)2  0.0001458  0.0000817   1.785  0.08870 .  
# ---
#   Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# 
# Residual standard error: 4.002 on 21 degrees of freedom
# Multiple R-squared:  0.4999,  Adjusted R-squared:  0.4523 
# F-statistic:  10.5 on 2 and 21 DF,  p-value: 0.0006919

Коэффициент квадратичного члена равен 0.0001458, близок к нулю, но статистически значимо отличается от нуля на уровне 0,1 (значение p = 0.08870). Поэтому модель как бы чувствует, что есть провал.

Вы можете построить подгонку следующим образом:

# plot the model
ggplot(d1, aes(x_groups, scores)) + 
  geom_point() +
  geom_smooth(formula = y ~ poly(x, 2, raw = TRUE),
              method = "lm")

Вы можете видеть это как сглаженную версию реального паттерна (1-й график).

введите здесь описание изображения

person AntoniosK    schedule 09.12.2017