Пуассонова GLM для нецелых чисел - R - PullRequest
1 голос
/ 19 февраля 2020

Я надеялся получить совет по GLM с семьей Пуассона.

У меня есть набор данных с количеством укусов, взятых каждым человеком за период времени. Поскольку наблюдаемые люди питались в течение разных периодов времени, когда я рассчитал частоту прикусов для каждого человека как количество укусов в минуту, я получил нецелые числа. Теперь, основываясь на том, что я прочитал до сих пор, я все еще должен быть в состоянии сделать GLM с семьей Пуассона. Тем не менее, я сталкиваюсь с ошибками, и я думаю, что это может быть потому, что R не нравится тот факт, что я использую нецелые числа. У кого-нибудь есть совет?

Example <- structure(list(Species = c("Fish1", "Fish2", "Fish3", "Fish4", 
"Fish5", "Fish6", "Fish7", "Fish1", "Fish2", "Fish3", "Fish4", 
"Fish5", "Fish6", "Fish7", "Fish1", "Fish2", "Fish3", "Fish4", 
"Fish5", "Fish6", "Fish7"), Site = c(1, 1, 1, 1, 1, 1, 1, 2, 
2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3), Bite_Rate = c(3.5, 7.5, 
0, 0, 2.45, 5.5, 6.5, 6.5, 7.5, 8.03, 32.1, 15.6, 18.2, 19.1, 
20.5, 20.5, 3.5, 5.7, 6.7, 23.2, 0)), class = c("spec_tbl_df", 
"tbl_df", "tbl", "data.frame"), row.names = c(NA, -21L), spec = structure(list(
   cols = list(Species = structure(list(), class = c("collector_character", 
   "collector")), Site = structure(list(), class = c("collector_double", 
   "collector")), Bite_Rate = structure(list(), class = c("collector_double", 
   "collector"))), default = structure(list(), class = c("collector_guess", 
   "collector")), skip = 1), class = "col_spec"))

str(Example) # check structure 
Example$Species<-as.factor(Example$Species) # set species as a factor 
str(Example) # check structure 
glm<-glm(Species~Bite_Rate, data=Example, family = poisson) # create the GLM

Сообщение об ошибке при запуске GLM:

Error in if (any(y < 0)) stop("negative values not allowed for the 'Poisson' family") : 
  missing value where TRUE/FALSE needed
In addition: Warning message:
In Ops.factor(y, 0) : ‘<’ not meaningful for factors 

На самом деле у меня нет никаких отрицательных значений, что немного сбивает меня с толку. Любой совет будет высоко ценится!

РЕДАКТИРОВАТЬ: Основываясь на комментариях, я обновил мой пример данных, так что он имеет количество укусов и время наблюдения в секундах


Example <- structure(list(Species = structure(c(1L, 2L, 3L, 4L, 5L, 6L, 
7L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 1L, 2L, 3L, 4L, 5L, 6L, 7L), .Label = c("Fish1", 
"Fish2", "Fish3", "Fish4", "Fish5", "Fish6", "Fish7"), class = "factor"), 
    Site = c(1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 3, 3, 3, 
    3, 3, 3, 3, 3), Bites = c(0, 10, 18, 17, 6, 0, 1, 0, 19, 
    12, 7, 3, 5, 1, 5, 0, 10, 18, 17, 7, 25), Observed_Seconds = c(50, 
    33, 47, 20, 17, 10, 14, 21, 48, 10, 50, 33, 47, 20, 17, 10, 
    14, 21, 48, 10, 90)), row.names = c(NA, -21L), spec = structure(list(
    cols = list(Species = structure(list(), class = c("collector_character", 
    "collector")), Site = structure(list(), class = c("collector_double", 
    "collector")), Bites = structure(list(), class = c("collector_double", 
    "collector")), Observed_Seconds = structure(list(), class = c("collector_double", 
    "collector"))), default = structure(list(), class = c("collector_guess", 
    "collector")), skip = 1), class = "col_spec"), class = c("spec_tbl_df", 
"tbl_df", "tbl", "data.frame"))

Спасибо вы!

1 Ответ

2 голосов
/ 19 февраля 2020

Вам необходимо получить секунды (знаменатель), с которыми вы нормализовались, и фактические укусы (количество).

Далее укажите минуты как смещение и обратите внимание, что ваша переменная ответа находится слева от ~:

fit = glm(Bites ~ Species,offset=log(Observed_Seconds),
family=poisson,data=Example)

Вы можете посмотреть сводку:

summary(fit)

   Coefficients:
             Estimate Std. Error z value Pr(>|z|)    
SpeciesFish1  -2.8679     0.4472  -6.413 1.42e-10 ***
SpeciesFish2  -1.1436     0.1857  -6.158 7.35e-10 ***
SpeciesFish3  -0.5738     0.1581  -3.629 0.000284 ***
SpeciesFish4  -0.7732     0.1543  -5.011 5.42e-07 ***
SpeciesFish5  -1.3269     0.1961  -6.766 1.33e-11 ***
SpeciesFish6  -1.7198     0.2887  -5.958 2.56e-09 ***
SpeciesFish7  -1.5244     0.1925  -7.921 2.35e-15 ***

Кажется очень важным, но не следует проверять, не являются ли данные чрезмерно рассредоточенными, а также учитываются ли другие факторы (например, Сайт):

fit_quasi = glm(Bites ~ Species + factor(Site),offset=log(Observed_Seconds),
          family=quasipoisson,data=Example)
summary(fit_quasi)

Coefficients:
              Estimate Std. Error t value Pr(>|t|)  
(Intercept)    -2.9754     1.0713  -2.777   0.0167 *
SpeciesFish2    1.8487     1.1434   1.617   0.1319  
SpeciesFish3    2.2731     1.1152   2.038   0.0642 .
SpeciesFish4    2.1246     1.1205   1.896   0.0823 .
SpeciesFish5    1.3533     1.1604   1.166   0.2662  
SpeciesFish6    1.2754     1.2658   1.008   0.3336  
SpeciesFish7    0.8922     1.1719   0.761   0.4612  
factor(Site)2  -0.2325     0.5132  -0.453   0.6587  
factor(Site)3   0.6118     0.4677   1.308   0.2154  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for quasipoisson family taken to be 5.521045)

Если оно следует пуассону, дисперсия будет около 1, но в этом случае он рассредоточен.

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...