r - GAM with mrf smooth - 错误(nb/polys 区域名称和数据区域名称不匹配

Question

在@GavinSimpson 撰写的 super 博客之后，我正在尝试拟合 2015 年波兰地方政府选举结果。 https://www.fromthebottomoftheheap.net/2017/10/19/first-steps-with-mrf-smooths/我将 xls 与 shp 数据连接到 6 位标识符(可能有前导 0)。我把它作为一个文本变量。编辑，我简化了标识符，现在使用从 1 到 nrow 的序列来简化我的问题。

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 代码让我这样做。