Pop-Con 1.0.5

Pop-Con：一种在现场频谱上可视化基因型分布的工具

简介

pop-con是一个工具，用于从snp和indel变体的变体调用格式（vcf）文件[1]中可视化站点频谱（sfs）的基因型剖面。 我们将基因型谱定义为给定位点上所有个体的基因型列表。 例如，"0/1,1/1,0/1,0/1"是对应于下面vcf线的基因型剖面：

Dst_b1v03_scaf00008 54329 . G A 1970.68 PASS AC=6;AF=0.500;AN=12;BaseQRankSum=-6.140e-01;ClippingRankSum=0.00;DP=144;ExcessHet=4.8280;FS=6.723;MLEAC=5;MLEAF=0.417;MQ=60.00;MQRankSum=0.00;QD=15.28;ReadPosRankSum=-8.460e-01;SOR=0.440 GT:AD:DP:GQ:PL 0/1:16,22:38:99:558,0,409 1/1:1,11:12:3:309,3,0 0/1:26,28:54:99:750,0,617 0/1:8,11:19:99:286,0,203

其中个体1、3和4为杂合子（"0/1"），个体3为该位点的替代等位基因（"1/1"）的纯合子。

pop-con接受一个vcf文件作为输入，并使用cyvcf2python包[2]读取和解析vcf文件 要通过cyvcf2读取和解析VCF文件，必须使用以下两个命令行之一进行压缩和索引：

bgzip input_vcf_file.vcf
tabix -p input_vcf_file.vcf.gz

或

bcftools view input_vcf_file.vcf -Oz -o input_vcf_file.vcf.gz
bcftools index input_vcf_file.vcf.gz

Pop-Con是在使用gatk[3]和read2snp[4]变量调用程序生成的VCF文件上开发的。必须为用于生成VCF文件的变量调用方提供选项--tool（gatk是默认变量调用方）。

--read和--marft筛选选项可用于测试不同的基因型筛选集：

• --read选项将读取覆盖率阈值列表作为输入。读取覆盖率低于阈值的基因型将被标记为lowcov
• --marft选项将次要等位基因读取频率阈值列表作为输入。少数等位基因读取频率低于阈值的基因型将被标记为lowmarf

对于每个过滤组合（"读取覆盖率"+"次要等位基因读取频率"），Pop-Con绘制2个图形（用于snp和indel）：

1. 所有个体都被基因分型并通过读取覆盖率和次要等位基因读取频率过滤的位点。
2. 所有的位点至少有一个个体被基因分型并通过了读取覆盖率和次要等位基因读取频率过滤。

对于图1，生成两个图：

• 用-max/--max-u profiles选项设置表示每个sfs-pic的x最典型基因型谱的观察比例的上图。
• 在hardy-weinberg平衡下，代表每个pic基因型分布的期望比例的较低的图。如果所观察到的基因型谱是上hardy-weinberg期望值的y倍，则基因型谱为红色；如果在平衡点（在y倍上和y倍下），则为蓝色；如果在平衡点（在y倍上和y倍下），则为白色。其中，y通过选项-fold/--hwe_fold_change设置

对于图2，只产生具有观察到的基因型分布比例的SFS图。此图可用于检测低阅读覆盖率、低次要等位基因阅读频率或无基因型个体的异常porpotion，但目前显示效果不理想…

Pop-Con在带有python2（版本==2.7）和python3（版本>=3.4）的Linux上受支持，当然在带有pypi python包管理器的Macintosh和Windows上也受支持，但未经测试。

安装

来源于Github

这假设您已经安装了python模块cyvcf2、numpy和matplotlib（cfrequirements.txt）。

git clone https://github.com/YoannAnselmetti/Pop-Con.git
cd Pop-Con
python ./Pop-Con

使用PIP（推荐）

sudo pip install Pop-Con

手动安装

在"setup.py"所在的弹出窗口文件夹中，运行：

sudo python setup.py install

