Correction: Statistical Dissection of Cyto-Nuclear Epistasis Subject to Genomic Imprinting in Line Crosses

Copyright: ß 2014 The PLOS ONE Staff. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.


Introduction
One of the central foci in biological study is to unravel the genetic secrets of complex traits of agricultural, evolutional and biomedical importance.Quantitative trait locus (QTL) mapping has been the major tool for this purpose over decades [1][2][3].In QTL mapping, the identified QTLs are chromosome segments harboring potential genetic variants that could give rise to phenotypical manifestation.Large successes have been witnessed in the past [4].However, there are still many phenomena that could not be explained by Mendelian genetics, leading to the new exploration of research focus on epigenetics [5].
Genomic imprinting, one of the major epigenetic phenomena, plays key roles in controlling embryonic growth and development [6,7].Let subscript letter M and F denote the parental origin of two alleles in a diploid organism, then a locus with two alleles A and a is thought to be imprinted if two heterozygotes A M a F and a M A F have different expressions [8].The malfunction of imprinted genes could potentially lead to abnormal characters such as cancers or other genetic disorders [9].
Genomic imprinting effect is considered as one type of parentof-origin effect due to allelic effect with specific parental origin.In contrast to this, maternal effect or cytoplasmic effect is also considered as one type of parent-of-origin effect in which the offspring's expression is influenced by maternal parent.For example, a mother's genotype, even if not transmitted to her offspring, may influence in utero conditions and increase risk and/ or interact with genetic predisposition for particular diseases among those offspring [10].For cytoplasm, it contains a wide variety of organelles such as mitochondria and chloroplast (in plant).Almost all the reactions of cellular metabolism take place in such an environment.It has been demonstrated that cytoplasm plays a central role in coordinating the activities of nuclear genetic materials [11][12][13][14][15]. Thus, the identification of cyto-nuclear interaction could shed novel insights into the genetic and epigenetic control of phenotypic variation.A number of empirical studies have documented the significant contribution of cytonuclear interaction to phenotypic variation in organisms such as wheat, rice, mice, yeast and Drosophila [16][17][18][19].
On the other hand, the existence of such parent-of-origin effects may lead to incorrect interpretation of the (marginal) effects of particular genes when performing genetic mapping studies, unless such effects are appropriately accounted for in the analysis [20].Statistical methods for dissecting genomic imprinting effect has been extensively studied in literature (e.g., [21][22][23][24][25][26]). Tang et al. [27] developed a model to evaluate cyto-nuclear interaction effect based on experimental crosses.However, how the two types of parent-of-origin effects, one in nuclear level and one in cytoplasmic level, interact to influence phenotypic variation is largely unknown due to the lack of proper statistical models.
In this work, we discuss potential scenarios of parent-of-origin effects, and present a statistical method to dissect the cyto-nuclear interaction effects subject to genomic imprinting.The developed framework is based on experimental crosses which can be realized through two different designs, reciprocal backcross and F 2 design.When an F 2 design is considered, sex-specific difference in recombination fractions, which is initially discovered in [28,29] and later observed in many species, is incorporated into our model to distinguish the genetic differences between two reciprocal heterozygote.Such information can be found in literature, such as the female-to-male recombination ratio of 1.6:1 for human [30], 1.4:1 for pig [31], 1.4:1 for dog [32], 1.25:1 for mouse [33] on average across the whole genome.A genome-wide scan for the identification of iQTL mapping cyto-nuclear epistasis is performed based on the interval mapping theory, and parameters are estimated based on the framework of maximum likelihood method implemented with the EM algorithm [34].Extensive simulations are conducted to evaluate the performance of our model under different scenarios, such as different sample sizes, different heritability levels, and different gene effects.The utility of the method is illustrated by applying it to a genome-wide scan of four traits in an F 2 family derived from two inbred mouse strains.

