JP6675164B2 - Mutation judgment method, mutation judgment program and recording medium - Google Patents

Description

  The present invention provides a mutation detection method executed by a computer to determine the presence or absence of a mutation in a target gene, a mutation detection program executed by a computer to determine the presence or absence of a mutation in the target gene, and the mutation detection program. The present invention relates to a recorded recording medium.

  In recent years, technology for determining the presence or absence of a gene mutation using a sequencer has been developed, and its usefulness as a tool for diagnosing a disease has been recognized. That is, a region including a known mutation in a specific disease-related gene is sequenced, and the presence or absence of the mutation is determined based on the obtained sequence data, whereby a disease can be diagnosed.

  An example of such a technique for determining the presence or absence of mutation is amplicon sequencing using a next-generation sequencer. In amplicon sequencing, first, an exon region of a target gene is subjected to PCR amplification, and the obtained PCR product is sequenced by a next-generation sequencer. Then, by mapping a read obtained by sequencing to a predetermined reference sequence, the presence or absence of a mutation can be evaluated and determined.

  Such a technique for determining the presence / absence of mutation has the advantage that the presence / absence of mutation in a large number of samples can be determined at low cost, but the problem that erroneous determination may occur due to sequence errors and the like. There is. In particular, as in the case of detecting mutations in ctDNA (circulating tumor DNA) derived from cancer cells in the bloodstream of a cancer patient, samples containing only a small amount of the target DNA or contaminations other than the target DNA are often used. This problem becomes remarkable when mutation is detected in a sample having a high background such as a sample containing the mutation.

  In contrast, Kukita et al. Have reported that they have established a method for detecting a mutation in the epidermal growth factor receptor (EGFR) gene, which is a cancer-related gene, in ctDNA (Non-Patent Document 1). ). In the method of Non-Patent Document 1, a threshold is calculated for each mutation location based on a statistical model estimated from sequence data of a normal sample having no mutation, and the sequencing result of the sample to be determined is compared with the threshold. To detect the mutation.

Kukita et al., PLoS One. 8 (11) e81468, 2013.

  In the technique described in Non-Patent Document 1, a threshold value is calculated in advance from sequence data of a normal sample, and thereafter, a sample to be determined is sequenced, and mutation is detected based on the result. Specifically, a statistical model is estimated as a probability distribution of the number of mutant reads detected in 100,000 reads from sequence data of a normal sample, and a threshold value is calculated from the statistical model. The sample is subjected to deep sequencing of 100,000 or more reads, the obtained result is converted into the number of mutant reads detected in 100,000 reads, and then compared with a threshold value.

  At this time, the number of reads obtained from the sample may differ for each sample and each region. For example, when sequencing is performed after simultaneously PCR-amplifying a plurality of regions for a sample, the amplification rate in each region may not be constant, and as a result, the number of reads obtained from the sample differs for each sample and for each region. Would.

  In the technique described in Non-Patent Document 1, a statistical model is estimated in advance as a probability distribution of the number of mutated reads detected in a fixed number of reads, and the sequencing result of a sample is read in a fixed number of reads. If the number of reads is different for each sample and each region, the judgment criteria may change for each sample and each mutation, and accurate judgment may not be possible. .

  In view of the above situation, a main object of the present invention is to provide a technique for judging the presence or absence of mutation with high accuracy even when the number of reads obtained from a sample is different for each sample and each region. is there.

  In order to solve the above problems, a mutation determination method according to one embodiment of the present invention is a mutation determination method that is performed by a computer to determine the presence or absence of a mutation in a target gene, wherein a polynucleotide derived from a specimen is From a plurality of reads obtained by sequencing, a corresponding read corresponding to a specific region of the target gene, and an extraction step of extracting a mutated read having a specific mutation in the specific region among the corresponding reads. A calculating step of calculating a threshold value for determining the presence or absence of a mutation at a predetermined significant level from the number of extracted corresponding reads, and, when the number of extracted mutant reads exceeds the calculated threshold value, Determining the presence of the mutation.

  According to the above configuration, every time a corresponding lead is extracted, a threshold is calculated based on the number of extracted corresponding leads, and the presence or absence of a mutation at a predetermined significance level is determined using the actual number of corresponding leads. The threshold value for determination can be calculated more appropriately. Thus, even if the number of reads obtained from a sample differs for each sample and each region, the presence or absence of a mutation can be determined with high accuracy.

  In the mutation determination method according to one aspect of the present invention, the computer stores statistical model information indicating a statistical model associated with the specific mutation, and refers to the statistical model information in the calculation step. Then, a threshold value for determining the presence or absence of a significant difference of a predetermined significance level from the number of corresponding leads may be calculated.

  According to the configuration, by referring to the statistical model information stored in the computer, using the statistical model associated with the specific mutation, the threshold for determining the presence or absence of the mutation at a predetermined significance level It can be suitably calculated.

  In the mutation determination method according to one aspect of the present invention, the statistical model information may include a type of a statistical model associated with the specific mutation.

  According to the above configuration, since the statistical model information includes the type of the statistical model (for example, Poisson distribution, negative binomial distribution, and the like), it is more preferable to use a probability distribution model suitable for each mutation. Can be calculated.

  In the mutation determination method according to one aspect of the present invention, the statistical model information may further include normal sample data relating to the specific mutation.

  According to the above configuration, since the statistical model information includes the normal sample data, the threshold value for determining the presence or absence of the mutation at a predetermined significance level can be preferably calculated.

In the mutation determination method according to one aspect of the present invention, the computer is stored in association with the specific mutation, minimum value information indicating whether to use the minimum value of the threshold,
In the calculation step, when the minimum value information indicates that the minimum value of the threshold is to be used, the minimum value of the threshold calculated from the number of corresponding leads is determined by referring to a predetermined statistical model. The threshold may be calculated.

  According to the above configuration, it is possible to prevent the threshold from becoming excessively small. Thereby, the reliability of the determination result can be improved.

  In the mutation determination method according to an aspect of the present invention, the computer stores a plurality of reference sequences, and the extracting step maps the plurality of reads to each of the plurality of reference sequences. It is preferred to include.

  According to the above configuration, by mapping the plurality of reads to each of the plurality of reference sequences, a corresponding read corresponding to a specific region of the target gene from the plurality of reads, Mutant reads having the above mutations can be suitably extracted.

  In the mutation determination method according to one aspect of the present invention, when the specific mutation is a single nucleotide substitution mutation, the reference sequence includes a specific reference sequence that is a sequence of the specific region, and the extraction step In the method, the read mapped to the specific reference sequence is extracted as the corresponding read, and it is determined whether or not each of the extracted corresponding reads has the specific mutation. The read may be extracted as the mutant read.

  According to the configuration, when the specific mutation is a single-base substitution mutation, a corresponding lead and a mutated lead can be suitably extracted.

  In the mutation determination method according to one aspect of the present invention, when the specific mutation is at least one mutation of insertion and deletion, the reference sequence is a specific reference sequence that is a sequence of the specific region. A detection reference sequence which is a sequence in which at least one of insertion and deletion has occurred in the specific reference sequence. In the extraction step, the read and the detection reference sequence mapped to the specific reference sequence are included. The mapped read may be extracted as the corresponding read, and the read mapped to the detection reference sequence may be extracted as the mutant read.

  According to the configuration, when the specific mutation is a mutation consisting of at least one of an insertion and a deletion, a corresponding lead and a mutated lead can be suitably extracted.

  In the mutation determination method according to one aspect of the present invention, the specific mutation is a mutation consisting of at least one of a predetermined base insertion and a predetermined base deletion, and the detection reference sequence A plurality of detection reference sequences corresponding to mutually different mutations, and in the extraction step, among the detection reference sequences, the read mapped to the detection reference sequence corresponding to the specific mutation is read. It may be extracted as the mutant read.

  According to the above configuration, when the specific mutation is a mutation consisting of at least one of insertion of a predetermined base and deletion of a predetermined base, in other words, the specific mutation is a specific mutation. In the case of a type of insertion and / or deletion mutation, the corresponding lead and the mutated lead can be suitably extracted. As a result, insertion and / or deletion mutations can be detected for each type as compared with the prior art in which it was only possible to determine whether or not insertion and / or deletion mutations exist. Can be obtained. In addition, since a threshold for determining the presence or absence of an insertion and / or deletion mutation is calculated for each type, more accurate determination can be performed.

  In the mutation determination method according to an aspect of the present invention, in the extraction step, among the mutation reads mapped to a single detection reference sequence, a mismatch between the detection reference sequence and the detection reference sequence When the ratio of reads exceeds a predetermined ratio, the mismatched read may be excluded from the corresponding read and the mutant read.

  According to the above configuration, when the specific mutation is at least one mutation of insertion and deletion, a mismatched read having a sequence mismatch with a single detection reference sequence is regarded as a low-reliability read. By removing it, the reliability of the determination can be improved.

  In particular, the determination can be performed efficiently by excluding only the case where the mismatched leads occupy the entirety exceeds the predetermined ratio and have a large influence on the determination result.

  In the mutation determination method according to one aspect of the present invention, in the extraction step, the reads mapped to a plurality of the reference sequences may be excluded from the corresponding reads and the mutation reads.

  According to the above configuration, the reliability of the determination can be improved by removing the read mapped to the plurality of reference sequences as a low-reliability read.

  The mutation determination method according to one aspect of the present invention may further include a warning step of outputting a warning when at least one of the leads is excluded from the corresponding leads in the extraction step.

  According to the above configuration, when a lead with low reliability is removed, a warning can be output to allow the user to recognize the reliability of the determination.

  In the mutation determination method according to one aspect of the present invention, the target gene may be an EGFR gene, and the plurality of reference sequences may include each of the nucleotide sequences shown in SEQ ID NOs: 1 to 101. .

  According to the above configuration, a mutation in the EGFR gene can be suitably detected. In particular, the reference sequences include detection reference sequences (SEQ ID NOS: 5 to 13, 15 to 18, 20 to 33, 35 to 101) corresponding to almost all types of known insertion and / or deletion mutations. In addition, since a detection reference sequence (SEQ ID NOs: 14, 19, and 34) is included in consideration of a read error, insertion and / or deletion mutation can be determined with high accuracy.

  The mutation determination method according to one aspect of the present invention may further include, before the extraction step, an update step of updating the reference sequence based on information obtained from an external database.

  According to the above configuration, since a reference sequence based on the latest knowledge from an external database can be used, a desired mutation can be suitably detected.

  In the mutation determination method according to one aspect of the present invention, in the extraction step, the lead shorter than a predetermined length may be removed before extracting the corresponding lead.

  According to the above configuration, the reliability of the determination can be improved by removing the short leads.

  In the mutation determination method according to one aspect of the present invention, the polynucleotide may be a PCR product, and the sequencing may be a sequencing in which a plurality of reads are obtained by overlappingly reading the same region. Good.

  According to the above configuration, the determination of the mutation can be suitably performed.

  The mutation determination program according to one embodiment of the present invention is a mutation determination program executed by a computer to determine the presence or absence of a mutation in a target gene, the plurality of reads being obtained by sequencing a polynucleotide derived from a specimen. From, a corresponding read corresponding to a specific region of the target gene, and, among the corresponding reads, an extraction step of extracting a mutated read having a specific mutation in the specific region, and from the number of extracted corresponding reads A calculation step of calculating a threshold for determining the presence or absence of a mutation at a predetermined significance level, and a determination step of determining that the specific mutation is present when the number of extracted mutation reads exceeds the calculated threshold. , Run.

  According to the above configuration, the same effects as those of the mutation determination method according to one embodiment of the present invention can be obtained.

  Further, a computer-readable recording medium that records the mutation determination program according to one embodiment of the present invention is also included in the scope of the present invention.

  According to the present invention, the presence or absence of mutation can be determined with high accuracy even when the number of reads obtained from a sample differs for each sample and each region.

