A variant filtering method based on multi-sample whole exome sequencing

Abstract

The present invention provides a variation filtering method based on multi-sample whole exon sequencing, and belongs to the technical field of second-generation sequencing variation analysis. The present invention includes the steps of constructing a variation blacklist, filtering low-quality and high-frequency variations, performing variation filtering according to the blacklist, variation annotation, and sample splitting. The present invention can filter out frequently occurring false positive sites and polymorphic variation sites by running a bcftools filtering command without relying on other databases and variation annotation software, making full use of the quality and frequency information of the sample data itself, and significantly reducing the number of invalid sites that need to be analyzed in the next step of analysis. The present invention reduces the consumption of computer resources, shortens the running time, and avoids the occurrence of false negative pathogenic variation sites while filtering out benign and high-frequency variations.

Description

Technical Field

[0001] The present invention belongs to the technical field of second-generation sequencing variation analysis, and specifically relates to a variation filtering method based on multi-sample whole exome sequencing. Background Art

[0002] With the development of high-throughput sequencing technology, the number of bases that can be produced by a single chip in a single sequencing can exceed 1T. Based on 10G bases of a single whole exome sequencing sample, it can accommodate more than 100 samples. Whole Exome Sequencing (WES) refers to a technology based on high-throughput sequencing after DNA capture and enrichment of the exon region of the human genome.

[0003] The original WES sequence undergoes bioinformatics analysis including sequence preprocessing, sequence alignment, variant identification, information annotation and variant filtering. Genetic analysts then comprehensively consider information such as population frequency, variant quality, the extent of the variant's impact on the protein, phenotypic associations and case reports to find 1-2 variant sites that cause the corresponding disease from millions of variants. The above is the current common WES data analysis solution.

[0004] The main filtering methods for WES variant sites currently include the following steps: 1) download the GnomAD database and construct a population frequency file that meets the requirements of the annotation software; 2) download the ClinVar database and construct a variant pathogenicity file that meets the requirements of the annotation software; 3) use variant annotation software, such as VEP (Ensembl Variant Effect Predictor), to perform database annotation on multiple vcf files output by variant identification, and obtain the population frequency and pathogenicity information of each variant site; 4) according to the ACMG guidelines, construct quality and population frequency thresholds and pathogenicity rules, filter out low-quality, high-frequency and benign variants, and obtain variant result files that can be used for medical genetic analysis.

[0005] In summary, the problems of existing mutation filtering solutions are as follows:

[0006] (1) Variant annotation is a process that consumes a lot of computer resources and time. Using variant annotation software to annotate the population frequency and pathogenicity information of each sample's vcf file will result in the same variant sites between samples being annotated multiple times, resulting in a large waste of computer resources and time.

[0007] (2) The quality of the second-generation sequencing technology varies between sequencing batches. The filtering scheme based on a single sample makes it difficult to filter out false positive sites in a specific genomic region caused by library construction, probe capture or sequencing in a single batch.

[0008] (3) Failure to set a whitelist during the filtering process may result in missing possible pathogenic variants. Summary of the invention

[0009] In view of the problems existing in the prior art, the present invention provides a variation filtering method based on multi-sample whole exome sequencing. The present invention constructs a blacklist and a whitelist of variation sites, combines the variation quality and variation detection frequency of the samples to be analyzed, and performs joint filtering on the variation sites of multiple samples.

[0010] The present invention first provides a variation filtering method based on multi-sample whole exome sequencing, which comprises the following steps:

[0011] 1) Select the variants with the CLINSIG fields of Pathogenic, Likely_pathogenic and Uncertain_significance in the variant classification database ClinVar as the whitelist;

[0012] The variants with a population frequency greater than 0.02 in the GnomAD database were selected, and the union of the variants with the CLNREVSTAT field of practice_guideline / reviewed_by_expert_panel and the CLINSIG field of Benign / Likely_benign in the ClinVar database was taken, and the difference between the whitelist and the whitelist was taken to construct the variant blacklist;

[0013] 2) Select quality threshold and frequency threshold to filter low-quality and high-frequency variants;

[0014] 3) Perform mutation filtering based on the blacklist;

[0015] 4) Annotate the results obtained in step 3);

[0016] 5) Split the result obtained in step 4) by sample to obtain the annotated result for each sample.

[0017] As a preferred embodiment of the present invention, the step 1) specifically comprises the following steps:

[0018] 1.1 Download the WES vcf file provided by the GnomAD database, extract the variant sites whose FILTER field is PASS and whose population frequency is greater than 0.02, and generate a "temporary blacklist 1";

