Sequence Variation and Haplotype Structure at the Human HFE Locus
Christopher Toomajian, Martin Kreitman

Abstract

The HFE locus encodes an HLA class-I-type protein important in iron regulation and segregates replacement mutations that give rise to the most common form of genetic hemochromatosis. The high frequency of one disease-associated mutation, C282Y, and the nature of this disease have led some to suggest a selective advantage for this mutation. To investigate the context in which this mutation arose and gain a better understanding of HFE genetic variation, we surveyed nucleotide variability in 11.2 kb encompassing the HFE locus and experimentally determined haplotypes. We fully resequenced 60 chromosomes of African, Asian, or European ancestry as well as one chimpanzee, revealing 41 variable sites and a nucleotide diversity of 0.08%. This indicates that linkage to the HLA region has not substantially increased the level of HFE variation. Although several haplotypes are shared between populations, one haplotype predominates in Asia but is nearly absent elsewhere, causing higher than average genetic differentiation among the three major populations. Our samples show evidence of intragenic recombination, so the scarcity of recombination events within the C282Y allele class is consistent with selection increasing the frequency of a young allele. Otherwise, the pattern of variability in this region does not clearly indicate the action of positive selection at this or linked loci.

HFE was the first gene to be associated with hereditary hemochromatosis, a recessive disease common in many populations of European descent and characterized by iron overload (Federet al. 1996). Recently, rare alleles in additional genes have been identified that are associated with the hemochromatosis phenotype (Roettoet al. 1999; Camaschellaet al. 2000; Njajouet al. 2001; Katoet al. 2001) and lead to distinct disorders. However, the majority of hereditary hemochromatosis cases in Europe are due to changes in the HFE gene (Beutleret al. 1996; Federet al. 1996; Carellaet al. 1997). HFE plays an important role in regulating iron levels in the body (Federet al. 1998; Salter-Cidet al. 1999). Although much progress has been made in determining the function of HFE, questions concerning its specific mechanism of iron regulation still remain (e.g., Drakesmith and Townsend 2000). The potential for iron requirements or availability to change in novel environments makes HFE a possible target of local adaptation.

A feature of HFE relevant to population genetic inference is its chromosomal context. HFE is found telomeric to the human leukocyte antigen (HLA) class-I region on chromosome 6p21 in an area referred to as the extended class-I region. The HLA has been the focus of polymorphism studies for decades, since the most polymorphic loci in the human genome are found here. Initial studies revealed that balancing selection has acted at a number of HLA genes (e.g., Hughes and Nei 1988), and its effect can be seen in the increase in variation at neighboring loci (e.g., Grimsleyet al. 1998; Sattaet al. 1999). However, as diversity studies have broadened to include additional loci, the observation of many peaks and valleys of nucleotide diversity suggests that simple models of diversifying selection acting on a handful of exons cannot explain the full complexity of variation in this region (e.g., Gaudieriet al. 2000). Even though the HFE gene is ∼4 Mb away from the highly polymorphic HLA-A locus, the genetic distance between them is ∼1 cM (Malfroyet al. 1997). In fact, the HFE gene was first localized to chromosome 6 on the basis of the association between hemochromatosis and HLA-A3 (Simonet al. 1976). Another report notes an extended HLA haplotype common in Europeans, A1-B8, which maintains associations with microsatellites distal to HFE (Worwoodet al. 1997). These associations between HFE and HLA alleles raise the question of whether selection on HLA alleles has influenced the pattern and level of variability found at HFE within and between populations.

Studies of HFE variation have focused on two amino acid polymorphisms that were discovered when the gene was mapped by Feder et al. (1996). One mutation, C282Y, disrupts an intramolecular disulfide bridge and renders the protein nonfunctional (Federet al. 1997). This mutation, found at a frequency of up to 10% in European populations (Merryweather-Clarkeet al. 1997), is by far the major mutation that leads to hemochromatosis. The H63D mutation is generally found at a higher frequency in Caucasians and also appears to be associated with hemochromatosis, although its penetrance is low (Risch 1997). Several other rare hemochromatosis-associated HFE mutations have been described (Bartonet al. 1999; de Villierset al. 1999; Muraet al. 1999). However, not all cases of hemochromatosis can be explained by the known HFE mutations, leaving open the possibility of additional minor disease-associated mutations. Effort to characterize HFE variation has concentrated on Caucasian populations, and the full spectrum of variation at HFE has not been considered. Therefore, it is difficult to make inferences about the forces governing this variation.

Several groups have proposed that selection has favored the C282Y mutation, but a detailed knowledge of the linked variation around this site is necessary to independently test this hypothesis at the nucleotide level. Two lines of evidence have led to the hypothesis of selection. One is based on the function of HFE, with the selective advantage for C282Y possibly stemming from its potential to prevent iron deficiency. The second is based on the seemingly incongruous observation that C282Y appears extremely young but is a relatively high-frequency mutation in European populations. Models of the decay of linkage disequilibrium (LD) over time estimate a young age (<100 generations) for the C282Y mutation (Ajiokaet al. 1997; Thomaset al. 1998). In contrast, the expected age of an allele at 5% frequency in a population with an effective size of 104 is >6000 generations (Kimura and Ohta 1973). Selection for the young C282Y mutation or hitchhiking with a linked mutation may have allowed it to reach its present frequency. While experiments to test the selective advantage of C282Y on the basis of functional differences are difficult and will not necessarily shed light on the historical reason for C282Y's high frequency, the patterns of nucleotide variability and LD in and around the HFE region can provide evidence for the selective advantage of particular mutations or the pattern of hitchhiking. For example, if a rapid change in the frequencies of alleles under diversifying selection at classical HLA loci were responsible for C282Y's high frequency, one might expect that other allele classes at HFE would display a lack of variation similar to that of the C282Y class.

In this report we describe the nucleotide variation and haplotype structure in HFE for a worldwide population sample. We test whether the pattern of variability is consistent with an equilibrium neutral model of evolution. We compare the level of population variation and differentiation seen at HFE with similar studies of human genes. The partitioning of linked variation into haplotypes allows us to address the related questions of the origin of different alleles and the forces that have acted to produce their current global distribution and frequency.

MATERIALS AND METHODS