It is a schematic diagram of a computer etc. which perform the mutation judging method concerning one embodiment of the present invention. It is a figure showing the flow of the mutation judging method concerning one embodiment of the present invention. It is a figure showing the flow of SNV detection processing concerning one embodiment of the present invention. FIG. 4 is a diagram illustrating a flow of an InDel detection process according to an embodiment of the present invention. It is a figure showing the flow of the selection processing of the statistic selection model concerning one embodiment of the present invention. It is a figure showing the flow of the threshold value decision processing concerning one embodiment of the present invention. It is a figure showing the flow of the threshold value decision processing of the Poisson distribution concerning one embodiment of the present invention. It is a figure showing the flow of the threshold value decision processing of the negative binomial distribution concerning one embodiment of the present invention. FIG. 4 is a diagram schematically illustrating the number of template arrays and the status of mapping according to the embodiment of the present invention. FIG. 4 is a diagram showing an example of an improvement in mismapping due to an increase in the number of template sequences using a tumor sample of published data of a document according to an embodiment of the present invention.

  Embodiments of the present invention will be described below. Note that the present invention is not limited to this.

[Definition of terms, etc.]
As used herein, the term "polynucleotide" can be referred to as "nucleic acid" or "nucleic acid molecule", and is intended to mean a polymer of nucleotides. Further, the “base sequence” can also be referred to as “nucleic acid sequence” or “nucleotide sequence”, and unless otherwise specified, intends a sequence of deoxyribonucleotides or a sequence of ribonucleotides. It also includes single-stranded or double-stranded polynucleotides.

  As used herein, the term “specimen” can be referred to as “sample,” which is used synonymously with a specimen or preparation in the art, and as a source of biological material (eg, an individual, a body fluid, a cell line, a tissue culture or Any preparation obtained from a tissue section is contemplated.

  In the present specification, the “polynucleotide derived from a specimen” is also referred to as a “polynucleotide sequence sample” and intends a polynucleotide obtained by processing a specimen. The treatment for obtaining the polynucleotide includes a treatment known in the art, such as a nucleic acid extraction treatment, a nucleic acid amplification treatment, and a reverse transcription treatment.

  As used herein, “read”, “read sequence” or “sequence read” refers to a polynucleotide sequence obtained by sequencing, and refers to a sequence output from a sequencer. As used herein, a read corresponding to a specific region of a gene refers to a lead having a sequence containing the same or an allowable range of mutation as the specific region of the gene. In the present specification, a read corresponding to a specific mutation of a gene is a read corresponding to a region containing the specific mutation of the gene, and refers to a read having a sequence containing the specific mutation.

  As used herein, “InDel (Insertion and / or Deletion)” means a mutation containing an insertion, a deletion, or both an insertion and a deletion. It may be described as "insertion and / or deletion mutation".

  In the present specification, “SNV (Single Nucleotide Variant)” means a single nucleotide substitution mutation.

  As used herein, "depth" refers to the total number of reads obtained by sequencing the same sequence or region.

  As used herein, a "reference" or "reference sequence" is used to determine which region on a gene a read corresponds to and / or which mutation on a gene the read corresponds to. This is the sequence to which the read is mapped. If the mutation is an SNV, the reference sequence can be a wild-type subsequence of the gene. When the mutation is InDel, the reference sequence may include a wild-type partial sequence of a gene and a sequence in which a known or unknown InDel is generated from the wild-type partial sequence. In the present specification, a sequence in which a known or unknown InDel is generated from a wild-type partial sequence of a gene is also referred to as “template” or “template sequence”.

  In this specification, "RE (read error rate)" means a read error rate, and indicates a ratio of the number of reads in which a base read error has occurred to the total number of sequence reads obtained by sequencing. Here, the reading error is an error that is read as a base different from the original sequence. When the sample is a PCR product, it includes both an error generated during sequencing and an error generated during PCR. I do.

  As used herein, "alignment" intends a state of the sequence identified as an order match of one or more nucleotides in the polynucleotide sequence of the sample relative to a reference sequence, which is a known sequence. Such alignment can be performed by a known computer algorithm.

  In this specification, the “statistical model” is also referred to as a “distribution model” or a “probability distribution model”, and when used in this specification, refers to a statistical model representing a probability distribution of a reading error rate.

  As used herein, “parameter” means a variable that characterizes a quantitative data set and / or a numerical relationship between the quantitative data sets, and is described as, for example, a variable applied to a distribution model. .

  As used herein, the term “threshold” is calculated using a distribution model described in the following embodiments, and in the determination of a mutation in a polynucleotide sequence of a sample, an arbitrary numerical value as the number of reads that acts as a cutoff value. Means The method of setting a threshold according to the present invention is described in detail in the following (Step 6. Setting of threshold) and in Examples and the like.

  As used herein, “mapping” means a process for aligning each read with a target reference sequence.

  As used herein, "positive" means that the read sequence has a true mutation to be determined.

  In the present specification, “false positive” means that the read sequence is determined to have a mutation even though it does not have a true mutation to be determined.

  As used herein, "negative" means that the lead sequence does not have the mutation of interest.

  As used herein, "subject" refers to human subjects as well as non-human subjects, such as mammals, invertebrates, vertebrates, fungi, yeast, bacteria, viruses, plants, and the like. Although the examples herein relate to human subjects, the concepts of the present invention are applicable to genomes from any non-human organism, such as an animal or plant, and are useful in fields such as medicine, veterinary medicine and animal science. is there.

  As used herein, "diagnose" or "diagnosis" refers to the identification made by a physician of a disease, disorder or condition based on signs and symptoms of a patient.

The “test” in the present specification refers to a test that does not require identification (diagnosis) by a physician, including a test subject, including a human being, for the presence or absence of a predisposition or the onset of a disease.
[Mutation determination method and mutation determination program]
The mutation determination method according to one embodiment of the present invention is a mutation determination method performed by a computer to determine the presence or absence of a mutation in a target gene, and can be performed by, for example, the computer 10 illustrated in FIG.

  As shown in FIG. 1, the computer 10 includes a processing unit 11, a storage unit 12, and a communication unit 13. The computer 10 executes the mutation determination program 100 stored in the storage unit 12 to execute the program. The mutation determination method according to the embodiment is performed. The mutation determination program 100 includes an instruction group 101, a reference sequence 102, minimum value information 103, and statistical model information 110 for executing each step of the mutation determination method. The statistical model information 110 includes a statistical model type 111 and normal sample data 112.

  Further, the computer 10 is connected to the sequencer 20, and receives a result of sequencing from the sequencer 20. In one embodiment, the computer 10 analyzes the sequencing data using IonPGM (registered trademark), and the mutation detection method according to the present embodiment uses IonPGM (registered trademark). It is executed for the run after the analysis of the sequencing data is completed. Further, the communication unit 13 communicates with the external server 30 via a network.

<About the mutation to be determined>
The types of mutations that can be determined by the mutation determination method according to the present embodiment are InDel and SNV. Therefore, the polynucleotide sequence of the sample includes those having any of the mutations of InDel, SNV, and a combination thereof. That is, the polynucleotide sequence of the sample may be a sequence having any one of insertion mutation, deletion mutation and single nucleotide substitution mutation, and may be a sequence having two or more mutations among these mutations at the same time. There may be. When two or more mutations are present, such as having both SNV and InDel mutations, each mutation may be individually detected by executing an independent detection flow for each mutation. In addition, for example, a sequence having both InDel and SNV, a sequence having a deletion after insertion, and the like are also targets of the determination in the present embodiment. The single nucleotide substitution mutation also includes a silent mutation. The mutation type to be subjected to the mutation determination method according to the present embodiment includes, for example, InDel, and a further example is Deletion. In the present embodiment, a desired number of those having a desired sequence can be set as a reference sequence, so that mismatches in mapping can be reduced, and the mutation of InDel is more preferable than a conventionally known technique. Can be detected.

  The length of the deletion region in the deletion mutation to be determined by the mutation determination method according to the present embodiment is not particularly limited and may be any length, but is, for example, 5 bp to 30 bp, and Is 10 bp to 30 bp. When the deletion region is 5 bp or more or 10 bp or more, the difference from the wild-type reference sequence is sufficiently large, and the determination can be made with higher accuracy.

  The length of the lower limit of the length of the insertion region in the insertion mutation to be determined by the mutation determination method according to the present embodiment is, for example, 5 bp or more, or 10 bp or more. When the insertion region is 5 bp or more or 10 bp or more, the difference from the wild-type reference sequence is sufficiently large, and the determination can be made with higher accuracy. Also, the upper limit length needs to be sufficiently shorter than the length of one sequence read.

<Sequencing technology and sequencing equipment>
The mutation determination method according to the present embodiment preferably targets sequencing data obtained by sequencing using the next-generation sequencing technology. That is, it is preferable to use a next-generation sequencer as the sequencer 20. The next-generation sequencer is a group of base sequence analyzers that have been developed in recent years, and greatly improved analysis by performing massively parallel processing of cloned amplified DNA templates or single DNA molecules in a flow cell. Have the ability. Specific examples include, but are not limited to, Ion PGM (registered trademark) of Life Technology Co., Ltd., and also include devices to be developed in the future.

  Further, the sequencing technology that can be used in the present embodiment may be a sequencing technology that acquires a plurality of reads by reading the same region twice (deep sequencing).

  Examples of sequencing techniques that can be used in this embodiment include ion semiconductor sequencing, pyrosequencing, sequencing-by-synthesis using a reversible dye terminator, sequencing-by-synthesis, and sequencing-by-synthesis. Sequencing techniques that can obtain a large number of reads per run based on sequencing principles other than the Sanger method, such as sequencing-by-ligation and sequencing by oligonucleotide probe ligation. .

  The sequencing primer used for sequencing is not particularly limited, and is appropriately set based on a sequence suitable for amplifying a target region. As for the reagent used for sequencing, a suitable reagent may be selected according to the sequencing technique and the sequencing device used.

  In one embodiment and the examples described below, the amplicon sequence is performed using Ion PGM (registered trademark) of Life Technology. In one embodiment, a sequencing primer for detecting a mutation in the EGFR gene is set in a region of an exon (exons 18, 19, 20, and 21) containing the mutation.

  The sequencing technique according to the present invention also includes a multiplex method in which sequences derived from a plurality of samples are simultaneously sequenced in one run. In the case of the multiplex method, sequence data can be identified for each sample during data analysis by adding a unique “barcode” sequence to each sample in order to distinguish data for each sample. The use of this multiplex method can dramatically reduce the time required for data acquisition even when analyzing a large number of samples.

<Sequence data>
All the read sequences obtained by sequencing using the above-described sequencing synthesizer 20 are used as sequencing data in the present embodiment. By the above-described sequencing, preferably 1,000 or more, more preferably 100,000 or more or more reads may be generated for each region to be determined. The detection sensitivity increases as the number of leads increases. The lead length may have a length depending on the platform used. For example, the sequence read length may be any length as long as it is compatible with the mapping software used, and may be, for example, 100 to 200 bp, and may include a paired-end read. This sequence data is mapped by aligning it with the reference sequence, and a mutation determination process is performed.

  In addition, as the sequence data according to the present invention, for example, data obtained in a FASTA format, a FASTQ format, or the like can be used.

<Determination of reference sequence>
The reference sequence 102 in the present embodiment is determined from a sequence including at least a region where a target mutation exists. One example of a reference sequence is the sequence of a specific region containing at least a mutation based on the sequence of interest. In particular, when the target genomic sequence contains a mutation in an exon region, the reference sequence may be limited to the exon region. The reference sequence used to detect SNV is a wild-type sequence (specific reference sequence) of the target gene, and the reference sequence used to detect InDel is a wild-type sequence (specific reference sequence) of the target gene. Reference sequence) and a sequence containing a known InDel mutation (reference sequence for detection). The specific reference sequence used to detect SNV may be a full-length wild-type genomic sequence of the gene of interest.

  The total length of the reference sequence is not particularly limited, and a sequence having a desired length can be used depending on the purpose. However, the time required to complete the analysis becomes longer depending on the length of the reference sequence.

  In the mutation determination method according to the present embodiment, the number of reference sequences used for one determination is not particularly limited, and is set according to the number of mutations to be detected.

  In one embodiment, the specific reference sequence may be generated for each exon containing the mutation to be determined.

  As an example of such a specific reference sequence, when the mutation determination method is to determine the presence or absence of a mutation in the EGFR gene, the specific reference sequence is a sequence of wild-type EGFR exon 18 shown in SEQ ID NO: 1. A specific reference sequence, a specific reference sequence which is the sequence of wild type EGFR exon 19 shown in SEQ ID NO: 2, a specific reference sequence which is a sequence of wild type EGFR exon 20 shown in SEQ ID NO: 3, and a wild type EGFR exon 21 which is shown in SEQ ID NO: 4 May be included.

  Further, in one embodiment, the reference sequence for detection may be prepared in advance by selecting only frequently-used mutants based on information on mutations obtained from a public database.

  As an example of such a detection reference sequence, when the mutation determination method is to determine the presence or absence of a mutation in the EGFR gene, the detection reference sequence includes SEQ ID NOs: 5 to 11, 13, 15 to 18, Reference sequences for detection shown in 22, 28, 31, 33, 35 and 39 may be included. These reference sequences for detection are obtained by reflecting, with respect to the specific reference sequence, which is the sequence of the wild-type EGFR exon 19 shown in SEQ ID NO: 2, the frequent InDel mutation obtained from the COSMIC database.

  In another embodiment, the sequence of the reference sequence for detection is a sequence that has been created in advance for almost all mutant types obtained based on information on mutations obtained from the latest public database. . By setting a large number of reference sequences corresponding to the target mutation in this way, it is possible to prevent mapping to an incorrect reference and detection as a false positive. As a result, the nucleic acid sequence contained in the sample can be mapped to a more appropriate base sequence corresponding to each variant.

  As an example of such a reference sequence for detection, when the mutation determination method is to determine the presence or absence of a mutation in the EGFR gene, the reference sequence for detection includes SEQ ID NOs: 5 to 13, 15 to 18, 20 to 20. 33, 35 to 101 can be included. These detection reference sequences reflect almost all InDel mutations obtained from the COSMIC database with respect to the specific reference sequence, which is the sequence of the wild-type EGFR exon 19 shown in SEQ ID NO: 2.

  In addition, due to the characteristics of the used sequencer, a specific read error (for example, GG is mistaken for G) may occur with high probability during sequencing. When there is a sequence in which such a read error can occur in the area to be mapped, a reference sequence for detection for responding to the read error may be further prepared.

  The detection reference array for responding to the read error is an array in which at least one read error has occurred with respect to the original array and is incorrect. Such a read error may be one that is known to occur with high probability in the sequencer used, regardless of whether it is an insertion mutation, a deletion mutation, or a single nucleotide substitution mutation. Further, the reference sequence for detection for coping with the read error does not need to be an array in which all possible read errors are reflected, and may be a sequence in which a part thereof is reflected.

  By preparing a detection reference sequence for such a read error, a read that should be positively determined as a mutation is not mapped to a detection reference corresponding to the target mutation because of including a read error. Can be prevented. Therefore, a more accurate determination can be made.

  As an example of such a reference sequence for detection, when the mutation determination method is to determine the presence or absence of a mutation in the EGFR gene, the reference sequence for detection includes (i) SEQ ID NOS: 5 to 13, 15 to 18 (Ii) In addition to the detection reference sequences shown in SEQ ID NOs: 13, 18, and 33, a read error in which GG is erroneously read as G is reflected in the detection reference sequences shown in SEQ ID NOS: 13, 18, and 33. In addition, reference sequences for detection corresponding to read errors shown in SEQ ID NOs: 14, 19, and 34 may be included. In this case, the reference sequence includes the reference sequences shown in SEQ ID NOs: 1 to 101.

  In one embodiment, the sequence information of the reference sequence is prepared in advance based on information on mutations obtained from a publicly known mutation information database based on a wild-type sequence obtained from a publicly available sequence information database. Is what is being done.

  As public sequence information databases, NCBI RefSeq (web page, http://www.ncbi.nlm.nih.gov/refseq/), NCBI GenBank (web page, http://www.ncbi.nlm.nih.gov) / genbank /), UCSC Genome Browser (web page, https://genome.ucsc.edu/) and the like. As publicly known mutation information databases, a COSMIC database (web page, http://www.sanger.ac.uk/genetics/CGP/cosmic/) and dbSNP (web page, http: //www.ncbi.nlm. nih.gov/SNP/). As the reference sequence, a database newly created based on other commercial databases, publicly known documents, experimental data, and the like may be used, or a plurality of these databases may be used in combination. In addition, the number of reference sequences, the length of each sequence, and the like are appropriately set according to the purpose by using sequences obtained from these databases and the like and detailed information on the sequences. For example, a reference sequence may be set for publicly known mutations in consideration of frequency information for each race or animal type. Publicly known mutation information databases having such information include HapMap Genome Browser release # 28 (web page, http://hapmap.ncbi.nlm.nih.gov/cgi-perl/gbrowse/hapmap28_B36/), Human Genetic Variation Browser (web page, http://www.genome.med.kyoto-u.ac.jp/SnpDB/index.html) and (1000 Genomes (web page, http://www.1000genomes.org/) From these databases, for example, Japanese mutation frequency information and the like can be obtained.

  In another embodiment, the reference sequence may be a sequence created based on information on a mutation obtained from the latest public sequence information database in the mutation determination method according to the present embodiment. That is, in one embodiment, in the mutation determination method according to the present embodiment, the communication unit 13 obtains information about mutation from the external server 30 such as a publicly known mutation information database and stores the specific reference stored in the storage unit 12. The method may include a step of creating a detection reference sequence (template sequence) in which the mutation is reflected in the sequence (wild-type sequence) and storing the reference sequence in the storage unit 12.

  For example, in one embodiment, when the mutation determination method is to determine the presence or absence of a mutation in the EGFR gene, the communication unit 13 sets the EGFR gene as an unacquired mutation from the COSMIC database. When the mutation “c.2235_2249del15” in the exon 19 is obtained, the communication unit 13 compares the wild-type sequence of the exon 19 (SEQ ID NO: 2) stored in the storage unit 12 with “c.2235_2249del15”. A new reference sequence for detection (SEQ ID NO: 8) is created by reflecting the mutation, and stored in the storage unit 12.

  In the case of having such a configuration, the nucleic acid sequence contained in the sample can be mapped to the latest mutant type, so that mismatches or erroneous determinations can be reduced, and more accurate mapping can be performed.

  In addition, the reference sequence is appropriately set so that two or more sequences having a certain mutation do not overlap (removal of the possibility of multi-hit). That is, it is sufficient to determine one reference sequence to correspond to one mutation. For example, different specimens are registered as different specimens on the database, but if they have the same mutation, one of them is removed.

<Data analysis flow>
Hereinafter, an example of the flow (data analysis flow) of the mutation determination method according to an embodiment of the present invention will be described in detail. Hereinafter, a flow for detecting a predetermined mutation in the EGFR gene using the EGFR gene as the target gene will be described, but the present invention is not limited to this.

(Step 1. Sequence read data input)
FIG. 2 is a diagram illustrating an example of a data analysis flow according to an embodiment of the present invention. First, the processing unit 11 inputs sequence read data to be used for the determination to the processing unit 11 (step S201). As the input data used for the analysis of the present invention, a sequence result file output after executing a known sequence analysis program can be used. That is, the processing unit 11 reads out the sequence result file and inputs the sequence read data to the processing unit 11.

  In one embodiment, a sequence result file (file name Unaligned BAM) output from IonPGM is used. Since this file is data before normal mapping is performed in IonPGM, it can be used as input data for executing the software of the present invention.

(Step 2. Lead removal)
As shown in step S202, before performing analysis, the processing unit 11 removes read data that is not suitable for use in analysis. A lead that is not suitable for use in analysis means, for example, a lead that is too short. For example, the processing unit 11 removes a read having a length of 70 bp or less from the reads included in the sequence read data input in step S201. Note that the present embodiment is not limited to this, and the processing unit 11 may remove a lead having a length of 80 bp or less or a lead having a length of 60 bp or less as another example. By removing a lead that is too short, a lead that is likely to cause a reading error or the like is removed, so that more accurate analysis can be performed.

(Step 3. SNV detection processing)
Subsequently, as illustrated in step S203, the processing unit 11 performs an SNV detection process. FIG. 3 is a diagram showing a detailed flow of an SNV detection process according to an embodiment of the present invention.

  As shown in FIG. 3, in the SNV detection flow, first, the processing unit 11 reads the sequence read from the reference sequence (the sequence of exons 18, 19, and 20 of wild-type EGFR) stored in the storage unit 12. To extract corresponding reads and mutant reads corresponding to specific regions (exon 18, exon 20, and exon 21) of the target gene (EGFR gene) (step S303, extraction step).

  Mapping can be performed using a known method and a computer program for mapping a read sequence to a reference sequence. In the present specification, “mapping a read sequence to a reference sequence” refers to comparing a read sequence with each reference sequence, and assuming that the read sequence is a reference sequence corresponding to the read sequence (for example, This means that a reference sequence with the highest mapping score is specified. In the present embodiment, the processing is performed by using the bwasw command of BWA 0.6.2 in an apparatus that executes a normal sequencing result analysis program of the IconPGM.

  Many known software or computer algorithms for aligning sequences are available for mapping. Such software includes, but is not limited to, BWA, Bowtie2 and BFAST. Although the algorithm is not particularly limited, Burrows Wheeler Transform (BWT), read hash, genome hash, merge sort, Smith-Waterman, and the like can be used.

  The processing unit 11 extracts, as a corresponding read, a read mapped to a reference sequence of a specific region (exon 18, exon 20, and exon 21) of the target gene (EGFR gene), and extracts the reference sequence and the sequence of the corresponding read. And the mutation in the corresponding read is detected. The storage unit 12 stores mutation information that defines a specific mutation to be detected in the present embodiment. The processing unit 11 compares the reference sequence with the sequence of the corresponding read mapped to the reference sequence. By comparing the obtained mutation with the mutation specified in the mutation information, it can be determined whether or not the corresponding read contains a specific mutation to be detected in the present embodiment. Thereby, the processing unit 11 can extract the corresponding read and the mutated read.

  Subsequently, in step S304, the processing unit 11 removes a mapping result with low reliability and an ambiguous mapping result. Specifically, the processing unit 11 removes the corresponding reads (and the mutated reads) mapped to the plurality of reference sequences from the corresponding reads (and the mutated reads) extracted in Step S303. Thereby, the reliability of the determination result can be improved. In the present embodiment, in a mapping result obtained by using the BWA, a read that hits a plurality of references is added to a BAM (sequence result file) as (mappingQuality = 0) and (QC fail). I have. The processing unit 11 may refer to such information.

  Subsequently, in step S305, the processing unit 11 sets the minimum read number required for mutation determination in the variable R. Such a value of the minimum number of reads can be defined in advance, for example, as a value such as 1000 or 10,000.

  Subsequently, the processing unit 11 repeats steps S306 and S314 while sequentially changing the mutations to be determined (mutation 1, mutation 2,..., Mutation N).

  First, in step S307, the processing unit 11 sets the variable X as normal sample data measured for each mutation, the variable M as a distribution model defined for each mutation, and the variable D as the total number of reads mapped to the target mutation location. (Number of corresponding reads) and the number of mutation reads mapped to the target mutation location are substituted for variable A.

  Next, in step S308, the processing unit 11 compares the variable D with the variable R, and if the variable D <the variable R (Yes), that is, the total number of reads mapped to the target mutation location in the variable D If (the number of corresponding reads) is less than the minimum number of reads necessary for mutation determination, it is determined that no detection (N / D) has been made in step S313, and the process proceeds to the next mutation determination. If variable D <variable R is not satisfied (No), the process proceeds to step S309.

  In step S309, the processing unit 11 assigns the data assigned to the variable X as “normal specimen data for the target mutation”, the data assigned to the variable M as the “distribution model of the target mutation”, and “observed A threshold value determination flow is executed using the value assigned to the variable D as the “total number of reads”, and a threshold value for determining the presence or absence of a mutation at a predetermined significance level is calculated (calculation step). Details of the threshold value determination flow will be described later.

  Subsequently, in step S310, the processing unit 11 compares the variable A with the threshold calculated in step S309, and when variable A ≧ the threshold (Yes), determines that the variable is positive (step S311; determination step). If the variable A is not equal to or larger than the threshold (No), it is determined to be negative (step S312).

  The above is performed for each mutation to be determined.

  In the case of a wild-type sequence having no mutation, it may be determined as a wild-type sequence. Further, when a plurality of different types of mutations are simultaneously detected as mapping candidates in the same sample, the mutations may be detected by executing the above-described detection flow independently for each mutant type. .

(Step 4. InDel detection process)
Subsequently, as illustrated in step S204, the processing unit 11 performs an InDel detection process. FIG. 4 is a diagram showing a detailed flow of the InDel detection process according to an embodiment of the present invention.

  As shown in FIG. 4, in the detection flow of InDel, first, the processing unit 11 reads the sequence read from the reference sequence (the sequence of the exon 19 of the wild-type EGFR and the detection target) stored in the storage unit 12. (The sequence of exon 19 of EGFR in which InDel is reflected) to extract corresponding reads and mutant reads (step S402, extraction step).

  The processing unit 11 extracts a read mapped to a reference sequence of a specific region (wild type or exon 19 having InDel) of a target gene (EGFR gene) as a corresponding read, and extracts a reference having a specific mutation (InDel). The corresponding reads mapped to the sequence are extracted as mutant reads.

  Subsequently, in step S403, the processing unit 11 removes the mapping result with low reliability and the ambiguous mapping result. Specifically, the processing unit 11 removes the corresponding reads (and the mutated reads) mapped to the plurality of reference sequences from the corresponding reads (and the mutated reads) extracted in Step S402. Thereby, the reliability of the determination result can be improved.

  In addition, in step S404, the processing unit 11 determines that, among the reads mapped to the reference sequence (template sequence) having InDel, 10% or more of the reads have a mismatch (mismatch) with the template. Next, the read having the mismatch (mismatched read) is removed from the corresponding read and the mutated read. Note that, in the IngPGM sequencing result analysis program, a read having a mismatch may cause the BAM to add information to the BAM as a QC fail, for example, and the processing unit 11 may refer to this. Note that the present embodiment is not limited to this, and the processing unit 11 has, as another example, the mismatch in the case where the read of the ratio of 5% to 20% or more has a mismatch. The leads may be removed. Thereby, the reliability of the determination result can be improved.

  In addition, for the removal of a lead having a mismatch, for example, a process of removing 10% or more of mismatched leads from data after removing 10% or more of mismatched leads is repeated until there is no mismatched lead. May be performed. As described above, by performing repeated removal, mismatched leads can be removed as much as possible.

  Note that a mismatch refers to the presence of one or more differences in the read sequence from the reference sequence. Here, the one or more differences are at least one base difference, at least one base deletion and / or insertion, and a combination thereof. These reads can cause false positives if there is a mismatch at or above the percentages described above. For example, reads containing new mutations not present in the reference sequence can be removed as mismatches.

  Subsequently, in step S405, the processing unit 11 sets the minimum read number required for the mutation determination in the variable R.

  Subsequently, the processing unit 11 repeats steps S406 and S414 while sequentially changing the mutations to be determined (mutation 1, mutation 2,..., Mutation N).

  First, in step S407, the processing unit 11 sets the variable X to normal sample data measured for each mutation, the variable M to a distribution model defined for each mutation, and the variable D to the total number of reads mapped to the target mutation location. (Number of corresponding reads) and the number of mutation reads mapped to the target mutation location are substituted for variable A.

  Next, in step S408, the processing unit 11 compares the variable D with the variable R, and when the variable D <the variable R (Yes), that is, the total number of reads mapped to the target mutation location in the variable D If (the number of corresponding reads) is less than the minimum number of reads required for mutation determination, it is determined that no detection (N / D) has been made in step S413, and the process proceeds to determination of the next mutation. If variable D <variable R is not satisfied (No), the process proceeds to step S409.

  In step S409, the processing unit 11 determines that the data substituted for the variable X as “normal specimen data for the target mutation” and the data substituted for the variable M as the “target mutation distribution model” A threshold determination flow is executed using the value substituted for the variable D as the “total number of reads”, and a threshold for determining the presence or absence of a mutation at a predetermined significance level is calculated (calculation step). Details of the threshold value determination flow will be described later.

  Subsequently, in step S410, the processing unit 11 compares the variable A with the threshold calculated in step S409, and when variable A ≧ the threshold (Yes), determines that the variable A is positive (step S411, determination step). If variable A is not equal to or larger than the threshold value (No), it is determined to be negative (step S412).

  The above is performed for each mutation to be determined. In the present embodiment, the InDel mutation to be determined is a specific type of InDel mutation (that is, a mutation consisting of at least one of insertion of a predetermined base and deletion of a predetermined base). The detection reference sequence includes a plurality of detection reference sequences corresponding to mutations different from each other. In step S402, the detection reference sequence corresponds to the specific type of InDel mutation to be determined among the detection reference sequences. The read mapped to the detection reference sequence to be detected is extracted as a mutant read. Thus, the presence or absence of a specific type of InDel mutation can be determined successfully. By this means, since the InDel mutation can be detected for each type as compared with the related art in which it was only possible to determine whether or not some InDel mutation exists, more detailed mutation information can be obtained. In addition, since a threshold value for determining the presence or absence of the InDel mutation is calculated for each type, more accurate determination can be performed.

  Further, in one embodiment, in parallel with the above-described InDel detection processing, in step S402, as in the related art, all the corresponding reads mapped to any of the detection references are extracted as mutant reads. It may be determined whether or not any InDel mutation exists. Accordingly, even when an InDel mutation that does not correspond to the detection reference exists, it is possible to detect the presence of some InDel mutation.

  In one embodiment, when the detection reference sequence includes the detection reference sequence for coping with the read error described above, in step S402, the detection corresponding to the mutation to be determined which does not include the read error is performed. Of the read mapped to the reference sequence for detection and the read mapped to the reference sequence for detection corresponding to the read error in which the read error is reflected on the reference sequence for detection, May be extracted as a mutation read corresponding to the mutation of Thereby, a read that is not mapped to the reference for detection corresponding to the mutation to be determined due to including a read error is also mapped to the reference sequence for detection to correspond to the read error, so It can be extracted as the corresponding mutant read. Thereby, erroneous determination can be prevented, and more accurate determination can be performed.

  Finally, the obtained mutation determination result is output. For example, the processing unit 11 may output a result to a Web browser in an HTML format. In addition, the processing unit 11 may output the determination result of the mutation as a text report.

<Distribution model selection work flow>
Before describing the threshold value determination flow in steps S309 and S409, a distribution model selection work flow for selecting a distribution model used to calculate a threshold value for each mutation will be described. In one embodiment, the distribution model selection work flow is performed in advance, and the result is stored in the statistical model information 110 of the mutation determination program 100.

  As a statistical model for setting a threshold, a model created in advance from normal sample data is applied. Here, the normal sample refers to a sample that has been proven not to have a mutation, and the normal sample data is an extraction step in the mutation determination method according to the present embodiment, using a normal sample as a determination target sample. (Steps S303, S402) indicates a set of the number of corresponding reads (total number of reads) and the number of mutated reads obtained by performing steps up to (Steps S303 and S402). The normal sample data includes data obtained from a plurality of normal samples. In the present embodiment, a statistical model that can predict background noise from such normal sample data is created in advance. Specifically, a distribution model to be used for determining a threshold is selected according to the flow shown in FIG.

  A procedure for creating a statistical model, selecting an appropriate statistical model, and criteria for selecting a statistical model will be described below. Hereinafter, a description will be given of a mode in which the computer 10 assists the determination of the worker, and finally the determination is made by the worker.

  The processing unit 11 obtains a set of the total number of reads and the number of mutated reads from the normal sample data for the target mutation (step S501). Subsequently, in order to assist comprehensive judgment, the processing unit 11 displays the relationship between the total number of reads and the number of mutated reads in the normal sample data on a display unit (not shown) as a scatter diagram and a histogram of the data (step). S502). Then, the worker visually confirms the state of the distribution in each model as a figure, thereby using the distribution as a judgment material for examining the validity (step S503).

  Subsequently, in order to further confirm the validity of the distribution model, the parameters of each model are estimated and the validity is examined. For example, the processing unit 11 first applies a Poisson distribution as a distribution model, and applies the total number of reads and the number of variant reads of normal sample data to a generalized linear model (GLM) of the distribution model. Then, the processing unit 11 calculates the PP plot, the KS (Kolmogorov-Smirnov (Kolmogorov-Smirnov)) test result, the deviance, and the AIC for the fitted model, and displays them on a display unit (not shown). The operator visually checks the displayed data to determine the validity (step S504). Next, the processing unit 11 applies a negative binomial distribution as a distribution model, and applies the same to the distribution model of the normal specimen data in the same manner as when the Poisson distribution is applied. After that, similarly to step S504, the PP plot, the KS test result, the deviation, and the AIC are displayed, and the operator visually confirms them (step S505).

  Next, the processing unit 11 checks whether or not each model has sparse data (step S506). Here, the sparse data refers to data in which the RER is mostly 0 and has a slightly smaller value, or is composed of all 0s.

  If it is determined in step S506 that the data is sparse (Yes), the “minimum value” model is used (step S507). Specifically, information indicating that the minimum value of the threshold is used is stored in the minimum value information 103 of the mutation determination program 100 in association with the mutation.

  If it is determined in step S506 that the data is not sparse data (No), the worker selects a suitable statistical model (step S508). At this time, (i) the better the shape of the PP plot is closer to a straight line, the better the (ii) the larger the p-value of the KS test, and (iii) the lower the AIC, the better the prediction. And (iv) overdispersion if the degree of residue deviation / degree of freedom greatly exceeds one.

  Then, the processing unit 11 stores the type of the statistical model selected by the operator in the statistical model type 111 of the mutation determination program 100, and stores the normal sample data in the normal sample data 112. In another embodiment, instead of storing the normal sample data, information indicating the statistical model created in step S504 or S505 may be stored in the statistical model information 110 of the mutation determination program 100. Further, information indicating that the minimum value of the threshold is not used may be stored in the minimum value information 103 of the mutation determination program 100 in association with the mutation.

<Correction of threshold>
In the present embodiment, a statistical model is created from a data set of analysis results of normal samples, and is used for calculating a threshold. The threshold may be corrected by excluding outliers from this data set on a fixed basis.

  As an example of a method of determining an outlier, first, the degree of conformity with the statistical model to be used is calculated for the data of the normal sample used for the threshold calculation. The degree of conformity with the statistical model is determined based on the p-value. The lower the p-value, the more out of the statistical model. The p-value is similarly calculated for all mutations. Subsequently, the p-values of all the mutations are arranged in the order of the numerical values, and the p-value as a reference for the outlier is determined. The outlier is, for example, p <0.02.

  Further, the statistical model may be corrected by feeding back data newly input in the mutation determination method as test data.

<Limit of detection (LOD)>
The detection limit used in the method of the present invention can be calculated from data of a normal sample and a positive sample. However, the setting of the detection limit may not be completely automated, and may be manually performed by those skilled in the art executing the program.

  In Non-Patent Document 1 described above, the threshold is determined in consideration of only an error of a normal sample. However, in order to actually set the detection limit, it is considered that not only data of a normal sample but also data of a positive sample are required. For the method of setting the detection limit, for example, determine the percentage of mutant reads such that the probability that the number of mutant reads in the sample to be determined falls below a threshold is equal to or lower than a certain level, and determine the minimum value as the detection limit. And so on. However, setting the detection limit is different from the setting of the threshold, and it is necessary to consider not only the possibility of erroneously determining a normal type as a mutant type but also the possibility of determining a mutant type as a normal type. .

As a method of calculating the LOD applied to the present invention, a known method (therascreen (registered trademark) KRAS RGQ PCR kit instruction manual (http://www.accessdata.fda.gov/cdrh_docs/pdf11/P110027c.pdf)) Can be used. For example, there is a method of calculating in advance by the following methods 1) to 5).
1) Dilution of a sample whose concentration of the target substance is known is performed in 5 to 10 steps.
2) For each dilution step, prepare samples so that 10-20 repetitions can be performed.
3) The target substance is detected, and the ratio of the number of positives to the number of repetitions in each dilution step is determined.
4) A function in which the concentration in each dilution step is an independent variable and the ratio of the number of positives to the number of repetitions is a dependent variable is determined.
5) From the function obtained in 4), obtain a dilution concentration such that the ratio of positives to the number of repetitions is 95%.

  Regarding the above 1), an example of the sample dilution step is as follows: 0.025%, 0.05%, 0.1%, 0.5%, 1%, and 100% when the sample stock solution is 100%. There are six stages. In the above 2), the number of repetitive experiments (N) is, for example, N = 3, 6, or 9. In the above 3), as a specific example, the ratio of the number of positive times is determined from a value such as 9 positives out of 12 repeated positives.

  Regarding 4), the function is obtained by, for example, logistic regression analysis. The model of the function is as follows:

