WO2022019331A1 - ウイルス変異予測装置、ウイルス変異予測方法、およびプログラム - Google Patents

ウイルス変異予測装置、ウイルス変異予測方法、およびプログラム Download PDF

Info

Publication number
WO2022019331A1
WO2022019331A1 PCT/JP2021/027331 JP2021027331W WO2022019331A1 WO 2022019331 A1 WO2022019331 A1 WO 2022019331A1 JP 2021027331 W JP2021027331 W JP 2021027331W WO 2022019331 A1 WO2022019331 A1 WO 2022019331A1
Authority
WO
WIPO (PCT)
Prior art keywords
mutation
virus
data
amino acid
sequence
Prior art date
Application number
PCT/JP2021/027331
Other languages
English (en)
French (fr)
Inventor
康悦 小笠原
Original Assignee
国立大学法人東北大学
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by 国立大学法人東北大学 filed Critical 国立大学法人東北大学
Priority to US18/017,039 priority Critical patent/US20230298700A1/en
Priority to DE112021003912.1T priority patent/DE112021003912T5/de
Priority to JP2022538042A priority patent/JPWO2022019331A1/ja
Publication of WO2022019331A1 publication Critical patent/WO2022019331A1/ja

Links

Images

Classifications

    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B20/00—ICT specially adapted for functional genomics or proteomics, e.g. genotype-phenotype associations
    • G16B20/20—Allele or variant detection, e.g. single nucleotide polymorphism [SNP] detection
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
    • G16B30/10—Sequence alignment; Homology search
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N20/00—Machine learning
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
    • G16B30/20—Sequence assembly
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B40/00—ICT specially adapted for biostatistics; ICT specially adapted for bioinformatics-related machine learning or data mining, e.g. knowledge discovery or pattern finding
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B40/00—ICT specially adapted for biostatistics; ICT specially adapted for bioinformatics-related machine learning or data mining, e.g. knowledge discovery or pattern finding
    • G16B40/20—Supervised data analysis
    • C—CHEMISTRY; METALLURGY
    • C12—BIOCHEMISTRY; BEER; SPIRITS; WINE; VINEGAR; MICROBIOLOGY; ENZYMOLOGY; MUTATION OR GENETIC ENGINEERING
    • C12Q—MEASURING OR TESTING PROCESSES INVOLVING ENZYMES, NUCLEIC ACIDS OR MICROORGANISMS; COMPOSITIONS OR TEST PAPERS THEREFOR; PROCESSES OF PREPARING SUCH COMPOSITIONS; CONDITION-RESPONSIVE CONTROL IN MICROBIOLOGICAL OR ENZYMOLOGICAL PROCESSES
    • C12Q1/00—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions
    • C12Q1/70—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions involving virus or bacteriophage

