入门

在比较早之前就看过一篇简书

[EEMS推测生境内生物的基因流 - 简书 (jianshu.com)](https://www.jianshu.com/p/9ef469ed7c98)

介绍过EEMS这款软件，点了收藏，直接搁置。最近才慢慢参悟该软件的使用，踩过很多坑，所以直接介绍一下EEMS。

乍练

首先上官网：https://github.com/dipetkov/eems
使用手册在Documentation中，有更详细的说明。

通晓

EMMS的安装比较麻烦，每一个模块都需要单独编译，这里只介绍我用到的三个模块。

#EMMS 主体可以直接git clone 下载 略
#首先是软件主命令 runeems_snps
#需要两款其他软件作为支持 
#Eigen进行线性代数计算
https://eigen.tuxfamily.org/index.php?title=Main_Page
#下载后解压 安装过程直接 粘贴下来：
Let's call this directory 'source_dir' (where this INSTALL file is).
Before starting, create another directory which we will call 'build_dir'.
Do:
  cd build_dir
  cmake source_dir
  make install
The "make install" step may require administrator privileges.
You can adjust the installation destination (the "prefix")
by passing the -DCMAKE_INSTALL_PREFIX=myprefix option to cmake, as is
explained in the message that cmake prints at the end.

#Boost进行随机数生成和生境几何形状
https://www.boost.org/
(不记得了 可以google 一下 how to install boost)

#两款软件装好后，找到/eems/runeems_snps/src Makefile
#修改依赖路径 make 就好 例：
EIGEN_INC = /data/01/user152/software/eigen-3.2.10/eigeneigen/include/eigen3/
BOOST_LIB = /data/01/user152/software/boost/lib
BOOST_INC = /data/01/user152/software/boost/include

#然后是bed2diffs模块的编译
#需要先安装  libplinkio (https://github.com/fadern/libplinkio)
#这个比较简单，略
#同样找到Makefile 修改依赖路径（在src里） 例：
PLINKIO = /data/01/user152/software/plinkio
EXE = bed2diffs_v1
OBJ = bed2diffs_v1.o data.o

CXXFLAGS = -I/data/01/user152/software/plinkio/include -O3 -Wall -Werror -fopenmp
LDFLAGS = -lplinkio

#最后是plot模块，跟着手册装就好
#遇到rgos装不上的情况，去跟管理员py一下，让root装一下就ok了。

至此，软件就搞定了。

小成

runeems_snp 需要三个主要输入文件：
*.diffs 遗传距离矩阵，可以通过VCF转换
*.coord 采样的坐标，经纬度就可
*.outer 栖息地坐标，同经纬度
现在我们一个个准备输入文件：

.diff

此处坑较多，务必注意

#首先要确定VCF中的个体是不是每个都有采样信息，没有的话需要删掉，外群也需要删掉。
#可以用vcftools 删除个体

#我在做的过程中做过很多次，只介绍跑通的这次的做法，要问我为什么这么做，别问，问就是这么做能跑通
1. 需要将VCF处理为下种模样

1. 第一列需要都改为1
2. POS 顺延
3. ID 为snp1 顺延

#之后运行PLINK
plink --vcf sample.vcf --allow-extra-chr --maf 0.05 --recode --out sample
plink --noweb --file sample --allow-extra-chr --make-bed --maf 0.05 --allele1234  --out sample

生成的smple.[bed/bim/fam]就是下一步的生成文件。

bed2diffs_v1 --bfile  ./sample

会生成第一个输入文件， sample.diffs

.coord

记录个体采样地理位置的文件：

'经度' \t'纬度'
建议经度在前，当然纬度在前也可以，不过后面需要注意这个的前后顺序

.outer

记录栖息地位置的文件
说人话就是
把所有样本包括进去的一个地理范围，类似下图

采样地点都在范围内。
可以用此网站准备改文件。

http://www.birdtheme.org/useful/v3tool.html

到此为止，最主要的三个文件也就准备好了，还有一类可选的文件。详细的可以阅读手册。

大成

输入文件准备好了就可以运行软件了。是一个类似配置文件的形式。

#params-run.ini
datapath = ./data/sample  # 输入文件前缀的路径
mcmcpath = ./out   #输出文件路径
#gridpath =   #刚提到的可选的输入文件前缀路径
nIndiv =48  #个体数
nSites =5343964  #位点数
nDemes =700 #grid 网格数700我觉得有点太慢了。。。。
diploid =true  #是否为二倍体
numMCMCIter =2000000
numBurnIter =1000000
numThinIter =9999

因为用到的是SNP数据，三个输入文件内并不能体现位点数，因此我用的是真实的位点数，可不可以瞎编我也不知道。
Demes 我理解的是最后网格数的多少，越多越慢，示例为300好像，我采用700，具体见下图。
剩下的三个建议就默认吧

#运行
runeems_snps --params params-run.ini

圆满

输出文件都有了，最重要的一点就是可视化。

###这里手册里介绍的很详细，不多加赘述
library("rEEMSplots")

eems_results <- file.path("out-1")  #输出文件路径
name_figures <- file.path("out")  #输出图片文件名前缀

eems.plots(mcmcpath = eems_results,
           plotpath = paste0(name_figures, "-output-PDFs"),  #输出为PDF
           longlat = TRUE,  # TRUE 指经度在前 纬度在后 FALSE反之
           out.png = FALSE, #输出为PDF
           #add.grid = TRUE,  #可选输入文件才有
           #col.grid = "gray90",  #可选输入文件才有
           #lwd.grid = 2,  #可选输入文件才有
           #add.outline = TRUE,  # 是否显示边界 （我觉得太丑了）
           #col.outline = "blue",
           #lwd.outline = 5,
           add.demes = TRUE,  #显示demes
           col.demes = "red",
           pch.demes = 5,
           min.cex.demes = 0.5,
           max.cex.demes = 1.5
           )

最后结果大致这样：

大道

判断结果好不好的重要标准是： 模拟的模型是否收敛
可视化后有一个图能用来判断。
doc里是这么解释的：

The eems.plots function in the rEEMSplots package generates a posterior trace plot. (a) This MCMC
chain obviously has not converged. (b) There is no indication that the MCMC chain has not converged. This is
not quite the same as there is evidence that the MCMC chain has converged. (c) To be confident that EEMS has
converged, we can run the MCMC sampler for more iterations. (d) Alternatively, to be confident that EEMS has
converged, we can run the MCMC sampler several times, each time starting from a different randomly initialized
parameter state.

说人话呢就是 最差也得是b图的样子。
为了达到收敛，有两种方法

1. 多跑几次，选最好的
2. 用上一次的结果作为下一次的输入，循环
我们就草率的选2哈！

#params-pre.ini
datapath = ./data/sample
mcmcpath = ./out
prevpath = ./out-1
#gridpath = 
nIndiv =48
nSites =5343964
nDemes =700
diploid =true
numMCMCIter =1000000
numBurnIter =0
numThinIter =9999

跑到能让自己舒服的结果就好。
感谢 软件作者 及 DumplingLucky 还有KZR