JP2010224815A - Retrieval algorithm for epistasis effect based on comprehensive genome wide snp information - Google Patents

Abstract

<P>PROBLEM TO BE SOLVED: To provide a method of identifying an SNP pair having synergistic epistasis effects within an actual time, even when main effects are not confirmed with respect to a binary phenotype by using genome side SNP data. <P>SOLUTION: The binary phenotype data of each specimen are input and stored (1). The genotype of M pieces of single nucleotide polymorphism (SNP) with respect to each specimen are input and stored (2). Genotype-categorized counting is calculated for every phenotype with respect to each SNP (3). The genotype-categorized counting is stored by determining the validity/invalidity of analysis continuation with respect to each SNP (4). Dominant/recessive type for the phenotype is determined and stored with respect to the SNP determined to be valid for analysis continuation (5). The epistasis effects are determined and stored based on a division table based on a division table based on the phenotype and the dominant/recessive type with respect to the two SNP determined to be valid for analysis continuation, and the SNP pair to be analyzed is changed according to an analytic procedure (6). The epistasis effects are verified by logistic regression analysis with respect to the SNP pair whose epistasis effects are determined (7). <P>COPYRIGHT: (C)2011,JPO&INPIT

Description

  The present invention provides data in which 1 million single nucleotide polymorphism (SNP) genotypes are observed on the genome wide for each subject having a label represented by a binary value as a phenotype. A method of rapidly identifying multiple SNPs that have an epistasis effect that synergistically affects the phenotype only when two SNPs are present simultaneously, although the individual SNPs do not affect the phenotype, and It relates to an identification program.

In various biological species, it is known that genes on the genome are involved in a phenotype that shows the biological characteristics of an individual. A single gene may act on the phenotype, but generally a plurality of genes may act on one phenotype. Epistasis has long been defined as a superordinate effect among non-additive gene effects as a gene action. Non-Patent Document 1: Bateson, Mendel's Principles of Heredity. See Cambridge University Press, Cambridge 1909. At present, it is regarded as “interaction between genes” and is an extremely important concept in clarifying the relationship between genes and phenotypes. Non-Patent Document 2: See Cordell. Hum Molecular Genet, Vol. 11, No. 20, 2463-2468, 2002. “Interaction” refers to the effect of individual genes acting on the phenotype independently. Is called synergistic epistasis, and if it is less than the effect of acting independently, it is called antagonistic epistasis. That is, when there is synergistic epistasis in a specific number of gene sets for a certain phenotype, the effect of the entire gene set is greater than the sum of the individual gene effects. It has been reported that as the genome becomes more complex, the effect of epistasis becomes synergistic, indicating the importance of considering epistasis in elucidating the relationship between the human genome and the phenotype. (See Non-Patent Document 3: Sanjuan and Elena, PNAS Vol.103, No.39, 14402-14405, 2006.)
The concept of epistasis has existed for a long time, but it may be difficult to actually search for the epistasis effect. An analysis method using a logistic model has also been proposed for analysis of a small number of gene sets. Non-Patent Document 4: See Cordell and Clayton. Am J Hum Genet, Vol. 70, No. 1, 124-141, 2002. There are penalized maximum likelihood estimation method and BSE method variable interval approach to detect epistasis effect when single gene has main effect on phenotype. See Non-Patent Document 5: Zhang, Shrinkage Estimation Method for Mapping Multiple Quantitative Trait Loci, Vol.33 No.10, Page.861-869, 2006. A method for analyzing the main effect at a locus and the interaction between two loci based on entropy has also been proposed. Non-Patent Document 6: See Dong et al. Eur J Hum Genet. Vol. 16, 229-235, 2008.

  However, when multiple genes are considered, each gene alone has no effect on the phenotype, but when multiple genes set exist and the effect appears only when there are multiple gene sets, the effect of each gene is searched and the result It is impossible to search for mutual effects of gene sets based on the above.

  Currently, a technology to comprehensively examine human genome mutations has been developed, and in order to identify human phenotypes, especially susceptibility to diseases and individual differences in side effects, single nucleotide polymorphism mutations ( Based on Single Nucleotide Polymorphizm), gene mutations have been examined genome-wide. Marchini et al. Pointed out that it was difficult to analyze genome data considering multiple loci as of 2005. Non-Patent Document 7: See Marchini et al. Nat. Genet. Vol. 37, No. 4, 413-417, 2005. In the latest technology (January 2009), about 900,000 SNPs have been examined for one patient. Until now, the concept of epistasis was “interaction between genes”, but when considering the effect of SNP, the concept of epistasis needs to be extended to “interaction between SNPs”. The smallest SNP set is the case of 2 SNPs, but there are about 500 billion combinations that take 2 SNPs from 900,000 SNPs, and in order to search for epistasis when there is no main effect as described above, About 500 billion streets must be examined.

Bateson, Mendel's Principles of Heredity.Cambridge University Press, Cambridge 1909. Cordell. Hum Molecular Genet, Vol. 11, No. 20, 2463-2468, 2002. Sanjuan and Elena, PNAS Vol.103, No.39, 14402-14405, 2006. Cordell and Clayton. Am J Hum Genet, Vol. 70, No. 1, 124-141, 2002. Zhang, Shrinkage Estimation Method for Mapping Multiple Quantitative Trait Loci, Vol.33 No.10, Page.861-869, 2006. Dong et al. Eur J Hum Genet. Vol.16, 229-235, 2008. Marchini et al. Nat. Genet. Vol. 37, No. 4, 413-417, 2005.

  The conventional method has the following problems.

  The problem with conventional methods is that until recently there has been no genome-wide technology for genotype detection, so even if the detection of the gene interaction effect between two specific genes or the number of genes to be examined increases. Since only the analysis of gene interaction effects on a limited number of genes, limited to a maximum of about 100, has been performed, the method for detecting gene interaction effects between two genes is an iterative calculation or complex search for solution search. Therefore, there is no method that can deal with a huge amount of gene data. In recent genome-wide analysis, about 900,000 SNP genotypes have been examined, and when examining epistasis, there are about 500 billion combinations that take 2 SNPs from 900,000 SNPs, and epistasis between two SNPs In order to explore, you have to investigate about 500 billion ways. Even if two SNPs are analyzed in 1 second, the calculation time of 15854 is required, which is practically impossible to analyze, and the conventional method has no main effect on genome-wide data. It is impossible to exhaustively search for epistasis between two SNPs.

  The object of the present invention is to complete the search for epistasis between two SNPs, even when there is no main effect on genotype data of about 1 million SNPs obtained by genome-wide analysis. To provide a high-speed identification method and a data analysis system.

  In the present invention, a total of M SNP genotype data observed from N specimens input via the input device (M is 500,000 or more) and phenotypic class data corresponding to each specimen, The result of determining superiority or inferiority based on class-specific genotype counts calculated from these data is stored in an internal storage device that can be accessed at high speed and referenced when necessary, eliminating unnecessary repeated calculations for the same SNP. ing.

  Furthermore, in the present invention, as a method for identifying an epistasis effect, a total of four types in a 2 × 2 contingency table configured by distinguishing genotypes by dominance and recessive for two SNP combinations by two phenotypes. Using summary numeric data, two odds ratio statistics that can be calculated for these four numeric data by a total of three computations are calculated and judged. In this way, by narrowing down information useful for epistasis identification and using statistics with extremely small amount of calculation, it is possible to construct a judgment method that can greatly reduce the calculation time, and comprehensive identification of epistasis effects Can be achieved in real time, and the object of the present invention can be achieved.

