实例教程 | 使用TDIP包对分类性状的缺失值进行插补

文摘   科学   2024-07-27 08:01   江苏  
TDIP包允许在系统发育树的帮助下或没有帮助的情况下插补性状数据集中的缺失值。通过使用这个包,用户可以模拟性状数据集,根据几种缺失机制包含缺失值,和/或应用各种插补方法。有关插补方法的细节可以查看上一篇推文(Method in Ecology and Evolution | 分类生物数据的基准插补方法)。本篇主要整理了TDIP包进行性状插补的基本工作流程。

  • 安装TDIP

所提供的代码将从gitHub安装TDIP。这需要用户在其计算机上安装devtools包。

# NADIA包是TDIP包的依赖之一,但这个包目前已从cran上移除了,需要自行下载,然后从本地安装# https://cran.r-project.org/src/contrib/Archive/NADIA/# 安装过程中可能还会提示缺少一些依赖,总之原则上提示缺少什么就安装什么# install.packages(c("mlr3", "mlr3pipelines", "paradox", "missMDA", "mlr3learners", "Amelia", "softImpute", "missRanger"))# 从github安装TDIP包devtools::install_github("Matgend/TDIP")
  • 数据

可以使用模拟数据或经验数据。要模拟性状数据集,可以使用函数data_simulator()。该函数需要2个参数:

1)由这些列组成的数据框:

nbr_traits:使用特定参数模拟的性状数量

class:性状类型(连续、non_eq_nominal(分类)和序数)

model:进化模型(BM1OU1ARDSYMERManual

states:离散性状的状态数(如果连续,则设为1

correlation:与其他性状相关或不相关的模拟性状组的相应指数

uncorr_traits:在“nbr_traits”中,它是不相关性状的数量

fraction_uncorr_traits:在“nbr_traits”中,不相关的性状的分数

lambdaPagel's lambda

kappaPagel's kappa

highCor:“manTrait”中定义的性状与模拟性状之间的相关率

manTrait:根据手动模型模拟的性状与相关性状的索引

2)包含出生率、死亡率和类群数的列表

函数data_simulator()将返回一个列表,其中包含模拟的性状数据、系统发育树、用于模拟的所有参数和输入。

# 加载TDIP包library(TDIP)
# 加载包中自带的数据集data(dataframe)
# 设置参数树,包含三个参数:0.4, 0.1和100param_tree <- list(0.4, 0.1, 100)
# 使用data_simulator函数生成模拟数据# param_tree参数指定参数树# dataframe参数使用加载的数据集simData <- data_simulator(param_tree = param_tree, dataframe = dataframe)
# 打印模拟数据的前几行# simData$FinalData 是模拟生成的最终数据集head(simData$FinalData)## I1.0/1 I2.0/1 I3.0/1 I4.0/1 F1.0/3 F2.0/3 F3.0/3 I1.0/4 I2.0/4 I3.0/4 F1.1/2 F2.1/2 F3.1/2## t23 2 1 2 1 0.203913671 -1.6171502 -0.3472887 2 2 0 1.13084411 0.388191354 1.974841## t24 2 1 2 1 0.569097928 -0.7370138 -0.1325044 2 2 0 0.26796591 -0.113133721 1.921079## t5 0 1 1 1 -0.003254024 0.5929791 -1.4368104 1 1 2 0.37160420 0.562254863 1.591107## t49 2 1 1 1 -0.308026821 0.6998594 0.2602554 0 2 0 -0.82149229 0.524967048 2.310056## t7 2 1 1 1 0.735442781 -0.1129450 -0.7334719 2 2 0 -0.03072779 -0.007397889 2.186394## t64 2 1 1 1 1.209052247 -0.8816450 -0.4864245 0 2 0 -1.01605308 0.287524388 1.256216
对于经验数据集,可以使用函数data_preprocessing()将所有字符列转换为因子,并创建一个包含数据集和系统发育树(如果提供)的列表。

# 加载simulatedData的数据集data(simulatedData)# 加载系统发育树数据集data(tree)
# 将simulatedData的行名设置为新的"Species"simulatedData$Species <- rownames(simulatedData)
# 使用data_preprocessing函数处理数据,不使用系统发育树# empTree = NULL 表示不使用系统发育信息# save = NULL 表示不保存处理结果empData <- data_preprocessing(simulatedData, empTree = NULL, save = NULL)
# 再次使用data_preprocessing函数处理数据,这次使用系统发育树empData <- data_preprocessing(simulatedData, empTree = tree, save = NULL)
  • 缺失值的生成

缺失值遵循缺失机制。在这个软件包中,用户可以根据4种缺失机制生成缺失值:

1)完全随机缺失(MCAR

2)随机缺失(MAR

3)非随机缺失(MNAR

4PhyloNa

# 设置缺失值比率为5%missingRate <- 0.05
# 生成完全随机缺失(MCAR)数据# 使用mcar_miss_meca函数在simData$FinalData的所有列中生成MCAR缺失值mcar_values <- mcar_miss_meca(missingRate = missingRate, ds = simData$FinalData, cols_mis = 1:ncol(simData$FinalData))
# 生成随机缺失(MAR)数据# 使用mar_miss_meca函数在simData$FinalData的第1列生成MAR缺失值,基于第3列的值mar_values <- mar_miss_meca(missingRate = missingRate, ds = simData$FinalData, cols_mis = 1, cols_ctrl = 3)
# 生成非随机缺失(MNAR)数据# 使用mnar_miss_meca函数在simData$FinalData的所有列中生成MNAR缺失值mnar_values <- mnar_miss_meca(missingRate = missingRate, ds = simData$FinalData, cols_mis = colnames(simData$FinalData))
# 生成基于系统发育关系的缺失(PhyloNa)数据# 使用phyloNa_miss_meca函数在simData$FinalData中生成基于系统发育树的缺失值phyloNa_values <- phyloNa_miss_meca(missingRate = missingRate, ds = simData$FinalData, tree = simData$TreeList$`0`)
函数na_insertion()可根据多个缺失机制同时模拟缺失值。

# 使用na_insertion函数生成缺失值missing_values <- na_insertion(  missingRate = missingRate,  # 缺失值比率  dataset = simData$FinalData,  # 完整数据集  missingTraits = ncol(simData$FinalData),  # 要生成缺失值的特征总数  MARTraits = 1,  # MAR机制中生成缺失值的变量  MARctrlTraits = 3,  # MAR机制中控制缺失值的变量  traitsNoNA = NULL,  # 不生成缺失值的特征  tree = simData$TreeList$`0`,  # 系统发育树  save = NULL  # 不保存结果到文件)
# 可视化缺失数据library(visdat)NaData <- missing_values$DataNaN
# 遍历NaData中的每个分区并绘制缺失数据图for(rdn in 1:length(NaData)){ partition <- NaData[[rdn]][[1]] # 获取当前分区 title <- names(partition) # 获取分区名称作为图表标题 plot(visdat::vis_miss(partition), main = title) # 使用vis_miss绘制缺失数据图}

  • 缺失值插补

缺失值的插补可根据几种策略采用几种插补方法。TDIP软件包中提供的方法如下:

1)系统发育插补方法

