具有MRF平滑功能的GAM-错误（nb /多边形区域名称与数据区域名称不匹配

Question

在@GavinSimpson撰写的超级博客之后，我正在尝试调整2015年波兰地方政府的选举结果。 https://www.fromthebottomoftheheap.net/2017/10/19/first-steps-with-mrf-smooths/ 我将xls上的shp数据与6位数字标识符（可能是前导0 s）结合在一起。我将其保留为文本变量。编辑，我简化了标识符，现在使用从1到n的序列来简化我的问题。

library(tidyverse)
library(sf)
library(mgcv)

# Read data
# From https://www.gis-support.pl/downloads/gminy.zip shp file

boroughs_shp <- st_read("../../_mapy/gminy.shp",options = "ENCODING=WINDOWS-1250",
                     stringsAsFactors = FALSE ) %>% 
  st_transform(crs = 4326)%>% 
  janitor::clean_names() %>% 
# st_simplify(preserveTopology = T, dTolerance = 0.01) %>% 
  mutate(teryt=str_sub(jpt_kod_je, 1, 6)) %>% 
  select(teryt, nazwa=jpt_nazwa, geometry)

# From https://parlament2015.pkw.gov.pl/wyniki_zb/2015-gl-lis-gm.zip data file
elections_xls <-
  readxl::read_excel("data/2015-gl-lis-gm.xls",
             trim_ws = T, col_names = T) %>% 
  janitor::clean_names() %>% 
  select(teryt, liczba_wyborcow, glosy_niewazne)

elections <-
  boroughs_shp %>% fortify() %>% 
  left_join(elections_xls, by = "teryt") %>% 
  arrange(teryt) %>%
  mutate(idx = seq.int(nrow(.)) %>% as.factor(), 
         teryt = as.factor(teryt)) 

# Neighbors

boroughs_nb <-spdep::poly2nb(elections, snap = 0.01, queen = F, row.names = elections$idx )
names(boroughs_nb) <- attr(boroughs_nb, "region.id")

# Model

ctrl <- gam.control(nthreads = 4) 
m1 <- gam(glosy_niewazne ~ s(idx, bs = 'mrf', xt = list(nb = boroughs_nb)), 
          data = elections,
          offset = log(liczba_wyborcow), # number of votes
          method = 'REML', 
          control = ctrl,
          family = betar()) 

这是错误消息：

    Error in smooth.construct.mrf.smooth.spec(object, dk$data, dk$knots) : 
  mismatch between nb/polys supplied area names and data area names
In addition: Warning message:
In if (all.equal(sort(a.name), sort(levels(k))) != TRUE) stop("mismatch between nb/polys supplied area names and data area names") :
  the condition has length > 1 and only the first element will be used

选举$ idx是一个因素。我正在使用它为boroughs_nb命名，以确保我具有相同数量的级别。我究竟做错了什么？

编辑：满足错误消息中提到的条件：

> all(sort(names(boroughs_nb)) == sort(levels(elections$idx)))
[1] TRUE

Answer 1

看来我解决了这个问题，也许不太意识到作为统计初学者它是如何做的。

首先，建模数据中不应出现任何一个 NA。有一个。之后，mcgv 似乎开始运行，但它花了很长时间（一刻钟），对我来说莫名其妙，只有当我将结数限制为k=50,结果不佳（更少或更多，并且没有返回任何结果）并发出警告时对结果持谨慎态度。然后我尝试删除offset=log(liczba_wyborcow)选民的偏移量，并将每 1000 票的无效票数设为我的预测变量。

elections <-
 boroughs_shp %>%  
 left_join(elections_xls, by = "teryt") %>% na.omit() %>% 
 arrange(teryt) %>% 
 mutate(idx = row_number() %>% as.factor()) %>% 
 mutate(void_ratio=round(glosy_niewazne/liczba_wyborcow,3)*1000)

既然它是一个计数，为什么不尝试将family = betar()gam 公式更改为poisson()- 仍然不是一个好的结果，然后更改为负二项式 family = nb() 现在我的公式看起来像

m1 <-
gam(
 void_ratio ~ s(
 idx,
 bs = 'mrf',
 k =500,
 xt = list(nb = boroughs_nb),
 fx = TRUE),
 data = elections_df,
 method = 'REML', 
 control = gam.control(nthreads = 4),
 family = nb()
)

现在看起来速度非常快，并且返回有效结果，没有警告或错误。在配备 4 核 Intel Core I7 6820HQ @ 2.70GHZ 16GB Win10 的笔记本电脑上，构建模型现在只需 1-2 分钟。

简而言之，我所做的更改是：删除单个 NA，从公式中删除偏移量并使用负二项式分布。

这是我想要实现的结果，从左到右依次是实际无效选票率、模型平滑后的比率以及指示差异的残差。mcgv 代码让我可以做到这一点。