In one embodiment of the present invention, the main effect is not confirmed for a phenotype having a binary class from genome-wide single nucleotide polymorphism (SNP) genotype data of more than 500,000 sites using a computer. A data analysis system that comprehensively identifies pairs of SNPs that have synergistic interactions (epistasis effects),
(1) Input to input a total of M SNP genotype data (M is more than 500,000) observed from N specimens with two types of phenotypes and the phenotype class corresponding to each specimen. Means,
(2) storage means for storing a phenotype class of N specimens and a total of M genotype data inputted via the input means (1);
(3) The i-th SNP stored in the storage means (2) is subjected to statistical processing of the two phenotype classes and genotype data for the sample N persons, and the class-specific genotype count is calculated. An arithmetic means for performing screening as a preprocessing step for determining the suitability of the i-th SNP analysis continuation based on the calculated minor allele count,
(4) When the calculation means (3) determines that “analysis continuation is suitable” as the SNP to be analyzed, whether the i-th SNP is dominant or inferior to the phenotype based on the calculated count Storage means for determining by statistical means, and storing the determination result on the suitability of analysis continuity and the dominant type / recessive type in the internal storage device;
(5) For each of the two classes of phenotypes based on the determination result regarding the dominant type and the recessive type determined by the statistical means, the i-th and j-th (j ≠ i, i = 1, j = 2) 2x2 contingency table in which the dominant type and recessive type of each two SNPs are determined, and an index for determining epistasis is calculated for the created 2x2 contingency table. Computing means for determining the presence or absence of an epistasis effect,
(6) When it is determined by the calculation means (5) that “the epistasis effect is present”, the determination result of “the epistasis effect is present” for the two SNPs is stored, and the j th The SNP is changed to the j + 1-th SNP, the process returns to the step of the recording means (4), the dominant type / recessive type of the j + 1-th SNP is determined by statistical means, and the calculation means (3) When the calculation of the steps is repeated and j + 1 reaches M, the analysis means for selecting the i-th SNP as the i + 1-th and the j-th as the i + 2-th SNP,
(7) In the calculation means (5), when it is determined that “there is an epistasis effect”, calculation means for confirming the synergistic epistasis effect by a multivariate analysis means using logistic analysis analysis;
It is a data analysis system characterized by comprising.

Another aspect of the present invention is:
Genome-wide single nucleotide polymorphism (SNP) genotype data of more than 500,000 sites using a computer, synergistic interaction (epistasis) even if the main effect is not confirmed for a phenotype having a binary class A data analysis method for comprehensively identifying SNP pairs having an effect),
(1) Input to input a total of M SNP genotype data (M is more than 500,000) observed from N specimens with two types of phenotypes and the phenotype class corresponding to each specimen. Steps,
(2) a storage step of storing in the storage means the phenotype class of N specimens and the total M genotype data inputted through the input step (2);
(3) According to the storing step (2), statistical processing of the two phenotype classes and genotype data for the sample N is performed on the i-th SNP stored in the storage means, and the class-specific genotype A calculation step for calculating a separate count and performing screening as a pre-processing step for determining the suitability of i-th SNP analysis continuation based on the calculated separate count by minor allele,
(4) In the calculation step (3), if it is determined that “analysis continuation is suitable” as the analysis target SNP, whether the i-th SNP is dominant or inferior to the phenotype based on the calculated count A storage step of determining by statistical means, and storing the determination result regarding the suitability of continuation of analysis and the dominant / recessive type in an internal storage device;
(5) In step (4), the i-th and j-th (j ≠ i) for each of the two classes of phenotypes based on the determination result regarding the dominant type and the recessive type determined by the statistical means. Create a 2x2 contingency table in which the dominant and inferior types of two SNPs (i = 1, j = 2) are determined as initial values, and calculate an index for determining epistasis for the created 2x2 contingency table A calculation step for determining the presence or absence of an epistasis effect based on this index,
(6) When it is determined in the calculation step (5) that “the epistasis effect is present”, the determination result of “the epistasis effect is present” for the two SNPs is stored, and the j th The SNP of j + 1 is changed to the j + 1-th SNP, and the process returns to step (4). The dominant / recessive type of the j + 1-th SNP is determined by statistical means, and the calculation of step (3) is repeated, j When +1 reaches M, an analysis step of selecting the i-th SNP as the i + 1-th and the j-th as the i + 2-th SNP; (7) In the calculation step (5), “Epistasis” When it is determined that “there is an effect”, the data analysis method includes a calculation step of confirming a synergistic epistasis effect by a multivariate analysis unit using logistic analysis analysis.

  Another embodiment of the present invention is a program for causing a computer to execute the data analysis method according to the present invention. Specifically, a program source in which a program that causes a computer to execute a series of numerical computations of steps constituting the data analysis method according to the present invention is stored on a recording medium readable by the computer It has the form.

  The effect of the present invention is that it is possible to efficiently evaluate the possibility of having a “synergistic effect” due to the combination of two SNPs even for SNPs for which it is difficult to determine “the presence or absence of a main effect”. . In particular, even when the size of the “sample group” used for the evaluation (N = n1 + n2) is small, a candidate “SNP pair” having a possibility of having a “synergistic effect” can be effectively selected.

It is a flowchart of the data analysis method which shows the procedure in the data analysis method concerning this invention. It is a figure which shows the determination algorithm of the dominance and inferiority in each SNP used at step (4) in the data analysis method concerning this invention. It is a figure explaining the structure of the 2x2 contingency table according to two phenotype classes with respect to two SNP which considered the superiority / inferiority of each SNP created in step (4) in the data analysis method concerning this invention. . In the data analysis method concerning this invention, it is a figure explaining the "2x2 contingency table according to phenotype class" produced in consideration of the dominance and inferiority of each SNP in step (4). Configuration of 2 x 2 contingency tables by two phenotype classes for two SNPs taking into account the dominant / recessive type of each SNP (for dominant and dominant types) In the data analysis method concerning this invention, it is a figure explaining the "2x2 contingency table according to phenotype class" produced in consideration of the dominance and inferiority of each SNP in step (4). Configuration of 2 x 2 contingency tables by two phenotype classes for two SNPs considering the dominant and recessive types of each SNP 2 (for dominant and recessive types) In the data analysis method concerning this invention, it is a figure explaining the "2x2 contingency table according to phenotype class" produced in consideration of the dominance and inferiority of each SNP in step (4). Configuration of 2 x 2 contingency tables by two phenotype classes for two SNPs considering the dominant and recessive types of each SNP 3 (for recessive and dominant types) In the data analysis method concerning this invention, it is a figure explaining the "2x2 contingency table according to phenotype class" produced in consideration of the dominance and inferiority of each SNP in step (4). Configuration of 2 x 2 contingency tables by 2 phenotype classes for 2 SNPs considering the dominant and recessive types of each SNP 4 (in case of recessive and recessive types) It is a figure explaining embodiment of the data analysis method concerning this invention, and is a figure explaining concretely the example and algorithm of a synergistic epistasis when there is no main effect of each SNP. It is a figure explaining embodiment of the data analysis method concerning this invention, and the analysis result is concretely shown about the example of a result with high prediction capability in 211 sets extracted from the combination of about 125 billion SNP. FIG.

  The data analysis system of the present invention and a data analysis method that can be implemented by using the data analysis system will be described in detail.

  In the input step (1) of the data analysis method according to the present invention, two types of phenotypes for N specimens and a total of M SNP genotype data observed from each specimen are associated with each specimen. Enter.

  In the recording step (2), the phenotypes of the two classes for the N samples input in the input step (1) and the total M genotype data of M SNPs observed from each sample are It is stored in an accessible internal storage device.

In the calculation step (3),
Based on the two classes of phenotypes for the N samples stored in the storage device and the genotype data of a total of M SNPs observed from each sample, the samples for the i-th SNP Statistical processing of two phenotype classes and genotype data for N people is performed, and class-specific genotype counts are calculated by an arithmetic unit.

  Based on the minor allele count calculated by the arithmetic unit, screening is performed as a pre-processing step for determining the suitability of the i-th SNP analysis continuation.

  The genotype data of each SNP is classified into 3 types, 2 homozygotes and 1 heterozygote in the population, depending on the type of two bases derived from mother and father. Here, these two types of homozygotes are represented by AA and aa, and one type of heterozygote is represented by Aa.

  Furthermore, a11, a12, and a13 are the counts for genotypes AA, Aa, and aa in the phenotype 1 class, respectively, and a21, a22, and a23 are the genotypes AA, Aa, and aa in the phenotype 2 class, respectively. Count.