Here, X is the logarithm of the concentration in each dilution step (base 2), p is 1-the threshold of the percentage of positives in the number of repetitions, and β ₀ and β ₁ are coefficients obtained by logistic regression analysis, respectively. It is. For example, if the threshold value is to be 95%, p = 0.05.

  Here, as an example, the dilution concentration at which the ratio of positives to the number of repetitions is 95% is obtained from the following equation (2).

  Further, a step of estimating a 95% confidence interval of the LOD by a bootstrap method or the like and confirming that the right end of the confidence interval is below a specific value may be performed.

The steps of the bootstrap method are as follows 1) to 4).
1) Restore and extract the same number (n) of samples from the original data set (N = n).
2) Logistic regression analysis is performed on the new data set obtained by the restoration extraction to estimate a threshold density.
3) 1) and 2) are repeated B (= 1,000 to 10,000) times.
4) A confidence interval is calculated from the distribution of the B estimated values obtained in 3). For example, if it is a 95% confidence interval, it is calculated as (100 ± 95) / 2, and the interval is set between the 2.5 and 97.5 percentiles.

  If the results are 100% positive in all dilution steps, correct calculations cannot be made. Therefore, in this case, an exception is set so as to be described as “<0.025%”.

<Threshold determination flow>
Next, a threshold is set. Here, the total number of reads obtained by sequencing differs between samples. For example, even in the multiplex method described above, the total number of reads varies due to variations in the sequencing reaction between samples and errors in the amount of samples. The total number of reads also differs depending on the sequence length to be sequenced. In the conventional sequence data analysis method, in order to equalize such variation in the total number of reads, in the related art, the number of reads is normalized between samples. For example, in Kukita et. Al. (Non-Patent Document 1), in order to set and determine a threshold value, a value obtained by converting the actual number of leads (Depth) into 100,000 leads is used. However, the threshold value actually changes depending on the depth. The threshold value decreases as the value of the depth increases (that is, the number of reads increases). Therefore, depending on the value of the depth, the parameter of the distribution model changes and the distribution of leads changes, and the threshold also changes. In the method of the present invention, as measures against this, (1) setting the minimum required number of leads, and (2) calculating the threshold from the actual number of leads at the time of determination without fixing the threshold are executed. .

  FIG. 6 shows a flow of a threshold value determination process according to an embodiment of the present invention. The processing unit 11 calculates a threshold from the distribution model selected based on the normal sample data. If the normal sample data is sparse data (data that is almost 0), the processing unit 11 calculates the minimum value of the threshold value calculated from the predetermined distribution model (Poisson distribution) as the threshold value. I do.

  First, the processing unit 11 inputs the normal sample data of the target mutation, the distribution model of the target mutation, and the observed total number of reads to a work memory for calculating a threshold (step S701).

  Next, the processing unit 11 refers to the minimum value information 103 and determines whether the target mutation uses the minimum value of the above-described threshold or uses another distribution model (step S702). When the minimum value of the threshold value is used, the processing unit 11 outputs, as a threshold value, a value (minimum value of the threshold value) obtained from the Poisson distribution where λ = total number of observed leads (Depth) / 100,000 (step S711). ).

  As described above, in the present embodiment, in order to increase the reliability of data, (i) when normal sample data is sparse data (data that is almost 0) at the time of determination of a distribution model, The value information 103 indicates that the minimum value of the threshold value is used. (Ii) When determining the threshold value, the following Table 1 obtained from a Poisson distribution where λ = total number of observed leads (Depth) / 100,000 Is obtained as the minimum value of the threshold value, and the minimum value of the threshold value is used as the threshold value.

  On the other hand, when using the distribution model, the processing unit 11 refers to the statistical model type 111 and selects either the Poisson distribution or the negative binomial distribution (step S703). Using the selected model, the processing unit 11 estimates the number of prediction errors with respect to the observed total number of reads, and outputs a threshold value (steps S704 to S713).

  Specifically, when the Poisson distribution model is selected, the processing unit 11 applies the total number of reads and the number of variant reads of the normal specimen data to the generalized linear model (GLM) of the Poisson distribution (step S704). Next, the processing unit 11 estimates the number of prediction errors with respect to the total number of reads observed using the fitted model, and sets it to λ (step S705). Based on this value, the processing unit 11 performs a threshold value calculation process based on the Poisson distribution (Step S706).

  When the negative binomial distribution is selected, the processing unit 11 applies the total number of reads and the number of variant reads of the normal sample data to the generalized linear model (GLM) of the negative binomial distribution (step S707). The processing unit 11 estimates the size parameter estimated using the fitted model, and sets it to θ (step S708). The processing unit 11 then estimates the number of prediction errors with respect to the total number of reads observed using the fitted model, and sets it to μ (step S709). The processing unit 11 performs a negative binomial distribution threshold calculation process based on these values (step S710). After that, the processing unit 11 outputs the threshold value (Step S713).

  The following describes the threshold setting flow in the case of using each of the Poisson distribution and the negative binomial distribution in more detail.

