#r #statistics #regression #lm #anova
Вопрос:
У меня есть список почтовых индексов и количество смертей от covid на почтовый индекс в фрейме данных (не реальные цифры, просто примеры).:
| Город | Весь |
|---|---|
| Ричмонд | 552 |
| Las Vegas | 994 |
| Сан-Франциско | 388 |
Я хочу посмотреть, есть ли какая-либо связь между почтовым кодом и общим количеством смертей.
Я создал модель LM, используя функцию LM()
mod_zip <- lm(Total ~ City, data=zipcode)
Но когда я вызываю сводку(mod_zip) Я получаю NA за все, кроме столбца оценки.
| Коэффициенты | Оценивать | Ошибка Std. | значение t | Pr(> | т) |
|---|---|---|---|---|---|
| Город Ричмонд | 2851 | NA | NA | NA | NA |
| СитиЛасВегас | -2604 | NA | NA | NA | NA |
| CitySanFran | -966 | NA | NA | NA | NA |
Что я делаю не так?
Комментарии:
1. Неясно, что вы ищете, но я не думаю
lm, что это правильный подход. Похоже, у вас есть диапазон чисел смертей для разных городов, и вы хотите проверить, больше ли разница, чем вы ожидали бы, если бы смерти произошли независимо от города. Возможно, более уместен подход: en.wikipedia.org/wiki/Chi-squared_test2. У вас на самом деле есть только одно наблюдение за городом? Это не будет иметь большого значения для линейного момделя
3. Да, это только одно наблюдение/ряд на город. У меня возникли проблемы с ошибками DF. Данные довольно стандартные, поэтому я подумал, что что-то вроде ANOVA было бы более подходящим, но не подумал о дисперсии/DF. Отличная мысль. Я попробую это в следующем разделе.
Ответ №1:
lm превратит коэффициент в один-горячие столбцы, так что у вас есть параметр для каждого города, кроме одного, и глобальный перехват.
Затем (предполагая, не видя ваших данных) вы пытаетесь оценить n точек данных с n параметрами, что ему удается сделать, но у него недостаточно степеней свободы, чтобы получить стандартную ошибку.
Упрощенный пример для воспроизведения:
df <- data.frame(x = LETTERS, y = rnorm(26), stringsAsFactors = TRUE)
fit <- lm(y~x, data = df)
summary(fit)
Вы увидите перехват и параметры от B до Z (26 параметров для 26 наблюдений), таким образом, степени свободы равны 0, следовательно, стандартные ошибки и связанные с ними показатели не поддаются вычислению.
Комментарии:
1. Это правильно. Я получаю ошибку, когда я также пытался использовать функции ANOVA. Каков был бы наилучший статистический подход к получению данных? Я намерен также увеличить численность населения и вычислить процент населения, но я чувствую, что сначала мне нужно это выяснить.
2. У Anova та же проблема, что и у lm, потому что на самом деле она делает то же самое. Вам придется использовать различные условные обозначения или выполнить некоторые функции (агрегировать почтовые индексы на более высоком уровне, например, штата, или превратить их в широту и долготу, или долю голосов республиканцев на последних выборах или что-то в этом роде) таким образом, чтобы решить фактический вопрос, на который вы пытаетесь ответить.
Ответ №2:
Похоже, вы хотите проверить, является ли город важным фактором для прогнозирования смертности. Другими словами, ожидали бы вы увидеть наблюдаемый диапазон значений, если бы каждая смерть имела равные шансы произойти в любом городе?
Моя интуиция на этот счет заключалась бы в том, что, безусловно, должна быть разница, основанная на различиях в демографических показателях, правилах, нормах, показателях вакцинации и характере инфекционного заболевания, которое распространяется больше, если изначально инфицировано больше людей.
Если вы хотите подтвердить эту интуицию, вы можете использовать моделирование. Допустим, во всех городах был одинаковый базовый уровень риска 800, и все различия были полностью обусловлены случайностью.
set.seed(2021)
Same_risk = 800
Same_risk_deaths = rpois(100, Same_risk)
mean(Same_risk_deaths)
sd(Same_risk_deaths)
Здесь наблюдаемое среднее значение действительно близко к 800, со стандартным отклонением около 3% от среднего значения.
Если бы вместо этого у нас была ситуация, когда в некоторых городах по какой-либо комбинации причин были разные факторы риска (скажем, 600 или 1000), то мы могли бы увидеть одно и то же среднее значение около 800, но с гораздо более высоким стандартным отклонением около 25% от среднего значения.
Diff_risk = rep(c(600, 1000), 50)
Diff_risk_deaths = rpois(100, Diff_risk)
mean(Diff_risk_deaths)
sd(Diff_risk_deaths)
Я полагаю, что ваши данные не похожи на первое распределение и вместо этого гораздо более разнообразны.