DNA samples: A total of 30 samples (60 chromosomes) were chosen to represent ancestry from African, Asian, and European peoples. Identifiers in parentheses indicate the sample numbers from the Coriell Cell Repositories' National Institute of General Medical Sciences Human Genetic Mutant Cell Repository (http://arginine.umdnj.edu). Samples without these identifiers were either collected at the University of Chicago or provided by other labs. The 10 African samples include five Mbuti Pygmies (NA10492–NA10496) from the Ituri forest in northeast Zaire and one sample each of Kikuya (NA00522), Ghanaian (NA02064A), Zulu (NA02476),!Kung (NA03043), and Luo (NA03190A) descent. The 10 Asian samples include five individuals of Chinese descent (including NA11321–NA11323), two samples of Korean descent (including NA00726), and one sample each of Filipino (NA10798), Khmer (NA11373), and Vietnamese (NA03037) descent. The 10 European samples include four samples from the previous study of Edwards et al. (1988) that come from Utah and its surrounding states and six samples collected in Chicago from mixed European ancestry. The samples of Edwards et al. (1988) are members of hemochromatosis pedigrees that we selected on the basis of an apparently normal phenotype and the lack of the C282Y HFE mutation. All other samples were chosen either from healthy volunteers or from individuals with disorders that have no known association with variation in the HFE gene. In addition to the 10 European population samples, two unrelated samples of European descent (NA14620 and NA14621) that were homozygous for the HFE C282Y mutation were selected for sequencing. The Institutional Review Board of the University of Chicago approved this project.

PCR and sequencing: The region under study consists of bases 43,385–54,657 of the human hereditary hemochromatosis region (GenBank accession no. U91328; Laueret al. 1997); for convenience we change the base numbering to 1–11,273, respectively. A large number of primers were designed throughout this region from the reference sequence, which does not carry the C282Y or H63D mutation. The primers and conditions used for amplification and sequencing are available upon request. For 23 of the 30 samples, overlapping diploid PCR products of ∼1 kb were sequenced to determine the identity of nucleotides 31–11,244 (excluding external primer sequence). PCR products were used as templates for dRhodamine terminator cycle-sequencing reactions that were subsequently cleaned and run on an ABI 377 automated sequencer (Applied Biosystems, Foster City, CA). Chromatograms were imported into Sequencher v. 3.0 (Gene Codes, Ann Arbor, MI) for manual assembly of contigs and identification of polymorphic sites. Each base in the study was called, using at least single-fold coverage sequencing reads for each strand, except for a few small regions where sequence repeats made the reads in one direction of poor quality and bases were called using information primarily from one strand. Sequence from the HFE gene was also obtained from one common chimpanzee (Pan troglodytes) from DNA provided by Dr. D. H. Ledbetter. Most PCR and sequencing primers worked with the chimpanzee sample, and where gaps remained new PCR primers were designed.

For the remaining seven samples, a long-range PCR product of 11,273 bp was amplified with the Expand Long Template polymerase mix (Roche Molecular Biochemicals, Indianapolis), cloned with the Topo XL PCR cloning kit (Invitrogen, Carlsbad, CA), and the whole insert was sequenced directly. This method leads to the sequencing of PCR errors and may produce hybrid sequences derived from the maternal and paternal alleles. Therefore, we confirmed each difference from the reference sequence by either sequencing or DHPLC analysis (Transgenomic, Omaha) of smaller PCR products produced from genomic DNA.

Haplotype determination: For the 4 samples from the Edwards et al. (1988) study, pedigree analysis was performed to determine haplotypes. For the 19 other samples with more than one heterozygous site, the 11-kb long-range PCR product was amplified and cloned as described above. For at least two clones per sample, we directly sequenced small regions containing the heterozygous sites. Comparison of these sequences with the corresponding diploid sequence was used to determine if any clones were hybrids of the maternal and paternal alleles. This occurred for only one sample, and the sequencing of the heterozygous sites from a third clone resolved which clone was hybrid and which were true alleles. For the 7 samples that were cloned before being sequenced, initially one clone per sample was sequenced. To find the allelic clone for each of these samples, individuals were screened for heterozygosity at intermediate frequency polymorphisms by DHPLC. All 7 samples were found to be heterozygous at one or more sites, and a second clone for complete sequencing that carried the other base at a heterozygous site was chosen. To ensure that neither of the sequenced clones were hybrids, some sites that appeared homozygous and matched the published sequence after sequencing the two clones were tested for homozygosity by DHPLC. Analysis proved homozygosity for all sites tested in this way and demonstrated no further hybrid clones.

Data analysis: The program DnaSP, ver. 3.53 (Rozas and Rozas 1999) was used to estimate parameters and perform statistical analyses unless otherwise noted. Length polymorphisms in repetitive DNA (most frequently mono- or dinucleotide repeats) could not always be determined with certainty, are not reported here, and were ignored for these analyses. Also, length differences in mono- or dinucleotide runs between human and chimpanzee were not scored, and fixed differences that interrupt these repeat stretches have also been excluded. Only single-nucleotide polymorphisms have been included in summaries of nucleotide variability, and one singleton insertion/deletion (indel) polymorphism found in complex sequence (site 9681) has been excluded. For the computation of diversity estimates, the number of observed mutations is used in place of the number of segregating sites. The confidence interval for θ was calculated as described in Kreitman and Hudson (1991). Significance for Tajima's D (1989), Fu and Li's D (1993), and Fay and Wu's H (2000) tests was assessed by comparison to the output of neutral coalescent simulations of 103 random samples with identical sample size and polymorphism level as the observed data, assuming constant population size and no recombination (which makes the tests conservative). The H test was performed using a program provided by J. Fay. The program K-Estimator v5.5 (Comeron 1999) was used to test the difference between human-chimp divergence levels with 104 Monte Carlo simulation replicates. Noncoding regions conserved between human and mouse (see results) were compared to simulation results of random divergence between two sequences of the same length and G + C content as the HFE noncoding sequence (J. Comeron and M. Kreitman, unpublished results). Significance was assessed by measuring the longest region with at least 75% identity in each of 103 replicates. Hudson-Kreitman-Aguadé (HKA) tests were performed via coalescent simulation using the program of J. Hey (http://lifesci.rutgers.edu/~heylab).