<Setting of threshold using Poisson distribution>
FIG. 7 shows a flow of a Poisson distribution threshold determination process according to an embodiment of the present invention.

It is assumed that the number of prediction errors is λ, and a value obtained by rounding down a fraction is int (λ) (= N) (steps S802 and S803). The observed total number of leads (D) is input (step S801). When the number N of predicted errors is smaller than the total number D of observed leads (N <D) (step S804), the probability (p) that the number of errors is N or more when the Poisson distribution of the parameter λ is assumed is calculated. It is calculated (step S805). If p is equal to or smaller than 2 ⁻⁵ (step S806), a threshold (T) is calculated from N (steps S808 and S810). If the prediction error number N is equal to or greater than the total read number D (step S804), a threshold (T) is calculated from the total read number D (steps S809 and S810). If p is greater than 2 ⁻⁵ (step S806), the probability p is calculated using N + 1, and the flow after the calculation of the probability p repeats the above steps (step S807).

<Setting of threshold using negative binomial distribution>
FIG. 8 shows a flow of the threshold value determination processing of the negative binomial distribution according to the embodiment of the present invention.

It is assumed that the number of prediction errors is μ, the size parameter is θ, and a value obtained by rounding down a fraction is int (μ) (= N) (steps S902 and S903). The observed total number of leads (D) is input (step S901). When the number N of predicted errors is smaller than the total number D of observed leads (N <D) (step S904), the probability (p) that the number of errors is N or more when the Poisson distribution of the parameter λ is assumed is determined. It is calculated (step S905). If p is ^2-5 or less (step S906), a threshold (T) is calculated from N (steps S908 and S910). If the predicted error number N is equal to or larger than the total read number D, a threshold (T) is calculated from the total read number D (step S909). If p is greater than 2 ⁻⁵ (step S906), the probability p is calculated using N + 1, and the flow after the calculation of the probability p repeats the above steps (step S907).

  In this way, the threshold is set from the created model such that the predetermined significance level (in this embodiment, the number of false positives is less than 1 out of 50,000 samples).

  In other words, the threshold is set for each different mutation, each different number of reads, and each different sequence read sample.