Definitions

  • the present invention relates to a virus mutation predictor, a virus mutation prediction method, and a program.
  • the present application claims priority based on Japanese Patent Application No. 2020-12556 filed in Japan on July 22, 2020, the contents of which are incorporated herein by reference.
  • the virus is characterized by its inability to self-proliferate, and can propagate using other cells. That is, the virus utilizes various enzymes such as the host polymerase to help the growth. It is known that there are DNA virus and RNA virus in the virus.
  • the DNA virus uses the virus genomic DNA of the host RNA polymerase to synthesize messenger RNA, synthesizes a protein, and the virus proliferates. It is known that DNA viruses have fewer gene mutations than RNA viruses because they have a mechanism to correct DNA replication errors that occur during proliferation.
  • RNA virus undergoes many mutations and changes as the infection spreads, as typified by influenza. That is, RNA viruses have more gene mutations than DNA viruses.
  • coronaviruses such as the new coronavirus (SARS-CoV-2) and SARS are also RNA viruses, and mutations have been observed.
  • SARS-CoV-2 new coronavirus
  • coronavirus has an RNA calibrating enzyme in the viral genome, large-scale gene deletions, base substitutions over several bases, and mutations are unlikely to occur. Therefore, it is known that there are many point mutations in coronavirus.
  • a point mutation is a change due to a deletion, substitution, or insertion of a base.
  • RNA editing enzyme is involved in the point mutation of RNA virus.
  • the point mutation is caused by RNA editing enzymes such as ADARs and APOBECs.
  • Point mutations in RNA viruses have been presented with results suggesting the involvement of ADARs in particular.
  • the 5'side 2 bases are expressed as -2 and the 3'side 2 bases are expressed as +2 for the part of the surrounding base sequence when the mutation part by the RNA editing enzyme is set to 0.
  • -2 to +2 base sequence is characteristic (see, for example, Non-Patent Document 1).
  • HA hemagglutinin
  • RNA viruses such as the new coronavirus cause mutations.
  • the antibody and antigen tests used in the diagnosis made prior to the viral mutation become ineffective, and the therapeutic agent becomes ineffective.
  • Viral mutations have the problem that the position of the mutation on the genome and the replaced base can only be known after the mutation has occurred. In order to create an antibody test or an antigen test kit, it was necessary to identify the mutation site after the mutation occurred and then create a new protein to be used for the antibody test or the antigen test. Therefore, it takes a lot of time to make a diagnostic drug and a therapeutic drug corresponding to a new mutation.
  • the present invention has been made in view of the above problems, and provides a virus mutation prediction device, a virus mutation prediction method, and a program capable of predicting a virus mutation in advance before the mutation occurs. With the goal.
  • the present invention includes the following aspects.
  • C (cytosine) or G (guanine) is extracted from the acquisition unit for acquiring the gene sequence data of the virus genome and the acquired gene sequence data of the genome, and mutation from C or G to U (uracil) is performed.
  • C or G is changed to U, it is confirmed whether there is an amino acid mutation, and the sequence with the amino acid mutation is separated as a non-synonymous substitution.
  • a separation unit that separates non-synonymous sequences as synonymous substitutions
  • a learning unit that learns using the synonymous substitution sequence data as training data
  • a prediction unit that predicts mutations in the virus using the learned results.
  • the base sequence of the extracted context changes with the extraction unit that extracts the context in which the mutation from G to A, A to G, U to C, T to C occurs or occurs, the amino acid mutation occurs.
  • the separator that separates the sequence with the amino acid mutation as a non-synonymous substitution and the sequence without the amino acid mutation as a synonymous substitution, and the sequence data of the synonymous substitution as training data.
  • a virus mutation prediction device including a learning unit for learning and a prediction unit for predicting a mutation of the virus using the learned result.
  • the virus mutation prediction device further includes a sampling unit that selects a predetermined number from the synonymous substitutions, and the learning unit uses the sequence data of the synonymous substitutions selected by the sampling unit as learning data.
  • the virus mutation predictor is a feature amount in which two bases out of four types of RNA bases A (adenine), U, G, and C are selected and characterized, and is a feature amount during learning. Further includes a feature amount addition selection unit for adding the feature amount used in the above, and the learning unit also uses the feature amount for the learning data.
  • the range of the context is -3 to +3 or more and -10 to +10 or less.
  • the virus is SARS-CoV-2.
  • the acquisition unit acquires the gene sequence data of the virus genome, and the extraction unit extracts C (cytosine) or G (guanine) from the acquired gene sequence data of the genome, and U from C or G. Extract the context in which the mutation to (uracil) occurs or occur, check if there is an amino acid mutation when the separator changes from C or G to U, and use the sequence with the amino acid mutation as a non-synonymous substitution. Separated, the sequence without the amino acid mutation is separated as a synonymous substitution, the learning unit learns using the sequence data of the synonymous substitution as the training data, and the prediction unit uses the learned result to describe the above.
  • a virus mutation prediction method that predicts virus mutations.
  • the acquisition unit acquires the gene sequence data of the genome of the virus, and the extraction unit acquires C (cytosine), G (guanine), A (adenin), U (uracil) from the acquired gene sequence data of the genome.
  • C cytosine
  • G guanine
  • A adenin
  • U uracil
  • T thymine
  • the context in which the mutation from G to A, A to G, U to C, or T to C occurs or occurs is extracted, and the nucleotide sequence of the extracted context changes in the separator. If so, it is confirmed whether or not there is an amino acid mutation, the sequence having the amino acid mutation is separated as a non-synonymous substitution, the sequence without the amino acid mutation is separated as a synonymous substitution, and the learning unit determines the sequence of the synonymous substitution.
  • a computer acquire the gene sequence data of the virus genome, extract C (cytosine) or G (guanine) from the acquired gene sequence data of the genome, and transfer from C or G to U (uracil).
  • C cytosine
  • G guanine
  • U uracil
  • a program that separates a sequence without a gene as a synonymous substitution trains the sequence data of the synonymous substitution as training data, and predicts a mutation of the virus using the learned result.
  • C cytosine
  • G guanine
  • A adenine
  • U uracil
  • T thymine
  • viral mutations can be predicted in advance before mutations occur.
  • SARS-CoV-2 virus [Outline of SARS-CoV-2 virus] Currently, vaccines, diagnostic methods, and therapeutic methods for SARS-CoV-2 are required. Vaccines and antibody tests are based on the SARS-CoV-2 protein (or gene sequence). According to genomic analysis, SARS-CoV-2 has several variants that fall into three types: A, B, and C. As a result, it is necessary to collect variants of SARS-CoV-2 for vaccine and antibody testing.
  • SARS-CoV-2 variants contain several gene mutations, but the effect of these mutations on infection is unknown. Mutations are introduced into the virus by self-renewal errors and cell-derived RNA editing enzymes. RNA editing enzymes are known to cause mutations in RNA viruses.
  • RNA editing enzymes such as adenosine deaminase (ADAR) that act on RNA, mRNA editing enzymes for apolypoprotein B, and catalytic polypeptides (APOBECs) have been studied in RNA virus infections.
  • ADAR is an enzyme that extracts an amino group from adenosine and converts it into inosine, and is a function that mainly acts on double-stranded RNA.
  • APOBECs a family of cytidine deaminase, are enzymes that extract amino groups from cytidine and convert them to uracil.
  • APOBECs have been reported to function using ssDNA as a substrate.
  • APOBEC1, APOBEC3A and APOBEC3G also recognize ssRNA as a substrate.
  • the place where the mutation can enter in the future and the base to be replaced are predicted. do. If a viral mutation can be predicted in advance, it will be possible to prepare a diagnostic agent or a therapeutic agent corresponding to the new mutation, and the diagnostic agent or the therapeutic agent can be applied immediately after the mutation occurs.
  • FIG. 1 is a diagram showing an example of the configuration of the virus mutation prediction device 1 according to the present embodiment.
  • the virus mutation prediction device 1 includes an acquisition unit 11, a storage unit 12, an extraction unit 13, a separation unit 14, a sampling unit 15, a feature amount addition selection unit 16, a learning unit 17, a prediction unit 18, and an output unit. 19 and an operation unit 20 are provided.
  • the virus mutation prediction device 1 acquires data from the DB (database) 2 via the network NW.
  • the virus mutation prediction device 1 learns the characteristics of the gene mutation from the acquired data and predicts the mutation.
  • the acquisition unit 11 is, for example, a wireless network circuit.
  • the acquisition unit 11 acquires data from DB2 (for example, GISAID (International Promotion Organization for Bird Influenza Information Sharing; https://www.gisaid.org/)) via the network NW.
  • the data are, for example, the gene sequence of the world genome of SARS-CoV-2 and are plural.
  • the storage unit 12 stores the acquired acquired genomic data of SARS-CoV-2.
  • the storage unit 12 stores information indicating whether or not the regularization parameter C has been mutated.
  • the storage unit 12 stores the confirmation result of confirming whether or not there is an amino acid mutation when the C (cytosine) or G (guanine) is changed to U (uracil).
  • the storage unit 12 stores algorithms, programs, threshold values, and the like necessary for learning and prediction.
  • Extraction unit 13 extracts C from the acquired SARS-CoV-2 genome.
  • the extraction unit 13 also extracts from the acquired SARS-CoV-2 genome the context in which the C or G to U mutation occurs or occurs.
  • the context is a set of sequences of several bases before and after the mutation site.
  • Separation unit 14 extracts the mutant portion from C or G to U of the acquired SARS-CoV-2 genomic data, and maps the extracted mutant portion onto one genome.
  • the separation unit 14 stores information indicating whether C or G has been mutated in the storage unit 12.
  • the separation unit 14 confirms whether or not there is an amino acid mutation, and stores the confirmation result in the storage unit 12.
  • the separation unit 14 confirms whether or not there is an amino acid mutation, separates the sequence having the amino acid mutation as a non-synonymous substitution, and separates the sequence without the amino acid mutation as a synonymous substitution. ..
  • the sampling unit 15 selects the first predetermined number of amino acid substitutions without amino acid substitutions (synonymous substitutions). In order to suppress noise, the sampling unit 15 selects a second predetermined number, which is smaller than the first predetermined number, as learning data from the selected first predetermined number.
  • the sampling process does not necessarily have to be performed. In this case, all synonymous substitutions may be used for the training data. Further, the sampling unit 15 may select a first predetermined number of data having no amino acid substitution (synonymous substitution) and use this as training data.
  • the feature amount addition selection unit 16 adds a feature amount (parameter).
  • the feature amount will be described later.
  • the feature amount is an amount in which two bases are selected and characterized from the four types of RNA bases A, U, G, and C.
  • the learning unit 17 uses the selected second predetermined number as learning data and the rest of the first predetermined number as test data.
  • the learning unit 17 performs learning using the feature amount and the learning data.
  • the learning unit 17 does not have to use the feature amount for learning.
  • the learning unit 17 learns using an algorithm such as a neural network, a support vector machine, reinforcement learning, or deep learning. In addition, learning may be performed using artificial intelligence (AI: Artificial Interigence).
  • AI Artificial Interigence
  • the prediction unit 18 predicts a point mutation using the learned result.
  • the output unit 19 displays information indicating the result predicted by the prediction unit 18 on the image display device 3.
  • the image display device 3 may be, for example, a tablet terminal or the like.
  • the operation unit 20 is, for example, a touch panel sensor, a mouse, or the like provided on the image display device 3.
  • the operation unit 20 detects the operation result operated by the user.
  • FIG. 2 is a diagram showing the distribution of point mutations in the SARS-CoV-2 genome.
  • the upper figure of FIG. 2 is a diagram (g1) showing the position of each gene of the full-length ssRNA.
  • the histogram g2 at the bottom of FIG. 2 shows the number of mutations at each position.
  • the vertical axis is the number of mutations and the horizontal axis is the number of bases (bp).
  • the average number of point mutations per 150 nucleotides (bins) was about 28, but it was observed that the frequency of point mutations was high in some places.
  • FIG. 3 is a diagram showing the number of point mutations for each gene.
  • the horizontal axis is the gene name and the vertical axis is the number of mutations.
  • ORF-1a and ORF-1b had many point mutations.
  • FIG. 4 is a diagram showing the point mutation rate per 100 bases of each gene.
  • the horizontal axis is the gene name and the vertical axis is the point mutation rate per 100 bases.
  • the frequency of point mutations was highest in the 5'-untranslated region (UTR) and 3'-UTR.
  • FIG. 5 is a diagram showing the results of examining the mutated nucleobase.
  • the horizontal axis is the number of substituted bases after a point mutation, and the vertical axis is a base (A (adenine), U, G (guanine), C).
  • A adenine
  • U adenine
  • G guanine
  • C a base
  • FIG. 6 is a diagram showing the results of investigating from which base each base is mutated.
  • the horizontal axis is the original base and the number of substituted bases at the time of each point mutation, and the vertical axis is from base to base.
  • C and G particularly C
  • G is predominantly mutated to A
  • A is predominantly mutated to G
  • U is dominated by C.
  • C to U and G to A are introduced by APOBEC
  • a to G and U to C are introduced by ADAR.
  • the mutation from C to U is also written as CtoU.
  • FIG. 7 is a diagram showing a mutation pattern of each gene.
  • FIG. 8 is a diagram showing the number of mutations obtained by dividing the number of point mutations in each gene by the gene length.
  • the horizontal axis is the gene name.
  • the vertical axis of FIG. 7 is the number of mutations.
  • the vertical axis of FIG. 8 is the number of mutations per 100 bases. From FIGS. 7 and 8, the mutation of CtoU was predominant, although there were some differences for each gene.
  • CtoU and GtoA are consistent with APOBEC
  • AtoG and CtoU are consistent with the mutations introduced by ADAR. Therefore, the inventors investigated the context of one base upstream and downstream for these four mutations.
  • FIG. 9 is a diagram showing the characteristics of the base sequences on both sides of the point mutation in CtoU.
  • FIG. 10 is a diagram showing the characteristics of the base sequences on both sides of the point mutation in GtoA.
  • FIG. 11 is a diagram showing the characteristics of the base sequences on both sides of the point mutation in AtoG.
  • FIG. 12 is a diagram showing the characteristics of the base sequences on both sides of the point mutation in UtoC.
  • the horizontal axis is the base name
  • the vertical axis is the ratio [%] of each of A, U, G, and C.
  • the graph on the left shows the base (-1) on the 5'side of the mutation site
  • the graph on the right shows the base (-1) on the 3'side of the mutation site.
  • the horizontal direction is the position of the context. Further, each numerical value is the number of each base AUGC at each position. As shown in FIG. 13, A and U were very large before and after C to be replaced. The reason for this is considered to be the bias that SARS-CoV-2 contains a large amount of A and U (A is 30%, U is 32%).
  • FIG. 14 is a diagram showing an increase / decrease [%] from the expected value corresponding to each base in all C contexts of the SARS-CoV-2 sequence.
  • the horizontal direction is the position of the context.
  • U was high at positions +2 and +1 and G was high at -1 (p ⁇ 10 ⁇ -3, fisher's exact test).
  • G was high at -1 (p ⁇ 10 ⁇ -3, fisher's exact test).
  • p ⁇ 0.01, fisher's exact test there was less C (p ⁇ 0.01, fisher's exact test).
  • FIG. 15 is a diagram showing the ratio of the context of all cytosine residues in the unmasked region of the reference sequence.
  • FIG. 16 is a flowchart of a learning procedure by the virus mutation prediction device 1 according to the present embodiment.
  • Step S1 The acquisition unit 11 acquires the genomic data of SARS-CoV-2 from DB2 (for example, GISAID).
  • DB2 for example, GISAID
  • the acquisition unit 11 stores the acquired genomic data of SARS-CoV-2 in the storage unit 12.
  • Step S2 The extraction unit 13 selects C or G from the acquired SARS-CoV-2 genome.
  • the extraction unit 13 also extracts the context g11 (FIG. 17) in which the C or G to U mutation occurs or occurs from the acquired SARS-CoV-2 genome.
  • FIG. 17 is an image diagram of mapping and mutation record.
  • the context is, for example, three ways (-2 to +2, -3 to +3, -10 to +10).
  • Step S3 The separation unit 14 extracts the mutant portion from C or G to U of the acquired genomic data of SARS-CoV-2, and maps the extracted mutant portion onto one genome (FIG. 17).
  • Step S4 The separation unit 14 stores information indicating whether C or G has been mutated in the storage unit 12 (FIG. 17). For example, the separation unit 14 stores the case of mutation from C or G to U as 1, and stores C or G as 0 as a numerical value.
  • Step S5 When the separation unit 14 changes from C or G to U, it confirms whether or not there is an amino acid mutation, and stores the confirmation result in the storage unit 12.
  • step S6 When the separation unit 14 determines that there is an amino acid mutation (step S5; YES), the separation unit 14 proceeds to the process of step S6.
  • step S5; NO When the separation unit 14 determines that there is no amino acid mutation (step S5; NO), the separation unit 14 proceeds to the process of step S7.
  • Step S6 The separation unit 14 determines that the substitution is non-synonymous, and this data is also used for learning.
  • Step S7 The separation unit 14 determines that the substitution is synonymous, and uses this data for learning. Mutations were confirmed at 675 sites out of about 1800 sites that were synonymously substituted. After the processing, the separating unit 14 proceeds to the processing of step S8.
  • Step S8 The sampling unit 15 selects 1000 pieces (500 with mutation, 500 without mutation) without amino acid substitution (synonymous substitution) (first random sampling). The sampling unit 15 performs this random selection 5 times and selects 1000 amino acid substitutions without amino acid substitutions (synonymous substitutions).
  • Step S9 Generally, in machine learning, the learning data is often set to 60 to 80%, so the sampling unit 15 selects 800 out of 1000 selected pieces as training data (second random). sampling). The sampling unit 15 makes random selections five times to select 800 pieces. The sampling unit 15 does not have to perform this process.
  • Step S10 The learning unit 17 uses the selected 800 pieces as learning data and the remaining 200 pieces as test data.
  • the learning unit 17 uses the learning data even if there is no mutation.
  • Step S11 Feature amount addition
  • the number of feature quantities is an example and is not limited to this.
  • the feature amount is a combination of two bases selected in the context as shown in FIG.
  • Step S12 The learning unit 17 performs learning using the feature amount and the learning data.
  • Step S13 The prediction unit 18 predicts a point mutation using the learned result. The forecast will be described later.
  • the context may be -3 to +3 or more and -10 to +10 or less.
  • -3 to +3 or more and -10 to +10 or less include -4 to +4, ..., -9 to +9.
  • FIG. 18 is a diagram showing an example of a combination of two positions when a synonymous substitution (no amino acid mutation) is used. For example, in “1_G 4_G" on the first line, 1_G indicates G at position +1 and 4_G indicates G at position +4. Further, “-2_T 1_G” on the second line indicates the context of TNCG.
  • FIG. 19 is a diagram showing an example of the selected top 30 feature quantities.
  • the hatching g21 represents an increase, and the hatching g22 represents a decrease.
  • the feature amount selected is not limited to 30.
  • FIG. 20 is a diagram showing an example of the relationship between the context and the score when no feature amount is added or selected.
  • FIG. 21 is a diagram showing an example of the relationship between the context and the score when the feature amount is added / selected.
  • the horizontal axis is the context ⁇ (-2, +2), (-3, +3), (-10, +10) ⁇
  • the vertical axis is the score.
  • the points are the regularization parameter C values in logistic regression, 0.0001, 0.001, 0.01, 0.1, 1.0, 10.0, 100.0, in order from the left. It is 1000.0.
  • the regularization parameter indicates that the larger the value, the easier it is to learn.
  • the score indicates the correct answer rate, and the variation indicates the robustness against the bias of the data.
  • the prediction unit 18 calculates and predicts by multiplying by a coefficient according to the rank of the top 30 by adding the feature amount (top 30).
  • the features (among the top 30), there are really important ones and noise.
  • FIG. 22 is a diagram showing the context and the average value of the scores for each regularization parameter when the feature amount is added / selected.
  • FIG. 23 is a diagram showing the context and the standard deviation of the score for each regularization parameter when the feature amount is added / selected. In the comparison between contexts -2 to +2 and 3 to +3 as shown in FIGS. 22 and 23, the scores from -3 to +3 are higher and the variation (standard deviation) is smaller. A high score indicates a high percentage of correct answers, and a small variation indicates that the obtained results are highly valid, and is considered to be practical.
  • the score was higher and the variation was smaller in -10 to +10 than in Context-3 to +3. Therefore, the context is better from -3 to +3 than from -2 to +2, and better from -10 to +10 than from -3 to +3. That is, the context of -10 to +10 was the best.
  • FIG. 24 is a flowchart of the mutation prediction processing procedure according to the embodiment.
  • the above-mentioned learning is performed in advance at the time of prediction.
  • Step S101 The prediction unit 18 calculates the score of the predicted result, and causes the image display device 3 to display the calculated score via the output unit 19. As a result, the image display device 3 displays a graph of the relationship between the context and the score, as shown in FIG. 25, for example.
  • FIG. 25 is a diagram showing an example of information displayed on the image display device 3 at the time of mutation prediction.
  • the operation unit 20 outputs the selection information selected by the user to the prediction unit 18.
  • Step S103 The prediction unit 18 performs statistical processing as shown in FIG. 26 with a predetermined algorithm (for example, logistic regression) for the regular parameters of the selected context.
  • FIG. 26 is a diagram showing an example of the result of logistic regression calculation.
  • the vertical axis of FIG. 26 is the score, and the straight line g42 is the threshold value for the presence or absence of mutation.
  • the prediction unit 18 displays a graph as shown in FIG. 26 on the image display device 3.
  • Step S104 The user looks at the displayed image (FIG. 26) and selects a point having a mutation, for example, point g43.
  • the operation unit 20 outputs the selection information selected by the user to the prediction unit 18.
  • Step S105 The prediction unit 18 maps the selected point to the position g44 on one SARS-CoV-2 genome as shown in FIG. 27, and displays the mapped image on the image display device 3.
  • FIG. 27 is a diagram showing mutation records and mutation predictions.
  • Step S106 When the prediction unit 18 detects that the extraction portion is selected by operating the operation unit 20 in the displayed image (FIG. 27), the prediction unit 18 displays where in FIG. 26 it corresponds (backcast). function). In addition, the prediction unit 18 may display all of FIGS. 25 to 27 in one screen, or at least one may be displayed and switched to be displayed. In this embodiment, for example, FIG. 26 and FIG. 27 can be bidirectionally selected and mapped in this way.
  • the processing procedure shown in FIG. 24 is an example and is not limited to this.
  • the virus genome is searched based on the characteristic sequences of several bases before and after the mutation of the virus gene, and the mutations (C or G to U) so far are used as teacher data for machine learning to predict the mutation. I did it.
  • this correct answer rate is the correct answer rate that includes not only mutations due to RNA editing enzymes but also mutations, and by distinguishing between mutations and mutations that are only RNA editing enzymes, the correct answer rate for prediction of mutations by RNA editing enzymes is It is easy to imagine that the rate is higher.
  • the AUC (Area Under the Curve) score was used as the correct answer rate. The calculation of the AUC score and the like will be described later.
  • the present embodiment it is an invention that enables the development of an ultra-early diagnostic kit. Further, according to the present embodiment, not only the diagnostic kit but also the effect of the vaccine, the effect of the viral antibody drug, and the authentication or cancellation of the immune passport can be determined. In addition, according to the present embodiment, it is possible to select a candidate for a therapeutic agent, so that ultra-early treatment is possible.
  • variant-1 (5'-AUUUAUGUUCUUUCCC-3'; at2946-2965 region in EPI_ISL_419308)
  • variant-2 (5'-AUUUAUUGUUCUUUUCUCUUUCCC-3';'; EPI_ISL_418420 regions 14392-14411)
  • variant-4 (5'-AAACCUUUGAGAGAGUU-3'; EPI_ISL_419846 regions 22946-22965).
  • a U-free sequence (5'-GACAGAGAGAGAACAAG-3') was used as a negative control to induce TLR7-mediated cytokine production.
  • ssRNA synthesized by Japan Genetic Research Institute Co., Ltd. (Sendai City, Miyagi Prefecture) was used.
  • the human monocytic leukemia cell line THP-1 is an RPMI-1640 medium supplemented with 10% FCS, 55 mM 2-mercaptoethanol, 100 mM non-essential amino acids (NEAAs), 1 mM pyruvate and 20 mM ml-1 penicillin and streptomycin. Maintained in.
  • FIG. 28 is a diagram showing a phylogenetic tree.
  • the inventors collected the gene sequence from GISAID based on the Wuhan type (W) reported in the early stage, and created the phylogenetic tree of FIG. 28.
  • FIG. 28 is a diagram showing a phylogenetic tree.
  • the four sequences consist of a first variant (variant-1, Japanese type), a second variant (variant-2, Georgia type), a third variant (variant-3, French type), and a fourth variant (variant-4,). It is derived from (Australian type).
  • W shows the original SARS-CoV-2 sequence reported in Wuhan.
  • FIG. 29 is a diagram showing the mutation sites of various mutants in the genome of the four selected mutants and the positions of the RNA sequences used in the pseudo-infection model.
  • the lateral direction is (bp)
  • the downward triangle is VtoU (V is all bases except U)
  • the upward triangle is UtoV
  • the squares indicate the sequence of ssRNA used for cell stimulation.
  • the number of U in the full-length ssRNA of each SARS-CoV-2 mutant is significantly increased as compared with the original isolated strain.
  • the frequency of point mutations for U was much higher than the frequency of U for A, G, or C.
  • the ability of full-length mutated ssRNA to induce inflammatory cytokines is much greater than in the original isolate.
  • FIG. 30 is a diagram showing the induction of TNF- ⁇ production by ssRNA.
  • ssRNA For the measurement of human TNF- ⁇ , cells were cultured in the presence of PMA (0.2 ng / ml, Sigma Aldrich, St. Louis, MO, USA) and DOTAP (10 ⁇ g, Roche Diagnostics, Mannheim, Germany). It was stimulated with 160 (pmol) ssRNA.
  • PMA 0.2 ng / ml
  • DOTAP DOTAP
  • W-1 indicates an early Wuhan type
  • variant-1 indicates a mutant type
  • FIG. 31 is a diagram showing the induction of IL-6 production by ssRNA.
  • ssRNA For the measurement of human IL-6, cells were cultured in the presence of PMA (50 ng / ml) and stimulated with DOTAP (15 ⁇ g) with 480 (pmol) ssRNA. IL-6 production was measured after 48 hours of stimulation.
  • the ssRNA sequence lacking the U residue did not upregulate the production of TNF- ⁇ , as shown in FIG.
  • the increase in U number induced by point mutation increased cytokine production of variant-1, 3 and 4 compared to stimulation with Wuhan-type reference ssRNA sequences.
  • the production of IL-6 was lower than that of TNF- ⁇ , but the same tendency was observed in the production of IL-6.
  • the acquisition unit acquires the gene sequence data of the virus genome.
  • the extraction unit extracts C (cytosine) or G (guanine) from the acquired genomic sequence data, and extracts the context in which the mutation from C or G to U (uracil) occurs or occurs.
  • C or G to U uracil
  • the base sequence of the extracted context changes from C or G to U, it is confirmed whether or not there is an amino acid mutation. Mutations caused by RNA editing enzymes are thought to occur with or without amino acid mutations because they act directly on genomic RNA to induce mutations.
  • the separation unit separates the sequence having the amino acid mutation as a non-synonymous substitution, and separates the sequence without the amino acid mutation as a synonymous substitution. Then, the learning unit learns by using the sequence data of the synonymous substitution as the learning data, and the prediction unit predicts the virus mutation by using the learned result.
  • FIG. 32 is a diagram showing an example of processing contents and an example of processing procedure of the analysis program according to the present embodiment.
  • the vertical direction is the main processing and the horizontal direction is the processing procedure.
  • the analysis program reads the file to be analyzed (step S211), sets the explanatory variables and objective functions (step S212), and defines the function for creating the feature quantity (step S212).
  • step S213) the range of the base sequence and the parameters for grid search are set (step S214).
  • the objective variable is the presence or absence of mutation, and the explanatory variables are the dummy base sequence and the base ratio.
  • the function for creating the feature amount takes a range of the base sequence (example: -3 to +3) as an argument, and the base ratio ("A", "G", "C", "T” for one record is the whole. % Is included) is a function to calculate.
  • the analysis program creates features (step S221) in the learning process (step S220), optimizes parameters by grid search (step S222), and executes cross-validation and learning of various models (step S223). ), Calculate the AUC scores of various models (step S224).
  • the ratio of the base is calculated using the function for creating the feature amount, and the base sequence is made into a dummy variable by using the function for dummyizing the variable specified in the argument.
  • the ACU score is the area of the part below the curve of the graph when the ROC (Receiver Operating Characteristic Curve) curve is created. For example, it takes a value from 0 to 1, and the closer the value is to 1, the better the discriminant ability. Indicates high.
  • the analysis program outputs the AUC scores of various models in the accuracy evaluation process (step S230) (step S231), and calculates the summary statistic of the AUC score (step S232).
  • the analysis program represents the coefficients of the regression equation in a histogram and plots them in a box plot (step S241), and plots the ROC curves of various models (step S242).
  • FIG. 33 is a diagram showing an example of hyperparameter values of each model optimized by performing a grid search for each range of the base sequence.
  • the model is LightGBM, a technique that combines logistic regression with decision trees and gradient boosting.
  • the analysis condition of FIG. 33 is that the number of cross-validations is five.
  • the range of the base sequence is -2 to +2, -3 to +3, -5 to +5, and -10 to +10.
  • the hyperparameters used for the verification of logistic regression are C: [0.0001, 0.001, 0.01, 0.1, 1,10,100,1000], and the hyperparameters used for the verification of LightGBM are.
  • Num_leaves [10,31,64], learning_late: [0.01,0.1,1].
  • the strength of regularization in logistic regression tends to increase as the range of the base sequence expands.
  • the hyperparameter "learning_rate" of LightGBM was constant at 0.01.
  • FIGS. 34 to 37 show the analysis results of A, C, G, T, A_percent, G_percent, C_percent, and T_percent.
  • A_percent, G_percent, C_percent, and T_percent are the ratio of bases per record.
  • 0 to 4 are coefficients for each crossing number, and for example, a bar graph of 0 is a coefficient at the time of the first cross-validation.
  • the horizontal axis is a parameter and the vertical axis is a rate.
  • the vertical axis is the rate.
  • FIGS. 34 to 36 are diagrams showing the coefficients of the regression equation in a histogram in the range of the base sequence ⁇ 10 to +10.
  • FIG. 37 is a diagram in which a histogram of the coefficients of the regression equation is plotted in a box plot in the range of the base sequence ⁇ 10 to +10.
  • the values of -6T, -2G, -2T, -1G, -1T, + 5A and the like were large.
  • the values of -2T and + 1G were large.
  • the values of -2T, -1G, and + 1G were large.
  • the values of -2T, -1G, and + 1G were large.
  • the values of -2T, -1G, -1T, + 1G and the like were large. It should be noted that such a correlation coefficient was used to visualize the weight for each base, which will be described later.
  • FIG. 38 is a diagram showing an outline and features of the compared learning models. As shown in FIG. 38, the models are logistic regression, SVM (Support Vector Machine), decision tree, random forest, XGBoost, and LightGBM.
  • FIG. 39 is a diagram showing an example of the result of analyzing the summary statistic of the AUC score for each model.
  • the correlation coefficients of the logistic regressions of the ranges -2 to +2, -3 to +3, -5 to +5, and -10 to 10 of the base sequence are used, and the third decimal place and the following are rounded down.
  • Image g101 of FIG. 39 shows ROC (ROC_xgt) of XDBost, ROC (ROC_tree) of a decision tree, and ROC (ROC_lgb) of LightGBM.
  • Image g102 shows ROC (ROC_svm) of SVM, ROC (ROC_tf) of random forest, and ROC (ROC_llr) of logistic regression.
  • ROC_xgt of XDBost
  • ROC_tree ROC
  • ROC_lgb ROC
  • Image g102 shows ROC (ROC_svm) of SVM, ROC (ROC_tf) of random forest, and
  • the scores are 55.4% for the base sequence range -10 to +10, 56.0% for the base sequence range -2 to +2, and 56.6% for the base sequence range -3 to +3.
  • the range of the base sequence -5 to +5 was 56.2%.
  • the overall logistic regression score was high. In other models, scores of about 52 to 57% were obtained.
  • FIG. 40 is a diagram showing an example of an AUC score before processing.
  • FIG. 41 is a diagram showing an example of the AUC score after processing.
  • the correlation coefficients of the logistic regressions of the ranges -2 to +2, -3 to +3, -5 to +5, and -10 to +10 of the base sequence are used to indicate the third decimal place and below. Truncate to show summary statistics.
  • the following variables that do not exist in the data before processing were deleted, and then the AUC score after processing was calculated.
  • the deleted variables are A_percent, G_percent, C_percent, and T_percent.
  • the AUC score before processing when logistic regression is used is about 51 to 54%, but the AUC score after processing is improved to about 56 to 57%.
  • FIGS. 42 to 43 show an example of the comparison result in the range of the base sequence -2 to +2.
  • FIG. 42 is a diagram showing the ROC curve of each model having the range of the base sequence -2 to +2; the number of cross-validations of the first time.
  • FIG. 43 is a diagram showing the ROC curve of each model having the range of base sequence -2 to +2; the number of cross-validations is the second.
  • the algorithms used are Logistic Regression, SVM, Decision Tree, Random Forest, XGBoost, and Light GBM shown in FIG. 38.
  • the line g201 is an XGBoost
  • the line g202 is a decision tree
  • the line g203 is a Light GBM
  • the line g205 is an SVM
  • the line g205 is a random forest
  • the line g206 is a logistic regression ROC curve.
  • the program that realizes the function of the virus mutation prediction device 1 so as to perform the above-mentioned analysis and the like has the following functions.
  • I. The first function that reads the file to be analyzed and deletes the "1" record that is not used in the analysis.
  • II. Execute the second function for calculating the base ratio, calculate the base ratio of the data read in I, and store it in a new variable.
  • III. Of the data read in I, the variables of the base sequence (for example, columns C to V of the file) are converted into dummy variables by using the third function.
  • a grid search is executed using the fourth function to optimize the parameters of various models (Fig. 33).
  • V. Perform 5-fold cross-validation using the fifth function.
  • the variables II and III are set as explanatory variables, and the presence or absence of mutation (for example, column B of the file) among the data read by I is set as the objective variable in the first method, and learning of each model is executed.
  • the first method performs machine learning by designating the test data to be classified as the first argument and the correct answer as a result of classifying into the second argument.
  • the AUC score of each model is calculated using the sixth function.
  • the summary statistic of the AUC score of each model is calculated by the second method for extracting statistical information (for example, FIGS. 38 to 43).
  • the third method is used to plot the coefficients of logistic regression (eg, FIGS. 34-36).
  • the third method is a method that outputs the confidence interval as an error bar with the average value of the given vector (array composed of numerical values) as the height.
  • X Using the third method, the coefficients are plotted in a box plot (eg, FIG. 37).
  • XI The ROC curve of each model is plotted using the plotting fourth method (eg, FIGS. 42-43).
  • the above-mentioned functions, functions, and methods of I to XI are examples, and the present invention is not limited to these.
  • FIG. 44 is a diagram showing an example of a method of dividing learning data by cross-validation five times. How to divide training data and test data is a very important issue. Therefore, in the present embodiment, the training data and the test data are divided as shown in FIG. 44, and the training data and the test data are exchanged for each intersection for learning.
  • FIG. 45 is a diagram for explaining a method of measuring generalization performance.
  • Stratified KFold was performed as shown in FIG. 45 as a method for measuring generalization performance.
  • the data is divided into training and testing while maintaining the distribution ratio.
  • the examples shown in FIGS. 44 and 45 are examples, and are not limited to these.
  • FIG. 46 is a boxplot based on the range of each base sequence and each learning model when mutating from G to U.
  • FIG. 47 is a boxplot based on the range of each base sequence and each learning model when mutating from G to A.
  • FIG. 48 is a boxplot based on the range of each base sequence and each learning model when mutating from A to G.
  • FIG. 49 is a boxplot based on the range of each base sequence and each learning model when mutating from U to C (or T (thymine) to C).
  • TtoC is expressed in DNA notation, but U to C in RNA notation.
  • xgb indicates XGBoost
  • Tree indicates a decision tree
  • Lab indicates LightGBM
  • Svm indicates SVM
  • rf indicates random forest
  • Lr indicates logistic regression.
  • the average value of the correct answer rate in the range of the base sequence -10 to +10 is 56.4% for XGBoost, 53.0% for the decision tree, and 50.0% for LightGBM. SVM was 51.4%, random forest was 54.0%, and logistic regression was 54.0%. As shown in FIG. 46, in the case of the mutation from G to U, the result of the combination of the base sequence range -10 to +10 and the model XGBoost was the best.
  • the average value of the correct answer rate in the range of -5 to +5 of the base sequence is 62.2% for XGBoost, 57.0% for decision tree, and 62.8 for LightGBM. %, SVM was 52.6%, Random Forest was 64.2%, and Logistic Regression was 60.2%.
  • the average value of the correct answer rate in the range of base sequence -10 to +10 is 60.6% for XGBoost, 56.6% for decision tree, 61.6% for Light GBM, 54.4% for SVM, and random forest. Was 64.2% and logistic regression was 59.8%.
  • FIG. 47 in the case of the mutation from G to A, the result of the combination of the base sequence range -10 to +10 or -5 to +5 and the model random forest was the best.
  • the average value of the correct answer rate in the range of base sequence -2 to +2 is 58.0% for XGBoost, 56.4% for decision tree, 60.2% for Light GBM, and so on. SVM was 48.8%, random forest was 57.2%, and logistic regression was 58.2%.
  • the result of the combination of the base sequence range -2 to +2 and the model Light GBM was the best.
  • the average value of the correct answer rate in the range of -5 to +5 of the base sequence is 61.0% for XGBoost, 62.4% for decision tree, and 64 for LightGBM. 0.0%, SVM was 55.0%, Random Forest was 62.4%, and Logistic Regression was 62.6%.
  • SVM was 55.0%
  • Random Forest was 62.4%
  • Logistic Regression was 62.6%.
  • the result of the combination of the base sequence range -5 to +5 and Light GBM was the best.
  • XGBoost decision tree
  • LightGBM decision tree
  • SVM random forest
  • logistic regression logistic regression
  • a program for realizing all or part of the functions of the virus mutation predictor 1 in the present invention is recorded on a computer-readable recording medium, and the program recorded on the recording medium is read into the computer system. By executing this, all or some of the treatments performed by the virus mutation prediction device 1 may be performed. Further, for machine learning, various learning methods such as a deep learning method may be used, or processing may be performed using artificial intelligence (AI: Artificial Interigence).
  • AI Artificial Interigence
  • the term "computer system” as used herein includes hardware such as an OS and peripheral devices. Further, the "computer system” shall also include a WWW system provided with a homepage providing environment (or display environment).
  • the "computer-readable recording medium” refers to a portable medium such as a flexible disk, a magneto-optical disk, a ROM, or a CD-ROM, and a storage device such as a hard disk built in a computer system.
  • a "computer-readable recording medium” is a volatile memory (RAM) inside a computer system that serves as a server or client when a program is transmitted via a network such as the Internet or a communication line such as a telephone line. In addition, it shall include those that hold the program for a certain period of time.
  • the above program may be transmitted from a computer system in which this program is stored in a storage device or the like to another computer system via a transmission medium or by a transmission wave in the transmission medium.
  • the "transmission medium” for transmitting a program refers to a medium having a function of transmitting information, such as a network (communication network) such as the Internet or a communication line (communication line) such as a telephone line.
  • the above program may be for realizing a part of the above-mentioned functions.
  • a so-called difference file (difference program) may be used, which can realize the above-mentioned function in combination with a program already recorded in the computer system.
  • Virus mutation prediction device 1 ... Virus mutation prediction device, 2 ... DB, 3 ... Image display device, 11 ... Acquisition unit, 12 ... Storage unit, 13 ... Extraction unit, 14 ... Separation unit, 15 ... Sampling unit, 16 ... Feature amount addition selection unit, 17 ... learning unit, 18 ... prediction unit, 19 ... output unit, 20 ... operation unit, A ... adenine, U ... uracil, G ... guanine, C ... cytosine, T ... thymine

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Medical Informatics (AREA)
  • Health & Medical Sciences (AREA)
  • Theoretical Computer Science (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Biophysics (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Biotechnology (AREA)
  • Evolutionary Biology (AREA)
  • General Health & Medical Sciences (AREA)
  • Software Systems (AREA)
  • Data Mining & Analysis (AREA)
  • Evolutionary Computation (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Artificial Intelligence (AREA)
  • Chemical & Material Sciences (AREA)
  • Analytical Chemistry (AREA)
  • Proteomics, Peptides & Aminoacids (AREA)
  • Bioethics (AREA)
  • Databases & Information Systems (AREA)
  • Epidemiology (AREA)
  • Public Health (AREA)
  • Computing Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Genetics & Genomics (AREA)
  • Molecular Biology (AREA)
  • Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)

Abstract

ウイルス変異予測装置は、ウイルスのゲノムの遺伝子配列データを取得する取得部と、取得したゲノムの遺伝子配列データからＣ（シトシン）またはＧ（グアニン）を抽出し、ＣまたはＧからＵ（ウラシル）への変異が起こるまたは起こったコンテキストを抽出する抽出部と、ＣまたはＧからＵに変わった時、アミノ酸変異があるかを確認し、アミノ酸変異がある配列を非同義置換として分離し、アミノ酸変異がない配列を同義置換であるとして分離する分離部と、同義置換の配列データを学習データに用いて学習する学習部と、学習された結果を用いて、ウイルスの変異を予測する予測部と、を備える。

Description

ウイルス変異予測装置、ウイルス変異予測方法、およびプログラム
　本発明は、ウイルス変異予測装置、ウイルス変異予測方法、およびプログラムに関する。
　本願は、２０２０年　７月２２日に、日本に出願された特願２０２０－１２５５６３号に基づき優先権を主張し、その内容をここに援用する。
　ウイルスは、自己増殖できないことが特徴であり、他細胞を利用して増殖することができる。すなわち、ウイルスは、宿主のポリメラーゼなどの種々の酵素を利用して増殖に役立てている。ウイルスには、ＤＮＡウイルスとＲＮＡウイルスが存在することが知られている。ＤＮＡウイルスは、ウイルスゲノムＤＮＡを、宿主のＲＮＡポリメラーゼを利用し、メッセンジャーＲＮＡを合成して、タンパク質を合成してウイルスは増殖する。ＤＮＡウイルスには、増殖の過程で生じたＤＮＡ複製のミスを修正する機構が備わっているので、ＲＮＡウイルスと比較すると遺伝子の変異が少ないことが知られている。
　ＲＮＡウイルスは、インフルエンザに代表されるように感染が伝播するにしたがって多くの変異が入りウイルスが変化していくことが知られている。つまり、ＲＮＡウイルスは、ＤＮＡウイルスに比べ遺伝子変異が多い。例えば、新型コロナウイルス（ＳＡＲＳ－ＣｏＶ－２）やＳＡＲＳ等のコロナウイルスもＲＮＡウイルスであり、変異が観察されている。しかし、コロナウイルスは、ＲＮＡ校正酵素をウイルスゲノム内に有しているため、大規模な遺伝子の欠失や数塩基にわたる塩基置換、変異は起こりにくい。そのため、コロナウイルスでは、点変異が多いことが知られている。なお、点変異とは、塩基の欠失、置換、挿入による変化である。
　ＲＮＡウイルスの点変異には、宿主のＲＮＡ編集酵素が関与していることが知られている。新型コロナウイルスの変異においては、ＲＮＡ編集酵素、ＡＤＡＲｓやＡＰＯＢＥＣｓなどにより点変異が引き起こされていることを示している。ＲＮＡウイルスの点変異では、特にＡＤＡＲｓの関与を示唆する結果が提示されている。加えて、ＲＮＡウイルスの点変異では、ＲＮＡ編集酵素による変異部分を０としたときの周囲の塩基配列の部分について５’側２塩基を－２と表記、３’側２塩基を＋２と表記すると、－２から＋２の塩基配列に特徴があることを示している（例えば、非特許文献１参照）。
　現在、ウイルスの変異予測については、インフルエンザウイルスについて予測が始まっており、ヘマグルチニン（ＨＡ）の構造を指標に変異予測がなされている。しかし、新型コロナウイルスなどＲＮＡ校正酵素をもつウイルスの変異の予測は、まだ行われていない。
Di Giorgio, S.,et al. Evidence for host-dependent RNA editing in the transcriptome of SARS-CoV-2. Science Advances: eabb5813, 2020.
　新型コロナウイルスなどのＲＮＡウイルスは変異を起こす。ウイルス変異が起こった場合は、ウイルス変異前に作成された診断に用いられる抗体検査や抗原検査が無効になる、および治療薬が無効になる。ウイルス変異は、ゲノム上の変異の位置や置換された塩基について、変異が起こってからしかわからないという問題がある。抗体検査や抗原検査キットを作成するためには、変異が起こったのち変異部位を特定してから、新たに抗体検査や抗原検査に利用するタンパク質を作成する必要があった。そのため、新たな変異に対応する診断薬や治療薬をつくるためには、多くの時間を要している。
　本発明は、上記の問題点に鑑みてなされたものであって、ウイルス変異を、変異が起こる前に事前に予測することができるウイルス変異予測装置、ウイルス変異予測方法、およびプログラムを提供することを目的とする。
　本発明は以下の態様を含む。
　［１］ウイルスのゲノムの遺伝子配列データを取得する取得部と、取得した前記ゲノムの遺伝子配列データからＣ（シトシン）またはＧ（グアニン）を抽出し、ＣまたはＧからＵ（ウラシル）への変異が起こるまたは起こったコンテキストを抽出する抽出部と、ＣまたはＧからＵに変わった時、アミノ酸変異があるかを確認し、前記アミノ酸変異がある配列を非同義置換として分離し、前記アミノ酸変異がない配列を同義置換であるとして分離する分離部と、前記同義置換の配列データを学習データに用いて学習する学習部と、学習された結果を用いて、前記ウイルスの変異を予測する予測部と、を備えるウイルス変異予測装置。
　［２］ウイルスのゲノムの遺伝子配列データを取得する取得部と、取得した前記ゲノムの遺伝子配列データからＣ（シトシン）、Ｇ（グアニン）、Ａ（アデニン）、Ｕ（ウラシル）またはＴ（チミン）を抽出し、ＧからＡ、ＡからＧ、ＵからＣ、ＴからＣへの変異が起こるまたは起こったコンテキストを抽出する抽出部と、抽出されたコンテキストの塩基配列が変化した場合、アミノ酸変異があるかを確認し、前記アミノ酸変異がある配列を非同義置換として分離し、前記アミノ酸変異がない配列を同義置換であるとして分離する分離部と、前記同義置換の配列データを学習データに用いて学習する学習部と、学習された結果を用いて、前記ウイルスの変異を予測する予測部と、を備えるウイルス変異予測装置。
　［３］前記ウイルス変異予測装置は、前記同義置換から所定数を選ぶサンプリング部、を更に備え、前記学習部は、前記サンプリング部によって選ばれた前記同義置換の配列データを学習データに用いる。
　［４］前記ウイルス変異予測装置は、ＲＮＡ塩基Ａ（アデニン）、Ｕ、Ｇ、Ｃの４種類のうち、２塩基が選ばれて特徴づけられた量である特徴量であって、学習の際に用いられる前記特徴量を追加する特徴量追加選択部、を更に備え、前記学習部は、前記特徴量も学習データに用いる。
　［５］前記ウイルス変異予測装置において、前記コンテキストの範囲は、－３から＋３以上、－１０から＋１０以下である。
　［６］前記ウイルス変異予測装置において、前記ウイルスは、ＳＡＲＳ－ＣｏＶ－２である。
　［７］取得部が、ウイルスのゲノムの遺伝子配列データを取得し、抽出部が、取得された前記ゲノムの遺伝子配列データからＣ（シトシン）またはＧ（グアニン）を抽出し、ＣまたはＧからＵ（ウラシル）への変異が起こるまたは起こったコンテキストを抽出し、分離部が、ＣまたはＧからＵに変わった時、アミノ酸変異があるかを確認し、前記アミノ酸変異がある配列を非同義置換として分離し、前記アミノ酸変異がない配列を同義置換であるとして分離し、学習部が、前記同義置換の配列データを学習データに用いて学習し、予測部が、学習された結果を用いて、前記ウイルスの変異を予測する、ウイルス変異予測方法。
　［８］取得部が、ウイルスのゲノムの遺伝子配列データを取得し、抽出部が、取得した前記ゲノムの遺伝子配列データからＣ（シトシン）、Ｇ（グアニン）、Ａ（アデニン）、Ｕ（ウラシル）またはＴ（チミン）を抽出し、ＧからＡ、ＡからＧ、ＵからＣ、ＴからＣへの変異が起こるまたは起こったコンテキストを抽出し、分離部が、抽出されたコンテキストの塩基配列が変化した場合、アミノ酸変異があるかを確認し、前記アミノ酸変異がある配列を非同義置換として分離し、前記アミノ酸変異がない配列を同義置換であるとして分離し、学習部が、前記同義置換の配列データを学習データに用いて学習し、予測部が、学習された結果を用いて、前記ウイルスの変異を予測する、ウイルス変異予測方法。
　［９］コンピュータに、ウイルスのゲノムの遺伝子配列データを取得させ、取得された前記ゲノムの遺伝子配列データからＣ（シトシン）またはＧ（グアニン）を抽出させ、ＣまたはＧからＵ（ウラシル）への変異が起こるまたは起こったコンテキストを抽出し、分離部が、ＣまたはＧからＵに変わった時、アミノ酸変異があるかを確認させ、前記アミノ酸変異がある配列を非同義置換として分離させ、アミノ酸変異がない配列を同義置換であるとして分離させ、前記同義置換の配列データを学習データに用いて学習させ、学習された結果を用いて、前記ウイルスの変異を予測させる、プログラム。
　［１０］コンピュータに、ウイルスのゲノムの遺伝子配列データを取得させ、
　取得された前記ゲノムの遺伝子配列データからＣ（シトシン）、Ｇ（グアニン）、Ａ（アデニン）、Ｕ（ウラシル）またはＴ（チミン）を抽出し、ＧからＡ、ＡからＧ、ＵからＣ、ＴからＣへの変異が起こるまたは起こったコンテキストを抽出させ、抽出されたコンテキストの塩基配列が変化した場合、アミノ酸変異があるかを確認し、前記アミノ酸変異がある配列を非同義置換として分離させ、前記アミノ酸変異がない配列を同義置換であるとして分離させ、前記同義置換の配列データを学習データに用いて学習させ、学習された結果を用いて、前記ウイルスの変異を予測させる、プログラム。
　本発明によれば、ウイルス変異を、変異が起こる前に事前に予測することができる。
実施形態に係るウイルス変異予測装置の構成の一例を示す図である。 ＳＡＲＳ－ＣｏＶ－２ゲノムにおける点突然変異の分布を示す図である。 遺伝子ごとの点突然変異の数を示す図である。 各遺伝子の１００塩基あたりの点突然変異率を示す図である。 変異した核酸塩基を調べた結果を示す図である。 それぞれの塩基がどの塩基から変異しているか調べた結果を示す図である。 各遺伝子の変異パターンを示す図である。 各遺伝子における点突然変異の数を遺伝子長で割った変異数を示す図である。 ＣｔｏＵにおける点変異の両側の塩基配列の特徴を示す図である。 ＧｔｏＡにおける点変異の両側の塩基配列の特徴を示す図である。 ＡｔｏＧにおける点変異の両側の塩基配列の特徴を示す図である。 ＵｔｏＣにおける点変異の両側の塩基配列の特徴を示す図である。 ＣからＵへの変異（ｎ＝２４０１）の上流下流３塩基ずつのコンテキストの特徴を示す図である。 ＳＡＲＳ－ＣｏＶ－２のシークエンスの全てのＣのコンテキストにおいて、それぞれの塩基に該当する期待値からの増減［％］を示す図である。 参照配列のアンマスク領域の全シトシン残基のコンテキストの比率を示す図である。 実施形態に係るウイルス変異予測装置による学習手順のフローチャートである。 マッピングと変異記録のイメージ図である。 同義置換（アミノ酸変異なし）を用いた場合の２つのポジションの組み合わせ例を示す図である。 選択された上位３０の特徴量の一例を示す図である。 特徴量追加無し・選択無しの場合のコンテキストとスコアの関係例を示す図である。 特徴量追加有り・選択有りの場合のコンテキストとスコアの関係例を示す図である。 特徴量追加有り・選択有りの場合のコンテキストと正則化パラメーター毎のスコアの平均値を示す図である。 特徴量追加有り・選択有りの場合のコンテキストと正則化パラメーター毎のスコアの標準偏差を示す図である。 実施形態に係る変異予測の処理手順のフローチャートである。 変異予測時に画像表示装置上に表示される情報の一例を示す図である。 ロジスティック回帰計算した結果例を示す図である。 変異記録と変異予測を示す図である。 系統樹を示す図である。 選んだ４つの変異型を各種変異型のゲノム上の変異部位と、疑似感染モデルの際に使用したＲＮＡ配列の位置を示す図である。 ｓｓＲＮＡによるＴＮＦ-α産生の誘導を示す図である。 ｓｓＲＮＡによるＩＬ-６産生の誘導を示す図である。 実施形態に係る解析プログラムの処理内容例と処理手順例を示す図である。 塩基配列の範囲ごとにグリッドサーチを行い、最適化した各モデルのハイパーパラメータの値の例を示す図である。 塩基配列の範囲－１０～＋１０；回帰式の係数をヒストグラムで表した図である。 塩基配列の範囲－１０～＋１０；回帰式の係数をヒストグラムで表した図である。 塩基配列の範囲－１０～＋１０；回帰式の係数をヒストグラムで表した図である。 塩基配列の範囲－１０～＋１０；回帰式の係数のヒストグラムを箱ひげ図でプロットした図である。 比較した学習モデルの概要と特徴を示す図である。 モデル別ＡＵＣスコアの要約統計量を解析した結果例を示す図である。 処理前のＡＵＣスコア例を示す図である。 処理後のＡＵＣスコア例を示す図である。 塩基配列の範囲－２～＋２；交差検証回数が１回目の各モデルのＲＯＣ曲線を示す図である。 塩基配列の範囲－２～＋２；交差検証回数が２回目の各モデルのＲＯＣ曲線を示す図である。 ５回の交差検証による学習データの分割方法例を示す図である。 汎化性能を測定する方法を説明するための図である。 ＧからＵへ変異する場合の各塩基配列の範囲と各学習モデルによる箱ひげ図である。 ＧからＡ（アデニン）へ変異する場合の各塩基配列の範囲と各学習モデルによる箱ひげ図である。 ＡからＧへ変異する場合の各塩基配列の範囲と各学習モデルによる箱ひげ図である。 ＵからＣ（Ｔ（チミン）からＣ）へ変異する場合の各塩基配列の範囲と各学習モデルによる箱ひげ図である。
　以下、本発明の実施の形態について図面を参照しながら説明する。なお、以下の実施形態では、対象ウイルスがＳＡＲＳ－ＣｏＶ－２の例を説明する。
　［ＳＡＲＳ－ＣｏＶ－２ウイルスの概略］
　現在、ＳＡＲＳ－ＣｏＶ－２に対するワクチン、診断法、治療法が求められている。ワクチンや抗体検査は、ＳＡＲＳ－ＣｏＶ－２のタンパク質（または遺伝子配列）をもとに作成される。ゲノム分析によれば、ＳＡＲＳ－ＣｏＶ－２には、Ａ、Ｂ、Ｃの３つのタイプに分類されるいくつかのバリアントがある。この結果、ワクチンや抗体検査のためには、ＳＡＲＳ－ＣｏＶ－２の変異型を収集する必要性がある。
　これらのＳＡＲＳ－ＣｏＶ－２バリアントに、いくつかの遺伝子変異が含まれているが、これらの変異が感染に及ぼす影響は不明である。ウイルスには、突然変異が自己複製時のエラーや細胞由来のＲＮＡ編集酵素によって、ウイルスに導入される。ＲＮＡ編集酵素は、ＲＮＡウイルスに変異を引き起こすことが知られている。
　ＲＮＡに作用するアデノシンデアミナーゼ（ＡＤＡＲ）などのＲＮＡ編集酵素やアポリポ蛋白質ＢのｍＲＮＡ編集酵素、触媒ポリペプチド（ＡＰＯＢＥＣｓ）は、ＲＮＡウイルス感染症において研究されてきた。ＡＤＡＲは、アデノシンからアミノ基を抽出してイノシンに変換する酵素であり、主に二本鎖ＲＮＡに作用する機能である。シチジンデアミナーゼの一族であるＡＰＯＢＥＣｓは、シチジンからアミノ基を抽出してウラシルに変換する酵素である。また、ＡＰＯＢＥＣｓは、ｓｓＤＮＡを基質として機能することが報告されている。さらに、ＡＰＯＢＥＣ１、ＡＰＯＢＥＣ３Ａ、ＡＰＯＢＥＣ３ＧもｓｓＲＮＡを基質として認識する。しかし、ＳＡＲＳ－ＣｏＶ－２変異体の変異が、宿主のＲＮＡ編集によって誘導されるか否かは、まだ不明である。
　そこで、本実施形態では、ＲＮＡ編集酵素に着目して、ウイルス遺伝子変異の前後数塩基の特徴的配列をもとにウイルスゲノムを検索することで、将来変異が入りうる場所、および置き換わる塩基を予測する。ウイルス変異を、事前に予測することができれば、新変異に対応する診断薬や治療薬を準備する時間ができ、変異が起こった後すぐに、診断薬や治療薬を適用できる。
　［ウイルスの点変異予測装置の構成例］
　図１は、本実施形態に係るウイルス変異予測装置１の構成の一例を示す図である。図１のように、ウイルス変異予測装置１は、取得部１１、記憶部１２、抽出部１３、分離部１４、サンプリング部１５、特徴量追加選択部１６、学習部１７、予測部１８、出力部１９、および操作部２０を備える。
　ウイルス変異予測装置１は、ＤＢ（データベース）２からネットワークＮＷを介してデータを取得する。ウイルス変異予測装置１は、取得したデータから遺伝子変異の特徴を学習させ、変異を予測する。
　取得部１１は、例えば無線ネットワーク回路である。取得部１１は、ＤＢ２（例えばＧＩＳＡＩＤ（鳥インフルエンザ情報共有の国際推進機構；https://www.gisaid.org/））からネットワークＮＷを介してデータを取得する。データは、例えばＳＡＲＳ－ＣｏＶ－２の世界のゲノムの遺伝子配列であり、複数である。
　記憶部１２は、取得された取得したＳＡＲＳ－ＣｏＶ－２のゲノムデータを記憶する。記憶部１２は、正則化パラメーターＣが変異を受けたか否かを示す情報を記憶する。記憶部１２は、Ｃ（シトシン）またはＧ（グアニン）からＵ（ウラシル）に変わった時、アミノ酸変異があるかが確認された確認結果を記憶する。記憶部１２は、学習、予測に必要なアルゴリズム、プログラム、閾値等を記憶する。
　抽出部１３は、取得したＳＡＲＳ－ＣｏＶ－２のゲノムから、Ｃを抽出する。抽出部１３は、取得したＳＡＲＳ－ＣｏＶ－２のゲノムから、ＣまたはＧからＵへの変異が起こる、または起こったコンテキストも抽出する。なお、コンテキストとは、変異部位の前後数塩基の配列のセットである。
　分離部１４は、取得したＳＡＲＳ－ＣｏＶ－２のゲノムデータのＣまたはＧからＵへの変異部分を抽出し、抽出した変異部分を１ゲノム上にマッピングする。分離部１４は、ＣまたはＧが変異を受けたか否かを示す情報を記憶部１２に記憶させる。分離部１４は、ＣまたはＧからＵに変わった時、アミノ酸変異があるかを確認し、確認結果を記憶部１２に記憶させる。分離部１４は、ＣまたはＧからＵに変わった時、アミノ酸変異があるかを確認し、アミノ酸変異がある配列を非同義置換として分離し、アミノ酸変異がない配列を同義置換であるとして分離する。
　サンプリング部１５は、アミノ酸置換のない（同義置換）のものを第１所定数選ぶ。サンプリング部１５は、ノイズを抑えるため、選択した第１所定数のうち、第１所定数より少ない第２所定数を学習データとして選ぶ。なお、サンプリング処理は、必ずしも行わなくてもよい。この場合は、同義置換全てを学習データに用いてもよい。また、サンプリング部１５は、アミノ酸置換のない（同義置換）のものを第１所定数選び、これを学習データとしてもよい。
　特徴量追加選択部１６は、特徴量（パラメーター）を追加する。なお、特徴量については後述する。例えば、特徴量は、ＲＮＡ塩基Ａ、Ｕ、Ｇ、Ｃの４種類のうち、２塩基が選ばれて特徴づけられた量である。
　学習部１７は、選ばれた第２所定数を学習データとし、第１所定数の残りをテストデータとする。学習部１７は、特徴量と学習データを用いて学習を行う。なお、学習部１７は、学習に特徴量を用いなくてもよい。なお、学習部１７は、例えば、ニューラルネットワーク、サポートベクターマシン、強化学習、ディープラーニング等のアルゴリズムを用いて学習する。また、人工知能（AI: Artificial Interigence）を用いて学習してもよい。
　予測部１８は、学習された結果を用いて、点変異を予測する。
　出力部１９は、予測部１８が予測した結果を示す情報を、画像表示装置３上に表示させる。なお、画像表示装置３は、例えばタブレット端末等であってもよい。
　操作部２０は、例えば、画像表示装置３上に設けられているタッチパネルセンサー、マウス等である。操作部２０は、利用者が操作した操作結果を検出する。
　［ＳＡＲＳ－ＣｏＶ－２の解析結果］
　ここで、発明者らが行ったＳＡＲＳ－ＣｏＶ－２の解析結果を説明する。発明者らは、ＧＩＳＡＩＤから収集したＳＡＲＳ－ＣｏＶ－２の世界のゲノム７８００遺伝子配列を網羅的に解析した。なお、収集の際、重複している配列、収集日が不明瞭な配列等は除外した。この結果、ＧＩＳＡＩＤから７８０４個の配列を取得した。
　まず、取得した配列に対して系統ネットワーク解析を行い系統樹作成した結果、５０００回以上の点突然変異の頻度が算出された。
　次に、これらの点突然変異の位置を解析した。図２は、ＳＡＲＳ－ＣｏＶ－２ゲノムにおける点突然変異の分布を示す図である。なお、図２の上の図は、全長ｓｓＲＮＡの各遺伝子の位置を示す図（ｇ１）である。図２の下のヒストグラムｇ２は、各位置での突然変異の数を示す。ヒストグラムｇ２において、縦軸は変異数であり、横軸は塩基数（ｂｐ）である。図２のように、１５０ヌクレオチド（ビン）あたりの点突然変異の平均数は約２８個であったが、いくつかの場所で点突然変異の頻度が高くなっていることが観察された。
　次に、各遺伝子における点突然変異の偏りをさらに解析するために、遺伝子ごとの点突然変異の数をカウントした。図３は、遺伝子ごとの点突然変異の数を示す図である。図３において、横軸は遺伝子名であり、縦軸は変異数である。図３のように、ＯＲＦ－１ａとＯＲＦ－１ｂは、点突然変異が多かった。
　しかし、図２のように、ＯＲＦ－１ａとＯＲＦ－１ｂは他の領域に比べて非常に長いため、より多くの突然変異が発生する可能性がある。このため、各遺伝子の１００塩基あたりの点突然変異率を推定した。図４は、各遺伝子の１００塩基あたりの点突然変異率を示す図である。図４において、横軸は遺伝子名であり、縦軸は１００塩基あたりの点突然変異率である。図４のように、遺伝子長で正規化した場合は、５’－非翻訳領域（ＵＴＲ）と３’－ＵＴＲで点突然変異の頻度が最も高かった。
　これらの結果は、点突然変異がＳＡＲＳ－ＣｏＶ－２変異体に存在することを示している。
　次に、発明者らは、遺伝子変異の可視化することで、遺伝子変異の特徴について解析した。
　図５は、変異した核酸塩基を調べた結果を示す図である。横軸は点突然変異後の置換塩基数であり、縦軸は塩基（Ａ（アデニン）、Ｕ、Ｇ（グアニン）、Ｃ）である。図５のように、Ｕへの変異が半分以上の割合を占めていることが分かった。
　図６は、それぞれの塩基がどの塩基から変異しているか調べた結果を示す図である。横軸は各点突然変異時の元の塩基と置換塩基数であり、縦軸は塩基から塩基である。この結果、Ｕへの変異は、Ｃ、Ｇ（特にＣ）からの変異が多いことが分かった。また、その他にもＧはＡ、ＡはＧ、ＵはＣへの変異が優位であることが分かった。ＣからＵとＧからＡはＡＰＯＢＥＣによって導入され、ＡからＧとＵからＣはＡＤＡＲによって導入されることが知られている。なお、実施形態では、例えばＣからＵの変異をＣｔｏＵとも書く。
　図７は、各遺伝子の変異パターンを示す図である。図８は、各遺伝子における点突然変異の数を遺伝子長で割った変異数を示す図である。図７と図８において、横軸は遺伝子名である。図７の縦軸は変異数である。図８の縦軸は１００ベースごとの変異数である。図７、図８より、遺伝子ごとに多少の違いはありながらも、ＣｔｏＵの変異が優位であった。
　さらに、図５～図８で観測された変異のうち、ＣｔｏＵとＧｔｏＡはＡＰＯＢＥＣ、ＡｔｏＧとＣｔｏＵは、ＡＤＡＲによって導入される変異と一致している。このため、発明者らは、この４つの変異について上流下流１塩基のコンテキストを調べた。
　図９は、ＣｔｏＵにおける点変異の両側の塩基配列の特徴を示す図である。図１０は、ＧｔｏＡにおける点変異の両側の塩基配列の特徴を示す図である。図１１は、ＡｔｏＧにおける点変異の両側の塩基配列の特徴を示す図である。図１２は、ＵｔｏＣにおける点変異の両側の塩基配列の特徴を示す図である。図９～図１２において、横軸は塩基名であり、縦軸はＡ、Ｕ、Ｇ、Ｃそれぞれの割合［％］である。また、図９～図１２において、左のグラフは突然変異部位の５’側の塩基（－１）を示し、右のグラフは突然変異部位の３’側の塩基（－１）を示す。
　図９のように、ＣｔｏＵの変異においては、５’側、３’側では、共にＡとＵが多く隣接していた。一方、図１０のように、ＧｔｏＡの変異において、５’側にはＡとＵが多く隣接し、３’側にはＧとＵが多く隣接していた。なお、これに相補的な配列は、［Ｃ／Ｕ］Ｃ、Ｃ［Ａ／Ｕ］である。
　図１１のように、ＡｔｏＧの変異においては、５’側ではＡとＵが多く隣接し、３’側ではＵとＧが多く隣接していた。また、図１２のように、ＵｔｏＣの変異においては、５’側にＵが多く隣接し、３’側にはＧとＵが多く隣接していた。なお、これに相補的な配列は、ＣＡ、［Ａ／Ｃ］Ｃである。
　次に、一番多く観測されたＣからＵへの変異（ｎ＝２４０１）の上流下流３塩基ずつのコンテキストについて、より詳細に検討した結果を説明する。図１３は、ＣからＵへの変異（ｎ＝２４０１）の上流下流３塩基ずつのコンテキストの特徴を示す図である。なお、０は突然変異部位を示し、負の数字と正の数字はそれぞれ上流側と下流側の部位を示す。図１３において、横方向はコンテキストの位置である。また、各数値は、各位置における塩基ＡＵＧＣそれぞれの数である。図１３のように、置換されるＣ前後にはＡとＵが非常に多くなっていた。この理由は、ＳＡＲＳ－ＣｏＶ－２にＡとＵが多く含まれているバイアスを受けていると考えられる（Ａが３０％、Ｕが３２％）。
　図１３のような特徴が得られたため、ＳＡＲＳ－ＣｏＶ－２のシークエンスの全てのＣのコンテキストについて調べ、期待値とし、それぞれの塩基について、該当する期待値からの増減［％］を調べた。図１４は、ＳＡＲＳ－ＣｏＶ－２のシークエンスの全てのＣのコンテキストにおいて、それぞれの塩基に該当する期待値からの増減［％］を示す図である。図１４において、横方向はコンテキストの位置である。図１４のように、位置＋２、＋１ではＵが高く、－１ではＧが高くなっていた（ｐ＜１０＾－３，　ｆｉｓｈｅｒ　ｅｘａｃｔ　ｔｅｓｔ）。一方、位置＋１では、Ｃが少なくなっていた（ｐ＜０．０１，ｆｉｓｈｅｒ　ｅｘａｃｔ　ｔｅｓｔ）。これは、変異を導入するＡＰＯＢＥＣの基質特異性を示唆している可能性がある。なお、Ｃ－ｔｏ－Ｕ変異におけるシトシン残基の上流側（－３）および下流側（＋３）のコンテキストの期待値は、参照配列のアンマスク領域の全シトシン残基のコンテキストの比率から算出した（図１５）。図１５は、参照配列のアンマスク領域の全シトシン残基のコンテキストの比率を示す図である。
　以上のような解析によって、以下の４つの遺伝子変異の特徴を発見した。
Ｉ．ウラシル（Ｕ）変異が多い
ＩＩ．シトシン（Ｃ）からウラシル（Ｕ）への変異が多い
ＩＩＩ．遺伝子変異にはＲＮＡ編集酵素が関与している
ＩＶ．ウラシル変異の前後１塩基から３塩基には特徴的な配列がある
　［学習手順］
　次に、ウイルス変異予測装置１による学習手順例を説明する。なお、実施形態では、ＳＡＲＳ－ＣｏＶ－２のゲノムを教師データとして使用した。図１６は、本実施形態に係るウイルス変異予測装置１による学習手順のフローチャートである。
　（ステップＳ１）取得部１１は、ＤＢ２（例えばＧＩＳＡＩＤ）からＳＡＲＳ－ＣｏＶ－２のゲノムデータを取得する。取得部１１は、取得したＳＡＲＳ－ＣｏＶ－２のゲノムデータを記憶部１２に記憶させる。
　（ステップＳ２）抽出部１３は、取得したＳＡＲＳ－ＣｏＶ－２のゲノムから、ＣまたはＧを選び出す。抽出部１３は、取得したＳＡＲＳ－ＣｏＶ－２のゲノムから、ＣまたはＧからＵへの変異が起こる、または起こったコンテキストｇ１１（図１７）も抽出する。図１７は、マッピングと変異記録のイメージ図である。なお、コンテキストは、例えば３通り（－２から＋２、－３から＋３、－１０から＋１０）である。
　（ステップＳ３）分離部１４は、取得したＳＡＲＳ－ＣｏＶ－２のゲノムデータのＣまたはＧからＵへの変異部分を抽出し、抽出した変異部分を１ゲノム上にマッピングする（図１７）。
　（ステップＳ４）分離部１４は、ＣまたはＧが変異を受けたか否かを示す情報を記憶部１２に記憶させる（図１７）。分離部１４は、例えば、ＣまたはＧからＵに変異した場合を１と記憶させ、ＣまたはＧのままを０として数値化して記憶させる。
　（ステップＳ５）分離部１４は、ＣまたはＧからＵに変わった時、アミノ酸変異があるかを確認し、確認結果を記憶部１２に記憶させる。分離部１４は、アミノ酸変異があると判別した場合（ステップＳ５；ＹＥＳ）、ステップＳ６の処理に進める。分離部１４は、アミノ酸変異がないと判別した場合（ステップＳ５；ＮＯ）、ステップＳ７の処理に進む。
　（ステップＳ６）分離部１４は、非同義置換であると判別し、このデータも学習に使用する。
　（ステップＳ７）分離部１４は、同義置換であると判別し、このデータを学習に使用する。なお、同義置換となる約１８００部位のうち６７５部位で変異を確認した。処理後、分離部１４は、ステップＳ８の処理に進める。
　（ステップＳ８）サンプリング部１５は、アミノ酸置換のない（同義置換）のものを１０００個（変異有５００、無５００）に選ぶ（第１ランダムサンプリング）。なお、サンプリング部１５は、このランダムな選択を、５回行ってアミノ酸置換のない（同義置換）のものを１０００個選択する。
　（ステップＳ９）一般的に機械学習の際には、学習データを６０から８０％とすることが多いため、サンプリング部１５は、選択した１０００個のうち８００個を学習データとして選ぶ（第２ランダムサンプリング）。なお、サンプリング部１５は、ランダムな選択を、５回行って８００個を選ぶ。なお、サンプリング部１５は、この処理は、行わなくてもよい。
　（ステップＳ１０）学習部１７は、選んだ８００個を学習データとし、残り２００個をテストデータとする。なお、学習部１７は、変異なしも学習データに使用する。
　（ステップＳ１１）特徴量追加選択部１６は、特徴量（パラメーター）を追加する。例えば、－１０から＋１０の塩基配列において、ＲＮＡ塩基は、Ａ、Ｕ、Ｇ、Ｃの４種類あり、２０塩基の配列があることから、特徴量は８０種類（＝４×２０）となる。このうち、２塩基を選んで特徴づけるため、８０の２乗の６４００種類あるが、組み合わせなので、特徴量はその半分の３２００通りになる。続けて、特徴量追加選択部１６は、３２００通りのパラメーターの中で、例えば上位３０を選び出す。なお、特徴量の個数は一例であり、これに限らない。特徴量追加選択部１６は、基準にｘ二乗検定を選択し、ＳｅｌｅｃｔＫＢｅｓｔ（ｃｈｉ２、Ｋ＝３０）を用いた。なお、この特徴量は、学習の際、スコア（ここでのスコアは正答率と同義）を向上させるために用いられる。なお、特徴量は、図１９のように、コンテキストの中で選ばれた２塩基の組み合わせである。
　（ステップＳ１２）学習部１７は、特徴量と学習データを用いて学習を行う。
　（ステップＳ１３）予測部１８は、学習された結果を用いて、点変異を予測する。なお、予測については後述する。
　なお、上述した例では、コンテキストが３通り（－２から＋２、－３から＋３、－１０から＋１０）の例を示したが、これに限らない。コンテキストは、－３から＋３以上、－１０から＋１０以下であればよい。なお、－３から＋３以上、－１０から＋１０以下とは、－４から＋４、・・・、－９から＋９を含む。
　図１８は、同義置換（アミノ酸変異なし）を用いた場合の２つのポジションの組み合わせ例を示す図である。例えば１行目の「１＿Ｇ　４＿Ｇ」は、１＿Ｇが位置＋１のＧを示し、４＿Ｇが位置＋４のＧを示す。また、２行目の「－２＿Ｔ　１＿Ｇ」は、ＴＮＣＧのコンテキストを示している。
　図１９は、選択された上位３０の特徴量の一例を示す図である。なお、ハッチングｇ２１は増加を表し、ハッチングｇ２２は減少を表している。なお、選択される特徴量は、３０個に限らない。
　［特徴有無によるスコアの比較］
　ここで、特徴量を追加しなかった場合と追加かした場合の学習結果のスコアの差を説明する。１０００部位のランダムサンプリングしたもののうち、８００部位を学習データとし、２００部位をテストデータとし、交差検証を行った（ｎ＝５）。結果を図２０と図２１に示す。
　図２０は、特徴量追加無し・選択無しの場合のコンテキストとスコアの関係例を示す図である。図２１は、特徴量追加有り・選択有りの場合のコンテキストとスコアの関係例を示す図である。図２０と図２１において、横軸はコンテキスト｛（－２、＋２）、（－３、＋３）、（－１０、＋１０）｝であり、縦軸はスコアである。また、各コンテキストにおいて点は、ロジスティック回帰における正則化パラメーターＣ値であり、左から順に０．０００１、０．００１、０．０１、０．１、１．０、１０．０、１００．０、１０００．０である。１つの点は交差検証のスコア（ｎ＝５）の平均値である。なお、正則化パラメーターは、値が大きいほど学習しやすいことを示している。また、スコアは正答率を示し、ばらつきはデータの偏りに対する頑健性を示す。
　特徴量追加無し・選択無しの場合は、図２０のようにコンテキストの範囲を増やしても学習結果のスコアが向上しなかった。これに対して、特徴量追加有り・選択有りの場合は、図２１の矢印ｇ３１のように、コンテキストの範囲を増やすと学習結果のスコアが向上した。このように、特徴量の追加、選択は、変異予測に有効であることが分かった。
　実施形態では、上述したように機械学習で変異を予測させるために、特徴量を加え８００個を学習させる。その時、予測部１８は、特徴量（上位３０）を加えたことで上位３０の順位に応じて、係数をかけて計算して予測している。特徴量（上位３０の中でも）の中でも、真に重要なものとノイズがある。
　Ｃ値を学習のしやすさと表現したが、その意味は、特徴量をもとに係数をかけて計算しているものの、その中にもノイズが含まれているので、Ｃ値で分類した（Ｃ値が小さければノイズなし、大きければノイズを含む）。
　例えば、Ｃ＝０．０００１であれば、ノイズを拾っていない学習のため、学習しきれていず、Ｃ＝１０００であれば、ノイズも拾って学習している、ということである。
　図２１では、適正のＣ値（－３から＋３、‐１０から＋１０とも）は、スコアが上がりきった最初の位置（Ｃ＝０．１か、１．０）が適正値と考えられる（真の特徴量の係数で計算できている可能性が高い）。なお、過学習は、ノイズを拾うので、ノイズをいかに少なくし、真に重要な特徴量を学ばせるかが重要となる。
　図２２は、特徴量追加有り・選択有りの場合のコンテキストと正則化パラメーター毎のスコアの平均値を示す図である。図２３は、特徴量追加有り・選択有りの場合のコンテキストと正則化パラメーター毎のスコアの標準偏差を示す図である。
　図２２、図２３のようにコンテキスト－２から＋２と、３から＋３との比較では、－３から＋３の方が、スコアが高いこと、かつ、ばらつき（標準偏差）が小さい。スコアが高いと正答率が高く、ばらつきが小さいと得られる結果の妥当性が高いことを示すため、実用的になっていると考えられる。
　さらに、コンテキスト－３から＋３よりも－１０から＋１０の方が、よりスコアが高く、ばらつきも小さかった。したがって、コンテキストは、－２から＋２よりも－３から＋３の方が良く、－３から＋３よりも－１０から＋１０が良い。すなわち、－１０から＋１０のコンテキストが一番良かった。
　［変異予測］
　次に、本実施形態における変異予測の例を説明する。図２４は、実施形態に係る変異予測の処理手順のフローチャートである。なお、予測に際して、上述した学習が予め行われている。
　（ステップＳ１０１）予測部１８は、予測した結果のスコアを算出し、算出したスコアを出力部１９を介して画像表示装置３に表示させる。この結果、画像表示装置３には、例えば図２５のような、コンテキストとスコアの関係のグラフが表示される。図２５は、変異予測時に画像表示装置３上に表示される情報の一例を示す図である。
　（ステップＳ１０２）利用者は、表示された画像（図２５）を見て、例えばコンテキスト－３から＋３のＣ＝０．１の領域ｇ４１を選択する。操作部２０は、利用者によって選択された選択情報を予測部１８に出力する。
　（ステップＳ１０３）予測部１８は、選択されたコンテキストの正則パラメーターに対して、所定のアルゴリズム（例えばロジスティック回帰）で図２６のような統計的処理を行う。図２６は、ロジスティック回帰計算した結果例を示す図である。図２６の縦軸はスコアであり、直線ｇ４２は変異有無の閾値である。予測部１８は、図２６のようなグラフを画像表示装置３上に表示させる。
　（ステップＳ１０４）利用者は、表示された画像（図２６）を見て、変異有りの点、例えば点ｇ４３を選択する。操作部２０は、利用者によって選択された選択情報を予測部１８に出力する。
　（ステップＳ１０５）予測部１８は、選択された点を、図２７のように１つのＳＡＲＳ－ＣｏＶ－２ゲノム上の位置ｇ４４にマッピングし、マッピングした画像を画像表示装置３上に表示させる。図２７は、変異記録と変異予測を示す図である。
　（ステップＳ１０６）予測部１８は、表示された画像（図２７）において抽出部分を操作部２０が操作されて選択されたことを検出した際、図２６のどこに相当するかを表示する（バックキャスト機能）。なお、予測部１８は、１つの画面内に、図２５～図２７の全てを表示してもよく、すくなくとも１つを表示して切り替えて表示するようにしてもよい。本実施形態では、このように、例えば図２６と図２７で双方向に選択、マッピングできる。
　なお、図２４に示した処理手順は一例であり、これに限らない。
　上述したように、新型コロナウイルスの世界のゲノム７８００遺伝子配列を網羅的に解析したところ、ウイルス遺伝子変異には特徴があることが判明した。その特徴は、１）ウラシル（Ｕ）変異が多いこと、２）シトシン（Ｃ）からウラシル（Ｕ）への変異が多いこと、３）遺伝子変異にはＲＮＡ編集酵素が関与していること、４）ウラシル変異の前後１塩基から３塩基には特徴的な配列があること、がわかった。また、コロナウイルスはＲＮＡ校正酵素をもつため、変異が点変異に限定され、かつＲＮＡ編集酵素による変異が顕在化していると考えられた。その結果、本実施形態では、ＲＮＡ編集酵素に着目して、ウイルス遺伝子変異の前後数塩基の特徴的配列をもとにウイルスゲノムを検索することで、将来変異が入りうる場所、および置き換わる塩基を予測することが可能となった。すなわち、本実施形態によれば、新型コロナウイルスの将来起こりうる変異を予測することができるようになった。
　本実施形態では、ウイルス遺伝子変異の前後数塩基の特徴的配列をもとにウイルスゲノムを検索し、これまでの変異（ＣまたはＧからＵ）を教師データとして、機械学習させ、変異を予測するようにした。
　この結果、本実施形態では、６０から７０％の精度正答率でウイルス変異を予測することが可能となった。ただし、この正答率は、ＲＮＡ編集酵素による変異のみならず突然変異も含んだ正答率であり、突然変異とＲＮＡ編集酵素のみの変異を区別することで、ＲＮＡ編集酵素による変異予測の正答率はより高率になっていることは容易に想像できる。なお、上記においてＡＵＣ（Area Under the Curve）スコアを正答率として用いた。ＡＵＣスコアの算出等については後述する。
　これにより、本実施形態によれば、ウイルス変異を、変異が起こる前に事前に予測できれば、ウイルス感染診断において、事前に診断キットを準備することができる。本実施形態によれば、超早期診断キットの開発を可能とする発明である。また、本実施形態によれば、診断キットのみならず、ワクチンの効果判定や、ウイルス抗体医薬の効果判定、免疫パスポートの認証や取り消しも可能とする。加えて、本実施形態によれば、治療薬の候補選択も可能となることから、超早期治療を可能とする。
　［検証結果］
　以下、上記の学習、予測について、検証した結果の一例を説明する。
　点変異によりウイルスゲノムのＵが増えることが明らかとなった。Ｕが増加することによる、炎症の増強が考えられるため、炎症性サイトカイン産生が変化するか否かを調べた。細胞刺激アッセイのために、４つの異なる配列、ＥＰＩ＿ＩＳＬ＿４１９３０８、ＥＰＩ＿ＩＳＬ＿４１５６４４、ＥＰＩ＿ＩＳＬ＿４１８４２０、およびＥＰＩ＿ＩＳＬ＿４１９８４６をＳＡＲＳ－ＣｏＶ－２バリアントから選択した。これらの変異配列は、それぞれ日本、ジョージア州、フランス、オーストラリアで検出されたものである。
　作業者は、４種類の変異体のそれぞれの一本鎖ＲＮＡ（ｓｓＲＮＡ）の全長の中から、Ｕへの変異が観察された１つの領域を抽出し、合成した。
　これらの異なるバリアントから得られたｓｓＲＮＡの配列は以下の通りであった。ｖａｒｉａｎｔ－１（５’－ＡＵＵＵＡＵＵＧＵＵＣＵＵＵＵＡＣＣＣ－３’；ａｔ２９４６－２９６５　ｒｅｇｉｏｎ　ｉｎ　ＥＰＩ＿ＩＳＬ＿４１９３０８）、ｖａｒｉａｎｔ－２（５’－ＡＵＵＵＡＵＵＧＵＵＣＵＵＵＵＵＣＵＵＵＵＡＣＣＣ－３’；ＥＰＩ＿ＩＳＬ＿４１５６４４の１１０４１～１１０６０領域）、バリアント－３（５’－ＵＵＵＣＵＡＣＡＧＵＧＵＣＣＣＡＣＵＵ－３’；ＥＰＩ＿ＩＳＬ＿４１８４２０の１４３９２～１４４１１領域）、およびバリアント－４（５’－ＡＡＡＣＣＵＵＵＧＡＧＡＧＡＧＵＵ－３’；ＥＰＩ＿ＩＳＬ＿４１９８４６の２２９４６～２２９６５領域）。
　変異したＳＡＲＳ－ＣｏＶ－２配列の対照として、参照配列（ＭＮ９０８９４７）と同じ領域を用いた。それぞれの異なる４種類の変異体に対応する参照配列は以下の通りであった。武漢－１（５’－ＡＵＧＵＡＡＵＧＵＵＣＵＣＣＣ－３’；ａｔ　３０２３－３０４２　ｒｅｇｉｏｎ）、武漢－２（５’－ＵＣＵＣＵＡＵＧＵＣＵＣＵＣＵＣＣＵＣＣＣ－３’；ａｔ　１１０６６－１１０８５　ｒｅｇｉｏｎ）、武漢－３（５’－ＵＣＵＣＵＡＵＣＡＧＵＣＣＣＵＣＣＣＵＣＣＵＣＵＣＵ－３’；ａｔ　１４３９０－１４４０９　ｒｅｇｉｏｎ、１１０６６－１１０８５領域）、武漢－３（５’－ＵＣＵＣＵＡＣＣＵＡＣＧＵＧＵＣＣＣＣＵＣＵ－３’；１４３９０－１４４０９領域）、および武漢－４（５’－ＡＡＡＣＣＣＵＡＣＵＵＵＧＵＡＧＡＧＡＧＵＡＵＡＵ－３’；２２９４６－２２９６５領域）。
　ＴＬＲ７媒介サイトカイン産生の誘導には、Ｕを含まない配列（５’－ＧＡＣＡＧＡＧＡＧＡＧＡＡＣＡＡＧ－３’）をネガティブコントロールとして用いた。検証には、株式会社日本遺伝子研究所（宮城県仙台市）が合成したｓｓＲＮＡを用いた。
　ヒト単球性白血病細胞株ＴＨＰ－１は、１０％ＦＣＳ、５５ｍＭ　２－メルカプトエタノール、１００ｍＭ　非必須アミノ酸（ＮＥＡＡｓ）、１ｍＭ　ピルビン酸および２０ｍＭ　ｍｌ－１の各ペニシリンおよびストレプトマイシンを添加したＲＰＭＩ－１６４０培地で維持した。
　４×１０＾５細胞を、９６ウェル平底プレートを用いて１５０μｌ　ＲＰＭＩで培養した。Ｙａｎ　Ｌｉらに準じて疑似感染モデルを行った。
　発明者らは、初期に報告された武漢型（Ｗ）をもとに、ＧＩＳＡＩＤより遺伝子配列を収集して、図２８の系統樹を作成した。図２８は、系統樹を示す図である。検証では、各ＲＮＡ配列の全長内でのＵへの点突然変異の頻度を調べるために、ＳＡＲＳ－ＣｏＶ－２変異体から以下の４つの異なる配列を選んだ。４つの配列は、第１バリアント（ｖａｒｉａｎｔ－１、日本型）、第２バリアント（ｖａｒｉａｎｔ－２、ジョージア型）、第３バリアント（ｖａｒｉａｎｔ－３、フランス型）、および第４バリアント（ｖａｒｉａｎｔ－４、オーストラリア型）に由来している。また、図２８において、Ｗは、武漢で報告されたオリジナルのＳＡＲＳ－ＣｏＶ－２配列を示す。
　図２９は、選んだ４つの変異型を各種変異型のゲノム上の変異部位と、疑似感染モデルの際に使用したＲＮＡ配列の位置を示す図である。図２９において、横方向は（ｂｐ）であり、下向き三角形はＶｔｏＵ（ＶはＵ以外のすべての塩基）であり、上向き三角形はＵｔｏＶであり、四角は細胞刺激に用いたｓｓＲＮＡの配列を示す。図２９のように、各ＳＡＲＳ－ＣｏＶ－２変異体の全長ｓｓＲＮＡ内のＵの数は、元の単離株と比較して有意に増加している。また、図２９のように、Ｕに対する点突然変異の頻度は、Ａ、Ｇ、またはＣに対するＵの頻度よりもはるかに高かった。このように、完全長の変異したｓｓＲＮＡの炎症性サイトカインを誘導する能力は、元の分離株よりもはるかに大きい。これらの結果は、ＳＡＲＳ－ＣｏＶ－２遺伝子変異が炎症性活性化の亢進を引き起こすメカニズムの一つである可能性を示唆している。
　これまでのいくつかの研究では、Ｕ－ｒｉｃｈ　ｓｓＲＮＡがＴＬＲ７シグナルを介して自然免疫細胞を刺激し、炎症性サイトカインを産生することが示されている。このように、点突然変異に起因する多数のＵ残基が、ヒトマクロファージによる炎症性サイトカインの誘導を促進しているのではないかという仮説を立てた。
　この仮説を検証するために、ＳＡＲＳ-ＣｏＶ-２変異体のＵリッチ領域で刺激されたヒト単球／マクロファージ細胞株ＴＨＰ-１におけるＴＮＦ-αおよびＩＬ（インターロイキン）-６の産生を解析した。図３０は、ｓｓＲＮＡによるＴＮＦ-α産生の誘導を示す図である。ヒトＴＮＦ－αの測定のために、細胞をＰＭＡ（０．２ｎｇ／ｍｌ、Ｓｉｇｍａ　Ａｌｄｒｉｃｈ、Ｓｔ．Ｌｏｕｉｓ、ＭＯ、ＵＳＡ）の存在下で培養し、ＤＯＴＡＰ（１０μｇ、Ｒｏｃｈｅ　Ｄｉａｇｎｏｓｔｉｃｓ、Ｍａｎｎｈｅｉｍ、Ｇｅｒｍａｎｙ）を用いて１６０（ｐｍｏｌ）のｓｓＲＮＡで刺激した。サイトカインの検出のため、ヒトＴＮＦ－αおよびＩＬ－６は、ＯｐｔＥＩＡセット（ＢＤ　Ｂｉｏｓｃｉｅｎｃｅ、Ｓａｎ　Ｄｉｅｇｏ、ＣＡ　ＵＳＡ）を用いて培養上清中で測定した。ＴＮＦ-αの産生は、１８時間の刺激後に測定した。
　なお、図３０と図３１において、例えばＷ－１は初期の武漢型を示し、ｖａｒｉａｎｔ－１は変異型を示す。
　図３１は、ｓｓＲＮＡによるＩＬ-６産生の誘導を示す図である。ヒトＩＬ-６の測定のために、細胞をＰＭＡ（５０ｎｇ／ｍｌ）の存在下で培養し、ＤＯＴＡＰ（１５μｇ）を用いて４８０（ｐｍｏｌ）のｓｓＲＮＡを用いて刺激した。ＩＬ-６の産生は、４８時間の刺激後に測定した。
　値は平均±ＳＤ（ｎ＝６）である。データは、類似の結果を有する２つの独立した実験の代表である。
　なお、フィッシャー厳密検定は、Ｐｙｔｈｏｎ　３ベースパッケージのｓｃｉｐｙ　１．４．１を用いて片側検定で行った。また、Ｍａｎｎ－Ｗｈｉｔｎｅｙ　Ｕ検定は、Ｐｒｉｓｍ　８　ソフトウェア（ＧｒａｐｈＰａｄ　Ｓｏｆｔｗａｒｅ、Ｓａｎ　Ｄｉｅｇｏ、ＣＡ）を用いて実施した。Ｐ＜０．０５の値は、有意性を示す。
　図３０のように、予想通り、Ｕ残基を欠いたｓｓＲＮＡ配列は、ＴＮＦ－αの産生をアップレギュレートしなかった。点突然変異によって誘導されたＵ数の増加は、武漢型の参照ｓｓＲＮＡ配列による刺激と比較して、ｖａｒｉａｎｔ－１、３、４のサイトカイン産生を増加させた。
　また、図３１のように、ＴＮＦ-αに比べてＩＬ-６の産生は低かったが、ＩＬ-６の産生には同様の傾向が観察された。これらの結果は、ＳＡＲＳ－ＣｏＶ－２ゲノム内のＵへの点変異が、ＴＮＦ－αやＩＬ－６などの炎症性サイトカインの産生増加を刺激する能力をもたらすことを示している。すなわち、Ｕ変異の予測は、炎症性サイトカイン産生の増強も予測できるため、患者の炎症の症状、重症化も判別が可能となる。
　そして、本実施形態では、取得部が、ウイルスのゲノムの遺伝子配列データを取得する。抽出部が、取得したゲノムの遺伝子配列データからＣ（シトシン）またはＧ（グアニン）を抽出し、ＣまたはＧからＵ（ウラシル）への変異が起こるまたは起こったコンテキストを抽出する。
　本実施形態では、上述したように、抽出したコンテキストの塩基配列がＣまたはＧからＵに変わった場合、アミノ酸変異があるかを確認する。ＲＮＡ編集酵素による変異は、ゲノムＲＮＡに直接作用して変異を誘導するためにアミノ酸変異の有無にかかわらずおこると考えられる。しかしその一方で、アミノ酸変異がある場合は、変異の要因にかかわらずウイルスの生存にかかわる変異があるため存在しないウイルス、存在しないゲノムデータもあるはずである。そのため、アミノ酸変異を含む変異データ自体が偏ったデータであると考えられる。従って、学習データするのには、アミノ酸変異のないデータを用いるのが妥当である。
　そこで、本実施形態では、分離部が、アミノ酸変異がある配列を非同義置換として分離し、前記アミノ酸変異がない配列を同義置換であるとして分離する。そして、学習部が、同義置換の配列データを学習データに用いて学習し、予測部が、学習された結果を用いて、ウイルスの変異を予測する。
　［解析プログラム］
　ここで、上述したウイルス変異予測装置１をソフトウエアプログラムである解析プログラムで実現した例を説明する。図３２は、本実施形態に係る解析プログラムの処理内容例と処理手順例を示す図である。図３２において、縦方向は主な処理であり、横方向は処理手順である。
　解析プログラムは、前処理（ステップＳ２１０）で、解析を行う対象のファイルを読み込み（ステップＳ２１１）、説明変数・目的関数の設定を行い（ステップＳ２１２）、特徴量作成用の関数を定義し（ステップＳ２１３）、塩基配列の範囲およびグリッドサーチ用パラメーターの設定を行う（ステップＳ２１４）。
　なお、目的変数とは変異の有無であり、説明変数はダミー化した塩基配列および塩基割合の２つである。特徴量作成用の関数は、例えば塩基配列の範囲（例：－３～＋３）を引数として、塩基割合（１レコードに対して「Ａ」「Ｇ」「Ｃ」「Ｔ」がそれぞれ全体の何％含まれているか）を算出する関数である。
　解析プログラムは、学習プロセス（ステップＳ２２０）で、特徴量の作成を行い（ステップＳ２２１）、グリッドサーチによるパラメーターの最適化を行い（ステップＳ２２２）、交差検証・各種モデルの学習を実行し（ステップＳ２２３）、各種モデルのＡＵＣスコアを算出する（ステップＳ２２４）。
　なお、特徴量の作成では、特徴量作成用の関数を用いて塩基の割合を算出し、引数に指定した変数をダミー化する関数を用いて塩基配列のダミー変数化を行う。また、ＡＣＵスコアは、ＲＯＣ（Receiver Operating Characteristic Curve）曲線を作成した時に、グラフの曲線より下の部分の面積であり、例えば０から１までの値をとり、値が１に近いほど判別能が高いことを示す。
　解析プログラムは、精度評価処理（ステップＳ２３０）で、各種モデルのＡＵＣスコアを出力し（ステップＳ２３１）、ＡＵＣスコアの要約統計量を算出する（ステップＳ２３２）。
　解析プログラムは、データ可視化処理（ステップＳ２４０）で、回帰式の係数をヒストグラムで表し箱ひげ図でプロットし（ステップＳ２４１）、各種モデルのＲＯＣ曲線をプロットする（ステップＳ２４２）。
　［モデルのハイパーパラメータの最適化の解析］
　次に、モデルのハイパーパラメータの最適化について解析した結果例を説明する。解析では、各モデルのハイパーパラメータについて、塩基配列の範囲ごとにグリッドサーチを行い、最適化した数値を算出した。
　図３３は、塩基配列の範囲ごとにグリッドサーチを行い、最適化した各モデルのハイパーパラメータの値の例を示す図である。図３２の例では、モデルはロジスティック回帰と、決定木と勾配ブースティングを組み合わせた手法であるＬｉｇｈｔＧＢＭである。図３３の解析条件は、交差検証の回数が５回である。塩基配列の範囲は、－２～＋２、－３～＋３、－５～＋５、および－１０～＋１０である。ロジスティック回帰の検証に使用したハイパーパラメータは、Ｃ：［０．０００１，０．００１，０．０１，０．１，１，１０，１００，１０００］であり、ＬｉｇｈｔＧＢＭの検証に使用したハイパーパラメータは、ｎｕｍ＿ｌｅａｖｅｓ：［１０，３１，６４］、ｌｅａｒｎｉｎｇ＿ｒａｔｅ：［０．０１，０．１，１］である。
　図３３のように、塩基配列の範囲が広がるにつれてロジスティック回帰における正則化の強さが増す傾向にある。また、ＬｉｇｈｔＧＢＭのハイパーパラメータ「ｌｅａｒｎｉｎｇ＿ｒａｔｅ」は０．０１で一定であった。
　［塩基配列の範囲別ロジスティック回帰の相関係数比較］
　次に、塩基配列の範囲－２～＋２、－３～＋３、－５～＋５および－１０～＋１０についてロジスティック回帰の相関係数を比較した結果の一例として塩基配列の範囲－１０～＋１０の結果を、図３４～図３７に示す。解析では、塩基配列の範囲別に各変数のロジスティック回帰の相関係数をプロットし、グループ間の比較を行った。なお、図３４～図３６では、Ａ、Ｃ，Ｇ、Ｔ、Ａ＿ｐｅｒｃｅｎｔ、Ｇ＿ｐｅｒｃｅｎｔ、Ｃ＿ｐｅｒｃｅｎｔ、およびＴ＿ｐｅｒｃｅｎｔの解析結果を示している。なお、Ａ＿ｐｅｒｃｅｎｔ、Ｇ＿ｐｅｒｃｅｎｔ、Ｃ＿ｐｅｒｃｅｎｔ、およびＴ＿ｐｅｒｃｅｎｔは、１レコードあたりの塩基の割合である。また、図３４～図３６において、０～４は交差回数毎の係数であり、例えば０の棒グラフが１回目の交差検証時の係数である。図３４～図３６において、横軸はパラメーター、縦軸は率である。図３７において、縦軸は率である。
　図３４～図３６は、塩基配列の範囲－１０～＋１０；回帰式の係数をヒストグラムで表した図である。図３７は、塩基配列の範囲－１０～＋１０；回帰式の係数のヒストグラムを箱ひげ図でプロットした図である。図３４～図３６のように、塩基配列の範囲－１０～＋１０では、－６Ｔ、－２Ｇ、－２Ｔ、－１Ｇ、－１Ｔ、＋５Ａ等の値が大きかった。
　また、塩基配列の範囲－２～２では、－２Ｔ、＋１Ｇの値が大きかった。塩基配列の範囲－３～３では、－２Ｔ、－１Ｇ、＋１Ｇの値が大きかった。塩基配列の範囲－５～５では、－２Ｔ、－１Ｇ、－１Ｔ、＋１Ｇ等の値が大きかった。
　なお、このような相関係数は、後述する塩基ごとの重みを可視化するために使用した。
　［モデル別ＡＵＣスコアの要約統計量］
　次に、モデル別ＡＵＣスコアの要約統計量を解析した結果例を説明する。解析では、学習アルゴリズムごとにＡＵＣスコアの要約統計量を算出した。
　図３８は、比較した学習モデルの概要と特徴を示す図である。図３８のように、モデルは、ロジスティック回帰、ＳＶＭ（Support Vector Machine）、決定木、ランダムフォレスト、ＸＧＢｏｏｓｔ、およびＬｉｇｈｔ　ＧＢＭである。
　図３９は、モデル別ＡＵＣスコアの要約統計量を解析した結果例を示す図である。なお、図３９では、塩基配列の範囲－２～＋２、－３～＋３、－５～＋５および－１０～１０それぞれのロジスティック回帰の相関係数を用いて、小数点第３位以下を切り捨てて要約統計量を示している。図３９の画像ｇ１０１は、ＸＤＢｏｏｓｔのＲＯＣ（ROC＿ｘｇｔ）、決定木のＲＯＣ（ROC_tree）、ＬｉｇｈｔＧＢＭのＲＯＣ（ROC_lgｂ)を示している。画像ｇ１０２は、ＳＶＭのＲＯＣ（ＲＯＣ＿ｓｖｍ）、ランダムフォレストのＲＯＣ（ＲＯＣ＿ｔｆ）、ロジスティック回帰のＲＯＣ（ＲＯＣ＿ｌｒ）を示している。図３９において、交差検証によってデータを５つに分割しているため、各要約統計量のｃｏｕｎｔが５になっている。ｍｅａｎはスコア（正答率と同義）であり、ｓｔｄは標準偏差であり、ｍｉｎは最小値で、ｍａｘは最大値である。
　図３９のように、スコアは、塩基配列の範囲－１０～＋１０が５５．４％、塩基配列の範囲－２～＋２が５６．０％、塩基配列の範囲－３～＋３が５６．６％、塩基配列の範囲－５～＋５が５６．２％であった。このように、全体的にはロジスティック回帰のスコアが高めであった。他のモデルでも、５２～５７％程度のスコアが得られた。
　次に、モデルとしてロジスティック回帰を用いた場合の処理前と処理後のＡＵＣスコアについて説明する。
　図４０は、処理前のＡＵＣスコア例を示す図である。図４１は、処理後のＡＵＣスコア例を示す図である。なお、図４０、図４１では、塩基配列の範囲－２～＋２、－３～＋３、－５～＋５および－１０～＋１０それぞれのロジスティック回帰の相関係数を用いて、小数点第３位以下を切り捨てて要約統計量を示している。また、図４０、図４１では、比較のため、処理前のデータに存在しない以下の変数を削除してから、処理後のＡＵＣスコアを算出した。なお、削除した変数（ハイパーパラメータ）は、Ａ＿ｐｅｒｃｅｎｔ，Ｇ＿ｐｅｒｃｅｎｔ，Ｃ＿ｐｅｒｃｅｎｔ，Ｔ＿ｐｅｒｃｅｎｔである。
　図４０、図４１のように、ロジスティック回帰を用いた場合の処理前のＡＵＣスコアは約５１～５４％程度であるが、処理後のＡＵＣスコアは約５６～５７％程度と向上している。
　［モデル別ＲＯＣ曲線］
　次に、モデル別ＲＯＣ曲線を用いて解析した結果例を説明する。解析では、塩基配列の範囲－２～＋２、－３～＋３、－５～＋５および－１０～＋１０それぞれに対して、学習アルゴリズム毎にＲＯＣ曲線をプロットし、モデル間の比較を行った。比較結果の一例として、塩基配列の範囲－２～＋２の比較結果例を図４２～図４３に示す。図４２は、塩基配列の範囲－２～＋２；交差検証回数が１回目の各モデルのＲＯＣ曲線を示す図である。図４３は、塩基配列の範囲－２～＋２；交差検証回数が２回目の各モデルのＲＯＣ曲線を示す図である。図４２～図４３において、横軸は誤答率（１．０＝１００％）、縦軸はスコア（１．０＝１００％）である。なお、用いたアルゴリズムは、図３８に示したロジスティック回帰(Logistic Regression)、ＳＶＭ、決定木（Decision Tree）、ランダムフォレスト(Random Forest)、ＸＧＢｏｏｓｔ、およびＬｉｇｈｔ　ＧＢＭである。なお、図４２～図４３において、線ｇ２０１はＸＧＢｏｏｓｔ、線ｇ２０２は決定木、線ｇ２０３はＬｉｇｈｔ　ＧＢＭ、線ｇ２０５はＳＶＭ、線ｇ２０５はランダムフォレスト、線ｇ２０６はロジスティック回帰それぞれのＲＯＣ曲線である。
　図４２～図４３、および塩基配列の範囲－３～＋３、－５～＋５および－１０～＋１０それぞれの結果から、いずれのモデルを用いても同等の結果が得られたが、他のモデルと比較するとＳＶＭのＲＯＣの面積の差異が塩基配列の範囲によって大きかった。
　［機械学習の実装例］
　上述した解析等を行えるように、ウイルス変異予測装置１の機能を実現したプログラムは、以下のような機能を備える。
　Ｉ．解析対象のファイル読み込み、解析で用いない”１”のレコードを削除する第１関数。
　II.塩基割合算出用の第２関数を実行し、Iで読み込んだデータの塩基割合を算出し、新たな変数に格納する。
　ＩＩＩ．Ｉで読み込んだデータのうち、塩基配列の変数（ファイルの例えばＣ～Ｖ列）を、第３関数を用いてダミー変数化する。
　ＩＶ．第４関数を用いてグリッドサーチを実行し、各種モデルのパラメーターを最適化する（図３３）。
　Ｖ．第５関数を用いて５分割交差検証を実行する。
　ＶＩ．ＩＩ、ＩＩＩの変数を説明変数、Ｉで読み込んだデータのうち変異の有無（ファイルの例えばＢ列）を目的変数として第１メソッドに設定し、各モデルの学習を実行する。なお、第１メソッドは、第一引数に分類対象のテストデータを、第二引数に分類した結果の正しい答えを指定することで、機械学習を行う。
　ＶＩＩ．ＶＩの学習結果をもとに、第６関数を用いて各モデルのＡＵＣスコアを算出する。
　ＶＩＩＩ．各モデルのＡＵＣスコアの要約統計量を、統計情報を抽出する第２メソッドで算出する（例えば図３８～図４３）。
　ＩＸ．第３メソッドを用いて、ロジスティック回帰の係数をプロットする（例えば図３４～図３６）。なお、第３メソッドは、与えられたベクトル（数値で構成される配列）の平均値を高さとして、信頼区間をエラーバーとして出力するメソッドである。
　Ｘ．　第３メソッドを用いて、係数を箱ひげ図でプロットする（例えば図３７）。
　ＸＩ．プロットする第４メソッドを用いて、各モデルのＲＯＣ曲線をプロットする（例えば図４２～図４３）。
　なお、上述したＩ～ＸＩの機能、関数、メソッドは一例であり、これに限らない。
　［学習データの分割方法、汎化性能を測定する方法］
　次に、学習データの分割方法、汎化性能を測定する方法を説明する。
　図４４は、５回の交差検証による学習データの分割方法例を示す図である。
　学習データやテストデータをどのように分割するのは非常に重要な問題である。このため、本実施形態では、図４４のようにトレーニングデータとテストデータとを分割し、交差毎にトレーニングデータとテストデータとを入れ替えて学習を行うようにした。
　図４５は、汎化性能を測定する方法を説明するための図である。
　本実施形態では、汎化性能を計測する方法として図４５のようにＳｔｒａｔｉｆｉｅｄＫＦｏｌｄを行った。この処理では、分布の比率を維持したままデータを訓練用とテスト用に分割する。
　なお、図４４、図４５に示した例は一例であり、これに限らない。
　［ＧｔｏＵ，ＧｔｏＡ，ＡｔｏＧ，ＵｔｏＣ］
　上述した例では、Ｃ（シトシン）またはＧ（グアニン）からＵ（ウラシル）への変異が起こるまたは起こったコンテキストを抽出する例を説明したが、これに限らない。以下に他の変異例に対する学習結果例を、図４６～図４９に示す。なお、この場合は、ＧからＵ、ＧからＡ、ＡからＧ、またはＵからＣ（またはＴ（チミン）からＣ）への変異が起こるまたは起こったコンテキストを抽出する。なお、学習や推定では、ＲＮＡ表記されているものについてＵ（ウラシル）を抽出し、ＤＮＡ表記されているものについてＴ（チミン）を抽出する。
　図４６は、ＧからＵへ変異する場合の各塩基配列の範囲と各学習モデルによる箱ひげ図である。図４７は、ＧからＡへ変異する場合の各塩基配列の範囲と各学習モデルによる箱ひげ図である。図４８は、ＡからＧへ変異する場合の各塩基配列の範囲と各学習モデルによる箱ひげ図である。図４９は、ＵからＣ（またはＴ（チミン）からＣ）へ変異する場合の各塩基配列の範囲と各学習モデルによる箱ひげ図である。なお、図４９では、ＤＮＡ表記でＴｔｏＣと表記しているが、ＲＮＡ表記ではＵからＣである。
　なお、以下の説明において、ｘｇｂはＸＧＢｏｏｓｔ、Ｔｒｅｅは決定木、ＬａｂはＬｉｇｈｔ　ＧＢＭ、ＳｖｍはＳＶＭ、ｒｆはランダムフォレスト、Ｌｒはロジスティック回帰を示す。
　ＧからＵへの変異の場合、例えば、塩基配列の範囲－１０～＋１０の正答率の平均値は、ＸＧＢｏｏｓｔが５６．４％、決定木が５３．０％、Ｌｉｇｈｔ　ＧＢＭが５０．０％、ＳＶＭが５１．４％、ランダムフォレストが５４．０％、ロジスティック回帰が５４．０％であった。
　図４６のように、ＧからＵへの変異の場合は、塩基配列の範囲－１０～＋１０、モデルＸＧＢｏｏｓｔの組み合わせの結果が最も良かった。
　また、ＧからＡへの変異の場合、例えば、塩基配列の範囲－５～＋５の正答率の平均値は、ＸＧＢｏｏｓｔが６２．２％、決定木が５７．０％、Ｌｉｇｈｔ　ＧＢＭが６２．８％、ＳＶＭが５２．６％、ランダムフォレストが６４．２％、ロジスティック回帰が６０．２％であった。また、塩基配列の範囲－１０～＋１０の正答率の平均値は、ＸＧＢｏｏｓｔが６０．６％、決定木が５６．６％、Ｌｉｇｈｔ　ＧＢＭが６１．６％、ＳＶＭが５４．４％、ランダムフォレストが６４．２％、ロジスティック回帰が５９．８％であった。
　図４７のように、ＧからＡへの変異の場合は、塩基配列の範囲－１０～＋１０または－５～＋５、モデルランダムフォレストの組み合わせの結果が最も良かった。
　ＡからＧへの変異の場合、例えば、塩基配列の範囲－２～＋２の正答率の平均値は、ＸＧＢｏｏｓｔが５８．０％、決定木が５６．４％、Ｌｉｇｈｔ　ＧＢＭが６０．２％、ＳＶＭが４８．８％、ランダムフォレストが５７．２％、ロジスティック回帰が５８．２％であった。
　図４８のように、ＡからＧへの変異の場合は、塩基配列の範囲－２～＋２、モデルＬｉｇｈｔ　ＧＢＭの組み合わせの結果が最も良かった。
　Ｕ（またはＴ）からＣへの変異の場合、例えば、塩基配列の範囲－５～＋５の正答率の平均値は、ＸＧＢｏｏｓｔが６１．０％、決定木が６２．４％、Ｌｉｇｈｔ　ＧＢＭが６４．０％、ＳＶＭが５５．０％、ランダムフォレストが６２．４％、ロジスティック回帰が６２．６％であった。
　図４９のように、Ｕ（またはＴ）からＣへの変異の場合は、塩基配列の範囲－５～＋５、Ｌｉｇｈｔ　ＧＢＭの組み合わせの結果が最も良かった。
　以上のように、本実施形態の手法を用いる場合は、ＸＧＢｏｏｓｔ、決定木、Ｌｉｇｈｔ　ＧＢＭ、ＳＶＭ、ランダムフォレスト、ロジスティック回帰を学習モデルとして用いることができる。この結果、本実施形態によれば、学習された結果を用いて、点変異を精度良く予測することでできる。
　また、本実施形態によれば、ＧからＵへ変異する場合に加え、ＧからＡに変異する場合、ＡからＧに変異する場合、およびＴからＣに変異する場合にも、本実施形態の手法を用いて、学習された結果を用いて点変異を予測することができる。
　なお、上述した説明と図におけるコンテキストの表記について説明する。
　本明細書中でコンテキストの表記は、変異部分を０とし、上流側をマイナス（－）、下流側をプラス（＋）で表記している。また、図や明細書の中では、プラス表記がある記載と無い記載（例えば“１＿Ｇ”と“＋１＿Ｇ”）があるが、同じコンテキストを表している。また、図と明細書において、例えば“１＿Ｇ”と“＋１＿Ｇ”と“１Ｇ”と“＋１Ｇ”のように数字とアルファベットとの間にアンダーバーがある表記と無い表記とが混在しているが、これらは同じコンテキストを表している。
　また、塩基配列の範囲については、例えば－２から＋２について、明細書と図で“－２－＋２”または“－２～＋２”と記載している。
　なお、本発明におけるウイルス変異予測装置１の機能の全てまたは一部を実現するためのプログラムをコンピュータ読み取り可能な記録媒体に記録して、この記録媒体に記録されたプログラムをコンピュータシステムに読み込ませ、実行することによりウイルス変異予測装置１が行う全ての処置または一部の処理を行ってもよい。また、機械学習には、ディープラーニング法など種々の学習法を用いてもよく、人工知能（AI: Artificial Interigence）を用いて処理を行ってもよい。なお、ここでいう「コンピュータシステム」とは、ＯＳや周辺機器等のハードウェアを含むものとする。また、「コンピュータシステム」は、ホームページ提供環境（あるいは表示環境）を備えたＷＷＷシステムも含むものとする。また、「コンピュータ読み取り可能な記録媒体」とは、フレキシブルディスク、光磁気ディスク、ＲＯＭ、ＣＤ－ＲＯＭ等の可搬媒体、コンピュータシステムに内蔵されるハードディスク等の記憶装置のことをいう。さらに「コンピュータ読み取り可能な記録媒体」とは、インターネット等のネットワークや電話回線等の通信回線を介してプログラムが送信された場合のサーバやクライアントとなるコンピュータシステム内部の揮発性メモリ（ＲＡＭ）のように、一定時間プログラムを保持しているものも含むものとする。
　また、上記プログラムは、このプログラムを記憶装置等に格納したコンピュータシステムから、伝送媒体を介して、あるいは、伝送媒体中の伝送波により他のコンピュータシステムに伝送されてもよい。ここで、プログラムを伝送する「伝送媒体」は、インターネット等のネットワーク（通信網）や電話回線等の通信回線（通信線）のように情報を伝送する機能を有する媒体のことをいう。また、上記プログラムは、前述した機能の一部を実現するためのものであってもよい。さらに、前述した機能をコンピュータシステムにすでに記録されているプログラムとの組み合わせで実現できるもの、いわゆる差分ファイル（差分プログラム）であってもよい。
　以上、本発明を実施するための形態について実施形態を用いて説明したが、本発明はこうした実施形態に何等限定されるものではなく、本発明の要旨を逸脱しない範囲内において種々の変形および置換を加えることができる。
　１…ウイルス変異予測装置、２…ＤＢ、３…画像表示装置、１１…取得部、１２…記憶部、１３…抽出部、１４…分離部、１５…サンプリング部、１６…特徴量追加選択部、１７…学習部、１８…予測部、１９…出力部、２０…操作部、Ａ…アデニン、Ｕ…ウラシル、Ｇ…グアニン、Ｃ…シトシン、Ｔ…チミン

Claims (10)

  1. 　ウイルスのゲノムの遺伝子配列データを取得する取得部と、
    　取得した前記ゲノムの遺伝子配列データからＣ（シトシン）またはＧ（グアニン）を抽出し、ＣまたはＧからＵ（ウラシル）への変異が起こるまたは起こったコンテキストを抽出する抽出部と、
    　ＣまたはＧからＵに変わった時、アミノ酸変異があるかを確認し、前記アミノ酸変異がある配列を非同義置換として分離し、前記アミノ酸変異がない配列を同義置換であるとして分離する分離部と、
    　前記同義置換の配列データを学習データに用いて学習する学習部と、
    　学習された結果を用いて、前記ウイルスの変異を予測する予測部と、
    　を備えるウイルス変異予測装置。
  2. 　ウイルスのゲノムの遺伝子配列データを取得する取得部と、
    　取得した前記ゲノムの遺伝子配列データからＣ（シトシン）、Ｇ（グアニン）、Ａ（アデニン）、Ｕ（ウラシル）またはＴ（チミン）を抽出し、ＧからＡ、ＡからＧ、ＵからＣ、ＴからＣへの変異が起こるまたは起こったコンテキストを抽出する抽出部と、
    　抽出されたコンテキストの塩基配列が変化した場合、アミノ酸変異があるかを確認し、前記アミノ酸変異がある配列を非同義置換として分離し、前記アミノ酸変異がない配列を同義置換であるとして分離する分離部と、
    　前記同義置換の配列データを学習データに用いて学習する学習部と、
    　学習された結果を用いて、前記ウイルスの変異を予測する予測部と、
    　を備えるウイルス変異予測装置。
  3. 　前記同義置換から所定数を選ぶサンプリング部、を更に備え、
    　前記学習部は、前記サンプリング部によって選ばれた前記同義置換の配列データを学習データに用いる、
    　請求項１または請求項２に記載のウイルス変異予測装置。
  4. 　ＲＮＡ塩基Ａ（アデニン）、Ｕ、Ｇ、Ｃの４種類のうち、２塩基が選ばれて特徴づけられた量である特徴量であって、学習の際に用いられる前記特徴量を追加する特徴量追加選択部、を更に備え、
    　前記学習部は、前記特徴量も学習データに用いる、
    　請求項１から請求項３のうちのいずれか１項に記載のウイルス変異予測装置。
  5. 　前記コンテキストの範囲は、－３から＋３以上、－１０から＋１０以下である、
    　請求項１から請求項４のいずれか１項に記載のウイルス変異予測装置。
  6. 　前記ウイルスは、ＳＡＲＳ－ＣｏＶ－２である、
    　請求項１から請求項５のいずれか１項に記載のウイルス変異予測装置。
  7. 　取得部が、ウイルスのゲノムの遺伝子配列データを取得し、
    　抽出部が、取得された前記ゲノムの遺伝子配列データからＣ（シトシン）またはＧ（グアニン）を抽出し、ＣまたはＧからＵ（ウラシル）への変異が起こるまたは起こったコンテキストを抽出し、
    　分離部が、ＣまたはＧからＵに変わった時、アミノ酸変異があるかを確認し、前記アミノ酸変異がある配列を非同義置換として分離し、前記アミノ酸変異がない配列を同義置換であるとして分離し、
    　学習部が、前記同義置換の配列データを学習データに用いて学習し、
    　予測部が、学習された結果を用いて、前記ウイルスの変異を予測する、
    　ウイルス変異予測方法。
  8. 　取得部が、ウイルスのゲノムの遺伝子配列データを取得し、
    　抽出部が、取得した前記ゲノムの遺伝子配列データからＣ（シトシン）、Ｇ（グアニン）、Ａ（アデニン）、Ｕ（ウラシル）またはＴ（チミン）を抽出し、ＧからＡ、ＡからＧ、ＵからＣ、ＴからＣへの変異が起こるまたは起こったコンテキストを抽出し、
    　分離部が、抽出されたコンテキストの塩基配列が変化した場合、アミノ酸変異があるかを確認し、前記アミノ酸変異がある配列を非同義置換として分離し、前記アミノ酸変異がない配列を同義置換であるとして分離し、
    　学習部が、前記同義置換の配列データを学習データに用いて学習し、
    　予測部が、学習された結果を用いて、前記ウイルスの変異を予測する、
    　ウイルス変異予測方法。
  9. 　コンピュータに、
    　ウイルスのゲノムの遺伝子配列データを取得させ、
    　取得された前記ゲノムの遺伝子配列データからＣ（シトシン）またはＧ（グアニン）を抽出させ、ＣまたはＧからＵ（ウラシル）への変異が起こるまたは起こったコンテキストを抽出し、
    　分離部が、ＣまたはＧからＵに変わった時、アミノ酸変異があるかを確認させ、前記アミノ酸変異がある配列を非同義置換として分離させ、アミノ酸変異がない配列を同義置換であるとして分離させ、
    　前記同義置換の配列データを学習データに用いて学習させ、
    　学習された結果を用いて、前記ウイルスの変異を予測させる、
    　プログラム。
  10. 　コンピュータに、
    　ウイルスのゲノムの遺伝子配列データを取得させ、
    　取得された前記ゲノムの遺伝子配列データからＣ（シトシン）、Ｇ（グアニン）、Ａ（アデニン）、Ｕ（ウラシル）またはＴ（チミン）を抽出し、ＧからＡ、ＡからＧ、ＵからＣ、ＴからＣへの変異が起こるまたは起こったコンテキストを抽出させ、
    　抽出されたコンテキストの塩基配列が変化した場合、アミノ酸変異があるかを確認し、前記アミノ酸変異がある配列を非同義置換として分離させ、前記アミノ酸変異がない配列を同義置換であるとして分離させ、
    　前記同義置換の配列データを学習データに用いて学習させ、
    　学習された結果を用いて、前記ウイルスの変異を予測させる、
    　プログラム。
PCT/JP2021/027331 2020-07-22 2021-07-21 ウイルス変異予測装置、ウイルス変異予測方法、およびプログラム WO2022019331A1 (ja)

Priority Applications (3)

Application Number Priority Date Filing Date Title
US18/017,039 US20230298700A1 (en) 2020-07-22 2021-07-21 Device for predicting mutation of virus, method for predicting mutation of virus, and program
DE112021003912.1T DE112021003912T5 (de) 2020-07-22 2021-07-21 Vorrichtung zum prognostizieren einer mutation eines virus, verfahren zum prognostizieren einer mutation eines virus, und programm
JP2022538042A JPWO2022019331A1 (ja) 2020-07-22 2021-07-21

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
JP2020-125563 2020-07-22
JP2020125563 2020-07-22

Publications (1)

Publication Number Publication Date
WO2022019331A1 true WO2022019331A1 (ja) 2022-01-27

Family

ID=79729156

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2021/027331 WO2022019331A1 (ja) 2020-07-22 2021-07-21 ウイルス変異予測装置、ウイルス変異予測方法、およびプログラム

Country Status (5)

Country Link
US (1) US20230298700A1 (ja)
JP (1) JPWO2022019331A1 (ja)
DE (1) DE112021003912T5 (ja)
TW (1) TW202217830A (ja)
WO (1) WO2022019331A1 (ja)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2024171375A1 (ja) * 2023-02-16 2024-08-22 富士通株式会社 情報処理プログラム，情報処理方法および情報処理装置

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2020125563A1 (zh) 2018-12-20 2020-06-25 厦门凯浦瑞运动器材有限公司 一种简易型多功能训练器

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
MATYÁŠEK ROMAN, KOVAŘÍK ALEŠ: "Mutation Patterns of Human SARS-CoV-2 and Bat RaTG13 Coronavirus Genomes Are Strongly Biased Towards C>U Transitions, Indicating Rapid Evolution in Their Hosts", GENES, vol. 11, no. 7, pages 1 - 13, XP055890032, DOI: 10.3390/genes11070761 *
SALAMA MOSTAFA A., HASSANIEN ABOUL ELLA, MOSTAFA AHMAD: "The prediction of virus mutation using neural networks and rough set techniques", EURASIP JOURNAL ON BIOINFORMATICS AND SYSTEMS BIOLOGY, vol. 2016, no. 1, 1 December 2016 (2016-12-01), pages 10, XP055890036, DOI: 10.1186/s13637-016-0042-0 *
WRIGHT ERIK S., LAKDAWALA SEEMA S., COOPER VAUGHN S.: "SARS-CoV-2 genome evolution exposes early human adaptations", BIORXIV, 26 May 2020 (2020-05-26), pages 1 - 17, XP055890025, Retrieved from the Internet <URL:https://www.biorxiv.org/content/10.1101/2020.05.26.117069v1.full.pdf> [retrieved on 20220210], DOI: 10.1101/2020.05.26.117069 *

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2024171375A1 (ja) * 2023-02-16 2024-08-22 富士通株式会社 情報処理プログラム，情報処理方法および情報処理装置

Also Published As

Publication number Publication date
DE112021003912T5 (de) 2023-07-13
TW202217830A (zh) 2022-05-01
JPWO2022019331A1 (ja) 2022-01-27
US20230298700A1 (en) 2023-09-21

Similar Documents

Publication Publication Date Title
CN112133365B (zh) 评估肿瘤微环境的基因集、评分模型及其应用
Gregori et al. Viral quasispecies complexity measures
Gan et al. Omicron Spike protein has a positive electrostatic surface that promotes ACE2 recognition and antibody escape
CN111785328A (zh) 基于门控循环单元神经网络的冠状病毒序列识别方法
Barrie et al. Elevated genetic risk for multiple sclerosis emerged in steppe pastoralist populations
WO2022019331A1 (ja) ウイルス変異予測装置、ウイルス変異予測方法、およびプログラム
CN114686591A (zh) 基于基因表达情况的肺鳞癌免疫治疗疗效预测模型及其构建方法和应用
CN113903395A (zh) 一种改进粒子群优化的bp神经网络拷贝数变异检测方法与系统
CN112037863B (zh) 一种早期nsclc预后预测系统
Kumar et al. Genetic affinities and adaptation of the South-west coast populations of India
Souza et al. Detecting clustered independent rare variant associations using genetic algorithms
CN116486913A (zh) 基于单细胞测序从头预测调控突变的系统、设备和介质
Barrie et al. Elevated genetic risk for multiple sclerosis originated in Steppe Pastoralist populations
CN115458054A (zh) 基于sigma+预测突变致病性分析方法、系统及设备
CN115691666A (zh) 基于sigma预测突变致病性分析方法、系统及设备
WO2018223185A1 (en) Methods of determining the likelihood of hepatitis b virus recrudescence
JP2024501141A (ja) 遺伝子データを分析するためのコンピュータ実施方法および装置
KR20220103819A (ko) 개인의 생물학적 상태를 예측하기 위한 시스템, 방법 및 유전자 시그니처
Morales-Arce et al. Inferring the distribution of fitness effects in patient-sampled and experimental virus populations: two case studies
CN112397200A (zh) 一种非综合征型唇腭裂遗传风险预测模型
Hua et al. Combining protein-protein interactions information with support vector machine to identify chronic obstructive pulmonary disease related genes
NL2030705B1 (en) Method for establishing comparative transcriptomics database of animal models of coronavirus infections
CN117743957B (zh) 一种基于机器学习的Th2A细胞的数据分选方法及相关设备
Ismaeel et al. Enhancement of a Novel Method for Mutational Disease Prediction using Bioinformatics Techniques and Backpropagation Algorithm
Wong et al. Hybrid feature selection to identify important single nucleotide polymorphism for asthma classification

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 21846790

Country of ref document: EP

Kind code of ref document: A1

ENP Entry into the national phase

Ref document number: 2022538042

Country of ref document: JP

Kind code of ref document: A

122 Ep: pct application non-entry in european phase

Ref document number: 21846790

Country of ref document: EP

Kind code of ref document: A1