Fst/Pi/Dxy的计算、基因岛屿与曼哈顿图(1) - 个人笔记
Posted
tags:
篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了Fst/Pi/Dxy的计算、基因岛屿与曼哈顿图(1) - 个人笔记相关的知识,希望对你有一定的参考价值。
参考技术A 待更新 2020.5.10-1.把Fst值较高的区域(比如前百分之1)当成岛屿(或者Dxy值,这里以Fst为例) (但为什么要这样做我并不明白)。
0. Fst/Pi/Dxy是什么
略
1.Fst值(还有Dxy,Pi)的计算
1.1 工具
https://github.com/simonhmartin/genomics_general 使用的是这里的脚本(也可以用vcftools)
wget https://github.com/simonhmartin/genomics_general/archive/refs/heads/master.zip解压即可
1.2 转换格式
./VCF_processing/parseVCF.py -i try.vcf -o try.geno 得到try.geno文件
.geno文件有这几种格式,对应的,更改下一步的 -f 参数 我这里是-f phased
1.3 计算
python popgenWindows.py -w 20000 -m 200 --windType coordinate -g try.geno -o try.csv -f phased -T 5 -p A A1,A2,A3,A4,A5,A6 -p B B1,B2,B3,B4,B5,B6,B7,B8,B9 --popsFile pop.file
-w是窗口值大小;-g是输入文件即1.2的输出文件;-o即输出文件;-m即使MinSites,这里-m 200的意思是:如果这个窗口里面的sites(SNP)低于200,就把这个窗口扔掉;--windType即窗口类型,coordinate 这个类型是指窗口的坐标恒定,但是里面的sites数量可能不同;;另一种窗口类型是 sites ,即每个窗口里面的sites数量相同,但是坐标不同;-T是线程数;-p则指定群体以及群体包含了哪些个体。 --popsFile这个文件的内容是两列,一列是个体名,一列是群体名。
得到try.csv
1.4 结果如下
1.5 画图(R)
对于Pi
待补充,正在学习
对于Fst和Dxy可以简单粗暴的把一个窗口当成一个SNP,手动加一列就行,然后拿去画图。
qqman、CMplot这两个包均可实现
https://www.jianshu.com/p/1edc44fe922f qqman
https://www.burning.net.cn/article/article-29 CMplot
但是我的疑问是在于一个窗口中有很多个Sites,如果Sites = SNP,那这样的画法明显不合适
遗传分化一些基本概念
群体遗传学中衡量群体间分化程度的指标有很多种,最常用的就是Fst指数。
Fst指数,由F统计量演变而来。F统计量(FIS,FIF,FST)主要有三种。
Fst是针对一对等位基因,如果基因座上存在复等位基因,则需要用Gst衡量,基因差异分化系数(gene differentiation coefficient,Gst)。
假定有s个地方群体,第k个地方群体相对大小为wk,第k个地方群体中第i个等位基因频率为qk(i),杂合体频率观察值为hk,那么整个群体中观察到的杂合体频率平均值HI,地方群体为理想群体的期望杂合体频率平均值HS,整个群体为理想群体的期望杂合体频率HT,分别为:
- FIS,是HI相对于HS减少量的比值,即地方群体的平均近交系数。
- FST,是HS相对于HT减少量的比值,即有亲缘关系地方群体间的平均近交系数。
- FIT,是HI相对于HT减少量的比值,即整个群体的平均近交系数。
三者在数量上的关系
从配子间亲缘关系角度分析,FST和FIT分别相当于地方群体和整个群体中携带的一对等位基因是同源的概率,而FST是从两个地方群体中任意抽取的两个配子是同源的概率。从两个地方群体中任意抽取的两个配子是同源的概率大,表明两个地方群体的遗传组成相似,分化程度低;反义,分化程度高。
FST取值范围[0,1],最大值为1,表明等位基因在各地方群体中固定,完全分化;最小值为0,意味着不同地方群体遗传结构完全一致,群体间没有分化。
Fst值分析遗传分化
Wright建议,实际研究中,FST为0~0.05:群体间遗传分化很小,可以不考虑;
FST为0.05~0.15,群体间存在中等程度的遗传分化;
FST为0.15~0.25,群体间遗传分化较大;
FST为0.25以上,群体间有很大的遗传分化。
以上是关于Fst/Pi/Dxy的计算、基因岛屿与曼哈顿图(1) - 个人笔记的主要内容,如果未能解决你的问题,请参考以下文章