<Alert for mutation judgment>
In one embodiment, for example, when there is a lead excluded from the corresponding leads in steps S304, S403, and S404, the processing unit 11 may output the determination result including an alert, for example (alert process). ). The alert is a display notifying that something is in an abnormal state. The processing unit 11 determines that many reads have the same base substitution when the presence of a mutation that does not exist in a known database is detected or due to a sequencing error or the like. Is output if it has In one aspect, an alert can be displayed when the gene sequence of interest is a cancer-related gene, for example, when there is a minor mutation that is not in the cancer gene database such as the COSMIC database. Further, from another viewpoint, the alert may be displayed when there is a mutation that is a known mutation but is not set as a determination target at the time of determination.

  More specifically, the conditions for displaying an alert are such that among the reads excluded in steps S304, S403, and S404, the degree of coincidence with a plurality of reference sequences is substantially the same, and it is determined that mapping to a plurality of reference sequences is possible. If a reference (InDel template) sequence is assigned to a plurality of read or identical references, or is assigned to a certain reference (InDel template) sequence, but there is a further mutation (SNV and / or InDel) to that reference sequence, There is a case where a judgment (False Positive) exists.

<Computer and program>
The processing unit 11 of the computer 10 includes a CPU (Central Processing Unit) that executes instructions of a program (including the mutation determination program 100) that is software for realizing each function, a RAM (Random Access Memory) that expands the program, and the like. Have. The storage unit 12 of the computer 10 includes a ROM (Read Only Memory) or a storage device (these are referred to as “recording media”) in which the program and various data are recorded so as to be readable by a computer (or a CPU). Then, the object of the present invention is achieved by the computer 10 (or the processing unit 11) reading and executing the program from the recording unit 12. As the recording medium, a “temporary tangible medium” such as a tape, a disk, a card, a semiconductor memory, or a programmable logic circuit can be used. Further, the program may be supplied to the computer via an arbitrary transmission medium (a communication network, a broadcast wave, or the like) capable of transmitting the program. Note that the present invention can also be realized in the form of a data signal embedded in a carrier wave, in which the program is embodied by electronic transmission.

  The mutation detection program according to the present embodiment is not particularly limited, but may be provided in the form of, for example, a plug-in, a standalone application, middleware, a library, or the like.

  Examples of the recording medium that stores the mutation detection program according to the present embodiment include a read-only memory (ROM), a floppy (registered trademark) disk, a magnetic medium such as a hard disk, an optical disk, a magneto-optical disk, a CD-ROM, and a CD-ROM. Optical media such as R, DVDROM, Blu-ray disc, and the like are included.

  Further, the computer 10 may include one or more graphic boards for processing and outputting graphic information on the display means. The above components can be suitably interconnected by a bus within computer 10. The computer 10 further includes a communication unit 13 for communication with the external server 20, and a suitable interface for communicating with general-purpose external components such as a monitor, a keyboard, a mouse, a printer, and a network.