Genetic designs
Consider a design initiated with two inbreed lines with two segregating alleles A and a.Let subscript letter M and F represent the parental origin of offspring alleles inherited from the mother and father, respectively.A complete dissection of the cyto-nuclear interaction subject to imprinting needs experimental designs that can distinguish the quantitative variation between two heterozygotes A M a F and a M A F , and also against the cytoplasmic effect.For this purpose, a reciprocal backcross or F 2 design is proposed so that variations and interactions can be fully introduced.
Figure 1 illustrates the reciprocal backcross design.Let the maternal parent carrying genotype AA (denoted as P2) has a cytoplasmic effect in contrast to the maternal line carrying genotype aa (denoted as P1) as the reference line.In the diagram, individuals that carry the maternal effect coming from the maternal parent with genotype AA are denoted by gray squares.Four possible backcrosses can be initiated as illustrated in Fig. 1.As shown in the figure, any backcross offsprings coming from the middle two designs carry the cytoplasmic maternal effect derived from the AA genotype.
Figure 2 shows the F 2 design.Denote the left one as design S 1 , and the right one as design S 2 .The cross of two F 1 's generates four possible allele-specific F 2 genotypes.Assuming there is a cytoplasmic effect, F 2 offsprings may show different phenotypes depending on the genotype of the maternal parental lines.For example, if maternal cytoplasmic effect exists, the offspring phenotypic value for AA genotype may be different depending on whether it comes from the S 1 or S 2 design.In the F 2 design, the two reciprocal heterozygotes A M a F and a M A F cannot be distinguished in general.Sex-specific recombination difference in male and female needs to be considered in order to distinguish the two (Cui et al. 2006) [25].