In the screening as a pre-processing step executed by the arithmetic device, whether or not the i-th SNP analysis is continued is determined by a11, a12, a13, a21, a22, a23 according to the following (I) to (IV) If any one of the above conditions is satisfied, it is determined that the analysis is not continued. SNPs that are determined to be “non-continuation of analysis” are excluded from the subsequent analysis.
(I) a11 + a12 ≤ 1 or a11 + a13 ≤ 1 or a12 + a13 ≤ 1 (Formula 1)
(II) a21 + a22 ≦ 1 or a21 + a23 ≦ 1 or a22 + a23 ≦ 1 (Formula 2)
(III) a11 = 0 and a23 = 0 (Formula 3)
(IV) a13 = 0 and a21 = 0 (Formula 4)
The conditions (I) and (II) specify that the number of specimens of two genotypes among three genotypes in each phenotype class is 0 or 1, and (I) and (II ) If any of the above conditions is met, when combined with other SNP genotypes, the number of specimens will be zero in two of the three other genotypes of other SNPs, This is a list of cases where the epistasis criteria to be described cannot be clearly satisfied.

  The condition (IV) is that the index values OR1 and OR2 used in “determination of dominant type or inferior type” described later are 0 when a13 = 0 and a21 = 0. This corresponds to a condition that makes the calculation impossible.

  On the other hand, the condition of (III) is that the index values OR1 and OR2 used in the “determination of dominant type or recessive type” described later are zero when a11 = 0 and a23 = 0. As a result, since OR1 = 0 and OR2 = 0, this corresponds to a condition in which reliable “dominance / inferiority determination” cannot be performed.

In the memory step (4), first,
In the calculation step (3), if the SNP to be analyzed is determined to be “suitable for continued analysis”, the i-th SNP is compared with the phenotype based on the calculated phenotypic class-specific genotype count. Whether it is dominant or recessive is determined by statistical means.

  The genotype data of each SNP is classified into 3 types, 2 homozygotes and 1 heterozygote in the population, depending on the type of two bases derived from mother and father. Here, these two types of homozygotes are represented by AA and aa, and one type of heterozygote is represented by Aa.

In the determination means, a11, a12, and a13 are respectively counted with genotypes AA, Aa, and aa in the phenotype 1 class, and a21, a22, and a23 are respectively genotypes AA, Aa, aa, in the phenotype 2 class. Let aa be the count. Here, OR1 and OR2 are defined as follows.
OR1 = (a11 + a12) x a23 / (a13 x (a21 + a22)) (Formula 14)
OR2 = a11 x (a22 + a23) / ((a12 + a13) x a21) (Formula 15)
If OR1 and OR2 are compared, and OR1 value is greater than or equal to OR2 (OR1 ≧ OR2), if you have alleles of genotypes AA and Aa, the first class phenotype (for example, with side effects) It is judged as dominant type (or type 1) in order to express that it is likely to occur, and when OR1 is less than OR2 (OR1 <OR2), it is judged as recessive type (type 2). The recessive type is calculated and stored for SNP from No. 1 to SNP of M.

In the calculation step (5),
Based on the “determination result regarding dominant type / recessive type” determined by the statistical means, for each of the two classes of phenotypes, the i-th and j-th (j ≠ i, i = 1 as an initial value) , j = 2), a 2 × 2 contingency table in which the dominance and inferiority of each two SNPs is determined, and the following indices for determining epistasis are calculated for the created 2 × 2 contingency table.
R1 = (x ₁₁ x ₂₂ ) / (x ₁₂ x ₂₁ ) ≧ w ₁ and R2 = (y ₁₁ y ₂₂ ) / (y ₁₂ y ₂₁ ) <1 / w ₂ (Formula 5)

That is, in the determination of “presence / absence of epistasis effect” in step (5),
About the combination of “i-th SNP and j-th SNP”, which is the object of determination,
Follow the procedure below to calculate R1 = (x ₁₁ x ₂₂ ) / (x ₁₂ x ₂₁ ) and R2 = (y ₁₁ y ₂₂ ) / (y ₁₂ y ₂₁ ) as indicators,
Based on the calculated index, the determination of “existence of epistasis effect”
Index: R1 = (x ₁₁ x ₂₂ ) / (x ₁₂ x ₂₁ ) and R2 = (y ₁₁ y ₂₂ ) / (y ₁₂ y ₂₁ )
R1 = (x ₁₁ x ₂₂ ) / (x ₁₂ x ₂₁ ) ≥ w ₁ and R2 = (y ₁₁ y ₂₂ ) / (y ₁₂ y ₂₁ ) ≤ 1 / w ₂ (Formula 5)
When the above (Formula 5) is satisfied,
It is determined that there is an epistasis effect.

here,
x ₁₁ i-th and j th number of samples of SNP are both predominant form phenotypically class 1,
x ₁₂ is the i-th and the SNP sample number of the j-th SNP is recessive form in predominant form,
x ₂₁ is the number of specimens where the i-th and SNP are recessive and the j-th SNP is dominant,
x ₂₂ is the number of samples in which the i-th and j-th SNPs are both recessive,
y ₁₁ is the number of specimens with phenotype class 2 and i-th and j-th SNPs both dominant
y ₁₂ is the i-th and the SNP sample number of the j-th SNP is recessive form in predominant form,
y ₂₁ is the number of specimens where the i th and SNP are recessive and the j th SNP is dominant,
y ₂₂ is the number of samples in which the i-th and j-th SNPs are both recessive.

In the above (Formula 5), x ₁₁ , x ₂₂ , x ₁₂ , x ₂₁ , y ₁₁ , y ₂₂ , y ₁₂ , y ₂₁ are calculated according to the following procedure.

In the above (Formula 5), x ₁₁ , x ₂₂ , x ₁₂ , and x ₂₁ are determined by the combination of the dominant type and the recessive type of the i-th SNP and the j-th SNP in the phenotype. It is a count.

Similarly, y ₁₁ , y ₂₂ , y ₁₂ , and y ₂₁ are counts determined by a combination of the dominant type and the recessive type of the i-th SNP and the j-th SNP in the phenotype.

The dominant type is a model in which the AA and Aa genotypes are related to the phenotype class 1, and are described as A1 = (AA, Aa) and A2 = (aa). The recessive type is a model in which the genotype of aa is related to the phenotype class 1, and is described as A1 = (AA) and A2 = (Aa, aa). The dominant type of the jth SNP is a model in which BB and Bb genotypes are related, and is described as B1 = (BB, Bb), B2 = (bb). The j-th recessive type of the SNP is a model in which bb genotype is related to phenotype class 1, and is described as B1 = (BB), B2 = (Bb, bb). (See Figure 3)
At this time, for specimens having a phenotype of class 1, c11 is a count of specimens having genotype AA of i-th SNP and genotype BB of j-th SNP (when AA and BB are included). Yes, c12 is the count of the specimen having AA and Bb, and c13 is the count of the specimen having AA and bb. Similarly, c21 is Aa and BB, c22 is Aa and Bb, c23 is Aa and bb, c31 is aa and BB, c32 is aa and Bb, and c33 is aa and bb. These counts satisfy the following formula:
c11 + c12 + c13 + c21 + c22 + c23 + c31 + c32 + c33 = n1 (Formula 6)
For specimens with phenotype class 2, d11 is the count of specimens with genotype AA of i-th SNP and genotype BB of j-th SNP (if it has AA and BB), d12 Is the count of specimens with AA and Bb, and d13 is the count of specimens with AA and bb.

