Method and electronic device for detecting somatic variation from a single sample
By calculating the parameters of the variation sites in a single sample and building feature vectors, and combining the somatic cell variant discrimination model, the problem of difficulty in detecting somatic cell variants in the prior art is solved, achieving high accuracy and low cost detection effects.
Patent Information
- Application Number
- CN202310520615.6
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2023-05-09
- Publication Date
- 2025-03-18
- Estimated Expiration
- 2043-05-09
AI Technical Summary
The prior art is difficult to accurately detect somatic variation in single samples, especially small panel data, and requires paired tissues for sequencing, which increases the cost and difficulty of data acquisition.
By analyzing the down-machine data of the detection sample, the germline variant parameters and somatic variation parameters of the mutation site were calculated, the mutation site characteristic vector was constructed, and the type of the mutation site was determined using the somatic variation discrimination model.
It realizes accurate detection of somatic variability in a single sample without paired tissue sequencing, reducing costs and data requirements, while improving detection accuracy.
Smart Images

Figure CN116769888B_ABST
Abstract
Description
Technical Field
[0001] The present invention relates to the field of bioinformatics, and in particular to a method and an electronic device for detecting somatic cell variation from a single sample. Background Art
[0002] In the field of precision tumor treatment, NGS sequencing technology is being used more and more widely. Currently, major public databases have accumulated sequencing data for multiple types of cancer. Based on the sequencing results of tumor patients, doctors can determine whether the patients are suitable for targeted drug treatment. In this process, distinguishing germline mutations from somatic mutations from the mutation results detected in patients is a key step.
[0003] The most common method to distinguish whether a variant site is a germline variant or a somatic variant is to sequence the tumor tissue and the paired normal sample at the same time, and filter the germline variants based on the paired samples. However, this will bring several problems: paired tissues are often not available in the clinic; sequencing paired tissues at the same time will incur additional costs; and paired sample data is often missing when retrospectively studying previous data. This requires us to create a method that can analyze variant sites based on separate tumor tissues and detect somatic variants therein.
[0004] At present, the filtering of germline variants is mainly based on database filtering and the SGZ algorithm published by Foundation Medicine. However, the current major public databases such as 1000Genome, ExAC, and dbSNP have not collected data for large cohorts of the Chinese population, and it is difficult to filter out all germline variants based on the database. There are still methods that filter based on the mutation frequency of the mutation site, which assumes that the frequency of germline variants will be close to 0.5 or 1, and the judgment is often distorted. The SGZ algorithm published by Foundation Medicine requires the establishment of a probability model to simultaneously estimate the purity, ploidy, and local copy number of the tumor, and is not applicable to small panel data. Summary of the invention
[0005] The main purpose of the present invention is to provide a method and an electronic device for detecting somatic cell variation from a single sample, so as to solve the problem in the prior art that it is difficult to accurately detect somatic cell variation in a single sample, especially in small panel data.
[0006] In order to achieve the above object, according to a first aspect of the present invention, a method for detecting somatic mutations from a single sample tumor tissue is provided, the method comprising:
[0007] (1) Analyze the offline data of the test samples, obtain the filtered variant sites, and calculate the posterior probability that the filtered variant sites are germline variants and somatic variants, which are recorded as germline variant parameters and somatic variant parameters;
[0008] (2) Using the background noise parameters, germline variation parameters, and somatic variation parameters of each filtered variation site, a mutation site feature vector is constructed;
[0009] (3) Using the somatic variation discrimination model, each filtered variation site is determined to be a somatic variation or a germline variation based on the mutation feature vector of each filtered variation site.
[0010] Furthermore, (1) includes: (1.1) obtaining the offline data of the test sample, comparing the offline data with the reference genome, and obtaining BAM file data, wherein the test sample is a single sample tumor tissue; (1.2) analyzing the BAM file data to obtain the initial mutation site, the depth of the initial mutation site, and the number of reads carrying the initial mutation site; (1.3) filtering the initial mutation site to obtain the filtered mutation site, the depth of the filtered mutation site, and the number of reads carrying the filtered mutation site; (1.4) calculating the germline mutation parameters and the somatic mutation parameters.
[0011] Furthermore, the calculation of germline variation parameters and somatic variation parameters includes:
[0012] a. dividing the genes of the reference sample and the test sample into segments, and equally dividing each gene of the reference sample and the test sample to obtain multiple windows; calculating the copy number baseline of each window in the reference sample, calculating the read segment ratio of each window in the test sample, and using the copy number baseline of each window to correct the read segment ratio of the corresponding window to obtain a correction value; using the correction value to calculate the gene copy number corresponding to each gene in the test sample; the reference sample includes a white blood cell sample;
[0013] b. For each gene in the test sample, the actual mutation probability of the heterozygous germline variation carried on the gene is statistically obtained; using the gene copy number, the tumor purity and the minor allele copy number are calculated by the grid search method;
[0014] c. Calculate the theoretical probability of germline variation and somatic variation of the test sample based on tumor purity, gene copy number and minor allele copy number; calculate the posterior probability of each filtered variation site being a germline variation or a somatic variation based on the theoretical probability of germline variation, the theoretical probability of somatic variation, the depth of the filtered variation site and the number of reads carrying the filtered variation site based on the binomial distribution formula, and record the posterior probability of the filtered variation site being a germline variation or a somatic variation as the germline variation parameter and the somatic variation parameter, respectively.
[0015] Furthermore, the data off the machine of the test sample includes large set data or small set data. The large set data refers to sequencing data of an area larger than 1Mb captured by hybridization of a high-throughput sequencing platform, and the small set data refers to sequencing data of an area between 100Kb and 1Mb captured by hybridization of a high-throughput sequencing platform.
[0016] Furthermore, the filtering of the initial variant sites in (1.3) includes: data self-filtering and / or database-based filtering; preferably, data self-filtering includes: filtering and removing initial variant sites with a variation frequency of <0.2%; filtering and removing initial variant sites with chain preference; filtering and removing initial variant sites with a supporting read number of <3; preferably, the chain preference includes that the forward read or reverse read supporting the initial variant site is 0; preferably, database-based filtering includes: database filtering of the initial variant sites using the population databases 1000Genome, ExAC and the SNP database dbSNP; preferably, filtering and removing initial variant sites with a population frequency >0.01 in the population databases 1000Genome and ExAC, and filtering and removing initial variant sites annotated as SAO=1 in the SNP database dbSNP.
[0017] Furthermore, the calculation of the gene copy number in a includes:
[0018] a1. Use more than 20 white blood cell samples to build a reference set, and count the ratio of the number of reads in each window of the reference set to the total number of reads in the sample sequencing as the copy number baseline for each window, represented by RDreference;
[0019] a2. Divide the hybrid capture area of the test sample into multiple windows, and count the ratio of the number of reads of the test sample in each window to the total number of reads of the offline data, represented by RDtest;
[0020] a3. Calculate the corrected value RDcor of each window after copy number baseline correction, RDcor=(RDtest-RDreference) / RDreference; a4. Remove the maximum and minimum RDcor values of all windows in a gene, and calculate the average value mean(RDcor) of the RDcor of the remaining windows, gene copy number=2×mean(RDcor);
[0021] Preferably, b includes: for each gene, counting the actual mutation frequency of heterozygous germline variation carried on the gene, expressed as AFactual; calculating the theoretical mutation frequency of the gene carrying germline variation, expressed as AFtheoretical, AFtheoretical = p×Vi / (p×cn+2×(1-p)), p is the tumor purity, cn is the gene copy number, and Vi is the minor allele copy number of the gene; the grid search method includes: substituting the tumor purity from 0.01 to 0.99, with a step length of 0.01, and the value of Vi is an integer from 1 to cn-1, into the square sum calculation formula ∑(AFtheoretical-AFactual) in sequence 2 Calculation is performed in the formula, when the result of the sum of squares calculation formula is the smallest, the tumor purity taken in the formula of the sum of squares calculation is the real tumor purity, and the minor allele copy number taken in the formula of the sum of squares calculation is the real minor allele copy number; in the formula of the sum of squares calculation, the square of each gene is calculated and the sum of squares is calculated, the tumor purity of different genes is the same, and the minor allele copy numbers of different genes are the same or different; preferably, the statistical method of AF actual includes: counting the polymorphic sites with a population frequency of 0.01-0.99 included in the population database 1000Genome, ExAC and the SNP database dbSNP, if the mutation site is a heterozygous mutation site, then the mutation frequency of the polymorphic site of the test sample is recorded as AF actual; the preferred population frequency is 0.05-0.95; preferably, the method for judging the heterozygous mutation site includes: if the mutation frequency of the polymorphic site is in the range of 20%-80%, the polymorphic site is considered to be a heterozygous mutation site; if the mutation frequency of the polymorphic site is <20% or >80%, the polymorphic site is considered not to be a heterozygous mutation site;
[0022] Preferably, the theoretical probability of germline variation in c is expressed as AFgermline, AFgermline=(p×Vi+1-p) / (p×cn+2×(1-p)); the theoretical probability of somatic variation is expressed as AFsomatic, AFsomatic=p×Vi / (p×cn+2×(1-p)); wherein p is the tumor purity, Vi is the number of copies of the minor allele of the gene, and cn is the number of copies of the gene; preferably, the posterior probability that the variant site is a germline variation is expressed as P(y|G;AFgermline), P(y|G;AFgermline)=Bin(f, n,AFgermline), which is the germline variation parameter; the posterior probability that the variant site is a somatic variation is expressed as P(y|G;AFsomatic), P(y|G;AFsomatic)=Bin(f, n,AFsomatic), which is the somatic variation parameter; wherein n is the depth of filtering the variant site, and f is the number of reads carrying the filtered variant site.
[0023] Furthermore, the background noise parameter is the posterior probability that the mutation site is background noise, calculated based on the mutation site in the reference sample; preferably, the mutation site feature vector in (2) is {background noise parameter, germline variation parameter, somatic variation parameter};
[0024] Preferably, the method for obtaining the background noise parameter in (2) includes: performing mutation detection on a white blood cell sample to obtain the initial mutation sites of the white blood cells, removing the germline mutation sites in the initial mutation sites, and constructing a background noise library using the mutation frequencies of the remaining mutation sites; counting the background noise frequency of the filtered mutation sites in the background noise library, and if the background noise frequency does not exist, the background noise frequency is recorded as 0; using the depth of the filtered mutation site, the number of reads carrying the filtered mutation site, and the background noise frequency, the posterior probability that the filtered mutation site is the background noise is calculated through the binomial distribution formula, and the posterior probability is recorded as the background noise parameter of the filtered mutation site;
[0025] Preferably, constructing a background noise library comprises: using more than 100 white blood cell samples, performing hybridization capture sequencing on the white blood cell samples to obtain second-generation sequencing data, and comparing with the reference genome to obtain a white blood cell BAM file; analyzing the white blood cell BAM file to obtain the initial mutation site of the white blood cells, and counting the mutation frequency of each white blood cell sample at all sites in the hybridization capture area, if the white blood cell sample does not have a mutation at the first site, the background noise frequency of the first site is recorded as 0; if the mutation frequency of the white blood cell sample at the second site is greater than 20%, it is considered that the second site carries a germline variation, and the background noise frequency of the second site is recorded as 0; for the remaining sites that do not belong to the first site or the second site among all the sites, the true mutation frequency of the white blood cell sample at the remaining sites is recorded as the background noise frequency, and for each single site in the remaining sites, the true mutation frequency of the single site is the median of the mutation frequencies of all white blood cell samples at the single site; the background noise frequencies of all sites in the hybridization capture area constitute the background noise library;
[0026] Preferably, the calculation of the background noise parameter includes: determining the background noise frequency recorded at the filtered variant site in the background noise library, using the depth of the filtered variant site and the number of reads carrying the filtered variant site to calculate the posterior probability when the mutation site is a false positive, the posterior probability is expressed as P(y|G; AFnoise), P(y|G; AFnoise)=Bin(f, n, AFnoise), wherein AFnoise is the background noise frequency recorded by white blood cells at the corresponding site in the background noise library, n is the depth of the filtered variant site, f is the number of reads carrying the filtered variant site, and P(y|G; AFnoise) is the background noise parameter.
[0027] Furthermore, the construction of the somatic variation discrimination model includes: obtaining sites known to be background noise, sites known to be germline variations, and sites known to be somatic variations, constructing feature vectors of known mutation sites, and using the feature vectors of known mutation sites and corresponding mutation site types to train a decision tree model to obtain a somatic variation discrimination model; preferably, the construction of the somatic variation discrimination model in (3) includes: using more than 100 sites known to be background noise, more than 100 sites known to be germline variations, and more than 100 sites known to be somatic variations to calculate the feature vectors of known mutation sites.
[0028] In order to achieve the above-mentioned purpose, according to a second aspect of the present invention, there is provided an electronic device for detecting somatic variation from a single sample, the electronic device comprising a variation parameter calculation unit, a feature vector construction unit and a variation discrimination output unit; the variation parameter calculation unit is used to obtain filtered variation sites using the offline data of the detection sample, and calculate germline variation parameters and somatic variation parameters; the feature vector construction unit is used to construct a mutation site feature vector according to background noise parameters, germline variation parameters and somatic variation parameters; the variation discrimination output unit is used to input the mutation site feature vector into a somatic variation discrimination model to determine whether each filtered variation site is a somatic variation or a germline variation.
[0029] Furthermore, the variation parameter calculation unit includes a data comparison unit, a variation site analysis unit, a variation site filtering unit, and a feature calculation unit; wherein the data comparison unit is used to obtain the offline data of the test sample, compare the offline data with the reference genome, and obtain the BAM file data; the variation site analysis unit is used to analyze the BAM file data to obtain the initial variation site, the depth of the initial variation site, and the number of reads carrying the initial variation site; the variation site filtering unit is used to filter the initial variation site to obtain the filtered variation site, the depth of the filtered variation site, and the number of reads carrying the filtered variation site; the feature calculation unit is used to calculate the gene copy number, tumor purity, minor allele copy number, germline variation theoretical probability, somatic variation theoretical probability, germline variation parameters, and somatic variation parameters.
[0030] Furthermore, the feature calculation unit includes a gene copy number calculation unit, a tumor purity and minor allele copy number calculation unit and a variation parameter calculation unit; the gene copy number calculation unit is used to calculate the gene copy number of each gene in the test sample according to the reference sample and the test sample; the tumor purity and minor allele copy number calculation unit is used to statistically obtain the actual mutation probability of the heterozygous germline variation carried on each gene, and calculate the tumor purity and minor allele copy number by using the grid search method and the gene copy number; the variation parameter calculation unit is used to calculate the theoretical probability of germline variation and the theoretical probability of somatic variation of the test sample according to the tumor purity, the gene copy number and the minor allele copy number, and filter the variation sites according to the theoretical probability of germline variation, the theoretical probability of somatic variation and the depth of the filtered variation sites. The degree and the number of reads carrying the filtered variant sites are used to calculate the germline variation parameters and the somatic variation parameters; preferably, the gene copy number calculation unit includes a window division unit, a copy number baseline calculation unit and a gene copy number correction unit; the window division unit is used to divide the genes of the reference sample and the test sample into segments, and equally divide each gene of the reference sample and the test sample to obtain multiple windows; the copy number baseline calculation unit is used to calculate the copy number baseline of each window in the reference sample; the gene copy number correction unit is used to calculate the read number ratio of each window in the test sample, use the copy number baseline of each window to correct the read number ratio of the corresponding window, obtain the correction value, and use the correction value to calculate the gene copy number corresponding to each gene in the test sample;
[0031] Preferably, the tumor purity and minor allele copy number calculation unit includes an actual mutation probability statistics unit and a grid search unit; the actual mutation probability statistics unit is used to statistically detect the actual mutation probability of heterozygous germline variation of each gene in the sample; the grid search unit is used to calculate ∑(AF theory-AF actual) by changing the tumor purity and minor allele copy number using a grid search method 2 The minimum value of , to obtain the true tumor purity and minor allele copy number;
[0032] Preferably, the variation parameter calculation unit includes a theoretical probability calculation unit and a binomial distribution calculation unit; the theoretical probability calculation unit is used to calculate the theoretical probability of germline variation and the theoretical probability of somatic variation of the test sample according to the tumor purity, gene copy number and minor allele copy number; the binomial distribution calculation unit is used to obtain the theoretical probability of germline variation, the theoretical probability of somatic variation, the depth of filtered variation sites and the number of reads carrying the filtered variation sites, and use the binomial distribution formula to calculate the germline variation parameters and the somatic variation parameters.
[0033] Furthermore, the variation site filtering unit includes a data self-filtering unit and / or a database filtering unit; the data self-filtering unit is used to filter and remove initial variation sites with a variation frequency of <0.2%, filter and remove initial variation sites with chain preference, and filter and remove initial variation sites with a supporting read number of <3; preferably, the chain preference includes a forward read or a reverse read supporting the initial variation site of 0; the database filtering unit is used to perform database filtering on the initial variation sites using the population databases 1000Genome, ExAC and the SNP database dbSNP; preferably, the initial variation sites with a population frequency >0.01 in the population databases 1000Genome and ExAC are filtered and removed, and the initial variation sites annotated as SAO=1 in the SNP database dbSNP are filtered and removed.
[0034] Furthermore, the feature vector construction unit includes a background noise parameter acquisition unit and a vector output unit; the background noise parameter acquisition unit is used to count the background noise frequency of the filtered variation site in the background noise library, and if the background noise frequency does not exist, the background noise frequency is recorded as 0; the background noise parameters of the filtered variation site are obtained by using the depth of the filtered variation site, the number of reads carrying the filtered variation site and the background noise frequency, through the binomial distribution formula; the vector output unit is used to obtain the background noise parameters, germline variation parameters and somatic variation parameters to construct the mutation site feature vector.
[0035] In order to achieve the above object, according to a third aspect of the present invention, a computer-readable storage medium is provided, the storage medium comprising a stored program, wherein when the program is executed, the device where the storage medium is located is controlled to execute the above method.
[0036] In order to achieve the above objective, according to a fourth aspect of the present invention, a processor is provided, the processor being used to run a program, wherein the program runs the above method.
[0037] By applying the technical solution of the present invention, the above method uses the purity of tumor tissue and the copy number of the gene to correct the variation frequency, which can make up for the defects of insufficient germline sites included in the database and the easy misjudgment based on the original variation frequency of the variation site, thereby improving the accuracy of detecting somatic cell variation in a single sample. This method also does not require the calculation of the ploidy of the tumor tissue, and has relatively low requirements for the size of the panel. And there is no need to sequence the paired tissue, which can reduce the requirements for samples and save costs. BRIEF DESCRIPTION OF THE DRAWINGS
[0038] The drawings constituting a part of the present application are used to provide a further understanding of the present invention. The exemplary embodiments of the present invention and their descriptions are used to explain the present invention and do not constitute an improper limitation of the present invention. In the drawings:
[0039] Figure 1 A flow chart of a method for detecting somatic variation from a single sample according to an embodiment of the present invention is shown.
[0040] Figure 2 is a schematic diagram of an optional electronic device for detecting somatic variation from a single sample according to an embodiment of the present invention.
[0041] Figure 3 is a schematic diagram of an optional feature calculation unit 104 according to an embodiment of the present invention.
[0042] Figure 4 It is a hardware structure block diagram of an optional method for detecting somatic variation from a single sample according to an embodiment of the present invention. DETAILED DESCRIPTION
[0043] It should be noted that, in the absence of conflict, the embodiments and features in the embodiments of the present application can be combined with each other. The present invention will be described in detail below in conjunction with the embodiments.
[0044] Terminology explanation:
[0045] Somatic mutation: Somatic mutation, also known as acquired mutation, refers to mutation that occurs in somatic cells other than sex cells. It will not cause genetic changes in offspring, but can cause changes in the genetic structure of certain contemporary cells. It is an acquired mutation during the growth and development process or under the influence of environmental factors, and usually only some cells in the body carry it.
[0046] Germline mutation: Germline mutation, also known as germ cell mutation, is a mutation that originates from reproductive cells such as sperm or eggs, so usually all cells in the body carry this mutation.
[0047] Tumor purity: Cancer cells in a tumor sample will be mixed with a certain unknown proportion of normal cells. We call the proportion of cancer cells in a tumor sample tumor purity.
[0048] As mentioned in the background technology, the filtering of germline variations in the prior art is mainly based on the variation frequency of the database or the variation site. Filtering germline variations based on paired tissue samples will bring additional sequencing costs, and in many cases, the data of paired tissues may not be available at all. Moreover, the current major public databases such as 1000Genome, ExAC, and dbSNP have not collected data for large cohorts of the Chinese population. The records of germline variations in the Chinese population are not complete, and it is difficult to filter out all germline variations based on the database. The method of filtering based on the variation frequency of the variation site is often distorted. Based on the SGZ algorithm published by Foundation Medicine, it is necessary to infer the purity, ploidy, and local copy number of the tumor tissue, and the panel captured by hybridization needs to be greater than a certain threshold to ensure the accuracy of the results. It is not suitable for small panel data.
[0049] Therefore, in this application, the inventors tried to explore a new detection method, based on the existing technology, corrected the mutation frequency for tumor purity and copy number variation, and developed a new method for somatic cell mutation detection. Therefore, a series of protection schemes of this application were proposed.
[0050] In a first typical embodiment of the present application, a method for detecting somatic mutations from a single sample tumor tissue is provided, the method comprising:
[0051] (1) Analyze the offline data of the test samples, obtain the filtered variant sites, calculate the posterior probability that the filtered variant sites are germline variants, recorded as the germline variant parameter, and calculate the posterior probability that the filtered variant sites are somatic variants, recorded as the somatic variant parameter;
[0052] (2) Using the background noise parameters, germline variation parameters, and somatic variation parameters of each filtered variation site, a mutation site feature vector is constructed;
[0053] (3) Using the somatic variation discrimination model, each filtered variation site is determined to be a somatic variation or a germline variation based on the mutation feature vector of each filtered variation site.
[0054] Figure 1 is a flow chart of a method for detecting somatic variation from a single sample according to an embodiment of the present invention. Figure 1 As shown, the method includes the following steps.
[0055] First, the offline data of the test samples are analyzed to obtain the filtered variant sites, and the posterior probability that the filtered variant sites are germline variants is calculated, which is recorded as the germline variant parameter. The posterior probability that the filtered variant sites are somatic variants is calculated, which is recorded as the somatic variant parameter (S10).
[0056] Secondly, the background noise parameters, germline variation parameters and somatic variation parameters of each filtered variation site are used to construct a mutation site feature vector (S20). The mutation site feature vector {background noise parameter, germline variation parameter, somatic variation parameter}, wherein the background noise parameter is based on each site in the reference sample (including but not limited to the white blood cell sample), and the posterior probability that the corresponding site is background noise is calculated, and the posterior probability is the background noise parameter corresponding to the site. The background noise parameter is a parameter obtained by calculating the information of the background noise library and the filtered variation site. This background noise library and the background noise frequency of each site in the background noise library are fixed parameters, and there is no need to recalculate each time a somatic variation is detected.
[0057] Finally, the somatic variation discrimination model is used to determine whether each filtered variation site is a somatic variation or a germline variation based on the mutation feature vector of each filtered variation site (S30). The somatic variation discrimination model is a model obtained by training using known variation sites, and does not need to be rebuilt every time a somatic variation is detected. A feature vector of a known mutation site is constructed using sites known to be background noise, sites known to be germline variations, and sites known to be somatic variations. The feature vector and the determined mutation site type are used to train a decision tree model to obtain a somatic variation discrimination model that can determine whether a mutation site is a somatic variation.
[0058] In a preferred embodiment, (1) includes:
[0059] (1.1) Obtain the offline data of the test sample, compare the offline data with the reference genome, and obtain the BAM file data. The test sample is a single sample tumor tissue;
[0060] (1.2) Analyze the BAM file data to obtain the initial mutation site, the depth of the initial mutation site, and the number of reads carrying the initial mutation site;
[0061] (1.3) Filtering the initial variant sites to obtain the filtered variant sites, the depth of the filtered variant sites, and the number of reads carrying the filtered variant sites;
[0062] (1.4) Calculate germline variation parameters and somatic variation parameters.
[0063] In the above (1.1), the sequencing data of the test sample is obtained, and the data is compared with the reference genome to obtain the BAM file data. In the BAM data, by comparing with the reference genome, the data can be located to the sequencing position that it truly reflects, so that the subsequent analysis of the variant site can be accurately located to the specific base. After the comparison, the sequencing fragments with overlapping sequencing sites can be analyzed, thereby improving the accuracy of the sequencing data through the information of multiple sequencing fragments. The reference genome used in this application is the human reference genome, including but not limited to the genome HG19 or the genome HG38.
[0064] In the above (1.2), based on the BAM data obtained above, the BAM data is analyzed to obtain the initial mutation site, the depth of the initial mutation site and the number of reads carrying the initial mutation site.
[0065] In the above (1.3), the initial variant sites are filtered to obtain the filtered variant sites, and accordingly the depth (n) of the filtered variant sites and the number of reads carrying the filtered variant sites (f) can be obtained. The filtering method includes but is not limited to filtering using a database or filtering the data itself. Using a database for filtering can filter out the germline variant sites with numerical values in the prior art, reduce the amount of subsequent calculations required, and improve the calculation speed and accuracy of the method.
[0066] In the above (1.4), the germline variation parameters and somatic variation parameters are calculated using methods in the prior art.
[0067] In a preferred embodiment, calculating the germline variation parameter and the somatic variation parameter comprises:
[0068] a. dividing the genes of the reference sample and the test sample into segments, and equally dividing each gene of the reference sample and the test sample to obtain multiple windows; calculating the copy number baseline of each window in the reference sample, calculating the read segment ratio of each window in the test sample, and using the copy number baseline of each window to correct the read segment ratio of the corresponding window to obtain a correction value; using the correction value to calculate the gene copy number corresponding to each gene in the test sample; the reference sample includes a white blood cell sample;
[0069] b. For each gene in the test sample, the actual mutation probability of the heterozygous germline variation carried on the gene is statistically obtained; using the gene copy number, the tumor purity and the minor allele copy number are calculated by the grid search method;
[0070] c. Calculate the theoretical probability of germline variation and somatic variation of the test sample based on tumor purity, gene copy number and minor allele copy number; calculate the posterior probability of each filtered variation site being a germline variation or a somatic variation based on the theoretical probability of germline variation, the theoretical probability of somatic variation, the depth of the filtered variation site and the number of reads carrying the filtered variation site based on the binomial distribution formula, and record the posterior probability of the filtered variation site being a germline variation or a somatic variation as the germline variation parameter and the somatic variation parameter, respectively.
[0071] The above a includes: selecting a white blood cell sample as a reference sample, performing hybridization capture sequencing, and obtaining reference sample data. The reference sample data is divided into different window (bin) areas, and the copy number baseline RDreference of each window is calculated; the test sample is also divided into bin areas, and the read segment ratio RDtest of each bin in the test sample is calculated. RDreference is used to correct RDtest to obtain the correction value RDcor, RDcor=(RDtest-RDreference) / RDreference, and the correction value RDcor is used to calculate the gene copy number cn of each gene in the test sample. The calculation of the copy number includes: removing the maximum and minimum values of the RDcor values of all bins in a gene after baseline correction, and calculating the average value mean(RDcor) of the remaining values. The copy number cn of the gene=2*mean(RDcor);
[0072] The above b includes: fitting the tumor purity, including: constructing a probability model based on the actual mutation probability of the gene copy number and the heterozygous germline variation carried on the gene, and obtaining the tumor purity and the minor allele copy number through grid search; Germline variation refers to the variation already carried during the human embryonic development period, which is the variation carried by the cells in the body. Heterozygous germline variation refers to the same site on two homologous chromosomes, only one of which carries the variation, and the other does not carry the variation.
[0073] The above c includes: according to the tumor purity, gene copy number and minor allele copy number obtained above, the theoretical probability of germline variation AFgermline and the theoretical probability of somatic variation AFsomatic of the test sample are calculated, AFgermline=(p×Vi+1-p) / (p×cn+2×(1-p)), AFsomatic=p×Vi / (p×cn+2×(1-p)), where p is the tumor purity, cn is the copy number of the gene, and Vi is the minor allele copy number of the gene; then the obtained theoretical probability of germline variation AFgermline, theoretical probability of somatic variation AFsomatic, the depth n of the filtered variation site and the number of reads f carrying the filtered variation site are used to calculate the germline variation parameter P(y|G; AFgermline) and the somatic variation parameter P(y|G; AFsomatic) based on the binomial distribution formula, P(y|G; AFgermline)=Bin(f, n, AFgermline), P(y|G; AFsomatic)=Bin(f, n, AFsomatic).
[0074] The calculation formula for the binomial distribution is as follows:
[0075] .
[0076] In the calculation of germline variation parameters or somatic variation parameters, y represents frequency, G represents the current situation is a germline variation or somatic variation scenario, and Bin() represents the calculation using binomial distribution. Combined with the above binomial distribution calculation formula, when the variation site is a germline variation, the probability of a depth of n and a number of reads carrying the variation is f is C n f ×AFgermline f -(1-AFgermline) n-f .
[0077] In a preferred embodiment, the offline data of the test sample includes large set data or small set data, the large set data is sequencing data of a region larger than 1Mb captured by hybridization of a high-throughput sequencing platform, and the small set data is sequencing data of a region between 100Kb and 1Mb captured by hybridization of a high-throughput sequencing platform.
[0078] The existing algorithms for detecting somatic mutations all rely on large amounts of samples to determine whether a mutation is a somatic mutation. This makes it difficult to detect somatic mutations in small panel data, especially in sequencing data with hybridization capture regions between 100kb and 1Mb.
[0079] In a preferred embodiment, the filtering of the initial variant sites in (1.3) includes: data filtering itself and / or database-based filtering;
[0080] Preferably, data self-filtering includes: filtering and removing initial mutation sites with mutation frequency <0.2%; filtering and removing initial mutation sites with chain preference; filtering and removing initial mutation sites with supporting read number <3; preferably, chain preference includes that the forward read or reverse read supporting the initial mutation site is 0; preferably, database-based filtering includes: database filtering of initial mutation sites using population databases 1000Genome, ExAC and SNP database dbSNP; preferably, filtering and removing initial mutation sites with population frequency >0.01 in population databases 1000Genome and ExAC, and filtering and removing initial mutation sites annotated as SAO=1 in SNP database dbSNP.
[0081] The annotation "SAO=1" above indicates a known germline variation. Using the above database filtering method, using the data of known germline variations in the prior art, the germline variation sites in the single sample offline data are eliminated to prevent such germline variation from being erroneously calculated as somatic variation in the subsequent process, which can improve the accuracy of subsequent judgments. Eliminating here can also reduce the amount of subsequent calculations on these known sites and improve the efficiency of the method.
[0082] In a preferred embodiment, the calculation of the gene copy number in a includes:
[0083] a1. Use more than 20 white blood cell samples to build a reference set, and count the ratio of the number of reads in each window of the reference set to the total number of reads in the sample sequencing as the copy number baseline for each window, represented by RDreference;
[0084] a2. Divide the hybrid capture area of the test sample into multiple windows, and count the ratio of the number of reads of the test sample in each window to the total number of reads of the offline data, represented by RDtest;
[0085] a3. Calculate the correction value RDcor of each window after copy number baseline correction, RDcor=(RDtest-RDreference) / RDreference;
[0086] a4. Remove the maximum and minimum RDcor values of all windows in a gene, and calculate the average value of RDcor of the remaining windows, mean(RDcor). The gene copy number = 2×mean(RDcor);
[0087] Using the above method, the copy number baseline obtained by calculating the reference set is used to correct the copy number of the test sample, reducing the impact of sequencing bias caused by the difference in the sequence itself on the analysis results, thereby obtaining a more reliable gene copy number. Since the white blood cell samples and test samples used as the reference set are both human samples, the gene information in the two samples is slightly different, so the division method of multiple bins in each gene of the two samples can be consistent, thereby achieving the correction of the test sample through the reference set.
[0088] Preferably, b includes: for each gene, counting the actual mutation frequency of heterozygous germline variation carried on the gene, expressed as AFactual; calculating the theoretical mutation frequency of the gene carrying germline variation, expressed as AFtheoretical, AFtheoretical = p×Vi / (p×cn+2×(1-p)), p is the tumor purity, cn is the gene copy number, and Vi is the minor allele copy number of the gene;
[0089] The grid search method includes: the tumor purity is from 0.01 to 0.99, with a step length of 0.01, and the value of Vi is an integer from 1 to cn-1, which is substituted into the square sum calculation formula ∑ (AF theory-AF actual) 2 Calculation is performed in , when the result of the square sum calculation formula is the smallest, the tumor purity taken in the square sum calculation formula is the real tumor purity, and the minor allele copy number taken in the square sum calculation formula is the real minor allele copy number; in the square sum calculation formula, the square of each gene is calculated and the square sum is calculated, the tumor purity of different genes is the same, and the minor allele copy number of different genes is the same or different;
[0090] Preferably, the statistical method of AF actual includes: counting the polymorphic sites with a population frequency of 0.01-0.99 included in the population database 1000Genome, ExAC and the SNP database dbSNP, if the mutation site is a heterozygous mutation site, then the mutation frequency of the polymorphic site of the test sample is recorded as AF actual; preferably, the population frequency is 0.05-0.95;
[0091] Preferably, the method for determining a heterozygous mutation site includes: if the mutation frequency of the polymorphic site is within the range of 20%-80%, the polymorphic site is considered to be a heterozygous mutation site; if the mutation frequency of the polymorphic site is <20% or >80%, the polymorphic site is considered not to be a heterozygous mutation site.
[0092] The statistical method of the above AF actual includes: through the polymorphic sites covered by the detected gene with a population frequency between 1%-99% (preferably 5%-95%) in the database (such as 1000Genome, ExAC and SNP), it is considered to be the known germline variation sites carried by the gene. According to the offline data, the actual mutation frequency of these sites is analyzed and calculated. The above AF actual is the probability obtained for each germline variation site.
[0093] The actual mutation frequency of germline variation reflects the actual law. The theoretical mutation probability of the gene carrying germline variation is calculated using the tumor purity, gene copy number and minor allele copy number obtained above. ∑(AF theory - AF actual) is calculated by grid search method. 2 The tumor purity and minor allele copy number corresponding to the minimum result of . The tumor purity between different genes is the same, and the minor allele copy number can be the same or different. The smaller the sum of squares, the closer it is to the overall objective law. Therefore, the tumor purity and minor allele copy number corresponding to the minimum sum of squares can be considered the most accurate value.
[0094] Preferably, the theoretical probability of germline variation in c is expressed as AFgermline, AFgermline=(p×Vi+1-p) / (p×cn+2×(1-p)); the theoretical probability of somatic variation is expressed as AFsomatic, AFsomatic=p×Vi / (p×cn+2×(1-p)); wherein p is the tumor purity, Vi is the number of copies of the minor allele of the gene, and cn is the number of copies of the gene where the mutation site is located;
[0095] Preferably, the posterior probability that the variant site is a germline variant is expressed as P(y|G; AFgermline), P(y|G; AFgermline)=Bin(f, n, AFgermline), which is the germline variant parameter; the posterior probability that the variant site is a somatic variant is expressed as P(y|G; AFsomatic), P(y|G; AFsomatic)=Bin(f, n, AFsomatic), which is the somatic variant parameter; wherein n is the depth of the filtered variant site, and f is the number of reads carrying the filtered variant site.
[0096] The tumor purity, minor allele copy number, and gene copy number obtained above are used to calculate the theoretical probability of germline variation and the theoretical probability of somatic variation of the gene. The binomial distribution and the sequencing depth of the variant site obtained previously and the number of reads carrying the filtered variant site are further used to calculate the somatic variation parameters and germline variation parameters of the mutation site.
[0097] In a preferred embodiment, the background noise parameter is the posterior probability that the mutation site is background noise, calculated based on the mutation site in the reference sample;
[0098] Preferably, the mutation site feature vector in (2) is {background noise parameter, germline variation parameter, somatic variation parameter};
[0099] Preferably, the method for obtaining the background noise parameter in (2) includes: performing mutation detection on a white blood cell sample to obtain the initial mutation sites of the white blood cells, removing the germline mutation sites in the initial mutation sites, and constructing a background noise library using the mutation frequencies of the remaining mutation sites; counting the background noise frequency of the filtered mutation sites in the background noise library, and if the background noise frequency does not exist, the background noise frequency is recorded as 0; using the depth of the filtered mutation site, the number of reads carrying the filtered mutation site, and the background noise frequency, the posterior probability that the filtered mutation site is the background noise is calculated through the binomial distribution formula, and the posterior probability is recorded as the background noise parameter of the filtered mutation site;
[0100] Preferably, constructing a background noise library comprises: using more than 100 white blood cell samples, performing hybridization capture sequencing on the white blood cell samples to obtain second-generation sequencing data, and comparing with the reference genome to obtain a white blood cell BAM file; analyzing the white blood cell BAM file to obtain the initial mutation site of the white blood cells, and counting the mutation frequency of each white blood cell sample at all sites in the hybridization capture area, if the white blood cell sample does not have a mutation at the first site, the background noise frequency of the first site is recorded as 0; if the mutation frequency of the white blood cell sample at the second site is greater than 20%, it is considered that the second site carries a germline variation, and the background noise frequency of the second site is recorded as 0; for the remaining sites that do not belong to the first site or the second site among all the sites, the true mutation frequency of the white blood cell sample at the remaining sites is recorded as the background noise frequency, and for each single site in the remaining sites, the true mutation frequency of the single site is the median of the mutation frequencies of all white blood cell samples at the single site; the background noise frequencies of all sites in the hybridization capture area constitute the background noise library;
[0101] The above-mentioned first site and second site do not refer to a specific site. In the present application, sites without mutation in the white blood cell sample belong to the first site, and sites with mutation channels greater than 20% in the white blood cell sample belong to the second site.
[0102] Preferably, the calculation of the background noise parameter includes: determining the background noise frequency recorded at the filtered variant site in the background noise library, using the depth of the filtered variant site and the number of reads carrying the filtered variant site to calculate the posterior probability when the mutation site is a false positive, the posterior probability is expressed as P(y|G; AFnoise), P(y|G; AFnoise)=Bin(f, n, AFnoise), wherein AFnoise is the background noise frequency recorded by white blood cells at the corresponding site in the background noise library, n is the depth of the filtered variant site, f is the number of reads carrying the filtered variant site, and P(y|G; AFnoise) is the background noise parameter.
[0103] Since there are no somatic mutations in leukocytes, after removing germline mutations in leukocytes, the remaining mutations belong to background noise. Therefore, the background noise library is constructed using leukocyte samples to reflect the probability of sequencing errors at each site in high-throughput sequencing. The construction of the background noise library includes: hybridization capture sequencing of more than 100 leukocyte samples, and after obtaining the second-generation sequencing data, the sequencing data of each leukocyte sample is compared with the reference genome, and the leukocyte BAM files are obtained respectively. The leukocyte BAM files are analyzed to obtain the leukocyte initial mutation sites, and the background noise frequencies of all sites covered by the panel are counted. If there is no mutation at the leukocyte initial mutation site in the leukocyte sample, the background noise frequency of the site is recorded as 0; if the mutation frequency of the leukocyte sample at the leukocyte initial mutation site is greater than 20%, it is considered that the leukocyte initial mutation site carries germline mutations, and the background noise frequency of the site is recorded as 0; the true mutation frequency is used as the background noise frequency for the remaining sites, and the background noise frequency of the background noise library at a certain site is recorded as the median of the mutation frequency of all samples at the site, and the background noise frequencies of all sites constitute the background noise library. The background noise frequencies of all sites in leukocytes are summarized to obtain a background noise library containing sites and their corresponding background noise frequency information, which can reflect the probability of background noise such as high-throughput test errors in each site, and can more objectively and truly calibrate the test samples according to the differences of each site to prevent the influence of background noise on the detection of somatic mutations. The leukocyte sample panel for constructing the background noise library is the same as the panel of the sample to be tested.
[0104] In a preferred embodiment, the construction of the somatic variation discrimination model includes: obtaining sites known to be background noise, sites known to be germline variations, and sites known to be somatic variations, constructing feature vectors of known mutation sites, and training a decision tree model using the feature vectors of known mutation sites and corresponding mutation site types to obtain a somatic variation discrimination model;
[0105] Preferably, the construction of the somatic variation discrimination model in (3) includes: using more than 100 sites known to be background noise, more than 100 sites known to be germline variations, and more than 100 sites known to be somatic variations to calculate the feature vector of known mutation sites.
[0106] In a second typical embodiment of the present application, an electronic device for detecting somatic variation from a single sample is provided, the electronic device comprising a variation parameter calculation unit 10, a feature vector construction unit 20 and a variation discrimination output unit 30; the variation parameter calculation unit 10 is used to obtain filtered variation sites using the offline data of the detection sample, and calculate germline variation parameters and somatic variation parameters; the feature vector construction unit 20 is used to construct a mutation site feature vector according to background noise parameters, germline variation parameters and somatic variation parameters; the variation discrimination output unit 30 is used to input the mutation site feature vector into a somatic variation discrimination model to determine whether each filtered variation site is a somatic variation or a germline variation.
[0107] Figure 2 is a schematic diagram of an optional electronic device for detecting somatic cell variation from a single sample according to an embodiment of the present invention, such as Figure 2 As shown, the device includes a variation parameter calculation unit 10, a feature vector construction unit 20 and a variation discrimination output unit 30.
[0108] Among them, the variation parameter calculation unit 10 is used to analyze the offline data of the test sample, obtain the filtered variation site, and calculate the posterior probability that the filtered variation site is a germline variation and a somatic variation, which are recorded as the germline variation parameter and the somatic variation parameter respectively.
[0109] The feature vector construction unit 20 is used to obtain germline variation parameters and somatic variation parameters, and obtain or call built-in or externally stored background noise parameters to construct a mutation site feature vector.
[0110] The variation discrimination output unit 30 is used to obtain a mutation site feature vector, input the mutation site feature vector into a built-in or externally stored somatic variation discrimination model, and use the discrimination model to determine whether each filtered variation site is a somatic variation or a germline variation.
[0111] In a preferred embodiment, the variation parameter calculation unit 10 includes a data comparison unit 101, a variation site analysis unit 102, a variation site filtering unit 103 and a feature calculation unit 104; wherein the data comparison unit 101 is used to obtain the offline data of the test sample, compare the offline data with the reference genome, and obtain the BAM file data; the variation site analysis unit 102 is used to analyze the BAM file data to obtain the initial variation site, the depth of the initial variation site and the number of reads carrying the initial variation site; the variation site filtering unit 103 is used to filter the initial variation site to obtain the filtered variation site, the depth of the filtered variation site and the number of reads carrying the filtered variation site; the feature calculation unit 104 is used to calculate the gene copy number, tumor purity, minor allele copy number, germline variation theoretical probability, somatic variation theoretical probability, germline variation parameters and somatic variation parameters.
[0112] The data comparison unit 101 is used to obtain the offline data and sequencing depth of a single sample, compare the offline data with the reference genome, and output BAM file data. The data comparison unit 101 can be built-in or associated with the reference genome and comparison algorithm in the prior art.
[0113] The variant site analysis unit 102 is used to obtain BAM file data, analyze the BAM file data, obtain the initial variant site, and use the algorithm in the prior art to calculate the depth of the initial variant site and the number of reads carrying the initial variant site.
[0114] The variant site filtering unit 103 is used to obtain the initial variant site, filter the initial variant site, obtain the filtered variant site, and screen the depth of the filtered variant site and the number of reads carrying the filtered variant site from the depth of the initial variant site and the number of reads carrying the initial variant site obtained above. The variant site filtering unit 103 can be built-in or associated with databases such as 1000Genome, ExAC, SNP database dbSNP in the prior art, and has a built-in filtering algorithm to achieve filtering based on the database and / or based on the data itself.
[0115] The feature calculation unit 104 is used to obtain the BAM file data obtained in the above units, filter the variation sites, the depth of the filtered variation sites and the number of reads carrying the filtered variation sites, and calculate the gene copy number, tumor purity, minor allele copy number, germline variation theoretical probability, somatic variation theoretical probability, germline variation parameters and somatic variation parameters.
[0116] In a preferred embodiment, the feature calculation unit 104 includes a gene copy number calculation unit 1041, a tumor purity and minor allele copy number calculation unit 1042 and a variation parameter calculation unit 1043; the gene copy number calculation unit 1041 is used to calculate the gene copy number of each gene in the test sample based on the reference sample and the test sample; the tumor purity and minor allele copy number calculation unit 1042 is used to statistically obtain the actual mutation probability of the heterozygous germline variation carried on each gene, and calculate the tumor purity and minor allele copy number using the grid search method and the gene copy number; the variation parameter calculation unit 1043 is used to calculate the theoretical probability of germline variation and the theoretical probability of somatic variation of the test sample based on the tumor purity, the gene copy number and the minor allele copy number, and calculate the germline variation parameters and the somatic variation parameters based on the theoretical probability of germline variation, the theoretical probability of somatic variation, the depth of the filtered variation site and the number of reads carrying the filtered variation site.
[0117] Figure 3 is a schematic diagram of an optional feature calculation unit 104 according to an embodiment of the present invention, such as Figure 3 As shown, the feature calculation unit 104 includes a gene copy number calculation unit 1041, a tumor purity and minor allele copy number calculation unit 1042 and a variation parameter calculation unit 1043.
[0118] Among them, the gene copy number calculation unit 1041 is used to obtain the BAM file data of the test sample, and has the BAM file data of the reference sample built in, and each gene in the BAM file data of the test sample and the reference sample is equally divided into multiple windows (bins); the ratio of the number of reads in each bin of the reference sample to the total number of reads sequenced by the reference sample is calculated, that is, the copy number baseline in each window; the ratio of the number of reads in each bin of the test sample to the total number of reads sequenced by the test sample is calculated, that is, the ratio of the number of reads in each window; the copy number baseline in each window is used to correct the ratio of the number of reads in the corresponding window to obtain the correction value of each window; and the correction value is used to calculate the gene copy number corresponding to each gene in the test sample.
[0119] Tumor purity and minor allele copy number calculation unit 1042, used to obtain the gene copy number, and calculate the theoretical mutation frequency of the gene carrying germline variation, and has built-in actual mutation probability of heterozygous germline variation carried on the gene, and uses grid search method to calculate tumor purity and minor allele copy number of the gene;
[0120] The variation parameter calculation unit 1043 is used to obtain tumor purity, gene copy number, minor allele copy number, sequencing depth of variation sites and number of reads carrying filtered variation sites, and calculate the theoretical probability of somatic variation and the theoretical probability of germline variation, and then calculate the germline variation parameters and somatic variation parameters using binomial distribution.
[0121] Preferably, the gene copy number calculation unit 1041 includes a window division unit 10411, a copy number baseline calculation unit 10412 and a gene copy number correction unit 10413; wherein the window division unit 10411 is used to divide the genes of the reference sample and the test sample into segments, and equally divide each gene of the reference sample and the test sample to obtain multiple windows; the copy number baseline calculation unit 10412 is used to calculate the copy number baseline of each window in the reference sample; the gene copy number correction unit 10413 is used to calculate the read segment ratio of each window in the test sample, use the copy number baseline of each window to correct the read segment ratio of the corresponding window, obtain the correction value, and use the correction value to calculate the gene copy number corresponding to each gene in the test sample.
[0122] Preferably, the tumor purity and minor allele copy number calculation unit 1042 includes an actual mutation probability statistics unit 10421 and a grid search unit 10422; wherein the actual mutation probability statistics unit 10421 is used to statistically detect the actual mutation probability of the heterozygous germline variation of each gene in the sample; the grid search unit 10422 is used to calculate ∑(AF theory-AF actual) by changing the tumor purity and the minor allele copy number using a grid search method. 2 The minimum value of , to obtain the true tumor purity and minor allele copy number.
[0123] The above-mentioned actual mutation probability statistics unit 10421 can be built-in with the actual mutation probability of different genes in the known database, and can also match and obtain the actual mutation probability in the known data after obtaining the filtered mutation sites. The above-mentioned grid search unit 10422 has a built-in grid search algorithm in the prior art, which searches for and obtains the actual changes in tumor purity and minor allele copy number by changing the tumor purity and minor allele copy number. The built-in grid search algorithm includes but is not limited to changing the tumor purity from 0.01 to 0.99, with a step length of 0.01, and the value of Vi is an integer from 1 to cn-1, and is substituted into the square sum calculation formula ∑ (AF theory-AF actual) in sequence 2Calculation is performed in , when the result of the sum of squares calculation formula is the smallest, the tumor purity taken in the sum of squares calculation formula is the true tumor purity, and the number of minor allele copies taken in the sum of squares calculation formula is the true number of minor allele copies; in the sum of squares calculation formula, the square of each gene is calculated and the sum of squares is calculated, the tumor purity of different genes is the same, and the number of minor allele copies of different genes is the same or different.
[0124] Preferably, the variation parameter calculation unit 1043 includes a theoretical probability calculation unit 10431 and a binomial distribution calculation unit 10432; the theoretical probability calculation unit 10431 is used to calculate the theoretical probability of germline variation and the theoretical probability of somatic variation of the test sample according to the tumor purity, gene copy number and minor allele copy number; the binomial distribution calculation unit 10432 is used to obtain the theoretical probability of germline variation, the theoretical probability of somatic variation, the depth of filtered variation sites and the number of reads carrying the filtered variation sites, and use the binomial distribution formula to calculate the germline variation parameters and the somatic variation parameters.
[0125] In a preferred embodiment, the variation site filtering unit 103 includes a data self-filtering unit 1031 and / or a database filtering unit 1032; the data self-filtering unit 1031 is used to filter and remove initial variation sites with a variation frequency of <0.2%, filter and remove initial variation sites with chain preference, and filter and remove initial variation sites with a supporting read number of <3; preferably, the chain preference includes the forward read or reverse read supporting the initial variation site being 0; the database filtering unit 1032 is used to perform database filtering on the initial variation sites using the population databases 1000Genome, ExAC and the SNP database dbSNP; preferably, the initial variation sites with a population frequency >0.01 in the population databases 1000Genome and ExAC are filtered and removed, and the initial variation sites annotated as SAO=1 in the SNP database dbSNP are filtered and removed.
[0126] In a preferred embodiment, the feature vector construction unit 20 includes a background noise parameter acquisition unit 201 and a vector output unit 202; wherein, the background noise parameter acquisition unit 201 is used to count the background noise frequency of the filtered variation site in the background noise library, and if the background noise frequency does not exist, the background noise frequency is recorded as 0; the background noise parameter of the filtered variation site is obtained by using the depth of the filtered variation site, the number of reads carrying the filtered variation site and the background noise frequency, through the binomial distribution formula; the vector output unit 202 is used to obtain the background noise parameter, the germline variation parameter and the somatic variation parameter to construct the mutation site feature vector.
[0127] In a third typical embodiment of the present application, a computer-readable storage medium is provided, which includes a stored program, wherein when the program is executed, the device where the storage medium is located is controlled to execute the above method.
[0128] In a fourth typical implementation of the present application, a processor is provided, the processor being configured to run a program, wherein the program runs the above method.
[0129] It should be noted that, for the above-mentioned method embodiments, for the sake of simplicity, they are all described as a series of action combinations, but those skilled in the art should know that the present invention is not limited by the order of the actions described, because according to the present invention, certain steps can be performed in other orders or simultaneously. Secondly, those skilled in the art should also know that the embodiments described in the specification are all preferred embodiments, and the actions involved are not necessarily required by the present invention.
[0130] It can be known from the description of the above implementation methods that those skilled in the art can clearly understand that the present application can be implemented by means of software plus hardware devices such as detection devices. Based on such an understanding, the data processing part of the technical solution of the present application can be embodied in the form of a software product, and the computer software product can be stored in a storage medium such as ROM / RAM, a magnetic disk, an optical disk, etc., including a number of instructions for a computer device (which can be a personal computer, a server, or a network device, etc.) to execute the methods of each embodiment of the present application or some parts of the embodiments.
[0131] The present application can be used in many general or special computing system environments or configurations, such as personal computers, server computers, handheld or portable devices, tablet devices, multiprocessor systems, microprocessor-based systems, set-top boxes, programmable consumer electronic devices, network PCs, minicomputers, mainframe computers, distributed computing environments including any of the above systems or devices, etc.
[0132] The method provided in this application can be executed in a terminal, a computer terminal or a similar computing device. Taking running on a terminal as an example, Figure 4 FIG. 1 is a hardware structure block diagram of a terminal for a method of detecting somatic cell variation from a single sample according to an embodiment of the present invention. Figure 4 As shown, the terminal may include one or more ( Figure 4 Only one is shown in the figure) processor A2 (processor A2 may include but is not limited to a processing device such as a microprocessor MCU or a programmable logic device FPGA) and a memory A4 for storing data. Optionally, the terminal may also include a transmission device A6 and an input / output device A8 for communication functions. It can be understood by those skilled in the art that Figure 4The structure shown is for illustration only and does not limit the structure of the above terminal. Figure 4 More or fewer components as shown, or with Figure 4 Different configurations are shown.
[0133] The memory A4 can be used to store computer programs, for example, software programs and modules of application software, such as computer programs corresponding to methods such as read segment splicing, clustering, and consistency processing in the embodiments of the present invention. The processor A2 executes various functional applications and data processing by running the computer programs stored in the memory A4, that is, implementing the above method. The memory A4 may include a high-speed random access memory, and may also include a non-volatile memory, such as one or more magnetic storage devices, flash memory, or other non-volatile solid-state memory. In some instances, the memory A4 may further include a memory remotely arranged relative to the processor A2, and these remote memories may be connected to the terminal via a network. Examples of the above-mentioned network include, but are not limited to, the Internet, an intranet, a local area network, a mobile communication network, and combinations thereof.
[0134] The transmission device A6 is used to receive or send data via a network. The specific example of the above network may include a wireless network provided by a communication provider of the terminal. In one example, the transmission device A6 includes a network adapter (Network Interface Controller, referred to as NIC), which can be connected to other network devices through a base station so as to communicate with the Internet. In one example, the transmission device A6 can be a radio frequency (Radio Frequency, referred to as RF) module, which is used to communicate with the Internet wirelessly.
[0135] Obviously, those skilled in the art should understand that some modules or steps of the present application described above can be implemented in a general computing device, they can be concentrated on a single computing device, or distributed on a network composed of multiple computing devices, and optionally, they can be implemented with program codes executable by the computing device, so that they can be stored in a storage device and executed by the computing device, or they can be made into individual integrated circuit modules, or multiple modules or steps therein can be made into a single integrated circuit module for implementation. Thus, the present application is not limited to any specific combination of hardware and software.
[0136] In the present application, for the parameters involved in the above method, including but not limited to the depth n of the filtered variant site, the number of reads f carrying the filtered variant site, the copy number baseline RDreference, the read number ratio RDtest, the correction value RDcor, the mean value mean(RDcor), the copy number cn of the gene, the theoretical probability of germline variation AFgermline, the theoretical probability of somatic variation AFsomatic, the minor allele copy number Vi, the germline variation parameter P(y|G; AFgermline) and the somatic variation parameter P(y|G; AFsomatic), the character expressions of the parameters can be flexibly replaced, and the specific description method does not affect the substantive content of the parameters. We aim to show all the features and applications of the technology as much as possible and to describe it accurately in all aspects. Therefore, even if we use multiple languages and forms to express it, there is no need to change the amount of information or intrinsic characteristics contained in the parameters themselves. The scope of the protected rights covers all the specified values, intervals, ratios, etc. mentioned in these description methods, and is not limited to any deterministic representation method. Rather than being limited to a certain symbolic representation method.
[0137] The beneficial effects of the present application will be further explained in detail below in conjunction with specific embodiments.
[0138] Example
[0139] The above method for detecting somatic variants from a single sample:
[0140] 1. Obtain the sequence information of the target gene of the sample to be tested.
[0141] 2. Align the sequence obtained by sequencing with the human genome HG19 to obtain the original BAM file. Use Picard and Samtools to remove the duplicate sequences and sequences aligned to multiple positions to obtain the final BAM file.
[0142] 3. Analyze the final BAM file (variation detection) to obtain the initial variant sites, annotate the initial variant sites using the population database 1000Genome, ExAC and SNP database dbSNP, and then perform the following filtering:
[0143] (1) Filter out sites with a population frequency > 0.01 in the population databases 1000Genome and ExAC;
[0144] (2) Filter out variants whose variant sites are annotated as SAO=1 in the SNP database dbSNP.
[0145] 4. Calculation of background noise parameters:
[0146] (1) Perform mutation detection based on white blood cell samples to obtain initial mutation sites, remove germline mutation sites, and record the background noise frequency of the sites where the germline mutation sites are located as 0, the background noise of sites that do not carry mutations as 0, and the background noise of the remaining points as mutation frequencies. The background noise frequencies of all sites covered by the comprehensive panel are used to construct a background noise library;
[0147] (2) Count the background noise frequency of the initial mutation site of the sample to be tested in the background noise library (if it does not exist, the frequency is recorded as 0), and use the depth of the mutation site, the number of reads carrying the mutation, and the background noise frequency to calculate the posterior probability that the mutation site of the sample to be tested is the background noise through the binomial distribution formula. This probability is recorded as the background noise parameter of the initial mutation site;
[0148] 5. Calculate the copy number of genes covered by the panel based on the BAM file data:
[0149] (1) A reference set was constructed using more than 20 white blood cell samples. The ratio of the number of reads in each bin of the reference set samples to the total number of reads in the sample sequencing was counted as the copy number baseline for each bin, represented by RDreference.
[0150] (2) Divide the hybridization capture area of the sample to be tested into bins, and count the ratio of the number of reads of the test sample in each bin to the total number of reads of the sample sequencing, expressed as RDtest. Calculate the baseline-corrected value of each bin based on the baseline data of each bin: RDcor = (RDtest-RDreference) / RDreference; Calculate the copy number of the gene: remove the maximum and minimum values of the baseline-corrected RDcor values of all bins in a gene, and calculate the average value of the remaining values mean(RDcor). The copy number of the gene cn = 2*mean(RDcor);
[0151] 6. Calculation of purity of tumor tissue and copy number of minor alleles:
[0152] (1) Count the actual mutation frequency of heterozygous germline variants carried by each gene, and calculate the theoretical mutation frequency of the gene carrying germline variants AFtheoretical = p*Vi / (p*cn+2*(1-p)), where p is the tumor purity, cn is the copy number of the gene, and Vi is the copy number of the minor allele of the gene;
[0153] (2) Count the actual mutation frequency of heterozygous germline variants carried by each gene, and calculate the theoretical mutation frequency of the gene carrying germline variants: AFtheoretical = p*Vi / (p*cn+2*(1-p)), where p is the tumor purity, cn is the number of copies of the gene, and Vi is the number of copies of the minor allele of the gene;
[0154] (3) The sum of the squares of the difference between the theoretical and actual mutation frequencies of heterozygous germline mutations carried by the gene ∑(AF theory - AF actual) 2 , the possible values of Vi range from 1 to cn-1, and the tumor purity ranges from 0.01 to 0.99, with a step length of 0.01. The possible values of Vi from 1 to cn-1 are substituted in sequence. When the sum of the squares of the difference between the theoretical and actual mutation frequencies of heterozygous germline mutations is the smallest, the tumor purity and minor allele copy number taken are the tumor purity and minor allele copy number required;
[0155] 7. Calculate the probability of the variant site being assumed to be a germline variant and a somatic variant:
[0156] When the variant site is a germline variant, the theoretical frequency of the variant site is: AFgermline=(p×Vi+1-p) / (p×cn+2×(1-p));
[0157] When the variant site is a somatic variant, the theoretical frequency of the variant site is: AFsomatic = p × Vi / (p × cn + 2 × (1-p));
[0158] Wherein, AFgermline is the theoretical probability of germline variation, AFsomatic is the theoretical probability of somatic variation, p is the tumor purity, Vi is the copy number of the minor allele of the tumor, and cn is the copy number of the gene where the mutation site is located.
[0159] 8. Calculate the probability based on the binomial distribution and determine the source of the variant site:
[0160] Assuming that the depth of the variant site is n, and the number of reads carrying the mutation is f, when the variant site is a germline variant, the probability of measuring the current situation is P(y|G; AF germline ) = Bin(f, n, AFgermline); when the mutation site is a somatic mutation, the probability of measuring the current situation is P(y|G; AF somatic ) = Bin(f, n, AFsomatic), which are recorded as germline variation parameters and somatic variation parameters respectively. The background noise parameter, germline variation parameter and somatic variation parameter of each filtered variation site are used to construct the mutation site feature vector {background noise parameter, germline variation parameter, somatic variation parameter}.
[0161] 9. Determination of germline and somatic mutations
[0162] The above mutation site feature vector is input into the trained somatic cell variation discrimination model to obtain the specific variation type of the mutation site.
[0163] Specific applications include: collecting 20 pairs of tumor tissue-leukocyte samples for library construction and sequencing, using sequencing data to perform somatic variation detection on single tumor samples and on paired tumor tissue-leukocyte samples, and comparing the results of variation detection in the two processes to determine the performance of somatic variation detection based on single samples.
[0164] Specific applications:
[0165] (1) Using a single sample from Youxun Lung, the panel size was approximately 92 kb, and hybridization capture and on-machine sequencing were performed on 20 pairs of tumor tissue-leukocyte paired samples to obtain off-machine data;
[0166] (2) Align the sequenced sequence with the human genome HG19 to obtain the original BAM file, and use Picard and Samtools to remove the duplicate sequences and sequences aligned to multiple positions to obtain the final BAM file;
[0167] (3) Somatic variation detection was performed using single samples of 20 tumor tissues, with filtering performed using only the database, filtering using the SGZ algorithm published by Foundation Medicine, filtering using the method disclosed in the prior art (CN202210605149.7), and filtering using the methods 1-9 in the above embodiments. In addition, somatic variation detection of tumor tissue-leukocyte paired samples was used as the gold standard to evaluate the accuracy of the two somatic variation detection methods.
[0168] (4) A total of 84 variant sites were detected. The comparison results of the methods are shown in Table 1.
[0169] Table 1
[0170] Filter method Somatic mutation detection accuracy Database filtering (68 / 84)81% This application method (82 / 84)97.6% SGZ algorithm (78 / 84)92.9% Prior art methods Method failure
[0171] From the above description, it can be seen that the above-mentioned embodiments of the present invention achieve the following technical effects: The present invention proposes a method for somatic variation detection based on a single sample of tumor tissue, which can make up for the defects of insufficient germline sites included in the database and easy misjudgment based on the original variation frequency of the variation site. This method does not require sequencing of paired tissues, can reduce the requirements for samples, and save costs at the same time, to achieve the same purpose of somatic variation detection. Germline variations are filtered in combination with database variations and variation site frequencies to make somatic variation detection more accurate. The variation frequency is corrected according to the purity of the tumor tissue and the number of copies of the gene to improve the accuracy of subsequent detection. It helps researchers to conduct retrospective studies on previous research data.
[0172] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. For those skilled in the art, the present invention may have various modifications and variations. Any modification, equivalent replacement, improvement, etc. made within the spirit and principle of the present invention shall be included in the protection scope of the present invention.
Claims
1. A method for detecting somatic cell variation from a single sample tumor tissue for non-disease diagnosis purposes, characterized in that: The method comprises: (1) Analyze the offline data of the test sample to obtain the filtered variant sites, and calculate the posterior probability that the filtered variant sites are germline variants and somatic variants, which are recorded as germline variant parameters and somatic variant parameters; (2) constructing a mutation site feature vector using the background noise parameter, the germline variation parameter, and the somatic variation parameter of each of the filtered variation sites; (3) inputting the mutation site feature vector into a somatic variation discrimination model to determine whether each of the filtered variation sites is a somatic variation or a germline variation; The (1) includes: (1.1) Obtaining the offline data of the test sample, comparing the offline data with the reference genome, and obtaining BAM file data, wherein the test sample is the single sample tumor tissue; (1.2) analyzing the BAM file data to obtain the initial mutation site, the depth of the initial mutation site, and the number of reads carrying the initial mutation site; (1.3) filtering the initial variant sites to obtain the filtered variant sites, the depth of the filtered variant sites, and the number of reads carrying the filtered variant sites; (1.4) Calculating the germline variation parameter and the somatic variation parameter, including: a. dividing the gene segments of the reference sample and the test sample, and obtaining multiple windows by equally dividing each gene of the reference sample and the test sample; Calculate and obtain the copy number baseline of each window in the reference sample, Calculate the ratio of the number of read segments in each window in the test sample, Using the copy number baseline of each window, correcting the read number ratio of the corresponding window to obtain a correction value; The correction value is used to calculate the gene copy number corresponding to each gene in the test sample; The reference sample includes a white blood cell sample; b. For each gene in the test sample, statistically obtaining the actual mutation probability of the heterozygous germline variation carried on the gene, expressed as AF actual; The gene copy number is used to calculate the theoretical mutation frequency of the gene carrying germline variation, which is expressed as AF theory, wherein AF theory = p×Vi / (p×cn+2×(1-p)), where p is tumor purity, cn is the gene copy number, and Vi is the minor allele copy number of the gene; Using the grid search method, ∑(AF theory - AF actual) was calculated by varying tumor purity and minor allele copy number. 2 The minimum value of , to obtain the true tumor purity and minor allele copy number; c. Calculate the theoretical probability of germline variation and the theoretical probability of somatic variation of the test sample according to the gene copy number and the true tumor purity and the minor allele copy number calculated in b; the theoretical probability of germline variation is represented by AFgermline, AFgermline=(p×Vi+1-p) / (p×cn+2×(1-p)); the theoretical probability of somatic variation is represented by AFsomatic, AFsomatic=p×Vi / (p×cn+2×(1-p)); wherein p is the true tumor purity, Vi is the true minor allele copy number of the gene, and cn is the gene copy number; According to the theoretical probability of germline variation, the theoretical probability of somatic variation, the depth of the filtered variation site and the number of reads carrying the filtered variation site, the posterior probability when each filtered variation site is a germline variation or a somatic variation is calculated based on the binomial distribution formula, and the posterior probability when the filtered variation site is a germline variation or a somatic variation is recorded as a germline variation parameter and a somatic variation parameter, respectively; The posterior probability when the variant site is a germline variant is expressed as P(y|G; AFgermline), P(y|G; AFgermline)=Bin(f, n, AFgermline), which is the germline variant parameter; the posterior probability when the variant site is a somatic variant is expressed as P(y|G; AFsomatic), P(y|G; AFsomatic)=Bin(f, n, AFsomatic), which is the somatic variant parameter; wherein n is the depth of the filtered variant site, and f is the number of reads carrying the filtered variant site; In (2), the method for obtaining the background noise parameter includes: Perform mutation detection on the white blood cell sample to obtain the initial mutation site of the white blood cell, remove the germline mutation site in the initial mutation site, and use the mutation frequency of the remaining mutation sites to build a background noise library; count the background noise frequency of the filtered mutation site in the background noise library, and if it does not exist, the background noise frequency is recorded as 0; use the depth of the filtered mutation site, the number of reads carrying the filtered mutation site and the background noise frequency, and calculate the posterior probability when the filtered mutation site is the background noise through the binomial distribution formula, and record the posterior probability as the background noise parameter of the filtered mutation site; The posterior probability when the filtered variant site is background noise is expressed as P(y|G; AFnoise), P(y|G; AFnoise)=Bin(f, n, AFnoise), wherein AFnoise is the background noise frequency recorded by the white blood cells at the corresponding site in the background noise library, n is the depth of the filtered variant site, f is the number of reads carrying the filtered variant site, and P(y|G; AFnoise) is the background noise parameter; In (3), the construction of the somatic variation discrimination model includes: obtaining sites known to be background noise, sites known to be germline variations, and sites known to be somatic variations, constructing feature vectors of known mutation sites, and using the feature vectors of known mutation sites and corresponding mutation site types to train a decision tree model to obtain the somatic variation discrimination model.
2. The method according to claim 1, characterized in that The data of the test sample includes large set data or small set data. The large set data is sequencing data of a region larger than 1Mb captured by hybridization of a high-throughput sequencing platform, and the small set data is sequencing data of a region between 100Kb and 1Mb captured by hybridization of the high-throughput sequencing platform.
3. The method according to claim 1, characterized in that The filtering of the initial variant sites in (1.3) includes: data self-filtering and / or database-based filtering.
4. The method according to claim 3, characterized in that The data filtering itself includes: Filter and remove the initial variant sites with a variant frequency of <0.2%; Filtering to remove the initial variant sites with chain preference; The initial variant sites with supporting read numbers <3 were filtered out.
5. The method according to claim 4, characterized in that The strand preference includes a forward read or a reverse read supporting the initial variant site being 0.
6. The method according to claim 3, characterized in that: The database-based filtering includes: using the population database 1000Genome, ExAC and the SNP database dbSNP to perform database filtering on the initial variant sites.
7. The method according to claim 6, characterized in that The database-based filtering includes: filtering and removing the initial variation sites with a population frequency > 0.01 in the population database 1000Genome and the ExAC, and filtering and removing the initial variation sites annotated as SAO=1 in the SNP database dbSNP.
8. The method according to claim 1, characterized in that The calculation of the gene copy number in a includes: a1. Using more than 20 white blood cell samples to construct a reference set, the ratio of the number of reads in each window of the reference set to the total number of reads sequenced in the sample is counted as the copy number baseline of each window, which is represented by RDreference; a2. The hybrid capture region of the test sample is divided into a plurality of the windows, and the number of reads of the test sample in each window is counted as the ratio of the total number of reads of the data from the machine, represented by RDtest; a3. Calculating the correction value RDcor of each window after the copy number baseline correction, RDcor = (RDtest-RDreference) / RDreference; a4. Remove the maximum and minimum values of the RDcor values of all the windows in a gene, and calculate the average value mean(RDcor) of the RDcor of the remaining windows. The gene copy number = 2×mean(RDcor).
9. The method according to claim 1, characterized in that: In b, the grid search method includes: substituting the tumor purity from 0.01 to 0.99, with a step length of 0.01, and the value of Vi is an integer from 1 to cn-1, into the formula ∑ (AF theory-AF actual) in sequence: 2 Calculation is performed in the formula, when the result of the formula is the smallest, the tumor purity taken in the formula is the true tumor purity, and the minor allele copy number taken is the true minor allele copy number; the true tumor purity and minor allele copy number of each gene are calculated in the formula, the tumor purity of different genes is the same, and the minor allele copy numbers of different genes are the same or different.
10. The method according to claim 9, characterized in that The statistical method of AF actual includes: counting the polymorphic sites with a population frequency of 0.01-0.99 included in the population database 1000Genome, ExAC and the SNP database dbSNP, and if the mutation site is a heterozygous mutation site, the mutation frequency of the polymorphic site of the test sample is recorded as AF actual.
11. The method according to claim 10, characterized in that The population frequency is 0.05-0.
95.
12. The method according to claim 10 or 11, characterized in that: The method for determining the heterozygous mutation site includes: if the mutation frequency of the polymorphic site is within the range of 20%-80%, the polymorphic site is considered to be a heterozygous mutation site; if the mutation frequency of the polymorphic site is <20% or >80%, the polymorphic site is considered not to be a heterozygous mutation site.
13. The method according to claim 1, characterized in that The mutation site feature vector in (2) is {background noise parameter, germline variation parameter, somatic variation parameter}.
14. The method according to claim 1, characterized in that In (2), constructing the background noise library includes: Using more than 100 white blood cell samples, performing hybridization capture sequencing on the white blood cell samples to obtain second-generation sequencing data, and comparing with the reference genome to obtain a white blood cell BAM file; Analyze the leukocyte BAM file to obtain the initial mutation site of the leukocyte, and count the mutation frequency of each leukocyte sample at all sites in the hybrid capture area. If the leukocyte sample does not have a mutation at the first site, the background noise frequency of the first site is recorded as 0; If the mutation frequency of the leukocyte sample at the second site is greater than 20%, the second site is considered to carry a germline mutation, and the background noise frequency of the second site is recorded as 0; For the remaining sites among all the sites that do not belong to the first site or the second site, the true mutation frequency of the white blood cell sample at the remaining sites is recorded as the background noise frequency, and for each single site among the remaining sites, the true mutation frequency of the single site is the median of the mutation frequencies of all the white blood cell samples at the single site; The background noise frequencies of all sites in the hybrid capture area constitute the background noise library.
15. The method according to claim 1, characterized in that In (3), the construction of the somatic variation discrimination model includes: The known mutation site feature vector is constructed using more than 100 sites known to be background noise, more than 100 sites known to be germline mutations, and more than 100 sites known to be somatic mutations.
16. A computer device comprising a memory, a processor and a computer program stored in the memory, characterized in that: The processor executes the computer program to implement the steps of the method according to any one of claims 1 to 15.
17. A computer-readable storage medium having a computer program stored thereon, characterized in that: When the computer program is executed by a processor, the steps of the method according to any one of claims 1 to 15 are implemented.
18. A computer program product, comprising a computer program, characterized in that When the computer program is executed by a processor, the steps of the method according to any one of claims 1 to 15 are implemented.
Citation Information
Patent Citations
A Single-Sample Tumor Somatic Mutation Detection and TMB Detection Method Based on NGS Platform
CN114694750B
Method and device for detecting TMB
CN110689930A
Single-sample tumor somatic mutation discrimination and TMB detection method based on NGS platform
CN114694750A