[0019] 1.2 Download the vcf file provided by the ClinVar database, extract the sites whose CLNREVSTAT field is practice_guideline or reviewed_by_expert_panel and whose CLINSIG field is Benign or Likely_benign, and generate a "temporary blacklist 2";

[0020] 1.3 Take the union of the sites in "temporary blacklist 1" and "temporary blacklist 2" to generate "temporary blacklist 3";

[0021] 1.4 Extract the sites with the fields of Pathogenic, Likely_pathogenic and Uncertain_significance in the CLINSIG of the clinvar database to generate a “white list 4”;

[0022] 1.5 Take the difference between "temporary blacklist 3" and "whitelist 4" to generate "blacklist 5" in vcf format.

[0023] As a preferred embodiment of the present invention, step 2) is specifically:

[0024] Analyze the WES joint variant data of multiple samples. If the variant site meets any of the following conditions, it will be deleted, and the final retained result will be output as "Filter Result 1":

[0025] A. The maximum sequencing coverage depth DP of each sample is <10;

[0026] B. Total allele number AN>12 and variant frequency AF>0.5;

[0027] C. Total allele number AN>30 and variant frequency AF>0.2;

[0028] D. Total allele number AN>100 and variant frequency AF>0.1;

[0029] E. The number of samples with genotype GT = "het" is > 10, and the maximum value of the variation ratio VAF is < 0.2.

[0030] The present invention uses the information in the result vcf file of multi-sample joint identification, and runs a bcftools filtering command without relying on other databases and variation annotation software, so as to filter out frequently occurring false positive sites and polymorphic variation sites, make full use of the quality and frequency information of the sample data itself, and significantly reduce the number of invalid sites that need to be analyzed in the next step of analysis. The present invention pre-constructs a variation blacklist, realizes independence from variation annotation software, runs a bcftools filtering command once, can filter the variation sites of multiple samples, reduces the consumption of computer resources, shortens the running time, and avoids the occurrence of false negative pathogenic variation sites while filtering out benign and high-frequency variations.

[0031] In the process of constructing the variant blacklist, the present invention selects variants whose CLINSIG fields are Pathogenic, Likely_pathogenic and Uncertain_significance in the variant classification database ClinVar as the whitelist to prevent possible pathogenic variants from being missed in the blacklist filtering process. BRIEF DESCRIPTION OF THE DRAWINGS

[0032] Figure 1 A schematic diagram of a process for constructing a mutation blacklist according to the present invention;

[0033] Figure 2 The figure is a flow chart of the variation filtering method of the present invention. DETAILED DESCRIPTION

[0034] The present invention is further described below in conjunction with specific embodiments, but the protection scope of the present invention is not limited thereto.

[0035] The relevant terms involved in the present invention are explained as follows:

[0036] ACMG guidelines: specifically refers to the guidelines for the classification of human genome variations published by ACMG (The American College of Medical Genetics and Genomics) in the journal Genetics of Medicine in March 2015. For details, see reference 1

[0037] GnomAD: The full name of GnomAD is The Genome Aggregation Database. It is an academic alliance organization that collects and organizes various large-scale exome and whole genome sequencing data, and makes them available to the world for free.

[0038] ClinVar: An open genomic variation database to which every research institution can submit data. The submitted information will be reviewed and rated by a team of experts. The sites in the database are divided into 1 to 4 different star ratings based on the reliability of the annotation information. The higher the star rating, the higher the credibility.

[0039] vcf: The full name of the English name is The Variant Call Format, which is a standard file format used to store variant site information.

[0040] bcftools: A software for processing variant-related files such as vcf and vcf files.

[0041] VAF: Variant Allele Fraction, the value range is [0,1], and the calculation method is VAF = depth (ALT) / depth (ALT + REF), that is: variation ratio = the number of sequences that have mutated at the site / the total number of sequences covering the site. When it is a heterozygous variation, the closer the VAF is to 0.5, the higher the quality, and the closer it is to 0, the lower the quality; when it is a homozygous variation, the closer the VAF is to 1, the higher the quality, and the closer it is to 0, the lower the quality.

[0042] The present invention proposes a variation filtering method based on multi-sample whole exome sequencing, and its overall process includes the following steps:

[0043] 1) Construct a mutation blacklist;

[0044] 2) Select frequency threshold and quality threshold to filter low-quality and high-frequency variants;

[0045] 3) Perform mutation filtering based on the blacklist;