Similarly, d21 is Aa and BB, d22 is Aa and Bb, d23 is Aa and bb, d31 is aa and BB, d32 is aa and Bb, and d33 is aa and bb. These counts satisfy the following formula:
d11 + d12 + d13 + d21 + d22 + d23 + d31 + d32 + d33 = n2 (Formula 7)
Based on the determination results for dominant and recessive types,
Specifically, x ₁₁ , x ₂₂ , x ₁₂ , x ₂₁ , y ₁₁ , y ₂₂ , y ₁₂ , y ₂₁ are given in the following cases.
(i) The i-th SNP is the dominant type and the j-th SNP is the dominant type (see Fig. 4)
x ₁₁ = c11 + c12 + c21 + c22, x ₁₂ = c13 + c23, x ₂₁ = c31 + c32, x ₂₂ = c33,
y ₁₁ = d11 + d12 + d21 + d22, y ₁₂ = d13 + d23, y ₂₁ = d31 + d32, y ₂ 2 = d33 (Formula 8)
(ii) The i-th SNP is the dominant type and the j-th SNP is the recessive type (see Fig. 5)
x ₁₁ = c11 + c21, x ₁₂ = c12 + c13 + c22 + c23, x ₂₁ = c31, x ₂₂ = c32 + c33,
y ₁₁ = d11 + d21, y ₁₂ = d12 + d13 + d22 + d23, y ₂₁ = d31, y ₂₂ = d32 + d33 (Formula 9)
(iii) The i-th SNP is recessive and the j-th SNP is dominant (see Fig. 6)
x ₁₁ = c11 + c12, x ₁₂ = c13, x ₂₁ = c21 + c22 + c31 + c32, x ₂₂ = c23 + c33,
y ₁₁ = d11 + d12, y ₁₂ = d13, y ₂₁ = d21 + d22 + d31 + d32, y ₂₂ = d23 + d33 (Formula 10)
(iv) The i-th SNP is recessive and the j-th SNP is recessive (see Fig. 7)
x ₁₁ = c11, x ₁₂ = c12 + c13, x ₂₁ = c21 + c31, x ₂₂ = c22 + c23 + c32 + c33
y ₁₁ = d11, y ₁₂ = d12 + d13, y ₂₁ = d21 + d31, y ₂₂ = d22 + d23 + d32 + d33 (Formula 11)
X ₁₁ , x ₂₂ , x ₁₂ , x ₂₁ , y ₁₁ , y ₂₂ , y ₁₂ , y ₂₁ given in any of the above (i) to (iv) based on the determination result regarding the dominant type or the recessive type From these, the indices: (x ₁₁ x ₂₂ ) / (x ₁₂ x ₂₁ ) and (y ₁₁ y ₂₂ ) / (y ₁₂ y ₂₁ ) are calculated.

Indicators: (x ₁₁ x ₂₂ ) / (x ₁₂ x ₂₁ ) and (y ₁₁ y ₂₂ ) / (y ₁₂ y ₂₁ ) are valid using two SNPs that consider phenotype, dominance and inferiority. used in log-linear model of 2x2x2 contingency table, the maximum likelihood of under the hypothesis that there is no 3 factor interactions _{(x 11, x 22, x} 12, x 21, y 11, y 22, y 12, y 21) The estimator (z ₁₁ , z ₂₂ , z ₁₂ , z ₂₁ , v ₁₁ , v ₂₂ , v ₁₂ , v ₂₁ ) is based on satisfying the following formula.
Log (z ₁₁ z ₂₂ ) / (z ₁₂ z ₂₁ )? log (v ₁₁ v ₂₂ ) / (v ₁₂ v ₂₁ ) = 0 (Formula 16)

In the index (Equation 5) for determining epistasis in the above means, the “selectable range” of w ₁ and w ₂ is given under the following conditions.
n1-3 ≤ w ₁ ≤ (n1 / 2-1) ² , (n2-3) ≤ w ₂ ≤ (n2 / 2-1) ² (Formula 12)
Here, w ₁ is the cross product ratio (x ₁₁ x ₂₂ ) / (x ₁₂ x ₂₁ ) for n1 specimens of class 1, and the minimum value of (x ₁₂ x ₂₁ ) is x ₁₂ = 1, x ₂₁ = 1. Under this condition, when x ₁₁ = 1 or x ₂₂ = 1, the minimum value of the cross product ratio (x ₁₁ x ₂₂ ) / (x ₁₂ x ₂₁ ) is a value of (n1-3).

Also, under the condition of x ₁₂ = 1 and x ₂₁ = 1 where (x ₁₂ x ₂₁ ) is the minimum value, x ₁₁ + x ₂₂ = n1-2, so the maximum that (x ₁₁ x ₂₂ ) can take Considering the maximum value of f (x ₁₁ , x ₂₂ ) = x ₁₁ x ₂₂ = x ₁₁ (n1-2-x ₁₁ ), {(n1 / 2-1) ² } is obtained. Range of w ₂ in the same manner can be calculated.

Of the above “selectable ranges” of w ₁ and w ₂ , w ₁ = n1-3 and w ₂ = n2-3 are the “slowest conditions” and w ₁ = (n1 / 2− 1) ² and w ₂ = (n2 / 2-1) ² correspond to “the most severe conditions”. In the present invention, when selecting a “relaxed condition”, w ₁ and w ₂ are selected within the following range, for example.

_{n1-3 ≦ w 1 ≦ (n1 +} √n1) -3, n2-3 ≦ w 2 ≦ (n2 + √n2) -3 ( Formula 17)
The “selectable range” of w ₁ and w ₂ is given by (Equation 12) above, but the denominator of the cross product ratio is not always the minimum, that is, x ₁₂ = 1 and x ₂₁ = 1 may not be satisfied. using the minimum value of w ₁ and a maximum value of 1 / w _2. That is, w ₁ = n1-3 and w ₂ = n2-3.

Here, since w ₁ ≧ 0 and w ₂ ≧ 0, the condition of the number of samples in each class of phenotype is given by n1 ≧ 3 and n2 ≧ 3. The theoretical lower limit of the number of samples (n1 + n2) in the sample group is (n1 + n2) ≧ 6. Considering the ratio of class 1 and class 2 in the population, the estimated values of class 1 and class 2 are given by n1 / (n1 + n2) and n2 / (n1 + n2), respectively. The lower limit of the number (n1 + n2) is (n1 + n2) ≧ 3 [1+ {n1 / (n1 + n2)} / {n2 / (n1 + n2)}] when n1> n2, and (n1 + n2) ≧ 3 [1+ {when n1 <n2. n2 / (n1 + n2)} / {n1 / (n1 + n2)}].
The method for determining that there is a possibility of “main effect” is as follows.

For SNP1, the estimated risk for “Class 1” is
The overall average risk r ₁ is (a11 + a12 + a13) / (a11 + a12 + a13 + a21 + a22 + a23),
The risk r _1AA in AA is a11 / (a11 + a21),
The risk r _1Aa in Aa is a12 / (a12 + a22),
The risk r _1aa in Aa is a13 / (a13 + a23).
In that case, for example, r _1aa > r _1Aa > r ₁ > r _1AA , that is, a relation of [a13 / (a13 + a23)]> [a12 / (a12 + a22)]> [a11 / (a11 + a21)] If there is, it is determined that there is a possibility of “main effect”.

Also, (r _1aa / r _1Aa )> (r _1Aa / r _1AA ), that is, [a13 / (a13 + a23)] / [a12 / (a12 + a22)]> [a12 / (a12 + a22)] / If there is a relationship of [a11 / (a11 + a21)], it is judged that there is a high possibility of a clear “main effect”.

In the analysis step (6),
When it is determined in the calculation step (5) that “there is an epistasis effect”, the results for the two SNPs are stored, and when moving to the analysis of the next SNP, the j-th SNP is changed to the j + 1-th SNP. Change to step (4), determine the dominant / inferiority of the j + 1-th SNP by statistical means, repeat the calculation in step (3), and if j + 1 reaches M, i The i th SNP is selected as the i th and the i th +2 th SNP is selected and analyzed.

In the calculation step (7),
When it is determined in the calculation step (5) that there is an epistasis effect, the synergistic epistasis effect is confirmed by multivariate analysis means using logistic analysis analysis.

