这个生态位分化分析整个的运行过程时间还是挺长的，所以想写个脚本直接在服务器上提交一个任务，最开始是想直接用save函数把运行结果保存下来。但是这个结果非常大，保存到本地以后有60多个G（暂时还没太搞懂这个数据里都有什么）。最后能用到的结果是
代码语言：javascript复制result$p.values
result$p.values%>%as.data.frame()%>%rownames_to_column()
这个是 p 值，可以转换为数据框保存
代码语言：javascript复制result$reps.overlap
result$reps.overlap %>%
rownames_to_column() %>%
write_tsv("reps.overlap.tsv")
这个是作图用到的数据
第一行是 结果图垂直虚线的位置
单独把这个数据保存成一个数据框然后输出应该就可以了，图用这个数据来做
其余的行是频率分布直方图的数据
代码语言：javascript复制result$d.plot
result$i.plot
这两个是可能会用到的结果图，暂时没有想明白单独保存图文件也会非常大，不知道为啥
下面运行一下这个R中的帮助文档，保存数据然后自己作图
帮助文档中的代码
代码语言：javascript复制library(ENMTools)
library(tidyverse)
cyreni <- iberolacerta.clade$species$cyreni
monticola <- iberolacerta.clade$species$monticola
cyreni$range <- background.raster.buffer(cyreni$presence.points, 100000, euro.worldclim)
monticola$range <- background.raster.buffer(monticola$presence.points, 100000, euro.worldclim)
identity.test(cyreni, monticola, env = euro.worldclim, type = "glm",
f = pres ~ bio1 bio12, nreps = 100) -> result
result$reps.overlap %>%
rownames_to_column() %>%
write_tsv("reps.overlap.tsv")
读取数据作图
代码语言：javascript复制read_tsv("reps.overlap.tsv") %>%
filter(rowname != "empirical") -> dat
p1<-ggplot(data=dat,aes(x=D))
geom_histogram(color="grey",fill="#009f73")
geom_vline(xintercept = 0.471,lty="dashed")
theme_bw(base_size = 15)
theme(panel.grid = element_blank())
annotate(geom = "text",x=0.6,y=50,label="P=0.0099",size=10)
p2<-ggplot(data=dat,aes(x=I))
geom_histogram(color="grey",fill="#0072b1")
geom_vline(xintercept = 0.471,lty="dashed")
theme_bw(base_size = 15)
theme(panel.grid = element_blank())
annotate(geom = "text",x=0.6,y=50,label="P=0.0099",size=10)
library(patchwork)
p1 p2
欢迎大家关注我的公众号
小明的数据分析笔记本
小明的数据分析笔记本 公众号 主要分享：1、R语言和python做数据分析和数据可视化的简单小例子；2、园艺植物相关转录组学、基因组学、群体遗传学文献阅读笔记；3、生物信息学入门学习资料及自己的学习笔记！