<Gene sequence for mutation determination>
The gene sequence to which the mutation determination method of the present invention is applied is not particularly limited, but one example is a sequence of a specific disease-related gene. Examples of the disease-related gene sequence include, but are not limited to, a cancer-related gene sequence and the like. Examples of the cancer-related gene include, but are not limited to, an epidermal growth factor receptor (EGFR) gene. EGFR encodes a transmembrane glycoprotein and is known to be involved in the growth and maintenance of cancer in living organisms. In particular, the expression of EGFR has been confirmed in non-small cell lung cancer and the like. The present invention can be applied to genetic diagnosis and monitoring by targeting such genes, for example.

<Polynucleotide sequence sample>
The polynucleotide sequence sample used to determine the mutation may be derived from a sequence present in a cell, or may be a cell-free polynucleotide sequence sample. In one embodiment, the polynucleotide sequence sample is a polynucleotide sequence sample extracted and prepared from a biological sample. The polynucleotide sequence sample of the sample of the present invention also contains, for example, circulating tumor DNA (ctDNA) in blood.

  In one embodiment of the present invention, the nucleic acid sample applied to the determination method according to the present invention is, for example, one prepared from a biological sample such as a tissue, a body fluid and cells, and preferably one prepared from a body fluid. is there. A biological sample is taken from a subject. When a polynucleotide sequence sample is obtained from a biological sample, the nucleic acid may be purified or separated from the biological sample using a conventionally known method. The biological fluid is not limited, but the biological sample includes, for example, a cell sample, a tissue sample, a body fluid sample, and the like, and among them, a body fluid sample is preferable. Examples of the body fluid sample include a blood sample, a lymph fluid sample, and a cerebrospinal fluid sample, and a blood sample is preferable, and a peripheral blood sample is particularly preferable. The peripheral blood sample can be easily collected, for example, by puncturing a fingertip or the like, so that the burden on the subject is small. Similarly, a biological sample can be taken from a biopsy, swab, smear, or the like. When a blood sample is used, it is preferable to use serum or plasma separated from the collected blood for preparing a polynucleotide sequence sample.

  Further, the biological sample may be stored by a method suitable for the type of the biological sample, such as cryopreservation, if necessary.

  The polynucleotide sequence sample may be obtained from a known nucleic acid sequence information source. The sample may be prepared from a plurality of different individuals, a normal individual, a sample collected at different stages of the disease in a particular individual, or a biological sample collected from an individual such as an individual having a pathological predisposition.

  Further, the polynucleotide sequence sample may be a PCR product amplified by PCR using any suitable primer set. The primer set used for PCR is a set of a forward primer and a reverse primer that can amplify a target DNA molecule or gene when performing PCR. The primer used for preparing the sample of the present invention is not particularly limited as long as it can specifically detect the target gene, and is, for example, an oligonucleotide consisting of 10 bp to 50 bp, preferably 18 bp to 24 bp. . The oligonucleotide base sequence is determined based on the sequence information of the target gene. Here, it may be a primer set designed at a position capable of amplifying an exon region containing a mutation in the genome sequence of the subject.

  A plasmid amplified and cloned in a system using Escherichia coli may be used as a sample. Further, the polynucleotide sequence sample according to the present invention may contain a certain amount of contamination in addition to the gene sequence to be sequenced.