The epistasis effect fast identification method in the absence of main effects based on comprehensive genome-wide SNP information according to the embodiment of the present invention will be described below with reference to the drawings. In the following, an example is given in which 60 patients with breast cancer are given preoperative chemotherapy with taxol alone to identify two SNPs that exhibit synergistic epistasis effects related to the occurrence of side effects of peripheral neuropathy. explain.

  In FIG. 1, a series of data analysis realized by the cooperation of computer hardware and a computer program as described above is illustrated in a flowchart format.

  In the following examples, the present inventors applied the data analysis method of the present invention to a gene search having a synergistic epistasis effect related to side effects after performing preoperative chemotherapy in a breast cancer patient with informed consent. The data analysis method of the present invention was confirmed to be effective. When statistical verification by logistic regression analysis was performed on the identified SNP data, an epistasis effect was observed. As a pioneering invention, the present invention can identify SNP pairs exhibiting a synergistic epistasis effect in real time even when there is no main effect even for SNP data of 1 million locations.

  Next, an example of the present invention will be described in detail with reference to the results. Such an example corresponds to an example of an embodiment of the present invention. The technical scope of the present invention is not limited to the specific modes exemplified in the examples.

  In this example, an SNP having an epistasis effect on the side effects of an anticancer drug was identified based on the data of cancer patients who had undergone preoperative chemotherapy for an anticancer drug with informed consent. The purpose of this study is to investigate genes that occur, and to identify the SNPs that have epistasis effects related to side effects. In recent years, human single nucleotide polymorphisms have been typed in large quantities. For example, in the Affymetrix SNP 6.0 array (registered trademark), allele-specific hybridization using a DNA chip is performed as a typing method. Typing of 906,600 SNPs covering the entire genome in 5 days per sample is possible.

1. Materials and Methods The specimens used for the analysis were 60 cases of breast cancer that were treated with preoperative chemotherapy (taxol 80 mg / m ² / q1w) of taxol alone with informed consent at the Cancer Research Society Breast Surgery.

  Side effects were examined for peripheral neuropathy (numbness). Eight people have peripheral neuropathy (CTC grade 2 or higher) and 52 have no peripheral neuropathy (CTC grade 0 or 1).

  The average call rate, which is the proportion of SNPs that were successfully typed, was 99.5%, and the coincidence rate of typing results with duplicate samples was 99.98%. 909622 SNP data were divided into groups with and without peripheral neuropathy, and the difference in allele frequency between the two groups was tested by Fisher's exact test using a 2 × 2 contingency table. SNPs (p ≦ 0.0001) that showed a correlation with CTC grade 2 or higher peripheral neuropathy were 33 SNPs and 17 genes. These results are analysis for each SNP and can be easily calculated by a conventional method.

  Next, up to 500,000 SNPs were analyzed as combinations of two SNPs using the data analysis method and data analysis system according to the present invention. There are 124999750,000 ways to select two SNPs from 500,000 SNPs.

  Of the 909622 SNP data, 409169 SNPs satisfied any one of the conditions (i) to (iv), and it was judged as “analysis continuation failure”. This made it possible to delete a large calculation load.

n1 = 8 and n2 = 52, and w ₁ = n1−3, w ₂ = n2 corresponding to “the mildest condition” among the “selectable range” of w ₁ and w ₂ described above. -3 was used.

Assuming W ₁ = (8-3) = 5 and W ₂ = (52-3) = 49 and considering all combinations to select two SNPs from 500,000 SNPs, 211 SNP pairs were selected It was done.

When it is determined that there is a possibility of “main effect”,
H1: [a13 / (a13 + a23)]> [a12 / (a12 + a22)]> [a11 / (a11 + a21)]
age,
When it is judged that there is a high possibility of a clear “main effect”
H2: [a13 / (a13 + a23)] / [a12 / (a12 + a22)]> [a12 / (a12 + a22)] / [a11 / (a11 + a21)]
And
Each of 211 SNP pairs is represented by SNP-A and SNP-B,
When H1 is established in SNP-A, AH1 = 1, otherwise AH1 = 0.
BH1 = 1 when H1 is established in SNP-B, and BH1 = 0 when not established.

  The results of combining AH1 and BH1 for 211 SNP pairs are summarized in Table 2 as shown in Table 1.

  As shown in Table 1, the number judged to be “mainly effective” was 3 pairs, only 1.42% of the total. In both cases, 84 pairs were judged not to have the “main effect”, which was 39.8% of the total. The number judged that only one of them may have “main effect” was 124 pairs, which was 58.7% of the total.

Therefore, out of 211 cases obtained in this analysis, 98.6% were pairs in which one SNP was judged not to have a “main effect” possibility.
In the conventional analysis method, it was difficult to identify 84 pairs (about 40%) that were judged not to have a “main effect”, indicating the usefulness of this method.

  Similarly, the results of combining AH1 and BH1 for 211 SNP pairs are summarized in Table 2 as shown in Table 2.

  As shown in Table 2, the number judged to have a clear “main effect” was 17 pairs, which was 3.31% of the total. In both cases, 57 pairs were judged to have no clear “main effect” possibility, which was 27.0% of the total. The number judged that only one of them may have a “main effect” was 137 pairs, which was 64.9% of the total.

  Therefore, out of 211 cases obtained in this analysis, 96.7% were pairs in which one SNP was judged not to be “mainly effective”.

  In the conventional analysis method, it was difficult to identify 57 pairs (27%), both of which were judged to have the “no main effect”, indicating the usefulness of this method.

  Of these 211 SNP pairs, the main effect is not recognized, the synergistic epistasis effect is recognized, and the SNP pair with high prediction ability is shown in FIG. In FIG. 9, SNP-A represents the 69146th SNP out of 909622, and SNP-B represents the 97440th SNP. SNP-A1 and SNP-B1, SNP-A1 and SNP-B2, NP-A2 and SNP-B1 combination is not high risk, but only with SNP-A2 and SNP-B2 combination Was found to have reached 0.8 or higher. This high-risk group accounts for about 10% of the Japanese population when calculated from the allele frequency, and the combination of SNP-A2 and SNP-B2 has a high probability of side effects, so the use of this drug should be carefully considered It shows the results that can contribute to the progress of personalized medicine in the future.

  As a result of examining all combinations for selecting two SNPs from 500,000 SNPs, 211 SNP pairs selected are shown in Tables 3-1 to 3-211. Each table of Table 3-1 to Table 3-211 includes the selected SNP number, the values of R1 and R2 determined to have an epistasis effect according to (Equation 5), and two genos of groups with and without side effects. A 3 × 3 table by type and two 2 × 2 tables of the group with and without side effects, which are summarized into dominant and recessive types, are described.

Example of verification by logistic regression analysis The usage data are the 69146th SNP and the 97440th SNP in Table 3.

  In the logistic regression analysis for verification, the interaction term SNP12 = SNP69146 / SNP97440 of 69146th SNP69146 and 97440th SNP97440 out of 909622 is used as a variable for the model, and the model with intercept and SNP12 as variables As a result of estimating the regression variables based on the maximum likelihood method, the intercept is -2.833 (standard error 0.59), the regression coefficient of the interaction term SNP12 is 4.442 (standard error 1.246, the lower and upper limits of the 95% confidence interval are (1.948: 6.937)), and showed a statistically significant result with a significance probability p = 0.00074. The fitness of the logistic regression model was also statistically significant at p = 0.00002.

  In logistic regression analysis for verification, the variables used in the model were the 69146th SNP69146 variable, the 97440th SNP97440 variable, and the interaction term SNP12, and these regression variables were estimated simultaneously based on the maximum likelihood method. As a result, the intercept is -2.079 (standard error 1.06), the regression coefficient representing the main effect of SNP69146 is -1.252 (standard error 1.470, the lower and upper limits of the 95% confidence interval are (-4.197: 6.937)) and significant There was no statistically significant difference with probability p = 0.398. The regression coefficient representing the main effect of SNP97440 is -0.629 (standard error 1.480, the lower and upper limits of the 95% confidence interval are (-3.594: 2.337)), and a statistically significant difference is observed with significance probability p = 0.672. I couldn't. The regression coefficient of the interaction term SNP12 was 5.570 (standard error 2.104, 95% confidence interval (1.355: 9.786)), and a statistically significant result was shown with a significance probability p = 0.011. The fitness of the logistic regression model was also statistically significant at p = 0.00024.

Furthermore, an example of the “analysis program” of the present invention is shown below. The following part is a program example corresponding to the process of steps (4) to (6), excluding the determination step of “appropriateness of continuation of analysis” as an example of a main part in the program of the present invention. The following program has made it possible to search for SNP pairs that have a synergistic epistasis effect even when there is no main effect on the side effects of SNP.

An example of “analysis program”:
program SNP
integer IG * 4
dimension IDAT (1000000,60), Adata (100000)
dimension IT (3,3,2), IS (60)
character Adata * 20
OPEN (UNIT = 1, FILE = 'D: \ ptxPNP # Aold.txt')
OPEN (UNIT = 2, FILE = 'D: \ ptxPNP-out.txt')
OPEN (UNIT = 3, FILE = 'D: \ ptxPNP-out-begin.txt')

