pyvcf 中带的一个工具 比其他工具用着好些 其他filter我很信不过~~  自己写的功能又很有限 所以转投vcf_filter.py啦

Filtering a VCF file based on some properties of interest is a common enough operation that PyVCF offers an extensible script. vcf_filter.pydoes the work of reading input, updating the metadata and filtering the records.

usage: vcf_filter.py [-h] [--no-short-circuit] [--no-filtered]
              [--output OUTPUT] [--local-script LOCAL_SCRIPT]
              input filter [filter_args] [filter [filter_args]] ...

Filter a VCF file

positional arguments:
  input                 File to process (use - for STDIN) (default: None)

optional arguments:
  -h, --help            Show this help message and exit. (default: False)
  --no-short-circuit    Do not stop filter processing on a site if any filter
                        is triggered (default: False)
  --output OUTPUT       Filename to output [STDOUT] (default: <open file
                        '<stdout>', mode 'w' at 0x2b0f9435c150>)
  --no-filtered         Output only sites passing the filters (default: False)  # 如果不加这个参数,不通过的位点也会被保留下来,只不过在filter列有相应说明!

--local-script LOCAL_SCRIPT
                        Python file in current working directory with the
                        filter classes (default: None)

sq:
  Filter low quailty sites

--site-quality SITE_QUALITY
                        Filter sites below this quality (default: 30)

dps:
  Threshold read depth per sample

--depth-per-sample DEPTH_PER_SAMPLE
                        Minimum required coverage in each sample (default: 5)

avg-dps:
  Threshold average read depth per sample (read_depth / sample_count)

--avg-depth-per-sample AVG_DEPTH_PER_SAMPLE
                        Minimum required average coverage per sample (default:
                        3)

eb:
  Filter sites that look like correlated sequencing errors. Some sequencing
  technologies, notably pyrosequencing, produce mutation hotspots where
  there is a constant level of noise, producing some reference and some
  heterozygote calls. This filter computes a Bayes Factor for each site by
  comparing the binomial likelihood of the observed allelic depths under: *
  A model with constant error equal to the MAF. * A model where each sample
  is the ploidy reported by the caller. The test value is the log of the
  bayes factor. Higher values are more likely to be errors. Note: this
  filter requires rpy2

--eblr EBLR           Filter sites above this error log odds ratio (default:
                        -10)

snp-only:
  Choose only SNP variants

mgq:
  Filters sites with only low quality variants. It is possible to have a
  high site quality with many low quality calls. This filter demands at
  least one call be above a threshold quality.

--genotype-quality GENOTYPE_QUALITY
                        Filter sites with no genotypes above this quality
                        (default: 50)

懒得翻译 自己看吧

注意,里面有默认值:

比如sq 默认30, dps 默认5, avg-dps 默认3, eblr 默认-10, mgq默认50.

但是经过测试, 不指定的话, 好像也不会针对这些项进行过滤, 所以还是指定的好!

下面举一些具体的例子:

vcf_filter.py testfilter.vcf sq 30 > testfilter30.vcf

会发现结果文件中,不满足要求的在filter列显示的是sq30, 满足过滤条件的显示的是PASS.

vcf_filter.py-LMLPHP

经观察sq参数是准确可靠地。

vcf_filter.py testfilter.vcf sq --site-quality 30 dps --depth-per-sample 10 > testfilter30_10.vcf

对sq 30 和 dps10 进行过滤:

发现虽然有的dp是小于10的但让然被PASS了,这很尴尬。。。由于我过滤的条件有两个,但是在filter那列,如果两个条件都不满足,只会显示其中一个条件不满足,而没有第二个条件不满足的信息并不会出现。

vcf_filter.py-LMLPHP

下面单独测试一下dps:

vcf_filter.py testfilter.vcf  dps --depth-per-sample 10 > testfilter10.vcf

vcf_filter.py-LMLPHP

说明dps确实有问题。更郁闷的是不仅有些小于10的不会被过滤掉,连大于10的竟然却被过滤掉了。上图中DP=58被过滤掉了。。。。呵呵了 dps压根不用就是了。

下面单独测试一下avg-dps:

vcf_filter.py testfilter.vcf avg-dps --avg-depth-per-sample 10 > testfilteravg_10.vcf

vcf_filter.py-LMLPHP

经测试avg-dps是准确的,其实这个参数正是我们需要的,只要控制了每个样的深度,总的深度也会得到控制,若avg-dps为10, 也就是说总深度我们控制到了30.静平也测试了一下, 发现avg-dps也准确:

这个脚本就avg-dps可以用, 但是不符合我的需要,因为如果三个样本的深度分别是20 5 5 这种情况也满足于avg-dps为10,。

by freemao

FAFU

[email protected]

04-26 14:21