[0046] 4) Annotate the results obtained in step 3);

[0047] 5) Split the result obtained in step 4) by sample to obtain the annotated result for each sample.

[0048] Reference Figure 1 In the embodiment of the present invention, the specific steps of constructing the mutation blacklist include:

[0049] 1. Download the WES vcf file provided by the GnomAD database, use the bcftools tool to extract the variant sites whose FILTER field is PASS and whose population frequency is greater than 0.02, and generate a "temporary blacklist 1";

[0050] 2. Download the vcf file provided by the clinvar database, and use the bcftools tool to extract the sites whose CLNREVSTAT (review status) field is practice_guideline or reviewed_by_expert_panel, and whose CLINSIG field is Benign or Likely_benign, to generate a "temporary blacklist 2";

[0051] 3. Take the union of the sites in "temporary blacklist 1" and "temporary blacklist 2" to generate "temporary blacklist 3";

[0052] 4. Extract the CLINSIG field (pathogenicity conclusion) in the clinvar database to generate a "white list 4" for the sites of Pathogenic, Likely_pathogenic and Uncertain_significance;

[0053] 5. Take the difference between "temporary blacklist 3" and "whitelist 4" to generate "blacklist 5" in vcf format.

[0054] In an optional embodiment of the present invention, the step 2) comprises the following steps:

[0055] 2.1 Select an appropriate frequency threshold: Select 500 retrospective positive samples, and select multiple sample subsets from them according to different sample number gradients of 6, 15, and 50. Set different variant detection frequency thresholds for each gradient, sort the number of variants retained in each gradient corresponding to each subset with a positive detection rate of 100% from small to large, and select the threshold near the first quartile as the frequency threshold of the corresponding gradient;

[0056] 2.2 Select an appropriate quality threshold: Set multiple VAF thresholds for 500 positive samples, sort the number of variants retained in each gradient of each subset corresponding to a positive detection rate of 100% from small to large, and select the threshold near the first quartile as the VAF threshold of the corresponding gradient;

[0057] 2.3 Use the -e parameter of the view function of the bcftools tool to jointly identify variant information based on multiple samples, filter out variants with quality lower than the threshold, coverage depth lower than the threshold, and detection frequency higher than the corresponding threshold in the samples analyzed this time, and output "filter result 1".

[0058] The present invention is further described below in conjunction with a specific embodiment. Figure 2 , the WES joint identification variation data of 100 samples were analyzed. The specific steps are:

[0059] 1. Filter low-quality and high-frequency variants:

[0060] Use the view function of the bcftools software to determine if the variant sites in the WES joint identification variant data "vcf file 1" of these 100 samples meet any of the following conditions, then delete them and output "filter result 1":

[0061] A. The maximum DP (sequencing coverage depth) of each sample is <10;

[0062] B. AN (total number of alleles)>12 and AF (variant frequency)>0.5;

[0063] C. AN (total number of alleles)>30 and AF (variant frequency)>0.2;

[0064] D. AN (total number of alleles)>100 and AF (variant frequency)>0.1;

[0065] E. The number of samples with GT (genotype) = "het" is > 10, and the maximum VAF (variance ratio) is < 0.2.

[0066] 2. Filter mutations based on mutation blacklist:

[0067] Use the isec function of the bcftools software to take the difference between the mutation list in "Filter Result 1" and the mutation list in "Blacklist 5" constructed in the previous step, and output "Filter Result 2"

[0068] 3. Variant information annotation:

[0069] The variant annotation software VEP (Ensembl Variant Effect Predictor) was used to annotate the population frequency, pathogenicity, gene, disease and other information of "Filter Result 2" and output "Annotation Result 1";

[0070] 4. Split by sample:

[0071] The -s parameter of the view function of the bcftools software was used to split the "annotation result 1" by sample, and the result file "annotation result 2" after annotation of each sample was obtained for manual screening and medical genetics analysis.

[0072] In order to evaluate the effect of the present invention, the present invention analyzes the WES joint identification variation data of 100 samples, adopts the existing technology solution and the present invention solution on the same Linux server, adopts 10 processes in parallel, and each process is set to run in a 10-thread mode, and successively analyzes the variation filtering link, and compares the running time of the two solutions and the difference in the list of variant sites after filtering. The results are shown in Table 1:

[0073] Table 1

[0074]

[0075]

[0076] As can be seen from Table 1, the present invention can reduce the consumption of computer resources, shorten the running time, and avoid the occurrence of false negative pathogenic variant sites while filtering out benign and high-frequency variants.