Statistical parameterization
For a particular cross, let y j (j~1, Á Á Á ,n) denote the phenotypic value of interest.Following Tang et al. [27], the one-QTL genetic model can be expressed as, where m is the overall mean; c is the cytoplasmic effect; a, d and i are the additive, dominance and imprinting effects of a QTL, respectively; and i ca , i cd and i ci are the cytoplasm by additive, cytoplasm by dominance and cytoplasm by imprinting interactions, respectively; x 1j is an indicator variable, denoting x 1j ~1 for the AA maternal cytoplasm and x 1j ~{1 for the aa maternal cytoplasm; x 2j ,x 3j , Á Á Á ,x 7j are other indicators of various effects describing the additive, dominance and the interaction effect between the cytoplasm and genetic variables.
For the S 1 design in the F 2 population initiated with cross AA|aa, the mean genotypic values of four possible genotypes formed by different allelic combination from the two F 1 parents can be expressed as, Similarly for the S 2 design initiated with cross aa|AA, the genetic model can be expressed as, For the reciprocal backcross design which consists of B 1 , B 2 , B 3 and B 4 , the indicator variables in Eq. (1) describing different QTL genotypes are defined in Table 1.
For simplicity, we will use matrix form to rewrite the models.Let us denote b~(m,c,a,d,i,i ca ,i cd ,i ci ) T , m~(m 1 ,m 2 ,m 3 ,m 4 ,m 5 ,m 6 ,m 7 ,m 8 ) T : Then the relationship between the eight genetic means and the eight parameters can be written as where for the reciprocal backcross design, and for the F 2 design.For the purpose of illustration, the following estimation and inference is demonstrated through the F 2 design.The same procedure applies to BC design too.

The mixture model and the likelihood
Statistical methods for QTL interval mapping based on a mixture model traced back to the work by Lander and Botstein [1].In the mixture model, each observation y is modeled as a weighted mixture of J (known and finite) components, and each component, which corresponds to a certain genotype category depending on the underlying genetic design, follows a certain distribution f j with weight p j .Conditional on the marker genotype M and unknown parameters w and g, the density of the observed y has the following expression where p~(p 1 ,:::,p J ) T refers to the mixture proportions which are constrained nonnegative and P J j~1 p j ~1; w~(w 1 ,:::,w J ) T is a vector for the component-specific parameters, with w j being specific to jth component; and g consists of parameters (i.e., residual variance) that are common to all components.
For the F 2 design we described above, there are four genotypes at each locus (A M A F , A M a F , a M A F , and a M a F ).The genotype of the QTL is generally unobservable, but can be inferred by using the two flanking markers' information.Given the flanking marker genotypes of the ith individual, the conditional probabilities p i ~(p AM AF Di ,p AM aF Di ,p aM AF Di ,p aM aF Di ) of the QTL genotype can be calculated.These conditional probabilities become the mixture proportions in the mixture model (7).Let us denote p i ~(p 1Di ,p 2Di ,p 3Di ,p 4Di ) to simplify the notation.These conditional probabilities are expressed in terms of sex-specific recombination rates in order to distinguish the two reciprocal heterozygotes.Table 1.QTL genotypes and corresponding genetic components under different backcross designs.Please refer to Cui et al. [25] for the conditional probabilities of QTL genotypes given marker genotypes in terms of sex-specific recombination fractions for an F 2 design.Assume the total number of F 2 offsprings for design S 1 and S 2 are n 1 and n 2 respectively, and let n~n 1 zn 2 .Phenotype data for a certain quantitative trait can be observed and recorded as a vector y~(y 1 ,:::,y n ), where y (1) ~(y 1 ,:::,y n1 ) are from design S 1 and y (2) ~(y n1z1 ,:::,y n ) are from design S 2 .Marker information can be reorganized as matrix M~(M 1 ; M 2 ) T , where the jth (j~1,:::,n 1 ) row of M (the 1 jth row of matrix M) contains all the marker information of the jth individual under design S 1 and the kth (k~1,:::,n 2 ) row of M 2 (the (n 1 zk)th row of matrix M) contains all the marker information of the kth individual under design S 2 .Based on the mixture model ( 7), and with the independence assumption, the joint likelihood function for the F 2 family with total n individuals, constructed by combining design S 1 and S 2 together, can be formulated as, where the unknown vector contains the QTL position, QTL effects and residual variance, and the density function f j (j~1,:::,8) is assumed to follow a normal distribution with mean m j and common variance s 2 .More specifically, parameter vector can be divided into two subsets, l and g , where l describes the location of QTL and g contains all the genetic parameters, including QTL-effects vector b~(m,c,a,d,i,i ca ,i cd ,i ci ) T and residual variance s 2 in our model.

Parameter Estimation
To estimate the unknown parameters ~( l , g ), several algorithms could be implemented, such as Expectation-Maximization (EM), Newton Raphson and Fisher Scoring.Among all these methods, EM algorithm, which was initially developed by Dempster et al. [34], is most commonly used in QTL mapping study.In this paper, EM algorithm is applied to obtain the maximum likelihood estimates (MLEs).This procedure involves differentiating the log-likelihood function with respect to each unknown parameter, letting the derivatives equal to zero, and solving the log-likelihood equation for the corresponding parameter.Please read Appendix S1 for the detailed derivations of parameters estimation and the algorithm.
For the QTL position which is unknown in the model, we did not estimate parameters l directly.As commonly treated in QTL mapping studies, we applied a grid search approach to estimate the putative QTL position via scanning the entire linkage genome by 1 or 2 cM increment flanked by two markers and did a hypothesis testing at each putative position.A likelihood ratio or LOD profile plot can be generated to graphically display the LR or LOD test statistic for a putative QTL at each testing position.The genomic position which corresponds to a peak in the profile plot is the MLE of the QTL location, given that the peak passes the genome-wide threshold determined by the permutation tests detailed below.Bootstrap methods can be applied to assess the confidence interval of the estimated position [35].

Hypothesis Test
Testing the overall QTL effect on the quantitative trait is the first step toward a complete dissection of different genetic contributions to the trait.Once the MLEs of the parameters are obtained, the presence of QTL responsible for the variation of the quantitative phenotype can be tested by using the following hypotheses, H 0 : a~d~i~i ca ~icd ~ici ~0 The test statistic for testing the above hypotheses is calculated as the log-likelihood ratio test statistic (LR) of the full model (H 1 ) over the reduced model (H 0 ): where and denote the MLEs of the unknown parameters under H 0 and H 1 , respectively.Since a genome-wide scan involves multiple correlated tests, we use the permutation tests proposed by Churchill and Doerge [36] to find the threshold value.
If there is a QTL, a number of other hypothesis tests can also be performed to test the property of the detected QTL.To test the imprinting effect, we can simply formulate the hypothesis as H 0 : i~0 vs H 1 : i=0 to assess the mean difference of the two reciprocal heterozygotes.To test the cytoplasmic maternal effect, the hypothesis can be stated as H 0 : c~0 vs H 1 : c=0.The epistatic effects of all interaction terms can also be tested as Similarly, additive and dominance effects can be tested as If specific interest is focused, for instance, the interaction of imprinting and maternal effect, the hypothesis can be formulated as H 0 : i ci ~0 vs H 1 : i ci =0.All the above tests can be done by applying the likelihood ratio test in which the test statistic asymptotically follows a x 2 distribution with degrees of freedom equal to the difference of the parameters under the null and the alternative hypotheses.For example, when testing H 0 : a~d~0, the LR test statistics is compared with the x 2 2 cutoff with 2 degrees of freedom.

Monte Carlo simulation
Monte Carlo simulations were performed to investigate the statistical behavior of our model.We simulated an F 2 population, with one half of the population coming from design S 1 and the other half coming from design S 2 .A genome with 100 cM long linkage group, composed of 6 equidistant markers, was constructed.The position of QTL was assumed to be located at 48 cM away from the first marker on the linkage group.The marker genotypes in the F 2 population were simulated by mimicking sexspecific recombination fractions in mice, i.e., r M ~1:25r P [33].A series of simulation study with different sample size (n~400 vs n~800) and different heritability levels (H 2 ~0:1,0:25,0:4) was conducted to examine the impacts of parameter spaces on parameter estimation and testing power.These simulation designs, which were aimed to give a better understanding of model performance under different situations, can provide biologists some empirical evidences to design their experiments.
In the simulation study, the residual variance was calculated under different heritability levels.For the F 2 design, the genetic variances of various terms can be calculated as follows: ci =2zci cd , and the broad sense heritability level can be expressed as For given genetic parameters and the heritability level, the residual variance can be calculated as s 2 e ~(s 2 Q zs 2 cQ )(1=H 2 {1){s 2 c , from which the phenotype data can be generated.
The MLEs of the QTL position and effect parameters, based on 200 simulation replicates under different heritability levels and sample sizes, are displayed in Table 2.The square root of mean squared error (RMSE) of parameter estimates are given in parenthesis to show the estimation accuracy.As we expected, the accuracy of parameter estimates increases as the sample size and heritability level increase.For instance, the RMSE of estimated QTL position decreases from 13.28 to 3.81, an 71% increase in accuracy when the sample size increases from 400 to 800 under a fixed heritability level 0.1.The other parameter estimates show the same pattern.If we increase the heritability level when the sample size is fixed, a clear reduction in RMSE can be observed.For example, with 400 samples, the RMSE of estimated QTL position decreases from 13.28 to 4.99, then to 2.77 as H 2 gradually increases from 0.1 to 0.4.From the decreasing RMSEs of the parameter estimates, we observed that simply increasing sample size is less efficient increasing heritability level in order to increase the precision of parameter estimation.Since high heritability corresponds to small environmental variability [37], reducing environmental variation is of more practically important than just simply increasing sample size.
Note that we did not list the estimation of the imprinting effect i and cyto-imprinting interaction i ci , which are not estimable under the F 2 design.The reason is that the imprinting direction cannot be inferred from the F 2 design [25].However, we can still conduct hypothesis test to infer the imprinting effect as well as its interaction with cytoplasmic effect.To further investigate the testing performance of cytoplasmic and imprinting effects, we introduced two proportions, namely g c and g i , where g c ~s2 c =s 2 G and g i ~s2 i =s 2 G .We can evaluate the test power under different cytoplasmic and imprinting effect sizes.Given all other genetic parameter values fixed (as shown in Table 2), simple algebra shows that the cytoplasmic effect c and imprinting effect i can be calculated for a given value of g c or g i , i.e., , and i~ffi Based on 1000 simulation runs under different heritabilities, sample sizes and variation proportions, the power of cytoplasmic effect test, imprinting effect test, interaction effects test and additive/dominance effects test are listed in Table 3.As we expect, the test power increases with the increasing of sample size and heritability level.For example, the cytoplasmic testing power increases from 0.663 to 0.905 as sample size increases from 400 to 800, a 36.5% increase in power for fixed g c ~0:1 and H 2 ~0:25.The same pattern is observed for the imprinting test.As the proportion of variance explained by the cytoplasmic and imprinting effect increases, the power increases accordingly.Noted that when both g c and g i are zeros, the testing power corresponds to the type I error rate for the corresponding factor.From the table we can see that the size of imprinting test is well controlled under different sample sizes and heritability levels.For the cytoplasmic effect, the size is a little inflated under n~400, but it gets close to the nominal 5% level as sample size increases to 800, especially under large heritability (e.g., H 2 ~0:4).In sum, the simulation evidences show that the model performs reasonably well in both parameter estimation and testing.

A case study
We applied the model to a published F 2 cross data set based on design S 1 and design S 2 aimed to find QTLs that contribute to variation in quantitative traits related to colitis severity in 1L-10deficient mice [38].The data contain 411 F 2 mice derived from inbred strains, where 203 mice are from design S 1 and 208 mice are from design S 2 .Ninety-one markers were obtained with an average length of ,15 cM spanning acroos the 19 autosome chromosomes.For more information about the data, the readers are referred to the original paper [38].
It has been reported that on average the female chromosome is 25% longer in genetic distance between homologous loci than the male in mice [38].The sex-specific recombination fraction, expressed as r M ~1:25r P , was reconstructed based on the marker information (see Cui et al. for details [25]).The method was applied to four phenotypes, cecum total score (CTS), spleen/body weight ratio (SBWR), mesenteric lymph node(MLN)/body weight ratio(MBWR), and secretory IgA (SIgA) level, where cecum total score was graded by using colitis-related criteria, including severity, hyperplasia, ulceration and the percentage of area involved.Box-Cox transformation was applied to all traits before fitting the Gaussian-mixture model.The genome-wide LOD profile plots for the four phenotypes are shown in Figure 3, where the solid blue curves correspond to the LOD values and the dashed red lines correspond to the 5% genome-wide threshold values out of 1000 permutations.The LOD score is calculated as log 10 LR, where LR is obtained from equation ( 9) to test the null hypothesis: H 0 : a~d~i~i ca ~icd ~ici ~0.
As shown in Figure 3, one QTL on chromosome 3 is detected for cecum total score trait, two QTLs on chromosomes 3 and one on chromosome 1 are detected for spleen/body weight ratio trait, three QTLs on chromosome 3 are detected for MLN/body weight ratio trait, and two QTLs located on chromosome 3 are detected for SIgA trait.The QTL located at 60.6 cM on chromosome 3 is common to three traits.The one located at 52.6 cM for SIgA trait is very close to it.It is highly possible that it is the same QTL that controls the four traits.Such a pleiotropic effect needs to be further evaluated.It should be mentioned that the QTL detected in the original paper for the four traits is located at 61.8 cM on chromosome 3 [38], which is 1.2 cM away from the one we found.Such a difference may arise because of the capitalization of sexspecific recombination rates and different models fitted.
In addition to the QTL identified in our analysis and the original paper, some other major QTLs which are not detected in Farmer et al. [14], such as those at 52.6 cM and 38.6 cM on chromosome 3, stand out in our model and therefore need further investigation.Almost all the QTLs on chromosome 3 are clustered together, whose local LOD profiles of four traits are shown in Figure 4.
The marker interval for each QTL is listed in Table 4, which also tabulates the p-values under four different tests for each estimated QTL using permutation tests.From the test results, it can be seen that most QTLs have strong additive and dominance genetic effect, except for the spleen/body weight ratio QTL located at D1Mit156+31.1 cM on chromosome 1.This QTL shows evidence of cytoplasmic, imprinting as well as cyto-nuclear interaction effects (p-values for the three tests are 0.0282, 0.0296 and 0.0073, respectively), but shows no sign of additive and dominance effect.In addition, the MLN/body weight ratio QTL located at D3Mit78+5.6 cM shows evidence of cytoplasmic effect.In summary, we identified one QTL on chromosome 1 with evidence of cyto-nuclear interaction effect and this QTL also shows evidence of cytoplasmic and imprinting effect.Further functional validation is needed to confirm the results.

Discussion
The cytoplasmic environment influences the expression of nuclear information in a very complicated way, which is still an unravel mystery to many organisms.For example, Burgess and Husband have demonstrated great maternal contributions to the fitness of mulberry hybrids [39].While it is an important parent-oforigin effect affecting offspring fitness, genomic imprinting, another source of parent-of-origin effect can also lead to phenotypic variation.Increasing evidence from cytoplasmic substation and cell fusion experiments suggests that weakness of hybrids may connect with the interactions between cytoplasm and nuclear [40,41], and the evidences about interaction between cytoplasm and imprinting have also been observed [42,43].As the source of genetic variation for many traits, these two types of parent-of-origin effects are often confounded, making it difficult to distinguish without proper statistical dissection.Although the role of cross-talk between the two sets of factors on phenotypic variation has been recognized, which genes are involved in the process and in what form they respond to the cytoplasmic changes are still unclear.
In this paper, we developed a statistical model to evaluate the cytoplasmic environment and nuclear gene interaction subject to imprinting effect within the framework of QTL interval mapping.: not all equal to 0; and 4) H 0 : a~d~0 vs H 1 : not all equal to 0, respectively.For a given g c , all other effect values are fixed as 0.8 except for c, which can be calculated in terms of g c and other parameters.The hypothesized c value is 0 for g c ~0, 0.461 for g c ~0:1 and 0.679 for g c ~0:2.Similarly, the value of i, which depends on imprinting effect variation proportion g i , is 0 for g i ~0, 0.680 for g i ~0:1 and 1.020 for g i ~0:2.doi:10.1371/journal.pone.0091702.t003The model that considers eight genetic factors which measure the degree of imprinting, cytoplasmic, additive, dominance effects as well as the interaction effects among them, provides a complete dissection of cyto-nuclear epistasis subject to imprinting effect.A number of hypothesis tests can be performed not only to assess major genetic effect(s) responsible for phenotypic variation, but also to find the statistical evidence for the existence of imprinting, cytoplasmic effect as well as the cyto-nuclear interactions.Simulation study showed relative good performance of the model under the F 2 design, in which parameters are estimated efficiently with modest heritability and sample size.Low heritability level (H 2 ~0:1) and small sample size (n~400) result in large mean squared error of parameter estimation.This result is valuable in practice as we need to be careful about the interpretation of genetic effects obtained in real data analysis when the proportion of variance explained by the QTL is small (i.e., low heritability).Although our model cannot estimate the imprinting effect (so the cyto-imprinting interaction effect) due to the nature of the F 2 design, existence test of imprinting (or cyto-imprinting interaction effect) can be achieved.Nevertheless, such imprinting estimation problem can be solved under the reciprocal backcross design illustrated in Figure 1.
In the real data analysis, one QTL located on chromosome 1 were found to have significant cytoplasmic, imprinting effect and cyto-nuclear interaction effects for spleen/body weight ratio.Other than that, no imprinting effect was found for all other QTLs, and only one on chromosome 3 that shows cytoplasmic effect for the MLN/body weigh ratio trait.It is worth mentioning that several QTL on chromosome 3 are detected by our model, but they are located relatively close to each other, as shown in Figure 3.For the sake of cautiousness, we reported all of them.However, these detected clustered QTLs may be due to the limitation of interval mapping, which can be overcome by fitting a composite interval mapping model as following, y j ~mzcx 1j zax 2j zdx 3j zix 4j zi ca x 5j zi cd x 6j zi ci x 7j z X K k~1 (a k x 2jk zd k x 3jk zi k x 4jk zi ca k x 5jk zi cd k x 6jk zi ci k x 7jk )ze j where x 2jk , Á Á Á ,x 7jk are corresponding variables for the kth marker, assuming total K markers are selected for controlling background genetic effect.Although more variables are introduced in the model, theoretically some dimension-reduction techniques such as LASSO, can be applied to implement the variable selection for each marker before fitting them into the final model [44].The composite interval mapping is know for its improved resolution in QTL detection.Regardless of the potential limitations mentioned above, the integration of cyto-nuclear interactions into the QTL mapping framework provides a testable platform with feasible experimental design for biologists to test the existence of cytoplasmic and imprinting effects, as well as the interactions of interested.The proposed model will have important biological implications with potentials to lift a corner of the great veil of the genetic system.

Figure 3 .
Figure 3.The LOD profiles of the four traits across the 19 chromosomes using the linkage map constructed from microsatellite markers [38].The genomic positions corresponding to the peak of the curve are the MLEs of the QTL locations.doi:10.1371/journal.pone.0091702.g003

Table 2 .
The MLEs of the QTL position and effect parameters based on 200 simulation replicates under different heritabilities and sample sizes.
The squared roots of the mean squared errors (RMSE) of the MLEs are given in parentheses.The locations of the QTL is described by the map distances (in cM) from the first marker of the linkage group.The hypothesized s 2 value is 3.81 for H 2 ~0:1, 2.04 for H 2 ~0:25 and 1.26 for H 2 ~0:4.doi:10.1371/journal.pone.0091702.t002