从不同长度列随机抽取指定元素的R语言实现问题

数据结构说明 #

现有一个维度为1042×64的tibble，每一列代表一个两栖动物科，每一行是对应科下的物种名称。除了包含1042个物种的科外，其他列均用NA补全至1042行。前5行2列的数据示例如下：

> amphilist[1:5,1:2]
A tibble: 5 × 2
  `Allophrynidae_(1_genus;_3_species)` `Alsodidae_(3_genera;_26_species)`
  <chr>                                <chr>                              
1 Allophryne_relicta                   Alsodes_australis                 
2 Allophryne_resplendens               Alsodes_barrioi                    
3 Allophryne_ruthveni                  Alsodes_cantillanensis            
4 NA                                   Alsodes_coppingeri                 
5 NA                                   Alsodes_gargola 

已完成的准备工作 #

已通过循环计算每个科的物种丰富度，并通过条件语句确定了每个科需要抽取的物种数量，代码及结果如下：

amphilist <- read_xlsx("amphilist.xlsx", col_names = TRUE)

families <- colnames(amphilist)
family_n <- ncol(amphilist)
spcR <- vector(length = family_n)

for(i in 1:length(families)) {
  families.i <- families[i]
  spcR[i] <- colSums(amphilist[,families.i] > 0, na.rm = TRUE)
}

species_no <- data.frame(families, spcR)
species_no$choose <- ifelse(species_no$spcR > 50, ceiling(species_no$spcR/10), 
                            ifelse(species_no$spcR >= 5 & species_no$spcR <= 50,
                                   5, species_no$spcR))

> species_no[1:3,]
                                        families spcR choose
1             Allophrynidae_(1_genus;_3_species)    3      3
2               Alsodidae_(3_genera;_26_species)   26      5
3 Alytidae_(2_subfamilies;_3_genera;_12_species)   12      5

当前问题 #

尝试随机抽取时，结果全为NA，即使无NA的列也出现该问题。错误代码如下：

choose_no <- species_no$choose
set.seed(43)
for(i in 1:length(families)) {
  families.i <- families[i]
  choose_no.i <- choose_no[i]
  rand_amphilist <- amphilist[sample(amphilist[,i], 
                                     size = choose_no.i), ]
}

解决方案 #

错误原因 #

原代码中sample(amphilist[,i], size = choose_no.i)是对列中的元素（物种名或NA）进行采样，随后用这些采样得到的字符串作为行索引去提取tibble的行，这完全不符合行索引的规则（行索引应为数字），因此导致结果全为NA。

正确代码1：循环实现 #

先提取每列的非NA物种，再从中随机抽取指定数量，最后将结果合并并补全NA至1042行：

library(dplyr)

# 初始化结果tibble
rand_amphilist <- tibble()

set.seed(43)
for(i in 1:length(families)) {
  # 获取当前列的所有非NA物种
  current_species <- na.omit(amphilist[[i]])
  # 按指定数量随机抽取
  sampled_species <- sample(current_species, size = choose_no[i], replace = FALSE)
  # 转换为带列名的tibble
  sampled_tbl <- tibble(!!families[i] := sampled_species)
  # 补全NA到1042行（若需要保持原数据行数）
  sampled_tbl <- bind_rows(sampled_tbl, tibble(!!families[i] := rep(NA, 1042 - nrow(sampled_tbl))))
  # 合并到结果中
  rand_amphilist <- if(i == 1) sampled_tbl else bind_cols(rand_amphilist, sampled_tbl)
}

正确代码2：tidyverse函数式实现 #

用purrr的map2_dfc函数，更简洁地完成批量处理：

library(purrr)
library(dplyr)

set.seed(43)
rand_amphilist <- map2_dfc(amphilist, choose_no, function(col, n) {
  # 提取非NA物种并采样
  sampled <- sample(na.omit(col), size = n, replace = FALSE)
  # 补全NA到1042行并转为带列名的tibble列
  c(sampled, rep(NA, 1042 - length(sampled))) %>%
    as_tibble_col(column_name = colnames(amphilist)[which(colnames(amphilist) == colnames(col))])
})

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