原理是筛选Fst、dxy均在前1%的区域,并且pi值岛屿区与背景区无差异(P < 0.05)
######run fst-dxy-no-pi pipeline######
1. run combine-plot.py 得到窗口对齐的文件
exp python3 combine-plot.py ASPI.list ESPI.list dxy.listfst.list > combine.list
2. run 01.1top-value.py 得到Fst dxy 前1% 的值
exp python3 01.1top-value.py combine.txt
3. run 02.call-spceia-island.py 得到 Fst dxy 均在top 1% 的文件 以及背景区文件
exp python3 02.call-spceia-island.py combine.txt 0.324281(dxy value) 0.543(fst value) top low
4. 对 top 文件 pi 值进行排序 (AS 高到底;ES 高到底)
window 运行 得到 ASsorttop ESsorttop
5. run 03.fst-dxy-filter.py 得到P < 0.05 的top文件所包含的行数
exp python3 03.fst-dxy-filter.py ESsorttop low 2 0.071
exp python3 03.fst-dxy-filter.py ASsorttop low 2 0.063
head -n 1062 ASsorttop > ASnopiTop
head -n 1053 ESsorttop > ESnopiTop
6. 两文件取交集
exp sort ASnopiTop ESnopiTop | uniq -d | sort -k 1,1 -k 2,2n > finaltop
7. run 04.new-back.py 得到新的背景区
exp python3 04.new-back.py finaltop combine.txt finallow
8. 准备画图输入数据
cut -f 3,4 finaltop > t.tmp
cut -f 3,4 finallow > f.tmp
paste f.tmp t.tmp > plot.tmp
cat plot.in| awk '{print$3"\t"$1"\t"$4"\t"$2}' > in.plot
9. 画图 [显著性在第5步可看]
exp Rscript 05.plot.R in.plot
10. 得到最后的窗口文件
cat finaltop | cut -f 1,2 | awk '{print$1"\t"$2-1"\t"$2+9999}' > region.tmp
bedtools merge -i region.tmp > spc_region.bed
网友评论