bayenv
是一种贝叶斯方法,通过检测等位基因频率和生态变量之间的相关性或地理区域之间的极端等位基因频率差异来识别适应当地环境的位点。该方法可以从一组基因组变异中估计种群之间等位基因频率的协方差,在单个 SNP 上进行测试。
由于基因流,有效群体大小等群体演化历史的不同,造成了跨种群的等位基因频率有一定的相关性,使适应性受选择位点的检测变得复杂。如果结合于环境的相关性,可能更能解释位点的受选择情况。
分析案例:
软件下载
wget https://bitbucket.org/tguenther/bayenv2_public/src/master/bayenv2
输入文件格式介绍
SNPSFILE
包含每个SNP 在不同群体的等位基因数量,每一个位点由两行表示,分别是两个等位基因的数量。 等位基因 1 和等位基因 2 的总和为该群体样本量*2(不包括缺失数据)。该文件由
Tab
分隔,可以使用PGDSpider[Lischer & Excoffier, Bioinformatics 2012]把VCF
文件换成ENV
格式。
0 0 0 0 0
42 24 16 22 32
2 0 0 1 1
40 24 16 21 31
SNPFILE
与SNPSFILE
输入文件格式相同,但是只包含单个 SNP 的数据,顺序相同。
ENVIRONFILE
是环境变量文件。 每个环境变量在群体之间进行标准化,即减去平均值,然后除以标准差。 每一行都是在各群体中的环境变量,Tab
分隔。
MATRIXFILE
基因型协方差矩阵。 可以通过从程序的输出中复制和粘贴最后的矩阵。
SAMPLEFILE
是包含每个种群的样本大小的文件。 样本大小与它们在协方差矩阵中出现的顺序相同,即它们出现在 SNPSFILE
中的顺序。
这里展示了在52个群体里(列)的两个环境变量(行)的标准化数值。
-0.071 0.649 -0.963 -1.421 -1.008 0.359 -1.170 0.549 0.789 0.432 0.744 0.114 0.421 0.527 0.387 0.900 0.917 1.503 0.917 1.157 1.157 1.693 0.761 1.553 0.733 0.761 0.365 -1.666 -1.159 -0.098 -0.974 -0.439 -0.450 -1.147 -0.841 0.549 -0.361 -0.026 -1.198 0.739 -0.606 1.308 1.302 0.833 -0.662 -1.979 -1.599 -0.651 0.644 -1.655 -1.343 -1.276
0.16423 0.89006 0.18691 0.28331 -0.24972 0.93543 0.15856 0.02813 -0.26107 -1.07763 -0.51624 -1.80914 -1.33281 0.19825 -0.31777 0.73696 0.72562 2.08656 0.72562 1.87675 1.87675 -2.60302 0.14722 -1.28177 1.97882 0.14722 2.22265 0.25496 -0.41417 -0.55027 0.22661 -0.86782 -0.34045 -0.41984 -0.08528 -0.26107 -1.03227 -0.68069 0.02246 -1.41787 -0.00022 -0.74307 -0.73739 -1.41787 -0.83379 -0.10796 0.14155 1.06585 0.89006 0.11886 1.00348 0.26630
参数介绍
Option Description
-i SNPFILE: 输入SNP文件 (always required)
-e ENVIRONFILE: 环境变量文件 (required for test mode)
-m MATRIXFILE: 矩阵文件 (output of Bayenv2.0, required for test mode)
-s SAMPLEFILE: 每个群体的样本量 (only required for pool mode)
-k NUMRUNS: 迭代次数 (always required)
-r SEED: 随机种子 (integer)
-p NUMPOPS: 群体个数 (always required)
-n NUMENVIRON: 环境变量的个数 (required for test mode)
-t: rest mode (not matrix estimation), 对每一个SNP计算 Z, BF or ρ.
-x: pool mode (to be used when the input data originates from pooled sequenc- ing of populations), calculate Z instead of BF (when in test mode, only one environmental variable)
-o OUTFILE: (optional) 输出文件的名字 (only available in test mode)
-f: write standardized allele frequencies X into a file
-c: calculate ρ in addition to BF
-X: calculate X T X
-z: calculate Z in test mode for unpooled data (only one environmental variable)
矩阵估计
./bayenv2 -i SNPSFILE -p NUMPOPS -k 100000 -r 63479 > matrix.out
这将每 500 次迭代将协方差矩阵的当前绘制输出到 matrix.out 中。
协方差矩阵中的行和列以与它们在等位基因计数文件中出现的相同的总体顺序出现。 协方差矩阵可以通过R中的image
命令进行可视化,cov2cor
函数可以用于将协方差矩阵转换为相关矩阵。
环境相关性估计
准备环境变量文件
估计 SNPFILE 的等位基因频率和 ENVIRONFILE 中的五个环境变量的贝叶斯因子。 输入文件是每一个SNP。
./bayenv2 -i SNPFILE -m MATRIXFILE -e ENVIRONFILE -p NUMPOPS -k 100000 -n 5 -t -r 429
如果没有请求特定的 OUTFILE (-o),则 SNP 的结果贝叶斯因子将附加到名为 bf 的文件中。
结果文件:
rs316 (X) (X) (X) (X)
第一列给出了 SNPFILE, 后续列给出了每个环境变量的估计贝叶斯因子, 输出会附加到文件中,因此针对不同 SNP 的运行将按顺序出现在文件中。
以下脚本可以为SNPSFILE的所有SNPs转换成每个SNPFILE,计算BFs
./calc bfs.sh SNPSFILE ENVIRONFILE MATRIXFILE NUMPOPS NUMITER NUMENVIRON
它将 SNPSFILE 的所有 SNP 拆分为单独的文件,然后运行 Bayenv2。
非参数检验
贝叶斯因子背后的线性模型可能不正确,或者异常值可能会误导模型。对于这种情况,使用-c参数,Bayenv2.0 计算标准化等位基因频率 (X),从中去除群体之间的协方差结构,计算 Spearman 等级相关系数 ρ。
./bayenv2 -i SNPFILE -m MATRIXFILE -e ENVIRONFILE -p NUMPOPS -k 100000 -n 1 -t -c -r 24542
此示例计算单个环境变量的贝叶斯因子、Spearman 的 ρ 和 Pearson 的相关系数。如果没有指定 OUTFILE (-o),结果将附加到名为 bf 的文件中。第一列给出了 SNPFILE 的名称,然后是一列显示贝叶斯因子、一列显示 Spearman 的 ρ 和一列显示 Pearson 的相关系数 r。如果使用了多个环境变量,则附加列(BF、ρ 和 rS)将附加到输出文件中。
参考:
- 官网
- manual
-
Robust identification of local adaptation from allele frequencies.
Genetics 2013 Sep;195(1):205-20. Günther T, Coop G. -
Using Environmental Correlations to Identify Loci Underlying Local Adaptation.
GENETICS August 1, 2010 vol. 185 no. 4 1411-1423. Graham Coop, David Witonsky, Anna Di Rienzo and Jonathan K. Pritchard.