Rphyloparspi_continuous_traits_phylo()

corHMMpi_categorical_traits_phylo()

2)机器学习

missForestmissForest_phylo()

kNNkNN_phylo()

预测均值匹配:mice_phylo()

3)深度学习

GAINgain_phylo()

Polytomous regressionmice_phylo()

上述每种方法都可以按照3种策略进行应用:

1)无系统发育信息

2)有作为特征向量的系统发育信息

3)“2-step”:将系统发育插补方法的输出结果与缺失值估算前的不完整数据合并作为输入。

函数missing_data_imputation()提供了一种通过多种插补方法和策略对缺失值进行插补的方法。每种插补方法也可以单独调用。

# 定义一个包含三种插补方法的字符向量methods <- c("mice_phylo", "missForest_phylo", "kNN_phylo")# 设置插补策略为"NP"(非系统发育)strategies <- "NP"# 设置方差分数为0.95,表示系统发育树特征向量解释的方差比例varfrac <- 0.95
# 调用missing_data_imputation函数进行缺失数据插补imputedData <- missing_data_imputation(ImputationApproachesNames = methods, # 指定使用的插补方法 data = missing_values$DataNaN$MCAR$`MCAR/13/0.05`, # 输入包含缺失值的数据集,这里使用完全随机缺失(MCAR)的数据 tree = simData$TreeList$`0`, # 输入系统发育树对象 strategies = strategies, # 指定插补策略 varfrac = varfrac, # 设置方差分数 save = NULL)# 设置save参数为NULL,表示不保存结果到文件imputedData
  • 硬投票

然后,函数hard_voting()提供了使用插补数据集应用硬投票分类器的可能性。

# 从imputedData对象中选择特定的插补数据集namesMethods <- c("MICE", "MissForest", "KNN")  # 定义插补方法的名称namesToselect <- paste(namesMethods, "NP", sep = "_")  # 创建要选择的数据集名称namesImputedData <- names(imputedData)  # 获取imputedData对象中所有数据集的名称datasets <- imputedData[which(namesImputedData %in% namesToselect)]  # 选择指定的数据集
# 执行硬投票hv <- hard_voting(listOfdataframe = datasets)  # 对选定的数据集应用硬投票方法hv
  • 误差计算

连续变量的误差计算是均方根误差(RMSE),而分类变量的误差计算是误分类条目的比例(PFC)。
# 计算包含由 missForest 插补的 MCAR 的数据的误差errors <- imputation_error(imputedData = imputedData$MissForest_NP, # 使用missForest方法插补后的数据集                           trueData = simData$FinalData, # 原始完整数据集                           missingData = missing_values$DataNaN$MCAR$`MCAR/13/0.05`,  # 包含MCAR缺失值的数据集                           imputationApproachesName = methods[2], # 插补方法名称                           dataframe = simData$Dataframe) # 用于模拟的数据框
errors## trait missForest_phylo## 1 I1.0/1 0.0000000## 2 I2.0/1 0.0000000## 3 I3.0/1 0.6000000## 4 I4.0/1 0.4000000## 5 F1.0/3 0.8053965## 6 F2.0/3 0.6036047## 7 F3.0/3 0.6121908## 8 I1.0/4 0.4000000## 9 I2.0/4 0.0000000## 10 I3.0/4 0.6000000## 11 F1.1/2 0.5367540## 12 F2.1/2 0.6819793## 13 F3.1/2 0.6687828
  • 参考来源

  1. Gendre, M., Hauffe, T., Pimiento, C., & Silvestro, D. (2024). Benchmarking imputation methods for categorical biological data. Methods in Ecology and Evolution, 00, 1–15. https://doi.org/10.1111/2041-210X.14339

Biodiversity Monitoring
生物多样性;监测保护;群落生态;生态统计;R语言;python。 主要分享一些前沿的文献和方法实例,更新看心情和时间。
 最新文章