2023-03-21 | Dsuite做D及f4-ratio统计
2023-03-20 本文已影响0人
千万别加香菜
Dsuite原理：
- D值（即ABBA统计量）和f4-ratio统计可以表示为适用于四个分类群的双等位基因SNP：P1,P2,P3,O，拓扑是 (((P1,P2),P3),O)。
- 其中外类群O携带祖先等位基因A，衍生等位基因用B表示。BBAA,ABBA,BABA分别代表四个分类群携带等位的三种模式。
- 在没有基因流的零假设下，由于具有相同频率的不完全谱系分类，预计P3与P1或P2共享衍生等位基因B的两种模式ABBA和BABA的频率相等，如果ABBA和BABA的频率有显著差异则代表在P3和P1或P2间存在基因渐渗。
- D=(nABBA-nBABA)/(nABBA+nBABA)；在外群对于祖先等位基因A是固定的（外群中B的频率为0）假设下，D统计量是等位基因模式计数的归一化差异。
- 如果外群中衍生等位基因B不为0，则Dsuite的D值是Patterson’s D，适用于无根的四分类群树。
基因流推断 —— Dsuite | 生信技工 (yanzhongsino.github.io)
Dsuite的使用：
/home/sll/software/Dsuite/Build/Dsuite
寻找系统中高版本的libstdc++.so: find / -name "libstdc++.so*"
激活环境：因为服务器的/lib64/libstdc++.so.6版本过老，所以这里使用我的conda下的版本
export LD_LIBRARY_PATH=/home/sll/miniconda3/pkgs/libstdcxx-ng-12.1.0-ha89aaad_16/lib:$LD_LIBRARY_PATH
准备两列的表格，第一列为样本ID，第二列为群体ID(需指定外群Outgroup，也就是树的根)，想屏蔽的个体用第二列xxx表示，注意群体ID不要包含“. - 空格”等字符，可以有下划线_
否则，Fbranch会报错
1、Dsuite Dtrios模块：
为所有可能的种群/物种三重组合计算D和f4-ratio统计量（ABBA-BABA）
建议分染色体计算，然后使用DtriosCombine 模块将各染色体结果合并
/home/sll/software/Dsuite/Build/Dsuite Dtrios sample-select.vcf d.txt -t sample.ML.tree.treeout -o sample
-t 物种树文件，可用treemix生成,根为outgroup，且m设为0，不考虑基因流
-o sample：指定输出文件前缀，默认是sets
-p 5：如果样品中包含pool-seq数据，-p用于设置最小深度，设置后从等位基因深度估计群体的等位基因频率。
D和f4-ratio结果包含在.Dmin.txt文件中
DtriosCombine 模块对Dtrios结果合并：
/home/sll/software/Dsuite/Build/Dsuite DtriosCombine -t sample.ML.tree.treeout -o sample_all DminFile1.txt DminFile2.txt DminFile3.txt
-o, --out-prefix=OUT_FILE_PREFIX 输出文件前缀，默认为 "out"
-n, --run-name 看不懂
-t , --tree=TREE_FILE.nwk 树文件
-s , --subset=start,length 只进行指定长度部分的合并
2、Dsuite Dinvestigate：
用于对感兴趣渗入组合的基因组区域的D值的计算，看哪些区域发生了渗入
/home/sll/software/Dsuite/Build/Dsuite Dinvestigate -w 50,25 INPUT_FILE.vcf.gz SETS.txt test_trios.txt
Outputs D, f_d, f_dM, and d_f in genomic windows
SETS.txt文件有两列 : SAMPLE_ID POPULATION_ID
test_trios.txt包含三个群体（除外群,外群在SETS.txt文件中已经指定）的名称:
POP1 POP2 POP3
There can be multiple lines and then the program generates multiple ouput files, named like POP1_POP2_POP3_localFstats_SIZE_STEP.txt
-h, --help display this help and exit
-w SIZE,STEP --window=SIZE,STEP (required)设置移动的窗口及步长大小 (default: 50,25)
-n, --run-name run-name will be included in the output file name
3、Dsuite Fbranch：
是一种启发式方法，执行f-branch计算，用于解释f4-ratio相关结果
/home/sll/software/Dsuite/Build/Dsuite Fbranch sample.ML.tree.treeout sample_tree.txt > fbranch.out
fbranch.out：f-branch统计量保存成矩阵格式
用dtools.py脚本绘制f-branch图
/home/sll/software/Dsuite/utils/dtools.py fbranch.out sample.ML.tree.treeout --outgroup Outgroup --use_distances --dpi 1200 --tree-label-size 30
–outgroup：指定外类群（与fbranch.out和species.newick一致，一般是Outgroup）
–use_distances：画树时使用newick文件里节点距离
–dpi：设置png分辨率，有些期刊投稿要求1200，800，600不等；最好高点。
–tree-label-size：设置树节点标签大小
结果展示：
Fbranch.png
- 真正的物种树作为数据模拟的输入文件显示在图的侧边。物种树在 y 轴以“展开”的形式进行展示，所以每一个分枝包括内部分枝都指向矩阵中对应的行和推断的 f-brach 统计值。
- 图的上方和左侧为群体/物种系统发育树，其中左侧为展开的群体/物种树。矩阵中色块颜色深浅表示渗入比例，颜色越深表示渗入比例越高，越浅渗入比例越低。
脚本自用：
Dsuite.sh
export LD_LIBRARY_PATH=/home/sll/miniconda3/pkgs/libstdcxx-ng-12.1.0-ha89aaad_16/lib:$LD_LIBRARY_PATH
Dsuite="/home/sll/software/Dsuite"
$Dsuite/Build/Dsuite Dtrios sample-select.vcf d.txt -t sample.ML.tree.treeout -o sample
$Dsuite/Build/Dsuite Fbranch sample.ML.tree.treeout sample_tree.txt > fbranch.out
$Dsuite/utils/dtools.py fbranch.out sample.ML.tree.treeout --outgroup Outgroup --use_distances --dpi 1200 --tree-label-size 30