[0077] The above is only a preferred specific embodiment of the present invention, but the protection scope of the present invention is not limited thereto. Any changes or substitutions that can be easily thought of by a person skilled in the art within the technical scope disclosed by the present invention should be included in the protection scope of the present invention. Therefore, the protection scope of the present invention should be based on the protection scope of the claims.

Claims

1. A variation filtering method based on multi-sample whole exome sequencing, characterized in that: The steps include: 1) Select the variants with the CLINSIG fields of Pathogenic, Likely_pathogenic and Uncertain_significance in the variant classification database ClinVar as the whitelist; The variants with a population frequency greater than 0.02 in the GnomAD database were selected, and the union of the variants with the CLNREVSTAT field of practice_guideline / reviewed_by_expert_panel and the CLINSIG field of Benign / Likely_benign in the ClinVar database was taken, and the difference between the whitelist and the whitelist was taken to construct the variant blacklist; The step 1) specifically includes the following steps: 1.1 Download the WES vcf file provided by the GnomAD database, extract the variant sites whose FILTER field is PASS and whose population frequency is greater than 0.02, and generate "temporary blacklist 1"; 1.2 Download the vcf file provided by the ClinVar database, extract the sites whose CLNREVSTAT field is practice_guideline or reviewed_by_expert_panel and whose CLINSIG field is Benign or Likely_benign, and generate "temporary blacklist 2"; 1.3 Take the union of the sites in "temporary blacklist 1" and "temporary blacklist 2" to generate "temporary blacklist 3"; 1.4 Extract the sites with the fields of CLINSIG in the clinvar database as Pathogenic, Likely_pathogenic and Uncertain_significance to generate "white list 4"; 1.5 Take the difference between "temporary blacklist 3" and "whitelist 4" to generate "blacklist 5" in vcf format; 2) Select quality threshold and frequency threshold to filter low-quality and high-frequency variants; The step 2) specifically includes: 2.1 Select N retrospective positive samples, select multiple sample subsets according to different sample number gradients, set different variant detection frequency thresholds for each gradient, sort the number of variants retained in each gradient corresponding to each subset with a positive detection rate of 100% from small to large, and select the threshold near the first quartile as the frequency threshold of the corresponding gradient; 2.2 Set multiple VAF thresholds for N positive samples, sort the number of variants retained in each gradient of each subset corresponding to a positive detection rate of 100% from small to large, and select the threshold near the first quartile as the VAF threshold of the corresponding gradient; 2.3 Use the -e parameter of the bcftools tool view function to analyze the WES joint variant data of multiple samples, filter out variants with variant quality lower than the threshold, coverage depth lower than the threshold, and detection frequency higher than the corresponding threshold in the samples analyzed this time, and output "filter result 1"; 3) Filter mutations based on blacklist; The step 3) is specifically as follows: using the -e parameter of the view function of the bcftools tool, filtering out the variants in "filter result 1" that have an intersection with "blacklist 5" to generate "filter result 2"; 4) Annotate the results obtained in step 3); The step 4) is specifically as follows: using the variation annotation software VEP to annotate the population frequency, pathogenicity, gene, and disease information of "filter result 2" and output "annotation result 1"; 5) Split the results obtained in step 4) by sample to obtain the annotated results for each sample; The step 5) is specifically as follows: using the bcftools tool, the "annotation result 1" is split according to the sample to obtain the result file "annotation result 2" after the annotation of each sample, which is used for subsequent manual screening and medical genetics analysis.

2. The variation filtering method based on multi-sample whole exome sequencing according to claim 1, characterized in that: In step 1.1, the bcftools tool is used to extract variant sites where the FILTER field is PASS and the population frequency is greater than 0.02; In step 1.2, the bcftools tool was used to extract sites whose CLNREVSTAT field was practice_guideline or reviewed_by_expert_panel and whose CLINSIG field was Benign or Likely_benign.

3. The variation filtering method based on multi-sample whole exome sequencing according to claim 1, characterized in that: The step 2) is specifically as follows: Analyze the WES joint variant data of multiple samples. If the variant site meets any of the following conditions, it will be deleted, and the final retained result will be output as "Filter Result 1": The maximum sequencing coverage depth DP of each sample is <10; The total number of alleles AN>12 and the variant frequency AF>0.5; The total number of alleles AN >30 and the variant frequency AF>0.2; The total number of alleles AN >100 and the variant frequency AF>0.1; The number of samples with genotype GT="het" is >10, and the maximum value of the variation ratio VAF is <0.2.

Images