0

我正在运行一个模型,它的响应变量是计数数据和非正常(非常右尾)有很多零。

我正在使用 glmmTMB 包并运行带有泊松分布的初始 glmm,其中 id(代码)和月份(montyear)作为随机因素。

m0 <- glmmTMB(deg~ SE_score + species + season + TL +
                (1|code) + (1|monthyear), family=poisson(), data=node_dat)

我使用了check_overdispersion()并且毫不奇怪地发现我在模型中存在过度分散。我使用负二项分布 (nbinom2) 运行了第二个模型,并使用 DHARMa 包绘制了 QQ 图和残差,但我仍然有很多分散。

m1 <- glmmTMB(deg~ SE_score + species + season + TL +
          (1|code) + (1|monthyear), family=nbinom2(), data=node_dat) 
res1 <- simulateResiduals(m1)
plot(res1)

在此处输入图像描述

我检查了零通货膨胀并且没有任何通货膨胀,并尝试了 nbinom1 系列,但在所有测试(KS、分散和异常值)中给出了更高的 AIC 和显着偏差。

我不确定下一个处理这个问题的方法是什么,也不确定为什么残差看起来很奇怪。

我的数据的一个子集如下

structure(list(station = structure(c(25L, 22L, 53L, 45L, 24L, 
65L, 12L, 64L, 27L, 28L, 2L, 15L, 46L, 53L, 10L, 16L, 11L, 59L, 
19L, 19L, 19L, 8L, 25L, 46L, 45L, 36L, 39L, 67L, 66L, 52L, 10L, 
10L, 3L, 67L, 41L, 12L), .Label = c("AR01", "BE01", "BE02", "BEUWM01", 
"BL01", "BL02", "GCB01", "GCB02", "GCB03", "GCB04", "NI01", "NI01b", 
"NI02", "NI03", "PB01", "PB02", "PB03", "PB04", "PB05", "PB06", 
"PB07", "PB08", "PB09", "PB10", "PB11", "PB12", "PB13", "PB14", 
"PB15", "PB16", "PB17", "PB18", "PB19", "PB20", "PB21", "PB22", 
"PB23", "PB24", "PB25", "PB26", "PB27", "PB28", "PB29", "PB30", 
"PB4G01", "PB4G02", "PBUWM01", "PBUWM02", "SA01", "SA02", "SA02b", 
"SA03", "SA04", "SA05", "SA06", "SA07", "SA11", "SA4G01", "SAUWM01", 
"SB01", "SB02/AR02", "SB03/AR05", "SB04/AR06", "VB01", "VB02", 
"VB03", "VB04"), class = "factor"), monthyear = structure(c(56L, 
27L, 52L, 57L, 40L, 57L, 7L, 45L, 51L, 46L, 3L, 46L, 4L, 11L, 
43L, 34L, 57L, 37L, 33L, 42L, 47L, 21L, 56L, 11L, 51L, 30L, 54L, 
18L, 57L, 55L, 12L, 26L, 27L, 7L, 2L, 25L), .Label = c("2014/01", 
"2014/02", "2014/03", "2014/04", "2014/05", "2014/06", "2014/07", 
"2014/08", "2014/09", "2014/10", "2014/11", "2014/12", "2015/01", 
"2015/02", "2015/03", "2015/04", "2015/05", "2015/06", "2015/07", 
"2015/08", "2015/09", "2015/10", "2015/11", "2015/12", "2016/01", 
"2016/02", "2016/03", "2016/04", "2016/05", "2016/06", "2016/07", 
"2016/08", "2016/09", "2016/10", "2016/11", "2016/12", "2017/01", 
"2017/02", "2017/03", "2017/04", "2017/05", "2017/06", "2017/07", 
"2017/08", "2017/09", "2017/10", "2017/11", "2017/12", "2018/01", 
"2018/02", "2018/03", "2018/04", "2018/05", "2018/06", "2018/07", 
"2018/08", "2018/09", "2018/10", "2018/11", "2018/12"), class = "factor"), 
    deg = c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
    0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0), 
    gs = c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
    0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0), 
    btw = c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0), 
    code = structure(c(169L, 120L, 87L, 197L, 121L, 65L, 173L, 
    183L, 11L, 191L, 22L, 191L, 141L, 5L, 85L, 130L, 90L, 188L, 
    197L, 82L, 99L, 202L, 193L, 19L, 69L, 142L, 195L, 117L, 118L, 
    89L, 5L, 132L, 121L, 4L, 172L, 15L), .Label = c("2390", "13573", 
    "13574", "13575", "13576", "19318", "19319", "19321", "19322", 
    "19506", "19514", "19519", "19520", "19524", "25537", "25540", 
    "25541", "25543", "25546", "25549", "25552", "25553", "27583", 
    "27585", "27586", "27591", "27592", "27593", "27594", "27595", 
    "27596", "27597", "27600", "27601", "27605", "27607", "27608", 
    "27613", "27614", "27617", "27619", "27620", "27621", "27626", 
    "27627", "27629", "27630", "27631", "27632", "28608", "28611", 
    "28612", "28618", "28625", "28628", "28629", "28631", "28632", 
    "28633", "28638", "28641", "28644", "28662", "28672", "28674", 
    "52978", "54815", "54846", "54852", "54860", "54863", "54865", 
    "54866", "54868", "54877", "54882", "54883", "54884", "54886", 
    "54890", "54892", "54895", "54896", "54901", "54904", "54914", 
    "54919", "54920", "54922", "54925", "54931", "54932", "54938", 
    "54952", "54954", "54955", "54958", "54959", "54962", "59950", 
    "59953", "59954", "59955", "59957", "59958", "59959", "59960", 
    "59961", "59962", "59964", "59966", "59969", "59970", "59971", 
    "59972", "59973", "59975", "59976", "59979", "59981", "59988", 
    "2388", "12950", "12952", "12956", "12958", "12960", "12962", 
    "12964", "12966", "12968", "13577", "14203", "19320", "19523", 
    "25534", "25535", "25536", "25539", "25542", "25544", "25545", 
    "25547", "25548", "25550", "27584", "27588", "27589", "27590", 
    "27598", "27599", "27602", "27603", "27604", "27606", "27609", 
    "27610", "27611", "27615", "27616", "27618", "27622", "27624", 
    "27625", "28624", "28627", "28637", "28639", "28642", "28660", 
    "28670", "34176", "34177", "34178", "34179", "52975", "52977", 
    "54817", "54821", "54822", "54825", "54845", "54849", "54880", 
    "54887", "54889", "54893", "54898", "54899", "54905", "54911", 
    "54912", "54915", "54933", "54947", "54957", "54961", "59951", 
    "59963", "59968", "59978", "59991", "59992", "59993", "59994", 
    "59995"), class = "factor"), species = structure(c(2L, 1L, 
    1L, 2L, 1L, 1L, 2L, 2L, 1L, 2L, 1L, 2L, 2L, 1L, 1L, 2L, 1L, 
    2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 1L, 1L, 2L, 
    1L, 1L, 2L, 1L), .Label = c("Grey Reef Shark", "Silvertip Shark"
    ), class = "factor"), SE_score = c(0.13, 0.29, 0.08, 0.1, 
    0.14, 0.07, 0.18, 0.01, 0.12, 0.09, 0.09, 0.1, 0.18, 0.06, 
    0.08, 0.21, 0.09, 0.03, 0.16, 0.11, 0.15, 0.15, 0.13, 0.14, 
    0.22, 0.21, 0.08, 0.14, 0, 0.03, 0.1, 0.09, 0.09, 0, 0.31, 
    0.31), region = structure(c(6L, 6L, 7L, 6L, 6L, 9L, 5L, 9L, 
    6L, 6L, 2L, 6L, 6L, 7L, 4L, 6L, 5L, 7L, 6L, 6L, 6L, 4L, 6L, 
    6L, 6L, 6L, 6L, 9L, 9L, 7L, 4L, 4L, 2L, 9L, 6L, 5L), .Label = c("Acoustic Release", 
    "Benares", "Blenheim", "Grand Chagos Bank", "Nelsons Island", 
    "Peros Banhos", "Saloman", "Speakers Bank", "Victory Bank"
    ), class = "factor"), date = structure(c(1533078000, 1456790400, 
    1522537200, 1535756400, 1491001200, 1535756400, 1404169200, 
    1504220400, 1519862400, 1506812400, 1393632000, 1506812400, 
    1396306800, 1414800000, 1498863600, 1475276400, 1535756400, 
    1483228800, 1472684400, 1496271600, 1509494400, 1441062000, 
    1533078000, 1414800000, 1519862400, 1464735600, 1527807600, 
    1433113200, 1535756400, 1530399600, 1417392000, 1454284800, 
    1456790400, 1404169200, 1391212800, 1451606400), tzone = "", class = c("POSIXct", 
    "POSIXt")), month = c(8, 3, 4, 9, 4, 9, 7, 9, 3, 10, 3, 10, 
    4, 11, 7, 10, 9, 1, 9, 6, 11, 9, 8, 11, 3, 6, 6, 6, 9, 7, 
    12, 2, 3, 7, 2, 1), season = structure(c(1L, 2L, 1L, 1L, 
    1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    1L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 
    2L, 2L), .Label = c("dry.season", "wet.season"), class = "factor"), 
    sex = structure(c(3L, 3L, 2L, 3L, 2L, 2L, 3L, 3L, 3L, 2L, 
    3L, 2L, 3L, 2L, 2L, 3L, 2L, 3L, 3L, 3L, 2L, 2L, 2L, 3L, 3L, 
    2L, 2L, 2L, 3L, 2L, 2L, 3L, 2L, 2L, 2L, 2L), .Label = c("", 
    "F", "M", "U"), class = "factor"), TL = c(103, 131, 118, 
    86, 112, 123, 115, 137, 70, 104, 131, 104, 134, 100, 122, 
    117, 135, 151, 86, 106, 117, 132, 139, 125, 114, 108, 144, 
    114, 106, 132, 100, 119, 112, 94, 140, 88)), row.names = c(NA, 
-36L), na.action = structure(c(`107971` = 107971L), class = "omit"), class = "data.frame")
4

0 回答 0