The mutational relationships among the experimentally determined haplotypes were visualized by using the reduced median (RM) algorithm of the program Network 2.0 (http://www.fluxus-engineering.com; Bandeltet al. 1995). This algorithm is designed for use with nonrecombining DNA haplotypes, but it provides a convenient means of displaying haplotype relationships, especially for very similar haplotypes, when recombination is not very prevalent in a region. In constructing the network, the algorithm links haplotypes that differ at only one site and then assumes that mutational events proceed from a more frequent haplotype to a less frequent one to choose between equally parsimonious mutational paths linking more distantly related haplotypes. In cases where this assumption does not resolve alternate paths, the uncertainty is shown as a reticulation, or loop, in the network.

The program Arlequin 2.000 (Schneideret al. 2000) was used to perform an AMOVA analysis of the genetic differentiation among population samples. CHRM was estimated from coalescent simulations of a constant size population conditioned on the number of segregating sites with J. Wall's program hrmpg2, available on the Hudson lab homepage (http://home.uchicago.edu/~rhudson1/source/JWallCode.html). A total of 105 replicates are run for each tested value of C. Initially, 31 values of C were tested (0, 1, 2,..., 30) and then an additional 11–21 values were tested around the above maximum-likelihood estimate [e.g., 4.0, 4.1,..., 6.0 for an initial maximum-likelihood estimate (MLE) of 5]. ρCL (Hudson 2001) was estimated under a model without gene conversion with a program provided by R. Hudson. In assessing LD, we have assumed that at site 7633 a C to G mutation preceded a G to A mutation, while for site 3877 (also triallelic) the order of mutations is not clear, and we have excluded the site in populations segregating all three alleles. LD was calculated as D′ (Lewontin 1964) and significance was assessed by Fisher's exact test.

RESULTS

HFE diversity at the nucleotide level: A total of 41 variable sites were identified in the 11,214-bp region of the 60 chromosomes surveyed: 38 diallelic single-nucleotide polymorphisms (SNPs), two triallelic SNPs, and one diallelic single-base indel polymorphism. The two SNPs found triallelic in the pooled sample indicate that at least two mutations have occurred at these sites in the history of this sample. Of the 41 polymorphic sites, 8 segregate singletons, with the more frequent allele matching the chimpanzee sequence in each case. Only two SNPs were in exons: SNP 6724, which causes the H63D amino acid polymorphism (Federet al. 1996), and SNP 3470, which is a synonymous change. The bottom of Figure 1 indicates the location of all polymorphisms relative to the intron-exon structure of HFE.

Summary statistics describing the sequence diversity in the pooled and individual populations are presented in Table 1. Overall, average per-nucleotide expected heterozygosity, θw, for the total sample, estimated from the observed number of mutations (Watterson 1975), was 0.080% [0.040–0.148%, 95% confidence interval (C.I.)]. Nucleotide diversity (π), an estimate of θ based on the average pairwise sequence difference (Nei 1987), was nearly identical (0.084%). When the coding sequence is excluded, π increases by ∼8% to 0.091%. These estimates of θ assume an infinite-sites model, which is clearly violated since the data contain two SNPs that are triallelic. However, estimates of θ derived from finite-sites models lead to only negligible differences from the infinite-sites estimates for our data (Tajima 1996). Even when we allow the mutation rate to vary extensively among sites (with mutation rates following a γ distribution with parameter α= 0.1), the finite-sites estimates of θ are not substantially changed (π increases from 0.084 to 0.085%). Our estimate of θ for HFE is slightly lower than the average for fourfold degenerate coding sites in humans (0.11%; Li and Sadler 1991), but conforms to this value and the average estimated for similar gene regions in humans (0.081%; Przeworskiet al. 2000), for the whole genome (0.075%; International SNP Map Working Group 2001), and for chromosome 6 in particular (0.074%; International SNP Map Working Group 2001) on the basis of our confidence intervals.

Figure 1.

—Observed HFE haplotypes. Bases identical to the chimpanzee haplotype are marked with a dot. Haplotypes are numbered 1–19 to indicate relative frequency from highest to lowest and ordered so that closely related haplotypes are near each other. Population counts are shown at the right. Polymorphisms are grouped by different regions of the HFE locus, as indicated in the gene diagram at the bottom. The arrow marks the direction of transcription and its start site (the 5′ of the gene is at the right); solid regions are either flanking sequence or introns, striped regions are coding exons, and checked regions are the 3′ and 5′ UTRs. The locations of polymorphisms are shown below the gene diagram, with the location of C282Y (4762 in our notation; not observed in our population sample) indicated by a star.

Diversity in continental populations: When the coding region is excluded, π increases by nearly the same percentage (∼8%) for each population. As is commonly observed, the Africans have the highest nucleotide diversity and the largest number of population-specific SNPs, at 11. However, only 2 of these 11 are singletons, so that the value of θw, which can be greatly influenced by rare variants, does not rise above that of π. The nucleotide diversity for Europeans is only slightly lower than that for Africans, but they have many fewer population-specific SNPs than the Africans (4 vs. 11, Table 1). European HFE variability is not strictly a subset of that found in Africa, consistent with the previous finding of HFE polymorphisms with an apparent European origin (Fairbanks 2000). The Asian samples have the lowest nucleotide diversity, which is still 65% of the African value, indicating that all populations studied have substantial HFE variability. Asians also have relatively few population-specific polymorphisms (4, Table 1), but in this case, they are all singletons, contributing to a higher θw relative to π.

Allele frequency spectrum: Test statistics that utilize the frequency spectrum of alleles within a locus may detect departures from an equilibrium neutral model caused by demographic forces such as population growth, contraction, and subdivision or by the effect of diversifying or directional selection on linked sites. Tajima's (1989) D statistic compares the two estimates of θ described above, while Fu and Li's (1993) D statistic compares the number of singleton polymorphic sites with the estimate of θ based on segregating sites and incorporates outgroup information to infer derived alleles. Additionally, Fay and Wu's (2000) H statistic can detect departures in the frequency spectrum due to recent hitchhiking events. None of these three statistics, calculated for each population and for the pooled sample (Table 1), are significant at the 5% level, so the equilibrium model of neutral evolution cannot be rejected. However, it should be noted that the small size of the individual populations (n = 20) limits the power of these tests (Simonsenet al. 1995). The pattern observed for the different populations is informative, though, since negative values of Tajima's and Fu and Li's test statistics indicate an excess of low-frequency polymorphisms, which is expected under a model of human population growth from its ancestral size. Only the Asian population shows a slight excess of low-frequency polymorphisms. Among their variable sites, 17 have a minor allele frequency of 10% or less in Asians, with a comparison to chimpanzee sequence implying the rare alleles are derived. However, 8 of the remaining 10 derived SNP alleles have a frequency in Asians of >50%. It is possible that a selective sweep or hitchhiking has boosted these 8 SNP alleles to high frequency. This pattern in Asians may be consistent with drift and a stronger founder effect than is seen in the other populations. For the African, European, and pooled populations, these test statistics show no clear sign of population expansion when all observed SNP alleles are assumed neutral.

View this table:
TABLE 1

Diversity estimates and neutrality tests for HFE

Divergence from chimpanzee and haplotype variation: The SNP alleles in the total sample of 60 chromosomes were found to occur in 18 distinct haplotypes (19 if the singleton indel is included). These haplotypes are displayed in Figure 1 along with a haplotype composed of the ancestral state of each allele inferred from chimpanzee (P. troglodytes). The chimpanzee sequence for the complete region revealed 71 fixed differences from human, including 69 SNPs, a 2-bp indel, and a complex mutation involving a base change and a single-base indel at a neighboring site. Only one fixed difference was found in the coding region, and this was a synonymous change. The average number of nucleotide substitutions per site between human and chimpanzee is 0.690% (0.750% for noncoding sequence). A recent analysis of divergence levels between humans and chimpanzee reports an average distance of 1.03 ± 0.04% for introns, from 32 loci, with a combined length >41 kb once repetitive regions were removed (Chen and Li 2001). The same calculation for nonrepetitive sequence in HFE introns (4032 bp) is 0.623%, which is significantly lower than this average under a model where the mutation rate is constant among sites (P = 0.006). However, intron divergence values in 4 of the 32 loci from the Chen and Li study (including one region of 9556 bp) are lower than that for HFE, suggesting that the low divergence observed for HFE introns is not exceptionally low relative to the actual distribution of intron divergence values.

The low observed divergence might result from a high average degree of constraint for the whole HFE region, which is composed primarily of introns, untranslated regions (UTRs), and intergenic sequence. To address this question, we investigated divergence from the homologous sequence in mouse. Sánchez et al. (1998) reported conservation between human and rodent in the HFE promoter region, and we estimated conservation with the published mouse HFE sequence (GenBank AF007558) across our entire 11-kb region, using available global alignment and visualization tools (Mayoret al. 2000). Because the extent of divergence for noncoding sequence between human and mouse prevents the creation of a definitive global alignment, one cannot calculate a net divergence in the same way one would for the human/chimp comparison. Aside from the exons, two conspicuously conserved regions were found in noncoding sequence, one 74.6% identical over 114 bp in the 3′ UTR and one 74.5% identical over 188 bp in intron 5. These regions are significantly longer (P < 0.001) than what is expected for unconstrained regions of equal or higher identity due only to common ancestry and given the observed human/mouse KS for HFE (0.65). They provide strong evidence that some functional constraint exists outside the coding region that could contribute to the low divergence observed between human and chimp. The 3′ UTR conserved region contains two fixed differences between human and chimp while neither region contains polymorphisms in humans, although divergence with chimp and polymorphism within human throughout the whole region are too low to provide evidence for or against functional constraint of specific subregions.

View this table:
TABLE 2

Results of HKA tests

In addition to functional constraint, a low neutral mutation rate could contribute to the low divergence. Both mutation rates and sequence divergence are influenced by G + C content (Wolfeet al. 1989), but HFE has an intermediate G + C content (45.1%) so that low divergence caused by a low mutation rate is not expected. If the low divergence did reflect a relatively low neutral mutation rate, then we might expect to see a level of polymorphism lower than the average level observed. HKA tests (Hudsonet al. 1987) comparing the noncoding portion of the HFE region with the 10 “locus pairs” of Frisse et al. (2001) for each population are shown in Table 2. The locus pairs are noncoding and thought to evolve neutrally, providing a suitable reference for the HKA test. Although HFE has a polymorphism to divergence ratio higher than that of the pooled locus pair data, particularly in Asians and Europeans, at least 1 locus pair has a ratio higher than that of HFE in each population. None of the HKA tests are significant, so we cannot conclude that the ratio of polymorphism to divergence is significantly higher for HFE.

For every SNP in humans, one of the alleles was present in the homologous position in chimpanzee. In all but three cases, the more common SNP allele in the pooled sample corresponds to the inferred ancestral allele. For these exceptions (SNP alleles 11204C, 519A, and 7451A), the derived allele frequencies are 58, 67, and 72%, respectively. Inference of ancestral state based on only one outgroup can be incorrect, but at least 3 derived neutral alleles out of 42 are expected to have >50% frequency in a sample of this size. No haplotypes in humans have the same configuration at the 41 polymorphic sites as has the chimpanzee; that closest to this configuration is haplotype 9, found exclusively in Africans, which differs at 3 polymorphic sites and has an additional 71 fixed differences from the chimpanzee sequence.

Because the C282Y hemochromatosis mutation was not observed in the random population samples, we sequenced two Caucasians known to be C282Y homozygotes and observed that this mutation occurs on the haplotype 3 background. Haplotype 3 is found in all three continental populations and is the second most common haplotype in Europeans. No additional polymorphisms were discovered by sequencing the complete 11,214 bp of these two homozygotes (four chromosomes), consistent with the previous conclusion (e.g., Ajiokaet al. 1997) that samples carrying C282Y have a relatively recent common ancestry.

Number of haplotypes: Fu's (1997) Fs statistic compares the observed number of haplotypes in a sample to the number expected assuming an infinite-sites model of mutation under neutrality and no recombination and is useful for detecting population growth or hitchhiking. In each population as well as the pooled sample, no excess haplotype variation is observed at HFE, as seen by nonsignificant Fs values for all populations (see Table 1). In fact, in each case Fs is positive and therefore indicates a deficit of haplotypes given the observed level of nucleotide diversity. This finding is surprising, as both recurrent mutation and recombination have likely affected the samples, and both are expected to produce additional haplotypes. The deficit of haplotypes in this case suggests that either recombination and recurrent mutation have not increased the haplotype diversity of these samples greatly or other forces have kept the number of observed haplotypes low. However, Strobeck's statistic (1987), which tests for the opposite pattern of observing too few haplotypes, is also not significant for any population (P > 0.05). Again, no evidence for population expansion is apparent based on haplotype diversity in each population.

Haplotype network: Figure 2 displays an RM network of haplotypes constructed from the pooled samples. Most haplotypes are unique to one particular population, since the relatively small sample size of each population reduces the chance of including rare haplotypes in each population sample. But other features, such as a branch that leads to four haplotypes found exclusively in Africans and representing 40% of all African samples, suggest the population differentiation seen in this sample may be real. The branch leading to the chimpanzee haplotype (Figure 2, arrow) contains the fixed differences as well as site 11,204, which is polymorphic in humans. The presence of recombinant haplotypes complicates the inference of haplotype relationships and results in mutations that have occurred only once in the history of a sample to be displayed multiple times in the network. Of the 43 mutations inferred from the human sample (including the indel), 10 are found twice on the network. Excluding the loop, 5 mutations show up twice on the network, indicating either recurrent mutations or recombination events. The chance of recurrent mutation is appreciable, since 144 CpG sites are found in the region studied and two nucleotide sites segregate three alleles each. In fact, CpG sites have a mutation rate that is estimated from our data to be >15 times higher than that of other nucleotide sites, similar to results for the LPL gene (Templetonet al. 2000). A similarly high rate of mutation at these sites is apparent from the divergence data with chimp. Both triallelic SNPs are found in CpG sites, with 1 of the 2 inferred mutations at each site consistent with a transition from 5-methylcytosine to thymine. Of the diallelic sites that occur more than once in the network, 3 out of 10 are CpG sites but none are consistent with this type of transition. This makes recombination a more likely cause of the observed homoplasy. By pruning certain haplotypes from the network we can reduce the number of homoplasic sites and find evidence for which haplotypes are recombinant. Potential recombinant haplotypes include 17 and 19, both occurring only once and generating recurrent mutations. Haplotypes 1, 12, and 18 may all stem from one or more recombination events, as their exclusion will remove most instances of recurrent mutation. Potentially ancient recombination events that created haplotypes found at high frequency today make it difficult to produce a more accurate genealogy of HFE alleles. But the network facilitates the design of strategies to infer HFE haplotypes from samples by typing only a few polymorphisms.

Population subdivision: To investigate population differentiation at HFE, we describe the unequal distribution of variation among populations. In the sample, 19 (44%) derived SNP alleles were found in all three populations, while 11 were restricted to African populations, 5 to Asian populations, and 4 to European populations (a total of 20 or 47% restricted to one population). Also, European and Asian populations shared 4 SNP alleles that were not found in Africa, while Africans did not share any SNP alleles with only Europeans or Asians. This supports the shared ancestry of the Asian and European samples after their split from African populations, which is also apparent since Europeans and Africans share 1 haplotype in common and Europeans and Asians share 2. When alleles are grouped together as haplotypes, each haplotype is expected to have a more limited distribution. We see 14 haplotypes (74%) unique to single populations (7 in Africa, 3 in Asia, and 4 in Europe) while all three populations share only 2 haplotypes.

Figure 2.

—RM network of HFE haplotypes. Mutational relationships are indicated by lines linking the 19 unique haplotypes, represented in the network as circles. The size of each circle is proportional to the relative frequency of the haplotype in the total sample. Circle fill patterns show in which population(s) each haplotype was observed, with pie diagrams used for haplotypes found in multiple populations. Mutational differences between haplotypes are indicated on the branches of the network (amino acid replacement mutations are blocked, and mutations found on the network more than once are underlined). The arrow points to the node where the chimpanzee haplotype connects to this network.

Wright's (1931) FST statistic serves to quantify population differentiation by expressing the genetic variance among populations divided by the genetic variance of the total population. Using an AMOVA analysis based on polymorphic sites and in which our samples were split into three continental groups (Weir 1996), we estimate an FST value of 0.23 (P < 0.00001, by haplotype permutation). This value exceeds the average calculated for many other genes (Cavalli-Sforzaet al. 1994). The FST estimate based on treating the region as a single locus with many alleles corresponding to each different haplotype is 0.17 (P < 0.00001). For both calculations, continental populations are found to be significantly differentiated at HFE. Table 3 shows two estimates of pairwise population FST values. The comparison of African and Asian samples produces the highest values. The Asian samples are distinguished by the extremely high frequency (13/20) of haplotype 1, which is absent from Africans and at low frequency in Europeans. The African samples tend to carry haplotypes unique to Africa (15/ 20), and some of these haplotypes are not closely related to European and Asian haplotypes (see Figure 2). Both of these patterns contribute to the higher-than-average level of subdivision found at HFE.

Recombination and linkage disequilibrium: Since we have unambiguously determined the haplotypes at HFE, we are more confident in making conclusions about the evidence for recombination and LD in the region. When four gametic combinations are found in a sample of two-site haplotypes, this is an indication of recombination (crossing over between the two sites or gene conversion) or repeated mutation. Even with moderate mutation rate heterogeneity, one can assume the probability of repeated mutation at any one site is very small since this probability is no greater than that of a single mutation at the same site. A moderate level of recombination is consistent with the results of this four-gamete test, in which 25/528 site pairs (excluding singletons) have all four gametes found in the pooled sample (Table 4). A minimum number of recombination events (RM) can be inferred from the data to explain all instances of four gametes (Hudson and Kaplan 1985). The Asian and European populations show a number of site pairs with four gametes and, therefore, nonzero RM's (Table 1). Although no site pairs have four gametes in the African sample (RM = 0), recombination may still have occurred. In the pooled sample, one inferred crossingover event is between sites <1 kb apart (6567 and 7451). All instances of site pairs with four-gamete types in which the sites flank both sides of this interval (22/25) could be explained by this one recombination event.

View this table:
TABLE 3

Pairwise estimates of population subdivision (FST)

View this table:
TABLE 4

Significant pairwise linkage disequilibria and four-gamete site-pair counts

Since RM gives only a lower bound on the amount of recombination, we can use other methods to estimate the population recombination parameter C (= 4Nr, where N is the effective population size and r is the perlocus recombination rate per generation). Estimators of C that use patterns of sequence variation can avoid inaccuracies due to local variation in recombination rates found in estimates based on observed crossingover events between distant markers. CHRM, which uses the observed number of haplotypes and RM from data to estimate C under a model without gene conversion, performs well against other estimators of C (Wall 2000). It is ∼5 for HFE in different populations, although it is 0 in Africans, due to their lack of fourgamete site pairs (Table 1). Hudson (2001) has recently developed a composite-likelihood method for estimating C, ρCL, which performs better than another commonly used estimator of C (Hudson 1987) and that Frisse et al. (2001) have used to estimate crossing over and gene conversion rates in the human genome. When the gene conversion rate is held at 0, then this method's estimates of C are in good agreement with those produced by CHRM (Table 1). In the African sample, the likelihood of the ρCL estimate is not much greater than that of ρCL = 0, indicating little evidence for recombination and consistent with the CHRM value. This result is unusual, as the study of Frisse et al. (2001) found that this same estimate for the pooled data of 10 loci is much higher in Africans than in non-Africans. Their higher population recombination parameter estimate in Africans is also consistent with the finding that LD decays more rapidly in Africans than in non-Africans (Tishkoff et al. 1996, 1998, 2000; Kidd et al. 1998, 2000; Mateuet al. 2001; Reichet al. 2001).

A moderate level of LD is observed throughout the region, consistent with the estimates of the population recombination parameter. Due to the low frequency of most polymorphisms and our modest sample size, most pairwise LD comparisons (Table 4) do not have the power to detect significant LD at the 0.001 level (a value close to a global 0.05 level when multiple tests are considered). The pooled sample has the most power to detect LD, but pooling can also cause spurious LD due only to allele frequency differences between populations. D′ values for all alleles that are found at least twice in each of a pair of populations are highly correlated between populations (African vs. European, r = 0.999; African vs. Asian, r = 0.888; Asian vs. European, r = 0.718), with no cases of significant LD in opposite directions for pairs of populations. The high correlations are not surprising since most pairs of alleles are in complete LD in each population.

Figure 3 shows the location of site pairs in significant LD for the European population. LD appears evenly distributed throughout the region, with a minor concentration of significant LD at the 5′ end of the gene, particularly among two sites in the first intron and two sites in the 5′ flanking region, which have a perfect association (sites 9013, 10,047, 10,701, and 11,204). Plots from the other populations (not shown) reveal a similar, even distribution of LD, providing no evidence of a recombination hotspot within this region. The rate of decay of LD with distance varies depending on the measure of LD used and the minimum frequency of variants included in the analysis. We have plotted |D′| from the pooled population vs. distance for variants at 25% frequency or greater (Figure 4). This plot indicates that LD falls to one-half of its maximum value at ∼6 kb, a rate of LD decline that lies in between that of the ACE and LPL regions (Figure 1 in Przeworskiet al. 2000). Intermediate frequency variants tend to be the oldest ones in the total population, and therefore time has allowed more recombination events to break up their association. For lower-frequency variants, almost all site pairs have |D′| = 1. As expected, when a broader range of SNP frequencies is considered (10% and higher, see Figure 4), the distance at which LD reaches one-half of its maximum value substantially increases (>11 kb).

Figure 3.

—Significance of pairwise linkage disequilibrium for the European samples. Cell shading highlights site pairs with higher significance levels (as assessed by Fisher's exact tests without any correction for multiple testing). Singleton sites have been removed. Site 7633 is triallelic in Europeans; given our data we infer the order of the two mutations at this site to be C to G and then G to A. We therefore code all 7633A alleles as 7633G alleles with another mutation at the neighboring site, for the purpose of measuring LD.

Figure 4.

—Pairwise linkage disequilibrium as a function of physical distance. LD is measured by |D′|, which varies between 0 and 1, calculated for the pooled sample. Solid circles show |D′| only for pairs of polymorphic sites with both minor allele frequencies at 25% or greater. Solid triangles represent the average of these points in windows of 1 kb. Open triangles represent the average |D′| in windows of 1 kb when polymorphisms with minor allele frequency at 10% or greater are considered.

DISCUSSION

We have performed a full resequencing survey of nucleotide variation at the HFE locus using nonclinical samples from three major human population groups. The benefit of including noncoding variation in the study is that these polymorphisms are often more numerous and found at a broader range of frequencies, making them more informative in resolving which evolutionary forces have affected patterns of variation and governed the fate of alleles that alter the amount or function of a protein. Another significant aspect of our study is the experimental determination of haplotypes for an autosomal gene. These haplotypes provide a level of resolution greater than that of SNPs alone when drawing conclusions from genomic variability.

Summary of HFE variation: The results of any survey of population variation can be roughly separated into three categories: the level of variation, the frequency spectrum of that variation, and the haplotype structure of the variants. These categories are not independent, as they all reflect the underlying population history of a sequence of DNA, but they do capture different aspects of the data. Before this study, we had little information about the level of nucleotide variation at HFE. Protein polymorphisms seemed few and rare, but this could be due to high conservation imposed by the protein function and the slightly deleterious nature of most amino acid substitutions. The HFE gene is thought to lie in or near a region of low recombination, as Malfroy et al. (1997) have estimated a rate of 0.2 cM/Mb for the region between HLA and HFE. Despite this fact, we found that its level of polymorphism is about average for the genome. Both hitchhiking and background selection are expected to lower neutral variation in regions of low recombination and produce a positive correlation between levels of recombination and variation (Maynard Smith and Haigh 1974; Charlesworthet al. 1993). Evidence for this correlation in humans is debated, but the increasing number of polymorphism studies at individual loci and SNP data gathered from the Human Genome Project may resolve this debate (reviewed in Nachman 2001). If the correlation holds in humans, this suggests that either recombination in and directly flanking HFE is not severely reduced or a simple model of either genetic hitchhiking or background selection may not apply to this region. This second alternative might result from the nature of diversifying selection acting on the linked HLA region. This could produce a deeper genealogy for the HFE locus and raise its level of variation. When contrasted with the below-average divergence from chimpanzee, HFE variation appears slightly increased. This is not reflected in the HKA test and due to the stochastic nature of the mutation process may not have any real biological basis.

Test statistics reveal no major deviations from the equilibrium neutral expectation in the frequency spectrum of SNP alleles, although when the haplotype structure is considered, the Asian population shows an unusual pattern that is not seen in Africans or Europeans and may be consistent with a founder effect or hitchhiking. Thus, HFE provides no evidence for the long-term growth of the human population. However, on the basis of their study, Frisse et al. (2001) conclude that non-African populations are not in equilibrium. They summarize results that support a bottleneck in these populations but point out that more complicated demographic scenarios must be invoked to account for all of their non-African data. This appears true of the HFE data as well, as our non-Africans do not show the expected deficit of rare variants produced by a bottleneck.

HFE haplotype structure reveals some evidence for recombination, although fewer haplotypes than expected are observed on the basis of the number of variants. Przeworski and Wall (2001) have noted that estimates of C derived from sequence variation for human genes tend to be much higher than those based on experimentally measured crossing-over rates. This pattern is true for the HFE gene if we assume an effective population size of 104 and a crossing-over rate that is lower than the genomewide average (0.2 cM/Mb for the region between HLA and HFE; Malfroyet al. 1997). But, if the crossing-over rate in the HFE gene were 1 cM/ Mb (five times the estimated value), then the estimate of C would be 4.5, in good agreement with CHRM. Recent studies (Ardlieet al. 2001; Frisseet al. 2001; Przeworski and Wall 2001) find evidence that actual sequence data sets on an intragenic scale are more likely under a model of gene conversion and crossing over than under a model of crossing over alone. The data for HFE do not provide sufficient power to reject a model of crossing over without gene conversion, although gene conversion is known to have played a role in the allelic diversity of HLA genes (e.g., Zangenberget al. 1995). Our study detects the second most common HFE variant, H63D, at a moderate frequency in Europeans. Other reports have found this variant outside of Europe at a frequency consistent with gene flow from the Mediterranean region. We find this mutation on two different haplotypes, consistent with the conclusion that this allele has an origin much older than that of C282Y. The fact that HFE variation fits an equilibrium neutral model does not conclusively resolve the question of the history and potential fitness effects of HFE amino acid polymorphisms. Amino acid polymorphisms represent only a small proportion of HFE variation. Additionally, if the C282Y allele has been the target of positive selection, it is still far from fixation in any population, and therefore we expect the signature of this selection will be very subtle. Diversifying selection rapidly changing the frequency of alleles at linked HLA loci could have affected the frequency of several alleles at the HFE locus. But after investigating the pattern of HFE haplotype diversity, only haplotype 1 in the Asian populations (discussed below) indicates the lack of variation relative to the frequency of an allele class expected under this type of scenario, providing little evidence that the C282Y allele has increased in frequency due to this effect. An obvious conclusion based on the identity between the most common human HFE variant and the sampled chimpanzee protein is that no directional selection has been acting on the protein over the past few million years, although this does not exclude the possibility of much more recent or weaker selection on polymorphism.

The effect of sampling strategy: The structure of our population samples represents a compromise between sampling intensely from a few populations and broadly surveying many populations. For a number of analyses, samples are pooled on the basis of their continental origin as is frequently done for humans. This can have an effect on results when significant genetic subdivision exists among populations within a continent. Many genetic surveys are consistent with subdivision stronger in Africans than in Asians or Europeans (e.g., Jordeet al. 2000), although Zhao et al. (2000) suggest that the frequency spectrum of variation of a presumably neutral noncoding region on chromosome 22 indicates the opposite pattern. As implied by the study of Fu (1996) and demonstrated by Yu et al. (2001), population subdivision tends to reduce the proportion of low-frequency variants and therefore make Fu and Li's D positive. A positive DFL is found only in our African samples, while only the Africans have a significant FST value when it is calculated for individual populations. This indicates that the effect of population subdivision within continents may be limited to results from these samples. Additionally, a number of points indicate that our sampling scheme is not solely responsible for the observed patterns of variation within and between populations. First, for each continental group the sampling scheme is similar, with approximately one-half of the samples coming from one focal population and the other one-half pooled from several populations (or in the case of Europeans, from samples of highly mixed ancestry). Second, unusual patterns in our data are not due to the choice of the focal population for each continental group. In general, focal populations have slightly fewer SNPs and haplotypes than the other samples from the same continent, as would be expected given some amount of population differentiation within each continent. For example, the Asian samples have the lowest nucleotide diversity. This does not result from low diversity found only in the focal Chinese population, as the non-Chinese Asian samples have only a slight increase in nucleotide diversity (0.00058 vs. 0.00053 for the Chinese). The Asian sample is also characterized by low haplotype diversity, even given their low nucleotide diversity. The Chinese samples contain as many different haplotypes (four) as the non-Chinese Asian sample, and the high frequency of haplotype 1 is not unique to the Chinese sample (0.6 in Chinese vs. 0.7 in non-Chinese Asian samples). However, whether different sampling designs can lead to quite different conclusions about human diversity is a question worthy of further study.

Distribution of variation among continents: Relative to other species, the amount of genetic diversity that exists between different human populations is low. Cavalli-Sforza et al. (1994) have reported average values of FST for a large number of DNA polymorphisms (0.139) and non-DNA polymorphisms (0.119). In most cases, less differentiation is seen between populations within continents than between continents, consistent with simple isolation by distance. Wakeley (1999) has found evidence that, as might be expected, the level of migration in humans underwent an increase in the past. At the HFE locus, significant genetic differentiation is seen among continental populations, and FST is higher than the average values reported by Cavalli-Sforza et al. (1994). When FST for HFE is recalculated by including only the Mbuti, Chinese, and Utah samples it is actually higher (0.27 vs. 0.23), so that this high value is not an artifact of the diverse sample composition from each continent. However, these populations do not necessarily represent all of the variation within each continent, making a comparison with the results from a larger number of populations difficult.

The comparison of variability in different populations can provide evidence for local adaptive evolution but can be complicated by population history. The greater variability of populations from sub-Saharan Africa is seen in almost every study, and African samples have been contrasted to non-African samples to reveal differences in their population history. In a recent study of >300 genes (Stephenset al. 2001), the African-American sample was found to have >1.5 times the number of both SNPs and haplotypes than either the Asian or the Caucasian sample. Also, the African-American sample had >2 times the number of rare alleles than the other samples. Results from HFE show the most variation in African samples, but this difference is minor, with Asians and Europeans having nearly as many SNPs. Africans also have fewer rare SNPs compared to the other populations, in contrast to the Stephens study (population subdivision within our African samples could contribute to this result). Thus, the contrast between African and non-African samples is very subtle, becoming apparent only in the measure of nucleotide diversity since the African SNP alleles tend to be shifted toward intermediate frequencies. Due to the similar level of variability in the different populations and the finding of SNPs and haplotypes specific to each population, non-African variation is not just a subset of African variation.

When their haplotype structure is considered, the African samples are unusual. Although they do have the most haplotypes, this set of haplotypes is consistent with no recombination. Non-African populations show evidence of recombination, and this indicates that LD would decay over distance slower in the African population than in the other populations. Just the opposite is observed in other studies (Tishkoff et al. 1996, 1998, 2000; Kidd et al. 1998, 2000; Frisseet al. 2001; Mateuet al. 2001; Reichet al. 2001). Thus, a contrast in Africans between estimates of the effective population size based on variation vs. recombination data would suggest a departure from an equilibrium neutral model of evolution. This finding is remarkable, since similar departures are not often found in African populations (e.g., Frisseet al. 2001). Wall (2001) has reported that population structure and population bottlenecks can lead to underestimates of the population recombination parameter. Our finding could result if our sampling strategy has captured population structure in Africa more than previous studies. However, of the above-cited references, most share our strategy of broadly sampling multiple populations within continents.

The Asian samples provide a more striking contrast to the other populations. They have the lowest variation, due to both fewer observed haplotypes and the low frequency of many SNPs. Over 60% of SNP alleles in Asians are rare (at 10% frequency or lower). Remarkably, when the direction of mutation is inferred from the chimpanzee sequence, nearly 30% of Asian-derived SNP alleles are at a frequency >60%. This proportion of high-frequency-derived alleles could result from a past hitchhiking event, although the H statistic does not reach significance. The low number of haplotypes observed in Asians and the predominance of haplotype 1 (rarely found outside of Asia) suggest that this haplotype has risen to a high frequency rapidly since the split of the Asian populations from the other populations studied. The large number of low-frequency variants could mean that, had a selective sweep occurred, it happened long enough ago that many new segregating sites have since arisen and thus decreased the power of the H test (Przeworski 2002). However, when the pattern of haplotypes in all samples is considered, the haplotype diversity of Asia is not consistent with mutation arising solely from the haplotype 1 background, but instead suggests haplotype 1 never reached fixation in Asians or migration from other populations has introduced new haplotypes. The pattern of haplotype variation in Asians is also consistent with drift following a relatively strong founder event for this group.

The high frequency of haplotype 1 reveals the advantage of haplotype data over SNP sharing when assessing the similarity of populations. Haplotype 1 can explain a great deal of the population subdivision seen at HFE and can provide clues to how the interaction of selection with population history may have differed for individual populations. Haplotype 1 either originated within Asia and rose to a high frequency, spreading into Europe at low frequency, or occurred in a more ancient population but became prevalent in Asian lineages only through some sort of founder effect due to a bottleneck or hitchhiking. A broader sampling of populations could help resolve the history of this haplotype. Although not as informative as complete haplotype data, one SNP allele unique to haplotype 1 in our sample, SNP 4600G, has been genotyped in other populations. Rochette et al. (1999) have observed this allele at 8% in French, 20% in Sri Lankans, and 32% in Burmese. This report is consistent with our observation of haplotype 1 at low frequency in Europeans and also suggests an increasing gradient in this allele moving east in Eurasia. Beutler and West (1997) have also found a high frequency of this SNP allele in Asians (63%), although they observe it at nearly 13% in Europeans and at 15% in a small sample of African-Americans. These higher frequencies outside of Asia make the Asian-specific origin of SNP 4600G less likely. The Beutler and West study infers HFE haplotypes defined using SNP 4600 and two additional sites. In Asian samples 4600G is always (34/ 34) found on a haplotype consistent with our haplotype 1, while in the Caucasian samples, 4600G is inferred to appear on haplotype 1 in 7 of 9 cases. This pattern is consistent with a more recent increase in 4600G in Asia, since it has not had time to recombine onto other haplotypes, and supports the idea of a founder effect following a bottleneck or a hitchhiking event.

Conclusions: The discovery of these additional polymorphisms in the HFE region and the haplotypes they create may help identify other alleles that have an effect on the iron regulation phenotype. Several polymorphisms described in this study are found in the 3′ UTR of the messenger RNA and could conceivably affect mRNA stability or levels of protein translation. Regulatory regions that affect levels of transcription can be found in introns or flanking a gene, where most of our polymorphisms are found. Also, tightly linked regulatory polymorphisms could be in disequilibrium with our observed haplotypes. The effects of these polymorphisms in regulatory regions are likely to be quite subtle, but they could help explain the finding of hemochromatosis in individuals that carry only one copy of a hemochromatosis-associated allele.

Finally, knowing more about the levels of variation and recombination at HFE will help evaluate the peculiar pattern of high frequency and young age estimated for the C282Y allele. Evidence for intragenic recombination at HFE has been sparse. Our estimates of local recombination rates based on sequence variation allow a reevaluation of the high LD around the C282Y allele, providing support for its young age. Although this pattern seems consistent with a selective advantage for C282Y, until an appropriate population genetic test for positive selection is performed, alternative causes of the pattern such as drift or a population bottleneck or growth cannot be ruled out. We have begun to use the polymorphisms and haplotypes reported here to study the extent of LD between HFE alleles besides C282Y and markers in the several megabases around HFE. By analyzing LD found around other HFE alleles, as well as around alleles produced by neutral coalescent simulations for different demographic parameters, we can evaluate the claim of positive selection on the C282Y allele relative to other possible explanations of the allele's high frequency given its young apparent age.

View this table:
APPENDIX

References for previously described polymorphisms

Acknowledgments

We thank all those who agreed to donate DNA for this study and Ami Rice for help collecting samples. Additionally, we thank R. Ajioka and L. Jorde for providing European DNA samples from hemochromatosis pedigrees; D. Ledbetter for the chimpanzee sample; J. Fay and R. Hudson for providing computer programs; E. Stahl, M. Fullerton, and members of A. Di Rienzo's laboratory for helpful discussions; J. Comeron for help with the analysis of divergence from chimp and mouse; and J. Comeron, A. Di Rienzo, K. Dyer, M. Hamblin, and two anonymous reviewers for helpful comments on the manuscript. This work was supported by National Institutes of Health grant GM39355 to M.K. and a National Science Foundation Doctoral Dissertation Improvement Grant (DEB-0073297) to C.T. and M.K. C.T. was partially supported by a Howard Hughes Medical Institute predoctoral fellowship and by National Institutes of Health training grant T32 GM07197 (genetics and regulation).

Footnotes

  • The Pan troglodytes nucleotide sequence data from this article have been deposited with the EMBL/GenBank Data Libraries under accession no. AF447807.

  • Communicating editor: Y.-X. Fu

  • Received January 2, 2002.
  • Accepted May 3, 2002.

LITERATURE CITED

View Abstract