毎度おなじみ、僕です。
このブログが生きていることに驚いています。
1個体のゲノムから個体群動態を推定するPSMCという方法。
https://github.com/lh3/psmc
今ではいろんな発展系もありますが、やっぱり理論的にも使い方的にもシンプルでわかりやすい(相対的に)原作を使いたくなるときもあります。
日本語解説はこちらが詳しいです。
https://sites.google.com/site/hiromimatsumae/過去の情報/集団遺伝学
が、新しいsamtools系だとこの通りにはいかず、古いバージョンを併用する必要が出てきて面倒です。
特に公式のこの部分。
"""
samtools mpileup -C50 -uf ref.fa aln.bam | bcftools view -c - | vcfutils.pl vcf2fq -d 10 -D 100 | gzip > diploid.fq.gz
"""
これの新しいバージョンに対応したコマンドは以下です。
"""
bcftools mpileup -C50 -Ou -f ref.fa aln.bam | bcftools call -c - | vcfutils.pl vcf2fq -d 10 -D 100 | gzip > diploid.fq.gz
"""
bcftoolsで代用できました。
お試しあれ。