<Application to diagnosis>
Using the computer program and the mutation determination method of the present invention described above, it is possible to examine, predict, diagnose or monitor a disease, disorder or symptom. The results obtained by the computer program and the mutation determination method of the present invention can be used as a source for diagnosing a disease, disorder, or symptom performed by a physician, or predicting the response of a subject to treatment thereof.

  In addition, examples of a disease, disorder, or symptom that can be tested, predicted, or diagnosed using the present invention include, but are not limited to, cancer. In addition, the present invention can be applied to the investigation of diseases, disorders or symptoms of animals other than humans, for example, mammals, invertebrates, and vertebrates other than mammals. The method of the present invention can also be applied to, for example, prenatal genetic diagnosis and the like.

  Further, the computer program and the method for determining a mutation of the present invention require detection of any specific polynucleotide sequence, for example, detection of infectious disease, monitoring of viral load, viral genotyping, environmental testing, food testing, epidemiology, and forensic medicine. In various nucleic acid detection applications.

  The present invention is not limited to the embodiments described above, and various modifications are possible within the scope shown in the claims, and embodiments obtained by appropriately combining technical means disclosed in different embodiments. Is also included in the technical scope of the present invention.

(Reference Example 1)
When the target gene was the EFGR gene, the present inventors were able to obtain very good results using a large number of reference sequences as shown in FIGS. 9 and 10. That is, “106 types” in FIGS. 9 and 10 correspond to reference sequences corresponding to almost all mutations obtained from the COSMIC database, and “8 types” correspond to the reference sequence used in Patent Document 1. I do. The number of reference sequences was increased from the “eight versions”, and after the “23 versions” using references to frequently mutated mutations obtained from the COSMIC database, finally the “106 versions” were obtained. By setting the “106 types version”, mismatch and erroneous determination that may occur when the number of references is smaller than this can be prevented, and a favorable result can be obtained.

