Abstract
To identify putative relations between different genetic factors in the human genome in the development of common complex disease, we mapped the genetic data to an ensemble of spin chains and analysed the data as a quantum system. Each SNP is considered as a spin with three states corresponding to possible genotypes. The combined genotype represents a multispin state, described by the product of individual-spin states. Each person is characterized by a single genetic vector (GV) and individuals with identical GVs comprise the GV group. This consolidation of genotypes into GVs provides integration of multiple genetic variants for a single statistical test and excludes ambiguity of biological interpretation known for allele and haplotype associations. We analyzed two independent cohorts, with 2633 rheumatoid arthritis cases and 2108 healthy controls, and data for 6 SNPs from the HTR2A locus plus shared epitope allele. We found that GVs based on selected markers are highly informative and overlap for 98.3% of the healthy population between two cohorts. Interestingly, some of the GV groups contain either only controls or only cases, thus demonstrating extreme susceptibility or protection features. By using this new approach we confirmed previously detected univariate associations and demonstrated the most efficient selection of SNPs for combined analyses for functional studies.
1. Introduction
The amount of data from genetic studies is growing and its structure is becoming more complex due to the rapid development of new genotyping and sequencing techniques. However, current understanding of the correlation between found genotypes and observed phenotypes remains an obstacle for more efficient use of these data in medicine and biology. There are only a few methods for integration of genetic data in a reasonable way for statistical and functional analyses. An important feature of genetic data is that it is not a collection of random variations, but rather a system with intrinsic correlation between genetic variants. Investigations of study populations are often based on the frequencies of alleles rather than genotypes. However, these variants are distributed on two parallel carriers of information, paired sets of human chromosomes, and in most cases both alleles contribute to the phenotype. Traditionally, statistical analysis in genetics is based on several simplified models with limited options for integration of multiple variants in the same analysis. Although in recent years an exponential growth of different approaches to search for epistatic effects, or gene-gene interaction, in studies of complex traits is observed (see review [1]), a simultaneous consideration of multiple variants remains difficult due to “bottom-up” the design of the strategy of the search.
In this study we suggest and evaluate a genetic vector approach (GVA; briefly described in [2]) that would allow integration of available genotyping data in an unambiguous way, selection of the most representative combinations for the statistical analyses, and, finally, significant reduction of the number of variants that need to be inspected for functional studies. Our approach provides possibility to combine the advantage of joining several genetic markers in one combination (haplotype) with consideration to genotype (not allele) for each marker. Instead of ambiguous connection of allele or haplotype to certain phenotype in each individual, in GVA each individual has only a single genetic vector that corresponds to the phenotype. Thus the genotype-phenotype relation is unambiguous and it does not need a priori biological interpretation (dominance/recessiveness). GVA also permits expression of the traditional measures of the association for separate genotypes and, particularly, enables detection of gene-gene interaction.
2. Methods
2.1. Genetic Vectors for Association Studies
Due to the complexity of genetic data, multiple options for the analyses are available. However, the majority of genetic analyses are based on univariate tests for independent genetic markers and applications with simultaneous employment of several genetic markers are relatively rare. Based on linkage disequilibrium, haplotypes are considered to be a valuable tool [3], but this approach does not provide unambiguous assessment to phenotype, since individuals with heterozygotic state could be, at best, considered as an intermediate group or should be combined with homozygotic state in the absence of a biological hypothesis. We constructed an approach for genetic analyses that integrates a number of markers as a whole entity, regardless of the specific linkage on a particular chromosome.
The design of the model is based on integration of different genetic variants in a reasonable biological group by assigning each individual to a sequence of genotypes in relation to given genetic variants. Since the combinations available in the population are not random and reflect the population’s genetic structure, the number of these sequences of genotypes is significantly less than the number of individuals in a study population of reasonable size. The individuals identical in the order of genotypes represent a group, which may have differential representation of phenotypes. In a case-control design the frequencies of these groups could be statistically compared to identify those that are most likely to associate with a selected phenotype.
Imagine that we are interested in the degree of association to a certain phenotype of several single-nucleotide polymorphisms (SNPs), which may belong to the same, or to different chromosomes. A polymorphism is defined at the level of the whole study population, whereas for a particular person each SNP has a definite genotype value. Let us make the following convention: we choose and fix the order of the SNPs. Then each person will be characterized by the set of genotypes that can be thought of as a genetic vector (GV), characterizing a particular person. The number of genetic variations included in this vector is the length of the GV. The study population is then characterized by the set of GVs.
In an association study, each GV “contains” some number of controls and cases. Then the structure of the study population is described by the Table 1.
Thus, this is a simple way to avoid working with separate alleles, and even with genotypes; instead, the focus lies on their relevant combinations in the study population, that is, on the genetic vectors. One of the advantages worth mentioning is that if there exists some, say, genotype-genotype interactions, these will be taken into account automatically, since at this first stage the approach does not require any assumptions about the statistical nature of the GV’s constituents: this is simply one of the possible ways to represent the experimental data set. Since we now know the number of healthy controls and number of cases, as well as the total numbers of both, we can apply standard statistical machinery, testing the statistical hypotheses for whole GVs; particularly, in order to discover the “most promising” GVs, we can (and will) evaluate the odds ratios, the regions with the greatest fluctuations on the basis of, say, Fisher’s exact tests, and so on.
The GVA consists of several stages:(1)selection of genetic markers,(2)genotyping for two study populations,(3)assigning of genetic vector value to each individual,(4)sorting all individuals in the study population into groups, according to particular GVs and chosen phenotype parameters,(5)exclusion of GVs with zero number of individuals in any subgroup (subset),(6)computation of odds ratio and relative risk (OR/RR), pathogenic, and protective genetic vectors,(7)sorting GVs according to statistical confidence,(8)exclusion of GVs with low confidence,(9)analysis of the match between study populations: focus on replication and absence of opposite effect for the vector.
The mathematical details of the GVA itself are given below (Section 2.3–2.5).
2.2. Experimental Data Sets and GV Definitions
To demonstrate the capacity of the methodology we applied it to the analyses of data from two study populations: from the Swedish Epidemiological Investigation of Rheumatoid Arthritis (EIRA, 1820 cases and 947 controls) and from the North American Rheumatoid Arthritis Consortium (NARAC, 813 ACPA-positive cases and 1161 controls).
The EIRA study population is a population-based case-control study of incident cases of rheumatoid arthritis (RA) in which all patients fulfilled the American College of Rheumatology (ACR) 1987 criteria [4]. Controls were randomly selected from the Swedish national population registry, taking into consideration the patient’s age, sex, and residential area. More details about the EIRA study population have been described elsewhere [5].
The cases in the NARAC study population consisted of RA patients of self-reported white ancestry, who were randomly drawn from four different groups of patients, while controls were recruited from the New York Cancer Project [6, 7]. Both studies were conducted after obtaining approval from the Regional Ethics Committees and in accordance with the Declaration of Helsinki.
Below, we apply GVA for the analysis of EIRA’s and NARAC’s association to rheumatoid arthritis (RA) using the data for 6 SNPs from the HTR2A locus, rs6314, rs977003, rs1328674, rs2070037, rs6313, and rs6311 and for HLA-DRB1 shared epitope alleles. We defined SE alleles as any of HLA-DRB1 (not ), and ; previous study suggested very low frequency for non-SE HLA-DRB1 variants in Caucasian populations [8]. Conventional evaluation for HTR2A genetic markers was performed previously [9].
We defined the GV states as a product of SNP states, as shown in Table 2. In order to compare the GV structure between the EIRA and NARAC study populations and to perform statistical evaluation, we, first, assigned individuals to a particular GV group using individual genotypes and, second, normalized the number of people, for each GV (according to (10) below).
2.3. Mathematical Details of the Genetic Vector Approach (GVA)
Let us consider a population, consisting of controls ( = healthy) and ( = sick) cases. The number of genetic markers is a matter of preliminary research and will be assigned as . We choose here SNPs for simplicity, but the approach can be easily extended for multiallelic markers. There are no restrictions on the position of markers within chromosomes, and variants from different chromosomes could be used. In the whole study population, with the dizygotic human genome, each SNP can normally produce 3 different categorical values, or genotypes. However, these categories are exclusive and only one of the categories can be found in the DNA of an individual person. This is similar to one of the components (say, ) of spin , which can acquire three values, , and could be described by corresponding three states, . The description of the study population with a set of SNPs is similar to the description in physics of an ensemble of spin chains (see, e.g, [10]). Each of the spin chains contains spins and can be in different quantum states. The latter is described by the set of vectors:where represents other variables; the spin corresponds to the SNP, and its values correspond to the genotypes. Therefore, a combination of genotypes form a vector for each person with known genotypes and we refer to these combinations as genetic vectors (GVs).
Different statistical characteristics of the study population, such as the average frequency of people with a certain genotype and correlations between genotypes and phenotypes, can be obtained by, first, averaging certain spin-operator combinations on these multispin states and, second, averaging the result over the population of interest. However, it is more convenient for us to speculate not in terms of -operators, but, rather, in terms of the operators, defined as follows:
These variables can be interpreted as the operators of population numbers (PNs) of the -state in the cell of the genetic vector . The -operators are projection operators: . Notice that the projection operators are sometimes used for description of spin systems, since the spin operators can be easily expressed in terms of -operators. For example, -projecture of spin operator is ; that is, it is just the difference between the population-number operators of two states.
In genetics this multiple-spin state, describing the combination of genotypes, GV, and characterizing each person, could be obtained from the experimental genotyping. Then, one can compare all expected values of interest in groups via the expectation values of corresponding operators on the wave functions (GVs) for two ensembles, for example, controls and cases. We are dealing here with a static problem: all the states of “spin chains” are given; they do not depend on time; the total wave function of each of the ensembles = ctrl, case is a product of the wave functions of separate “chains,” . The expectation value of the population-number operator for a certain person is
Therefore, in general, the individual is characterized by numbers (, and , if other parameters are excluded) or by the matrix of PNs for each . We will refer to it as index matrix (IM): Thus, an individual is fully characterized by this matrix as described above. Each column, describing the filling of the “spin-one state”, can contain only one nonzero value equal to one and two others equal to zero. Each column represents a particular variation and, for a SNP, has three levels defined by nucleotides at the SNP position. The matrix for selected genotypes is a matter of experimental design and may reflect either a set of selected SNPs of interest or a sequence of all SNPs at a certain locus. For example, if we are interested in an ensemble of three SNP chains, with , , and , with , , and , where are corresponding genotypes, then some people may be in the state characterized by the genetic vector: Corresponding to this GV index matrix is
These IMs will be used in our further calculations. For statistical evaluation of differences between controls and cases we will develop two matrices of PNs: These, in fact, are the genotype-frequency maps for controls and cases for the study population of interest.
In the same fashion we can obtain the rate for people, who have a pair of certain genotypes and ; this rate can be calculated from the expression: The expression for triplets contains a product of three -operators and so on.
Each GV “contains” some number of controls and cases, whereas the structure of the study population is described by Table 1.
Notice that the set of GVs is orthonormalized and, therefore, can be considered as a basis set in -dimensional space, whereas each part of the study population, controls and cases, can be presented as a point in this space. If some other study population contains GVs that are not included in the set , the set should be complemented by these GVs in order to obtain a complete set for the extended population, including all populations in question.
2.4. Statistical Treatment of GVs
In order to avoid possible confusion it is important to remember that each GV describes the set of SNP values. If the population has a group of individuals with the same GV, we can use this fact in the statistical treatment of the population, describing the whole population in terms of these larger-scale variables, GVs. Then, each GV within a given study population is characterized by the number of healthy controls and the number of cases and we can introduce two frequencies: Comparing these, we can determine which of the GVs are more protective and which are more pathogenic. As a criterion for the separation of GV subgroups with protective and pathogenic behavior we choose here the odds ratio, defined by GV frequencies:
The difference between the standard definitions of OR and the one used here is discussed in the book [10, page 154]. Further, we have to choose some threshold , which reflects a chosen criterion of statistical confidence, and select all GVs that satisfy the conditionfor “pathogenic” GVs and for the “protective” ones, for which . The standard obvious threshold in (12) would be . However, statistical fluctuations in the vicinity of the threshold may make the choice of the groups, described by certain GVs, doubtful. In order to garantee that the selected GVs remain on a safe side, one can choose a harder criteria of the GV selection. Of course, the threshold for “protective” GVs not obligatory should be chosen equal to . Due to the different number of people in each GV, the threshold value can occur as -dependent. Let us denote the subset of GVs defined by the inequality (12) as (the details of this subset will be discussed in Sections 2.5 and 3.2). Then the combination of variants, responsible for a large OR, can be visualized by summing the IMs of the GVs that fulfil (12). This summation produces two matrices, which represent frequencies of genotypes in selected control and case groups. The cells in the matrix that contain the largest difference give us a possible indication of which of the genotypes, entering different GVs, are responsible for the large ORs displayed by the GVs.
2.5. The Choice of Thresholds
Let us now discuss a criterion for the choice of thresholds . In the first step, all individuals in the population are sorted in terms of corresponding genetic vectors (GV) of a length . Each GV contains of healthy controls (HCs) and cases, with a total of individuals in each GV. The whole population consists of controls and case individuals, with the total number of individuals in the population . Then for each GV we can write down the contingency Table 3 and the odds ratio based on this table. The full set of GVs can be sorted into two subsets, one with , which reflects the pathogenic combinations, and the rest with , corresponding to the protective combinations of genotypes.
At this stage we have to estimate the degree of confidence for each . In order to do this, we introduce the random variable , which describes a number of case individuals, selected randomly from persons, with the number kept fixed. Then the probability to pick up precisely this number for fixed values , , and is described by the hypergeometric distribution (as used in Fisher’s exact test):where are binomial coefficients. Constructing the complex random variable of interest, OR, for the chosen , we can calculate for each GV the average odds ratio and, say, variance, . This would be sufficient if we had normal distribution of the variable of interest, , but our case is far from this. One could also try to use the higher moments in order to build the confidence intervals for for each GV . However, due to asymmetry of the distribution for the odds ratio and the presence of GVs with a small total number of individuals in GV, this method is impractical. Instead, we build the distributions for the complex random variable , , and find the thresholds for the interval, where experiences the largest fluctuations (for the normal distribution this would be the interval , where is the standard deviation and is chosen accuracy) directly from the equation based on this distribution: We choose a 5% level of accuracy, that is, . Since the full number of steps in the discrete distribution for the total number of individuals in GV is restricted by , the best accuracy that can be achieved is also restricted by the magnitude of this step; the sum in (18) will consist of only one term with this minimal possible step. As normal, a decrease of also decreases the number of contributing GVs, decreasing the total statistical power.
The other way to estimate the degree of confidence for the experimentally found for each GV is to evaluate the sum of probabilities for all possible exceeding : The criteria are illustrated in Figure 1.
We use both criteria.
The GVs, which contain either only cases (, , i.e., “solely sick”) or only controls (, , i.e., “solely healthy”), require separate consideration. For simplicity we call them “zero GVs.”
3. Results and Discussion
3.1. GV Frequency Distributions
We have found that all individuals in two study populations, based on selected SNPs from HTR2A and HLA-DRB1 SE alleles for the chosen GV length (number of markers), can be ascribed to 161 GVs for EIRA and 163 GVs for NARAC; 131 of them are common GVs for both studies. The total basis of GVs for these both cohorts, therefore, consists of 193 GVs, while Figure 2 displays the results for common GVs only.
(a)
(b)
(c)
The rest of the GVs contain a relatively small number of individuals in each of the cohorts: 40 people of 2767 total in 30 GVs (1.45%) in EIRA and 36 people of 1974 total in 32 GVs (1.82%) in NARAC. Since each of these noncommon GVs contains just one or two people, it is obvious that these relatively rare GVs would be unlikely to give a significant contribution to the statistical characteristics of the cohorts. Thus, the vast majority of GVs in both cohorts are the same and it is reasonable to compare GV frequencies between the two study populations.
The data in Figure 2 demonstrate the high level of coherency of the GV distribution within each study population (EIRA and NARAC) and provides a rationale for statistical comparison between RA cases and controls. It is also evident that these two populations differ in the frequencies of specific GVs, which is likely to reflect a different genetic backgrounds or differences in the selection criteria. More specifically, as can be seen from Figure 2(c), frequencies of several GVs for the RA cases in the two studies differ significantly, with some GVs absent from NARAC, while relatively common in EIRA.
3.2. Stratification of the Common GVs and an Analysis of the Consistency
In order to facilitate identification and comparison between high-risk and protective GVs in the two study populations, we selected GVs which display or (see (18) in Section 2.5). However, the GVs with either or should be considered separately, since in this case OR cannot be formally defined. These subgroups of GVs with zero values will be described in Section 3.4. The distribution of GVs between the groups is shown in Table 4.
After exclusion of the GVs specific only for one study and those with zero frequency in any of subsets, we analysed the set of GVs with corresponding numbers of individuals, as shown in Table 5.
As seen from Table 5, despite the fact that the number of GVs decreased more than twice after exclusion of nonoverlapping GVs and GVs with zero values, the remaining GVs represent 76–88% of individuals. The data for the set of common GVs with both and are displayed in Figure 3.
(a)
(b)
In order to estimate quantitatively the degree of consistency of the data, we stratified the set of common GVs into four categories using OR, defined by (11):
The result of this is summarized in Table 6 and only GVs with consistency in direction of the effect (risk or protectivity) were selected for futher analyses; see Tables 7 and 8. At the next step, using (17) for the average OR for each GV and (18) at the level , we can construct the plot for the group of GVs: , and the “experimental” value of OR. Subsequently, for the group, where all , we used , . The results are displayed in Figure 4.
(a)
(b)
In this comparison, the NARAC group has a much greater size of ORs. However, some of the NARAC GVs happen to be inconsistent with the ones from EIRA. Therefore, only those GVs which exceed the upper 5% threshold in both cohorts were chosen for further analysis.
We first analyzed the risk groups, , represented in the left panels of Figure 4 and found that only 6 GVs, for which in both cohorts, give consistent results. Let us inspect the probabilities of occasionally finding an OR greater than in the selected GVs. The result is presented in Figure 5.
Thus, since these subgroups fulfill both criteria, namely, the GVs display OR beyond the fluctuation width near for chosen and show a low probability of finding ; the combination of variants corresponding to this subgroup GVs should be considered as most probable candidates for an association with RA.
Among the subgroup of “protective” genotype combinations, which are represented in the right panels of Figure 4, only the first eight GVs give consistent results in both the EIRA and NARAC subgroups. In contrast, the NARAC GVs from 9 to 15, which display , do not fulfil similar requirements in the EIRA subgroup. The values of for the consistent GVs in this subgroup are presented in Figure 6.
The probabilities shown in Figure 6 suggest that a protective role is played by the combinations of variants, presented by these GVs.
In summary, we found that there is a family of GVs in both study populations that follows the same direction of association, although the absolute values of the effects displayed in the NARAC study are more extreme than the ones from the EIRA study for both protective and risk GVs. We chose to use these GVs, with the aim of finding more specific effects from individual genetic variations for further optimization of analyses.
3.3. Candidate Combinations of Genotypes
We used all the available combinations from our genotyping data to differentiate GVs associated with disease. However, the functional consequences from these variations are likely to correspond only to a limited number of selected markers and most of the markers are not important for the association study. Therefore, for further studies it would be desirable to reduce the length of the GVs and to identify the most influential variants within GVs. By excluding unnecessary genetic variants from GVs, we could reconstruct shorter GVs, to facilitate inspection of the genotype content of the GVs in two contrast groups, chosen via two consequent statistical criteria.
In order to compare the results from the contrasting GVs of control and cases with different total numbers of individuals, we sum all IMs (see (4)) in each of categories 11 and 22 with the weights , where is number of individual cases (, or “sick”) or controls (, or “healthy”) in GV , belonging to the class ; here are the numbers of individuals in the selected group ; notice that here the normalization is different from the one used in (7):
As expected, our data demonstrate strong association of RA with SE alleles (see Table 9). At the same time, the SNP rs1328674 does not show an association with RA, in contrast to a previous report [9], while we detected an association of the SNPs rs6314, rs977003, rs2070037, rs6313, and rs6311 in the EIRA study. Since not all of these results agree with previously published evaluation at this locus, we compare it with independent data from NARAC; see Table 10.
As can be seen, the data for NARAC fully confirm the conclusions derived from the EIRA data. The numbers for NARAC in Table 10 differ from those for EIRA (see Table 9) only slightly; therefore, the tendencies revealed in both cohorts from the same sets of GVs can be considered to be cross-validated.
Thus, selection of shorter GVs for statistical evaluation demonstrated how the overall analysis could be optimized: , , , and remain indicative for the difference between groups, whereas , , and have no influence on association in the GVs of both EIRA and NARAC. These cropped GVs can be especially useful for optimizing the analysis of long GVs, since it can improve the statistical confidence of the results.
We identified significant heterogeneity of association of combinations of HTR2A variants, which depends on the absence or presence of a shared epitope allele in the GV. Several interesting observations have been found from application of the method to two RA study populations, EIRA and NARAC.
We found that the number of GVs that describe the study population for the chosen length of GV is relatively small: there are only 161 GVs required for description of the EIRA, with 2767 individuals, and 163 GVs for NARAC, consisting of 1974 people. A very high overlap between Swedish and North American study populations was detected in our analysis: 131 of GVs are common to both study populations.
Interestingly, some of the GVs contain only healthy controls and do not contain the case counterpart (subset of solely healthy persons), and vice versa; some of the GVs contain only cases and do not contain controls (solely ill). The statistical weight of these GVs, however, appeared to be small. Nevertheless, we have found similar behaviour in two independent cohorts and this observation deserves further investigation. Even after removing these “zero” GVs from consideration, the remaining set of common EIRA and NARAC GVs contains about 80% of the total number of individuals in each cohort.
Our literature search did not give much information about similar methods being used in genetic studies. The closest approach, haplotype analysis, is based on the classical genetic idea of transmission of the marker from parents to offspring and therefore the data in this analysis are arranged around the combination of genetic markers within a single chromosome. In an association study this approach may well serve in allelic mode, but it totally ignores the possible involvement of a second allele. In functional studies, selection of individuals by haplotypes is straightforward for homozygotic states but introduces multiple combinations in the case of heterozygosity and generates an extreme number of genotype groups in comparison with GVA. Experimental detection of haplotypes is a difficult procedure and the common statistical approach for assigning of haplotypes is never 100% unambiguous.
3.4. Genetic Vectors with Zero Values in One of the Subgroups, Experimental Data
In this section we illustrate previously defined “zero GV” detected in the EIRA and NARAC studies of RA. The number of individuals involved in each step of selection is decreased as shown in Table 11.
Thus, of 131 GVs common for EIRA and NARAC, 19 GVs happen to be zero GVs, that is, either and or, vice versa, and . They contain only about 2% of individuals from the study populations (62/2727 in EIRA and 43/1938 in NARAC). It is interesting to compare the contents of GVs from the “solely sick” (SSGV) and the “solely healthy” (SHGV) subgroups. Since the number of individuals belonging to the classes of interest is small, it would be nice to have a look at the pooled data. The inspection, however, shows that the part of the GVs that belong to the zero GV subgroup in, say, EIRA, often does not belong to this subgroup in NARAC, and vice versa. These GVs lose their status of zero GVs in the pooled data.
The insufficient number of individuals in SSGV and SHGV subgroups does not make it possible to draw statistically confident conclusions if we work with each study population separately. We can, however, investigate whether there exists some cross-validated tendency, that is, the one observed in both EIRA and NARAC.
This can be visualized via GV frequency weighted index matrices for zero GVs. Since we are interested only in the groups of zero GVs, we will perform averaging over these subgroups only. In other words, (13) should be modified as follows: where and are numbers of individuals in the GV belonging to the group SH or SS correspondingly; indices sh and ss mean “solely sick” and “solely healthy”; the numbers are total numbers of individuals in the SSGV and SHGV subgroups and and are index matrices of the GV from SH and SS subgroups. Then we can compare EIRA and NARAC results for case-control differences of index matrices, for zero GVs. The result is shown in Figure 7.
Full coherence in EIRA and NARAC data display is only SE (no SE is protective, while double SE strongly associates with RA, which is expected) and the SNP rs1328674. The positive value of the difference for frequencies of the genotype (rs1328674,CT) means the tendency to disease, whereas the negative one for the genotype (rs1328674,CC) indicates the tendency to protect. All other genotypes do not display consistency between the results for EIRA and NARAC. Thus, although the fact of zero GVs existence is interesting, we have to admit that, due to the absence of (i) the consistency for two study populations and (ii) the possibility of drawing statistically confident conclusions, analysis of zero GVs should be taken with caution when the numbers of observations are too low. In practice, it will need replication study in very big cohorts.
4. Conclusions
4.1. On GVA
Our approach, GVA, is the method for genetic analysis based on testing of association to a complex trait simultaneously for combinations of variants (multiple markers) inherent to individuals without any assumptions about their statistical relations. It is based on ascribing to each person a multiple-spin quantum state, GV, with subsequent description of the study population as an ensemble of spin chains with a given (determined by experimental data) set of allowed states. The full data set is sorted into subsets, each of which is described by a single GV and the numbers of cases and controls representing this GV. We successfully tested this approach for the set of experimental data for RA and validated the previously known association with shared epitope alleles.
GVA has several advantages:
First, since GVA works with the genotype combinations as a single entity, without decoupling of random variables, the issue of possible statistical independence (as well as interaction, influence of the linkage disequilibrium, etc.) of separate genotypes does not affect analyses.
Second, being formulated for genotypes, it is free of the uncertainties connected with the use of an allelic model or haplotypes: it is, therefore, not necessary to guess a mode of dominance. At the same time, the frequencies of separate genotypes, or their pairs, triples, and so forth, can be easily expressed in terms of GV variables, which means that other analyses could be easily performed by simplification of the GV algorithm.
Third, from the perspective of functional genetics and personalized medicine, this approach provides the opportunity to directly assign individuals belonging to the risk group of interest to a genetically determined subset defined by GV. Since GVA includes in the analysis only the genotype (phenotype) combinations which are met in the population study, the number of combinations that need to be analyzed is significantly decreased compared to the number of random combinations of genetic markers; therefore, the volume of necessary calculations is greatly reduced. The problem of multiple combinations of haplotypes is also resolved in GVA.
The probability distribution for OR of each GV, which we derived from the hypergeometric distribution for the random variable “number of cases in GV,” happens to be asymmetric and is not Gaussian. For this reason, the standard method for evaluating the confidence accuracy could not be used. Thus, we instead evaluated the thresholds of the OR interval where OR fluctuations exceed 5% for each OR separately. The statistics of each GV depends on three parameters: the total number of cases, the total number of controls in the study population, and the total number of individuals (cases plus controls), described by the GV in question. Only those of the GVs for which OR falls outside of the interval of large OR fluctuations have been selected for further analysis. The dependence of the OR statistics for a particular GV on the number of individuals “belonging” to this GV creates the obvious disadvantage of GVA: with an increase of the GV length (the number of SNPs included in the GV) the total number of GVs describing the study population approaches the number of individuals, since each person has unique DNA and, therefore, statistical analysis becomes invalid.
Regarding optimization of GVA and reducing the number of markers in GV (cropping the vector size), we suggested using index matrices, uniquely characterizing each GV. By using this approach a comparative analysis of the “genotype content” can be done for the sets of GVs of interest, identifying the sets of genotypes that influence the difference between the control and case groups.
It is worth noting that GVA does not require any assumptions about statistical coupling between variables. The study of interaction and synergetic effects between particular SNPs can also be performed in terms of GVA. Indeed, the answer to the question whether pairs, or triples, of some of SNPs produce a much stronger association to the phenotype of interest than is expected from independent factors can be obtained from the analysis of the correlation functions between SNPs. Examples of this type of analysis have been discussed in [2].
4.2. On RA in EIRA and NARAC
Our analysis reveals that the majority of GV groups, 76.7% (or 81.4% of individuals in all GVs, Table 6), are consistent for the effect between EIRA and NARAC. The degree of inconsistence could be possibly explained by low number of observation and randomness, but also by influence of different environmental factors or by different genetic architecture of the locus (between included SNPs and at flanking regions).
We confirm the well-known fact that shared epitope (SE) is strongly associated with RA. What seems to be significant here is that the whole set of genotypes displays quite large differences in frequencies in the presence and absence of SE. The SNPs rs2070037, rs6313, and rs6311 seem to play a more important role in RA formation than they were considered to do in a previous study [9].
List of Abbreviations
| EIRA: | Swedish Epidemiological Investigation of Rheumatoid Arthritis |
| GV: | Genetic vector |
| GVA: | Genetic vector approach |
| NARAC: | North American Rheumatoid Arthritis Consortium |
| OR: | Odds ratio |
| RA: | Rheumatoid arthritis |
| SNP: | Single-nucleotide polymorphism. |
Conflict of Interests
The authors declare that they have no conflict of interests.
Authors’ Contribution
Igor Sandalov and Leonid Padyukov are responsible for design of the study; Igor Sandalov is responsible for computations and model evaluation; Igor Sandalov and Leonid Padyukov are responsible for preparation of the paper.
Acknowledgments
This study was supported by the Swedish Research Council, Vinnova (COMBINE Grant), and Innovative Medicines Initiative Joint Undertaking, 115142-2 BTCURE. The authors would like to thank Janet Ahlberg for help with English editing of the paper.