如果在安装cyvcf2时遇到问题，或在执行弹出式菜单时收到以下错误消息：

Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/usr/local/lib/python3.6/site-packages/cyvcf2/__init__.py", line 1, in <module>
    from .cyvcf2 import (VCF, Variant, Writer, r_ as r_unphased, par_relatedness,
ImportError: /usr/local/lib/python3.6/site-packages/cyvcf2/cyvcf2.cpython-36m-x86_64-linux-gnu.so: undefined symbol: EVP_sha1

必须从源手动安装： https://github.com/brentp/cyvcf2/P>

用法

usage: Pop-Con [-h] -i VCF_FILE [-t VARIANT_CALLER]
               [-r READ_COVERAGE_THRESHOLD [READ_COVERAGE_THRESHOLD ...]]
               [-m MARFT [MARFT ...]] [-f VARIANT_CALLER_FILTERING]
               [-fmono MONOMORPH_FILTERING] [-p PREFIX] [-v VERBOSE]
               [-o OUTPUT_DIR] [-hf WRITE_HETEROZYGOSITY_FILE] [-sep SEP]
               [-max MAX_PROFILES] [-fold HWE_FOLD_CHANGE]

Pop-Con - A tool for Population genomic Conflicts detection.

optional arguments:
  -h, --help            show this help message and exit
  -i VCF_FILE, --input_file VCF_FILE
                        Variant Calling Format file containing variant calling
                        data.
  -t VARIANT_CALLER, --tool VARIANT_CALLER
                        Variant calling tool used to call variant. Values:
                        "read2snp" or "GATK". (Default: "GATK")
  -r READ_COVERAGE_THRESHOLD [READ_COVERAGE_THRESHOLD ...], --read READ_COVERAGE_THRESHOLD [READ_COVERAGE_THRESHOLD ...]
                        List of read coverage threshold filtering for each
                        genotype. Values range: [0,infinity[. (Default: 0)
  -m MARFT [MARFT ...], --marft MARFT [MARFT ...]
                        List of Minor Allele Read Frequency (MARF) threshold
                        for heterozygous genotypes filtering. Values range:
                        [0.0,0.5]. (Default: 0.0)
  -f VARIANT_CALLER_FILTERING, --variant_caller_filtering VARIANT_CALLER_FILTERING
                        Boolean to set if variant caller filtering is consider
                        or not. (if "True", sites with TAG column are not
                        empty are filtered out from analysis). (Default: True)
  -fmono MONOMORPH_FILTERING, --monomorph_filtering MONOMORPH_FILTERING
                        Boolean to filter out monomorph sites. ("True":
                        monomorph sites are filtered out from analysis. Reduce
                        drastically the execution time!!!). (Default: True)
  -p PREFIX, --prefix PREFIX
                        Experiment name (used as prefix for output files).
                        (Default: "exp1")
  -v VERBOSE, --verbose VERBOSE
                        Verbose intensity. (Default: 1)
  -o OUTPUT_DIR, --output OUTPUT_DIR
                        Output directory path. (Default: ./)
  -hf WRITE_HETEROZYGOSITY_FILE, --heterozygosity_file WRITE_HETEROZYGOSITY_FILE
                        Boolean to set if the heterozygosity files
                        (summarizing VCF file for each combination of read
                        coverage and MARF filtering) have to be written or
                        not. File used by script
                        "scripts/plot_SFS_geneotype_profiles.py" to plot SFS
                        without re-parsing VCF file. (Default: True)
  -sep SEP, --separator SEP
                        Separator used in genotype profiles. (Default: ",")
  -max MAX_PROFILES, --max_profiles MAX_PROFILES
                        Maximum number of genotype profiles displayed in SFS
                        plot. (Default: 10)
  -fold HWE_FOLD_CHANGE, --hwe_fold_change HWE_FOLD_CHANGE
                        Fold change value to define when an observed genotype
                        profile proportion is in excess/deficit compare to the
                        expected value under Hardy-Weinberg Equilibrium.
                        (Default: 2.0)

Source code and manual: http://github.com/YoannAnselmetti/Pop-Con

示例

下面是运行Pop-Con的最小命令行：

python Pop-Con -i input_vcf_file.vcf.gz

如果VCF文件包含snp和indel变量，则弹出窗口将在当前目录中创建两个目录"snp/"和"indel/"。

下面是目录"example/data/lineus_longissimus/"中VCF文件的示例：

Dst_b1v03_scaf00008 54329 . G A 1970.68 PASS AC=6;AF=0.500;AN=12;BaseQRankSum=-6.140e-01;ClippingRankSum=0.00;DP=144;ExcessHet=4.8280;FS=6.723;MLEAC=5;MLEAF=0.417;MQ=60.00;MQRankSum=0.00;QD=15.28;ReadPosRankSum=-8.460e-01;SOR=0.440 GT:AD:DP:GQ:PL 0/1:16,22:38:99:558,0,409 1/1:1,11:12:3:309,3,0 0/1:26,28:54:99:750,0,617 0/1:8,11:19:99:286,0,203

Dst_b1v03_scaf00008 54329 . G A 1970.68 PASS AC=6;AF=0.500;AN=12;BaseQRankSum=-6.140e-01;ClippingRankSum=0.00;DP=144;ExcessHet=4.8280;FS=6.723;MLEAC=5;MLEAF=0.417;MQ=60.00;MQRankSum=0.00;QD=15.28;ReadPosRankSum=-8.460e-01;SOR=0.440 GT:AD:DP:GQ:PL 0/1:16,22:38:99:558,0,409 1/1:1,11:12:3:309,3,0 0/1:26,28:54:99:750,0,617 0/1:8,11:19:99:286,0,203

Dst_b1v03_scaf00008 54329 . G A 1970.68 PASS AC=6;AF=0.500;AN=12;BaseQRankSum=-6.140e-01;ClippingRankSum=0.00;DP=144;ExcessHet=4.8280;FS=6.723;MLEAC=5;MLEAF=0.417;MQ=60.00;MQRankSum=0.00;QD=15.28;ReadPosRankSum=-8.460e-01;SOR=0.440 GT:AD:DP:GQ:PL 0/1:16,22:38:99:558,0,409 1/1:1,11:12:3:309,3,0 0/1:26,28:54:99:750,0,617 0/1:8,11:19:99:286,0,203

Dst_b1v03_scaf00008 54329 . G A 1970.68 PASS AC=6;AF=0.500;AN=12;BaseQRankSum=-6.140e-01;ClippingRankSum=0.00;DP=144;ExcessHet=4.8280;FS=6.723;MLEAC=5;MLEAF=0.417;MQ=60.00;MQRankSum=0.00;QD=15.28;ReadPosRankSum=-8.460e-01;SOR=0.440 GT:AD:DP:GQ:PL 0/1:16,22:38:99:558,0,409 1/1:1,11:12:3:309,3,0 0/1:26,28:54:99:750,0,617 0/1:8,11:19:99:286,0,203