(Reference Example 2)
When the target gene was the EFGR gene, the fitness of the Poisson distribution and the negative binomial distribution were tested for normal sample data of some mutations. There were both variants with good fit in the distribution.

  INDUSTRIAL APPLICABILITY The present invention can be used for examination, prediction or diagnosis of various gene-related diseases including cancer, or for prenatal genetic diagnosis. Further, it can be used in various nucleic acid detection applications in a field requiring detection of any specific polynucleotide sequence.

Claims (18)

A mutation determination method performed by a computer to determine the presence or absence of a mutation in the target gene,
From a plurality of reads obtained by sequencing of a polynucleotide derived from a specimen, a corresponding read corresponding to a specific region of the target gene, and a mutated read having a specific mutation in the specific region among the corresponding reads An extraction step of extracting
The threshold for determining the presence or absence of a mutation at a significance level of Jo Tokoro, a calculation step of calculating using the number of extracted said corresponding leads,
A determination step of determining that the specific mutation is present when the number of the extracted mutant reads exceeds the calculated threshold value. In the computer, statistical model information indicating a statistical model associated with the specific mutation is stored,
2. The mutation according to claim 1, wherein in the calculating step, a threshold value for determining the presence or absence of a significant difference of a predetermined significance level is calculated from the number of corresponding leads with reference to the statistical model information. 3. Judgment method.   3. The method according to claim 2, wherein the statistical model information includes a type of the statistical model associated with the specific mutation.   The mutation determination method according to claim 3, wherein the statistical model information further includes normal sample data relating to the specific mutation. In the computer, in association with the specific mutation, minimum value information indicating whether to use the minimum value of the threshold is stored,
In the calculation step, when the minimum value information indicates that the minimum value of the threshold is to be used, the minimum value of the threshold calculated from the number of corresponding leads is determined by referring to a predetermined statistical model. The method according to claim 1, wherein the threshold value is calculated. The computer stores a plurality of reference sequences,
The method according to any one of claims 1 to 5, wherein the extracting step includes mapping the plurality of reads to each of the plurality of reference sequences. The specific mutation is a single nucleotide substitution mutation,
The reference sequence includes a specific reference sequence that is a sequence of the specific region,
In the extracting step, the read mapped to the specific reference sequence is extracted as the corresponding read, and it is determined whether each of the extracted corresponding reads has the specific mutation, and the specific mutation is determined. 7. The mutation determination method according to claim 6, wherein the corresponding read having the mutation is extracted as the mutation read. The specific mutation is at least one mutation of insertion and deletion,
The reference sequence includes a specific reference sequence that is a sequence of the specific region, and a detection reference sequence that is a sequence in which at least one of insertion and deletion occurs in the specific reference sequence,
In the extracting step, the read mapped to the specific reference sequence and the read mapped to the detection reference sequence are extracted as the corresponding read, and the read mapped to the detection reference sequence is subjected to the mutation. The method according to claim 6, wherein the mutation is extracted as a lead. The specific mutation is a mutation consisting of at least one of insertion of a predetermined base and deletion of a predetermined base,
The detection reference sequence includes a plurality of detection reference sequences corresponding to different mutations from each other,
The mutation determination according to claim 8, wherein in the extraction step, among the detection reference sequences, the read mapped to the detection reference sequence corresponding to the specific mutation is extracted as the mutation read. Method.   In the extraction step, among the mutant reads mapped to the single detection reference sequence, when the ratio of mismatched reads having a sequence mismatch with the detection reference sequence exceeds a predetermined ratio 10. The method according to claim 8, wherein the mismatched read is excluded from the corresponding read and the mutant read.   The method according to any one of claims 6 to 10, wherein in the extraction step, the reads mapped to a plurality of the reference sequences are excluded from the corresponding reads and the mutation reads.   12. The mutation determination method according to claim 10, further comprising a warning step of outputting a warning when at least one of the leads is excluded from the corresponding leads in the extracting step. The target gene is an EGFR gene,
The method according to any one of claims 6 to 12, wherein the plurality of reference sequences include each of the base sequences shown in SEQ ID NOs: 1 to 101.   The mutation determination according to any one of claims 6 to 13, further comprising, before the extracting step, an updating step of updating the reference sequence based on information obtained from an external database. Method.   The method according to any one of claims 1 to 14, wherein, in the extracting step, the lead shorter than a predetermined length is removed before extracting the corresponding lead. The polynucleotide is a PCR product,
The mutation determination method according to any one of claims 1 to 15, wherein the sequencing is a sequence in which a plurality of the reads are obtained by reading the same region redundantly. A mutation determination program executed by a computer to determine the presence or absence of a mutation in the target gene,
From a plurality of reads obtained by sequencing of a polynucleotide derived from a specimen, a corresponding read corresponding to a specific region of the target gene, and a mutated read having a specific mutation in the specific region among the corresponding reads An extraction step of extracting
The threshold for determining the presence or absence of a mutation at a significance level of Jo Tokoro, a calculation step of calculating using the number of extracted said corresponding leads,
A determination step for determining that the specific mutation is present when the number of extracted mutation reads exceeds the calculated threshold value;   A computer-readable recording medium on which the mutation determination program according to claim 17 is recorded.

Images