write (3,300)
300 FORMAT (1H, 'Start')
CLOSE (UNIT = 3)

NN = 60
DO 500 I = 1,52
IS (I) = 0
500 CONTINUE
DO 505 I = 53,60
IS (I) = 1
505 CONTINUE

c ************************************************* *
read (1, *) Adata (IG)
read (1, *) Adata (IG)
read (1, *) n1, n2
IG = 1
5 continue @@
read (1, *, end = 99) Adata (IG), (IDAT (IG, K), K = 1, NN)
200 FORMAT (1H, F5.3,44I4, A20)

IG = IG + 1
c IF (IG.GT.500000) GOTO 99
GOTO 5
99 CONTINUE
IGEND = IG-1

DO 10 IG = 1, IGEND-1
DO 20 JG = IG + 1, IGEND

DO 40 kk = 1,2
DO 40 J = 1,3
DO 40 I = 1,3
IT (I, J, kk) = 0
40 CONTINUE

DO 30 K = 1, NN
IF (IDAT (IG, K) .EQ.-10.OR.IDAT (JG, K) .EQ.-10) GOTO 30
IF (IS (K) .EQ.0) GOTO 33
c AE (+): IS (K) = 1
IF (IDAT (IG, K) .EQ.0.and.IDAT (JG, K) .EQ.0) IT (1,1,1) = IT (1,1,1) +1
IF (IDAT (IG, K) .EQ.0.and.IDAT (JG, K) .EQ.1) IT (1,2,1) = IT (1,2,1) +1
IF (IDAT (IG, K) .EQ.0.and.IDAT (JG, K) .EQ.2) IT (1,3,1) = IT (1,3,1) +1
IF (IDAT (IG, K) .EQ.1.and.IDAT (JG, K) .EQ.0) IT (2,1,1) = IT (2,1,1) +1
IF (IDAT (IG, K) .EQ.1.and.IDAT (JG, K) .EQ.1) IT (2,2,1) = IT (2,2,1) +1
IF (IDAT (IG, K) .EQ.1.and.IDAT (JG, K) .EQ.2) IT (2,3,1) = IT (2,3,1) +1
IF (IDAT (IG, K) .EQ.2.and.IDAT (JG, K) .EQ.0) IT (3,1,1) = IT (3,1,1) +1
IF (IDAT (IG, K) .EQ.2.and.IDAT (JG, K) .EQ.1) IT (3,2,1) = IT (3,2,1) +1
IF (IDAT (IG, K) .EQ.2.and.IDAT (JG, K) .EQ.2) IT (3,3,1) = IT (3,3,1) +1
GOTO 30
c AE (-): IS (K) = 0
33 CONTINUE
IF (IDAT (IG, K) .EQ.0.and.IDAT (JG, K) .EQ.0) IT (1,1,2) = IT (1,1,2) +1
IF (IDAT (IG, K) .EQ.0.and.IDAT (JG, K) .EQ.1) IT (1,2,2) = IT (1,2,2) +1
IF (IDAT (IG, K) .EQ.0.and.IDAT (JG, K) .EQ.2) IT (1,3,2) = IT (1,3,2) +1
IF (IDAT (IG, K) .EQ.1.and.IDAT (JG, K) .EQ.0) IT (2,1,2) = IT (2,1,2) +1
IF (IDAT (IG, K) .EQ.1.and.IDAT (JG, K) .EQ.1) IT (2,2,2) = IT (2,2,2) +1
IF (IDAT (IG, K) .EQ.1.and.IDAT (JG, K) .EQ.2) IT (2,3,2) = IT (2,3,2) +1
IF (IDAT (IG, K) .EQ.2.and.IDAT (JG, K) .EQ.0) IT (3,1,2) = IT (3,1,2) +1
IF (IDAT (IG, K) .EQ.2.and.IDAT (JG, K) .EQ.1) IT (3,2,2) = IT (3,2,2) +1
IF (IDAT (IG, K) .EQ.2.and.IDAT (JG, K) .EQ.2) IT (3,3,2) = IT (3,3,2) +1
30 CONTINUE

ISNP1A = IT (1,1,1) + IT (1,2,1) + IT (1,3,1)
ISNP2A = IT (2,1,1) + IT (2,2,1) + IT (2,3,1)
ISNP3A = IT (3,1,1) + IT (3,2,1) + IT (3,3,1)

ISNP1B = IT (1,1,2) + IT (1,2,2) + IT (1,3,2)
ISNP2B = IT (2,1,2) + IT (2,2,2) + IT (2,3,2)
ISNP3B = IT (3,1,2) + IT (3,2,2) + IT (3,3,2)

JSNP1A = IT (1,1,1) + IT (2,1,1) + IT (3,1,1)
JSNP2A = IT (1,2,1) + IT (2,2,1) + IT (3,2,1)
JSNP3A = IT (1,3,1) + IT (2,3,1) + IT (3,3,1)
@
JSNP1B = IT (1,1,2) + IT (2,1,2) + IT (3,1,2)
JSNP2B = IT (1,2,2) + IT (2,2,2) + IT (3,2,2)
JSNP3B = IT (1,3,2) + IT (2,3,2) + IT (3,3,2)


C For SNP1
OR1 = FLOAT ((ISNP1A + ISNP2A) * ISNP3B) /
+ ((FLOAT (ISNP3A) +0.1) * (FLOAT (ISNP1B + ISNP2B) +0.1))
OR2 = FLOAT (ISNP1A * (ISNP2B + ISNP3B)) /
+ ((FLOAT (ISNP2A + ISNP3A) +0.1) * (FLOAT (ISNP1B) +0.1))

IF (OR1.GE.OR2) THEN
Itype = 1
ELSE
Itype = 2
ENDIF
35 CONTINUE
OR3 = FLOAT ((JSNP1A + JSNP2A) * JSNP3B) /
+ ((FLOAT (JSNP3A) +0.1) * (FLOAT (JSNP1B + JSNP2B) +0.1))
OR4 = FLOAT (JSNP1A * (JSNP2B + JSNP3B)) /
+ ((FLOAT (JSNP2A + JSNP3A) +0.1) * (FLOAT (JSNP1B) +0.1))

IF (OR3.GE.OR4) THEN
Jtype = 1
ELSE
Jtype = 2
ENDIF

37 CONTINUE
IF (Itype.EQ.1.and.Jtype.EQ.1) THEN
X11 = FLOAT (IT (1,1,1) + IT (1,2,1) + IT (2,1,1) + IT (2,2,1))
X12 = FLOAT (IT (1,3,1) + IT (2,3,1))
X21 = FLOAT (IT (3,1,1) + IT (3,2,1))
X22 = FLOAT (IT (3,3,1))
Y11 = FLOAT (IT (1,1,2) + IT (1,2,2) + IT (2,1,2) + IT (2,2,2))
Y12 = FLOAT (IT (1,3,2) + IT (2,3,2))
Y21 = FLOAT (IT (3,1,2) + IT (3,2,2))
Y22 = FLOAT (IT (3,3,2))
GOTO 77
ENDIF
IF (Itype.EQ.1.and.Jtype.EQ.2) THEN
X11 = FLOAT (IT (1,1,1) + IT (2,1,1))
X12 = FLOAT (IT (1,2,1) + IT (1,3,1) + IT (2,2,1) + IT (2,3,1))
X21 = FLOAT (IT (3,1,1))
X22 = FLOAT (IT (3,2,1) + IT (3,3,1))
Y11 = FLOAT (IT (1,1,2) + IT (2,1,2))
Y12 = FLOAT (IT (1,2,2) + IT (1,3,2) + IT (2,2,2) + IT (2,3,2))
Y21 = FLOAT (IT (3,1,2))
Y22 = FLOAT (IT (3,2,2) + IT (3,3,2))
GOTO 77
ENDIF

