本文为CSDN博主「生信菜鸟1号」的原创文章,原文链接:https://blog.csdn.net/weixin_45694863/article/details/126801831
SMC++是一个用于从全基因组序列数据中估计种群大小历史的程序,其可以进行多个样本的分析,v1.15.4 以后smc++不支持conda了,改为docker运行。
1.SMC++安装
SMC++可以作为Docker image直接运行,其安装方法也及其简单,就是运行一遍,当docker找不到SMC++image时自动下载,需要一定时间。
#简直不要太简单
sudo apt-get install -y python3-dev libgmp-dev libmpfr-dev libgsl0-dev #下载依赖库,我忘记用docker安装需不需要这些库了,但是下载了也没啥问题。
sudo docker run --rm -v $PWD:/mnt terhorst/smcpp:latest #下载
sudo docker images #看一下docker有哪些image
如下如,MSC++ image下载完毕
看一下具体参数
sudo docker run --rm -v $PWD:/mnt terhorst/smcpp:latest -h
有效种群大小主要用到vcf2msc,eastimate,plot
2.输入文件
SMC++的输入文件可以是vcf文件,非常棒。首先需要得到你分析的群体的vcf文件
vcftoosl --vcf xxx.vcf --recode --keep poplation.txt --out smc #population中时你的群体样本
bgzip smc.recode.vcf #压缩
tabix smc.recode.vcf.gz #为vcf构建tbi索引
之后就是格式转换,vcf转smc格式,用到了vcf2msc
for i in {1..19}
do
sudo docker run --rm -v $PWD:/mnt terhorst/smcpp:latest vcf2smc ./smc.recode.vcf.gz ./data/chr${i}.smc.gz $i S:S1,S2,S4,S5,S6,S7,S9
done
#输入vcf文件,输出smc.gz文件,$i为染色体,样本的排列是 population name:sample1,sample2,sample3
#值得注意的是mask文件,在MSMC2中需要,在MSC++也可以添加
3.分析
sudo docker run --rm -v $PWD:/mnt terhorst/smcpp:latest estimate --spline cubic --knots 15 --timepoints 1000 1000000 --cores 20 -o ./data/estimate-1/ 2e-8 ./data/*.smc.gz
#–spline 线条类型 cubic为平滑曲线吧
#–knots 节点,就是绘图曲线的点数
#–timepoints 时间范围,单位为多少代,如1000 1000000 为1000代到1000000代
#cores 计算核数
#2e-8为突变率
4.绘图
sudo docker run --rm -v $PWD:/mnt terhorst/smcpp:latest plot ./plot-IN-6e-9.pdf ./data/estimate-1/*.final.json -g 1 --ylim 0 100000000 -c
#.final.json为最终模型
#-c --csv 输出绘图文件,利用该文件在R中绘图,具体自己发挥
#-g 代数,如人为25,鼠为1
#–ylim Y轴范围
#–xlim X轴范围
总结:smc++不仅能多个文件分析,速度还快,还支持vcf文件,非常好。
stairway plot 分析
————————————————