使用R包Segregation计算县差异指数结果异常的问题排查

问题背景 #

我是中初级R用户，目标是计算美国所有县的Dissimilarity指数（D值），为减少运算量先运行了半数州的数据。但运行代码后，所有县的Dissimilarity估计值（D值）均为0.648，怀疑和手动创建COUNTY_FIPS变量缺少几何信息有关，希望获取修正代码或调整建议。

原代码 #

my_states <- c("AL", "AK", "AZ", "AR", "CA", "CO", "CT", "DE", "DC", "FL", "GA", "HI",
               "ID", "IL", "IN", "IA", "KS", "KY", "LA", "ME", "MD", "MA", "MI")

#my_states2 <- c ("MN", "MS", "MO", "MT", "NE", "NV", "NH", "NJ", "NM", "NY", 
               #"NC", "ND", "OH", "OK", "OR", "PA", "RI", "SC", "SD", "TN", "TX",
               #"UT", "VT", "VA", "WA", "WV", "WI", "WY") 

acs_data1 <- get_acs(
  geography = "tract",
  variables = c(
    white = "B03002_003",
    black = "B03002_004",
    asian = "B03002_006",
    hispanic = "B03002_012"), 
  state = my_states,
  geometry = TRUE,
  year = 2019
) 

seg_acs_data <- acs_data1 %>% 
  mutate(COUNTY_FIPS = substr(GEOID, 1, 5))

subsetseg <- seg_acs_data %>% filter(variable %in% c("white", "black"))
  
dissimilarity <- subsetseg %>% group_by(COUNTY_FIPS) %>%
  group_modify(~ 
                 dissimilarity(data = subsetseg,
                               group = "variable",
                               unit = "GEOID",
                               weight = "estimate"
                 )) 

原输出 #

# A tibble: 1,314 x 3
# Groups:   COUNTY_FIPS [1,314]
   COUNTY_FIPS stat    est
   <chr>       <chr> <dbl>
 1 01001       D     0.648
 2 01003       D     0.648
 3 01005       D     0.648
 4 01007       D     0.648
 5 01009       D     0.648
 6 01011       D     0.648
 7 01013       D     0.648
 8 01015       D     0.648
 9 01017       D     0.648
10 01019       D     0.648
# ... with 1,304 more rows

问题根源与修正方案 #

核心错误 #

问题和几何信息无关，根本原因是group_modify内部调用dissimilarity()时，传入了整个数据集subsetseg，而非分组后的子数据集.x。这导致每个县都在计算全量数据的D值，结果自然完全一致。

修正后的代码 #

只需修改group_modify中的data参数为.x：

dissimilarity <- subsetseg %>% 
  group_by(COUNTY_FIPS) %>%
  group_modify(~ 
                 dissimilarity(data = .x,  # 此处用分组后的子数据集.x替代subsetseg
                               group = "variable",
                               unit = "GEOID",
                               weight = "estimate"
                 ))

额外优化建议 #

1. 减少数据量：Dissimilarity指数不需要几何信息，且你只用到白人和黑人数据，可在get_acs阶段直接过滤变量并关闭几何数据加载：

acs_data1 <- get_acs(
  geography = "tract",
  variables = c(white = "B03002_003", black = "B03002_004"), # 仅保留需要的变量
  state = my_states,
  geometry = FALSE, # 无需几何数据，减少内存占用
  year = 2019
)

2. 验证结果：运行修正代码后，可通过filter(COUNTY_FIPS %in% c("01001", "06037"))查看不同县的D值是否存在差异，确认计算正常。

内容的提问来源于stack exchange，提问作者Richardmk