IF (Itype.EQ.2.and.Jtype.EQ.1) THEN
X11 = FLOAT (IT (1,1,1) + IT (1,2,1))
X12 = FLOAT (IT (1,3,1))
X21 = FLOAT (IT (2,1,1) + IT (2,2,1) + IT (3,1,1) + IT (3,2,1))
X22 = FLOAT (IT (2,3,1) + IT (3,3,1))
Y11 = FLOAT (IT (1,1,2) + IT (1,2,2))
Y12 = FLOAT (IT (1,3,2))
Y21 = FLOAT (IT (2,1,2) + IT (2,2,2) + IT (3,1,2) + IT (3,2,2))
Y22 = FLOAT (IT (2,3,2) + IT (3,3,2))
GOTO 77
ENDIF

IF (Itype.EQ.2.and.Jtype.EQ.2) THEN
X11 = FLOAT (IT (1,1,1))
X12 = FLOAT (IT (1,2,1) + IT (1,3,1))
X21 = FLOAT (IT (2,1,1) + IT (3,1,1))
X22 = FLOAT (IT (2,2,1) + IT (2,3,1) + IT (3,2,1) + IT (3,3,1))
Y11 = FLOAT (IT (1,1,2))
Y12 = FLOAT (IT (1,2,2) + IT (1,3,2))
Y21 = FLOAT (IT (2,1,2) + IT (3,1,2))
Y22 = FLOAT (IT (2,2,2) + IT (2,3,2) + IT (3,2,2) + IT (3,3,2))
ELSE
write (2, "Itype error exist")
ENDIF
77 CONTINUE

ICHECK = 0
RU = 50.0
RD = 1.0 / 50.0
X = FLOAT (ISNP1A + ISNP3A)
Y = FLOAT (ISNP1B + ISNP3B)

IF (X11 * X12 * X21 * X22.EQ.0.0) GO TO 20
IF (Y11 * Y12 * Y21 * Y22.EQ.0.0) GO TO 20

w1 = n1-3
w2 = n2-3
IF (X11 * X22 / (X12 * X21) .GE.w1
+ .and.Y11 * Y22 / (Y12 * Y21) .LT.w2)
+ ICHECK = 1

IF (ICHECK.EQ.1) NC = NC + 1
c
IF (ICHECK.EQ.1) THEN
WRITE (2,699) IG, JG
699 FORMAT (1H, 2I8)

WRITE (2,333) (X11 + 1.0) * (X22 + 1.0) / ((X12 + 1.0) * (X21 + 1.0)),
+ (Y11 + 1.0) * (Y22 + 1.0) / ((Y12 + 1.0) * (Y21 + 1.0))
333 FORMAT (1H, "Check", 2F8.3)

write (2,1000) OR1, OR2, Itype
write (2,1010) OR3, OR4, Jtype
1000 FORMAT (1H, 'OR, Itype', 2F6.1,2x, I3)
1010 FORMAT (1H, 'OR, Jtype', 2F6.1,2x, I3)

DO 70 I = 1,3
WRITE (2,700) (IT (I, J, 1), J = 1,3), (IT (I, J, 2), J = 1,3)
700 FORMAT (1H, 3I3, 3x, 3I3)
WRITE (2,698)
698 FORMAT (1H, '')
70 CONTINUE

WRITE (2,702) X11, X12, Y11, Y12
WRITE (2,702) X21, X22, Y21, Y22
702 FORMAT (1H, 2F3.0, 3x, 2F3.0)
WRITE (2,698)
END IF

20 CONTINUE
10 CONTINUE

write (2,15) IG
write (6,15) IG
15 FORMAT (1H, 'READ DATA =', I10)

write (2,25) NC
write (6,25) NC
25 FORMAT (1H, 'WRITE DATA =', I15)

CLOSE (UNIT = 1)
CLOSE (UNIT = 2)

stop
end

  As shown in the above examples, the data analysis method of the present invention identifies SNP combinations that have a high probability of side effects that cannot be detected by a single SNP, for example, considering the occurrence of side effects of a drug as a phenotype Therefore, the use of such SNP combinations should be carefully considered and can contribute to the development of personalized medicine in the future. In addition, there is a high possibility of giving knowledge to explore the mechanism of occurrence of side effects, and it has high applicability in that it can contribute to the progress of genomic science.

Claims (6)

