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.
经观察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那列,如果两个条件都不满足,只会显示其中一个条件不满足,而没有第二个条件不满足的信息并不会出现。
下面单独测试一下dps:
vcf_filter.py testfilter.vcf dps --depth-per-sample 10 > testfilter10.vcf
说明dps确实有问题。更郁闷的是不仅有些小于10的不会被过滤掉,连大于10的竟然却被过滤掉了。上图中DP=58被过滤掉了。。。。呵呵了 dps压根不用就是了。
下面单独测试一下avg-dps:
vcf_filter.py testfilter.vcf avg-dps --avg-depth-per-sample 10 > testfilteravg_10.vcf
经测试avg-dps是准确的,其实这个参数正是我们需要的,只要控制了每个样的深度,总的深度也会得到控制,若avg-dps为10, 也就是说总深度我们控制到了30.静平也测试了一下, 发现avg-dps也准确:
这个脚本就avg-dps可以用, 但是不符合我的需要,因为如果三个样本的深度分别是20 5 5 这种情况也满足于avg-dps为10,。
by freemao
FAFU