跟着NC论文学数据分析:R语言gradientForest包分析环境变量的重要性排序

论文

Genomic insights into local adaptation and future climate-induced vulnerability of a keystone forest tree in East Asia
论文中提供的代码
https://github.com/jingwanglab/Populus_genomic_prediction_climate_vulnerability/blob/2d314a9618a4e85d26030897c5b296dea5284e35/7-Genetic_offset/1GF.R#L3
没有找到论文中的数据
搜索了一下相关教程，找到了一个提供数据的教程
https://github.com/eriqande/merida-workshop-2022/tree/main
这个教程里写道这个分析还可以预测snp位点的重要性，有时间看看其中提到的论文
（这个github主页上还有好多其他教程）
这里分析需要准备两个数据
一个是环境数据，行是群体，列是环境变量，这里的群体是根据地理位置定的，在同一个区域内就是一个群体，这个群体的环境数据是一样的
另外一个数据是等位基因频率表，按照不用的群体算等位基因频率（这里没有搞懂是需要用突变基因型的等位基因频率还是用参考基因型的等位基因型频率）
利用vcf文件计算不同群体等位基因频率

如果有vcf文件按照不同的群体计算等位基因频率，使用plink2，准备一个群体分组文件

image.png
这个需要有表头
命令
~/biotools/plink2/plink2 –vcf ../rMVP/smoove_filtered.vcf –freq –pheno pheno02.txt –allow-extra-chr –allow-no-sex –loop-cats pop
最后的pop对应的群体分组文件里第二列的列名，这个根据自己的数据修改
输出结果是每个群体单独输出一个文件

image.png
计算的是突变的等位基因频率

image.png
多个群体自己写代码合并
gradiantForest R包安装

install.packages(“gradientForest”, repos=”http://R-Forge.R-project.org”)
分析代码

library(gradientForest)
library(tidyverse)

path<-“D:/Jupyter/IAAS/reference/LEA变异和环境因素关联/merida-workshop-2022-main/”
Gwifl<-read_csv(paste0(path,’data/6.0/wiflforest.allfreq.sample.csv’))
Ewifl<-read_csv(paste0(path,’data/6.0/wiflforest.env.csv’))

nSites<-21 lev <- floor(log2(nSites*0.368/2)) lev wiflforest<-gradientForest(cbind(Ewifl,Gwifl), predictor.vars=Ewifl %>% colnames(),
response.vars=Gwifl %>% colnames(),
ntree=10,
transform = NULL,
compact=T,
nbin=101,
maxLevel=lev,
trace=T)
环境变量的重要性排序

plot(wiflforest,plot.type=”Overall.Importance”)
这里怎么把数据提取出来自己作图，暂时没有搞明白

image.png
欢迎大家关注我的公众号
小明的数据分析笔记本