Genome-wide single nucleotide polymorphism (SNP) genotype data of more than 500,000 sites using a computer, synergistic interaction (epistasis) even if the main effect is not confirmed for a phenotype having a binary class A data analysis method for comprehensively identifying SNP pairs having an effect),
(1) Input to input a total of M SNP genotype data (M is more than 500,000) observed from N specimens with two types of phenotypes and the phenotype class corresponding to each specimen. Steps,
(2) a storage step of storing in the storage means the phenotype class of N specimens and the total M genotype data inputted through the input step (2);
(3) According to the storing step (2), statistical processing of the two phenotype classes and genotype data for the sample N is performed on the i-th SNP stored in the storage means, and the class-specific genotype A calculation step for calculating another count and performing screening as a pre-processing step for determining the suitability of i-th SNP analysis continuation based on the calculated minor allele count,
(4) In the calculation step (3), if it is determined that “analysis continuation is suitable” as the analysis target SNP, whether the i-th SNP is dominant or inferior to the phenotype based on the calculated count A storage step of determining by statistical means, and storing the determination result regarding the suitability of continuation of analysis and the dominant / recessive type in an internal storage device;
(5) In step (4), the i-th and j-th (j ≠ i) for each of the two classes of phenotypes based on the determination result regarding the dominant type and the recessive type determined by the statistical means. Create a 2x2 contingency table in which the dominant and inferior types of two SNPs (i = 1, j = 2) are determined as initial values, and calculate an index for determining epistasis for the created 2x2 contingency table A calculation step for determining the presence or absence of an epistasis effect based on this index,
(6) When it is determined in the calculation step (5) that “the epistasis effect is present”, the determination result of “the epistasis effect is present” for the two SNPs is stored, and the j th The SNP of j + 1 is changed to the j + 1-th SNP, and the process returns to step (4). The dominant / recessive type of the j + 1-th SNP is determined by statistical means, and the calculation of step (3) is repeated, When +1 reaches M, in the analysis step of selecting the i-th SNP as the i + 1-th and the j-th as the i + 2-th SNP, and (7) in the calculation step (5), the “epistasis” A data analysis method comprising: an operation step of confirming a synergistic epistasis effect by a multivariate analysis means using logistic analysis analysis when it is determined as “effective”. As described in step (3) of claim 1,
In the operation for screening as a pre-processing step for determining the suitability of continuing analysis of the i-th SNP, "suitability of continuing analysis of the i-th SNP" is determined according to the following procedure: 2. The data analysis method according to 1.
The genotype data of each SNP is classified into 3 types, 2 homozygotes and 1 heterozygote in the population, depending on the type of two bases derived from mother and father. Here, these two homozygotes are represented by AA, aa, and one heterozygote by Aa;
Furthermore, a11, a12, and a13 are the counts for genotypes AA, Aa, and aa in the phenotype 1 class where the total number of specimens is n1, respectively, and a21, a22, and a23 are the phenotypes where the total number of specimens is n2. The counts for genotypes AA, Aa, aa in two classes;
In screening as a pre-processing step executed by the arithmetic device, the determination of “appropriateness of continuing analysis of the i-th SNP” is as follows:
When a11, a12, a13, a21, a22, a23 satisfy any one of the following conditions (I) to (IV), it is determined as “analysis continuation failure”,
SNPs determined as “continuation failure” are excluded from the analysis after step (4) (I) a11 + a12 ≦ 1 or a11 + a13 ≦ 1 or a12 + a13 ≦ 1 (Formula 1)
(II) a21 + a22 ≦ 1 or a21 + a23 ≦ 1 or a22 + a23 ≦ 1 (Formula 2)
(III) a11 = 0 and a23 = 0 (Formula 3)
(IV) a13 = 0 and a21 = 0 (Formula 4) As described in step (5) of claim 1,
In determining whether or not there is an epistasis effect,
About the combination of “i-th SNP and j-th SNP”, which is the object of determination,
Follow the procedure below to calculate R1 = (x ₁₁ x ₂₂ ) / (x ₁₂ x ₂₁ ) and R2 = (y ₁₁ y ₂₂ ) / (y ₁₂ y ₂₁ ) as indicators,
Based on the calculated index, the determination of “existence of epistasis effect”
Index: R1 = (x ₁₁ x ₂₂ ) / (x ₁₂ x ₂₁ ) and R2 = (y ₁₁ y ₂₂ ) / (y ₁₂ y ₂₁ )
R1 = (x ₁₁ x ₂₂ ) / (x ₁₂ x ₂₁ ) ≥ w ₁ and R2 = (y ₁₁ y ₂₂ ) / (y ₁₂ y ₂₁ ) ≤ 1 / w ₂ (Formula 5)
When the above (Formula 5) is satisfied,
The data analysis method according to claim 2, wherein it is determined that “there is an epistasis effect”.
In the above (Formula 5), x ₁₁ , x ₂₂ , x ₁₂ , x ₂₁ , y ₁₁ , y ₂₂ , y ₁₂ , y ₂₁ are calculated according to the following procedure.
In the above (Formula 5), x ₁₁ , x ₂₂ , x ₁₂ , and x ₂₁ are determined by the combination of the dominant type and the recessive type of the i-th SNP and the j-th SNP in the phenotype. It is a count.
Similarly, y ₁₁ , y ₂₂ , y ₁₂ , and y ₂₁ are counts determined by a combination of the dominant type and the recessive type of the i-th SNP and the j-th SNP in the phenotype.
The dominant type is a model in which the AA and Aa genotypes are related to the phenotype class 1, and are described as A1 = (AA, Aa) and A2 = (aa). The recessive type is a model in which the genotype of aa is related to the phenotype class 1, and is described as A1 = (AA) and A2 = (Aa, aa). The dominant type of the jth SNP is a model in which BB and Bb genotypes are related, and is described as B1 = (BB, Bb), B2 = (bb). The j-th recessive type of the SNP is a model in which bb genotype is related to phenotype class 1, and is described as B1 = (BB), B2 = (Bb, bb).
At this time, for specimens having a phenotype of class 1, c11 is a count of specimens having genotype AA of i-th SNP and genotype BB of j-th SNP (when AA and BB are included). Yes, c12 is the count of the specimen having AA and Bb, and c13 is the count of the specimen having AA and bb. Similarly, c21 is Aa and BB, c22 is Aa and Bb, c23 is Aa and bb, c31 is aa and BB, c32 is aa and Bb, and c33 is aa and bb. These counts satisfy the following formula:
c11 + c12 + c13 + c21 + c22 + c23 + c31 + c32 + c33 = n1 (Formula 6)
For specimens with phenotype class 2, d11 is the count of specimens with genotype AA of i-th SNP and genotype BB of j-th SNP (if it has AA and BB), d12 Is the count of specimens with AA and Bb, and d13 is the count of specimens with AA and bb.
Similarly, d21 is Aa and BB, d22 is Aa and Bb, d23 is Aa and bb, d31 is aa and BB, d32 is aa and Bb, and d33 is aa and bb. These counts satisfy the following formula:
d11 + d12 + d13 + d21 + d22 + d23 + d31 + d32 + d33 = n2 (Formula 7)
Based on the determination results for dominant and recessive types,
Specifically, x ₁₁ , x ₂₂ , x ₁₂ , x ₂₁ , y ₁₁ , y ₂₂ , y ₁₂ , y ₂₁ are given in the following cases.
(i) The i-th SNP is the dominant type and the j-th SNP is the dominant type
x ₁₁ = c11 + c12 + c21 + c22, x ₁₂ = c13 + c23, x ₂₁ = c31 + c32, x ₂₂ = c33,
y ₁₁ = d11 + d12 + d21 + d22, y ₁₂ = d13 + d23, y ₂₁ = d31 + d32, y ₂₂ = d33 (Formula 8)
(ii) The i-th SNP is the dominant type and the j-th SNP is the recessive type
x ₁₁ = c11 + c21, x ₁₂ = c12 + c13 + c22 + c23, x ₂₁ = c31, x ₂₂ = c32 + c33,
y ₁₁ = d11 + d21, y ₁₂ = d12 + d13 + d22 + d23, y ₂₁ = d31, y ₂₂ = d32 + d33 (Formula 9)
(iii) The i-th SNP is recessive and the j-th SNP is dominant
x ₁₁ = c11 + c12, x ₁₂ = c13, x ₂₁ = c21 + c22 + c31 + c32, x ₂₂ = c23 + c33,
y ₁₁ = d11 + d12, y ₁₂ = d13, y ₂₁ = d21 + d22 + d31 + d32, y ₂₂ = d23 + d33 (Formula 10)
(iv) The i-th SNP is recessive and the j-th SNP is recessive
x ₁₁ = c11, x ₁₂ = c12 + c13, x ₂₁ = c21 + c31, x ₂₂ = c22 + c23 + c32 + c33
y ₁₁ = d11, y ₁₂ = d12 + d13, y ₂₁ = d21 + d31, y ₂₂ = d22 + d23 + d32 + d33 (Formula 11)

X ₁₁ , x ₂₂ , x ₁₂ , x ₂₁ , y ₁₁ , y ₂₂ , y ₁₂ , y ₂₁ given in any of the above (i) to (iv) based on the determination result regarding the dominant type or the recessive type From these, the indices: (x ₁₁ x ₂₂ ) / (x ₁₂ x ₂₁ ) and (y ₁₁ y ₂₂ ) / (y ₁₂ y ₂₁ ) are calculated.
Note that w ₁ and w ₂ described in (Equation 5) are specified in the following range.
n1-3 ≦ w ₁ ≦ (n1 / 2-1) ² , n2-3 ≦ w ₂ ≦ (n2 / 2-1) ² (Formula 12) Since w ₁ and w ₂ described in (Formula 5) of Claim 3 are specified within the range described in (Formula 12), w ₁ and w ₂ are specified under the most lenient conditions shown in (Formula 13). Do
w ₁ = n1-3, w ₂ = n2-3 (Formula 13)
The data analysis method according to claim 3. Genome-wide single nucleotide polymorphism (SNP) genotype data of more than 500,000 sites using a computer, synergistic interaction (epistasis) even if the main effect is not confirmed for a phenotype having a binary class A data analysis system for comprehensively identifying SNP pairs having effects),
(1) Input to input a total of M SNP genotype data (M is more than 500,000) observed from N specimens with two types of phenotypes and the phenotype class corresponding to each specimen. Means,
(2) storage means for storing a phenotype class of N specimens and a total of M genotype data inputted via the input means (1);
(3) The i-th SNP stored in the storage means (2) is subjected to statistical processing of the two phenotype classes and genotype data for the sample N persons, and the class-specific genotype count is calculated. An arithmetic means for performing screening as a preprocessing step for determining the suitability of the i-th SNP analysis continuation based on the calculated minor allele count,
(4) When the calculation means (3) determines that “analysis continuation is suitable” as the SNP to be analyzed, whether the i-th SNP is dominant or inferior to the phenotype based on the calculated count Storage means for determining by statistical means, and storing the determination result on the suitability of analysis continuity and the dominant type / recessive type in the internal storage device;
(5) For each of the two classes of phenotypes based on the determination result regarding the dominant type and the recessive type determined by the statistical means, the i-th and j-th (j ≠ i, i = 1, j = 2) 2x2 contingency table in which the dominant type and recessive type of each two SNPs are determined, and an index for determining epistasis is calculated for the created 2x2 contingency table. Computing means for determining the presence or absence of an epistasis effect,
(6) When it is determined by the calculation means (5) that “the epistasis effect is present”, the determination result of “the epistasis effect is present” for the two SNPs is stored, and the j th The SNP is changed to the j + 1-th SNP, the process returns to the step of the recording means (4), the dominant type / recessive type of the j + 1-th SNP is determined by statistical means, and the calculation means (3) When the calculation of the steps is repeated and j + 1 reaches M, the analysis means for selecting the i-th SNP as the i + 1-th and the j-th as the i + 2-th SNP,
(7) In the calculation means (5), when it is determined that “there is an epistasis effect”, calculation means for confirming the synergistic epistasis effect by a multivariate analysis means using logistic analysis analysis;
A data analysis system characterized by comprising:   The program which makes a computer perform the method as described in any one of Claim 1 to 4.