Abstract
Genotypefitness correlations (GFC) have previously been studied using allozyme markers and have often focused on shortterm processes such as recent inbreeding. Thus, models of GFC usually neglect marker mutation and only use heterozygosity as a genotypic index. Recently, GFC have also been reported (i) with DNA markers such as microsatellites, characterized by high mutation rates and specific mutational processes and (ii) using new individual genotypic indices assumed to be more precise than heterozygosity. The aim of this article is to evaluate the theoretical impact of marker mutation on GFC. We model GFC due to shortterm processes generated by the current breeding system (partial selfing) and to longterm processes generated by past population history (hybridization). Various mutation rates and mutation models corresponding to different kinds of molecular markers are considered. Heterozygosity is compared to other genotypic indices designed for specific marker types. Highly mutable markers (such as microsatellites) are particularly suitable for the detection of GFC that evolve in relation to shortterm processes, whereas GFC due to longterm processes are best observed with intermediate mutation rates. Irrespective of the marker type and population scenario, heterozygosity usually provides higher correlations than other genotypic indices under most biologically plausible conditions.
THE existence of correlations between individual genotypes at marker loci and fitnessrelated traits has caused much debate among evolutionary biologists. Such correlations were initially used as an argument in favor of selection acting on the maintenance of allozyme polymorphisms in the controversy that has historically opposed selectionists and neutralists (David 1998). Allozymes have been used for decades to detect correlations between multilocus heterozygosity (the number of heterozygous marker loci per individual) and fitnessrelated traits such as growth, viability, or physiological parameters (reviewed in Mitton and Grant 1984; David 1998). Such positive heterozygosityfitness correlations (HFC) have been reported for various organisms, including marine bivalves (Zouroset al. 1988), salmonid fishes (Learyet al. 1983), and pine trees (Ledig 1986). HFC have recently been reported using restriction fragment length polymorphism (RFLP) markers (Pogson and Fevolden 1998) and microsatellite markers (Bierne et al. 1998, 2000a).
The observation of significant HFC with noncoding DNA markers makes it clear that at least some of the correlations are not due to direct effects of the marker genes on the phenotype. Associative overdominance refers to any kind of HFC not due to a direct effect of marker genes, but to genetic associations between the markers and fitness genes (David 1998). The first kind of association is linkage desequilibrium due to genetic drift (correlation of allelic state within gametes). This has been identified as a possible cause of HFC in small populations when there is physical linkage between fitness genes and marker genes (Ohta and Kimura 1970; Pamilo and Palsson 1998). The second kind of association, identity disequilibrium due to variance in inbreeding (correlation of homozygosity between loci across the whole genome), has been identified as a major source of HFC in several theoretical models (Ohta and Cockerham 1974; Charlesworth 1991; Zouros 1993). Here we focus on this second process. Variance in inbreeding generates HFC because more inbred individuals are both more homozygous for their marker loci and less fit due to inbreeding depression. HFC may thus be a powerful tool to analyze the fitness consequences of inbreeding (Charlesworth and Charlesworth 1999; Pembertonet al. 1999). However, inbreeding itself may be caused by very different population processes. The most obvious is the mating system, which generates “shortterm” inbreeding, i.e., inbreeding caused by one or a few generations of consanguineous matings. This could explain HFC in large, partially selfing populations of pine trees (Ledig 1986; Bush and Smouse 1991). “Longterm” inbreeding, on the other hand, involves both recent coalescence events and coalescence events deeper in the pedigree. For example, when two isolated populations come into contact, hybrid offspring are more “outbred” than nonhybrid offspring. In this case, inbreeding may have built up during a long history of isolation of the parental source populations. Such a scenario was invoked to explain the HFC detected in the red deer population of the Isle of Rum (Coulsonet al. 1998) and in the harbor seal population breeding on Sable Island (Coltmanet al. 1998; Pembertonet al. 1999).
Previous models of associative overdominance implicitly neglect marker mutation, because they focus on shortterm inbreeding and on markers with low mutation rates (e.g, allozymes). Now that highly mutable markers (i.e., microsatellites) and longterm scenarios are included in HFC studies, a theoretical assessment of the importance of mutation is needed. This approach was implicitly followed by Coulson et al. (1998) who proposed the use of a new individual genotypic index, rather than heterozygosity, to account for the marker mutation properties (in their case, microsatellites). The model usually assumed for microsatellites is stepwise mutation (one repeat unit added or removed in each mutation event; Ohta and Kimura 1973; Valdèset al. 1993). This mutation model suggested the definition of an index, d^{2}, the squared difference in repeat units between the two alleles of an individual (Coulsonet al. 1998), whose distribution is closely related to the distribution of coalescence times under such a mutation model (Pritchard and Feldman 1996). Empirical studies have sometimes found d^{2} to correlate with fitness traits in samples where heterozygosity does not correlate significantly with these traits (Coltmanet al. 1998; Coulsonet al. 1998). It has thus been suggested that heterozygosity is suitable for detecting shortterm inbreeding, whereas d^{2} provides additional information about longterm inbreeding, due to the mixture of formerly isolated subpopulations (Pembertonet al. 1999; Coulsonet al. 1999; Marshall and Spalton 2000).
In this article we assess the significance of these arguments using a theoretical approach. The influence of marker mutation on genotypefitness correlations (GFC) due to inbreeding is investigated by comparing (i) different mutation models and mutation rates, (ii) different population inbreeding histories involving various time scales, and (iii) different genotypic indices related to the mutational processes of the marker.
THEORY
The model: The rationale of the model is as follows: Population history can generate variance in inbreeding among individuals, depending on individual pedigrees. The population can thus be partitioned into “inbreeding classes.” First, inbreeding levels are associated with genotypes at neutral markers. The latter can be obtained from identityinstate (IIS) relationships among marker alleles for a given marker mutation model and population scenario. Second, in the presence of inbreeding depression, the inbreeding level determines the value of the fitness trait. Thus, individual phenotype and genotype correlate through individual variation in the inbreeding level. For each mutation model, genotypic index, and population scenario assumed, the correlation coefficient ρ(X, W) between X (a given index) and the fitness trait W is derived. As
Mutation models and individual genotypic indices: In a second step we develop theoretical models that approximate mutational processes at DNA markers, such as microsatellites, RFLP markers, or neutral DNA sequences. Only the first two categories of markers have been used in GFC studies, but the third kind might also be used in the future.
Stepwise mutation model: Under a strict stepwise mutation model (SMM; Ohta and Kimura 1973) with mutation rate u, an allele with i repeat units is assumed to mutate only to the i − 1 or the i + 1 states, each with probability u/2 per generation. This model is classically taken as an approximation for mutation at microsatellite marker loci, for which u ranges from 10^{−6} to 10^{−2} per generation (Jarne and Lagoda 1996; Estoup and Angers 1998). Two indices are considered. The first is individual heterozygosity H, which takes the value 0 when the marker locus is homozygous and 1 when it is heterozygous. The second is d^{2}, the squared difference in repeat units between the two alleles of an individual (Coulsonet al. 1998).
Kalleles model: Under a model first formulated by Crow and Kimura (1970), there is a finite number K of possible allelic states, and each allele can mutate to any other at rate u/(K − 1). A twoallele model can be taken as an approximation for mutation at an RFLP locus, as only two alleles (“cut” and “uncut”) need to be distinguished, and for single nucleotide polymorphisms (SNPs; Kuhneret al. 2000). We assume the mutation rate of these DNA markers to be fairly low, even though no direct data are available. Estimates of the substitution rate per nucleotide site and per generation are typically of the order of 10^{−8} or less, depending on the organism under consideration (Li 1997; Drakeet al. 1998; Nachman and Crowell 2000). For a SNP site, u should thus be at most 10^{−8} per generation. Assuming that an RFLP polymorphism is due to nucleotide substitutions rather than to indels and that the RFLP locus is composed of few nucleotides, the mutation rate for an RFLP should be of the order of 10^{−8} or 10^{−7} per locus per generation. Although mutational processes affecting RFLP loci may be asymmetric, destroying any particular polymorphism more often than recreating it, we do not expect our model to be very sensitive to this asymmetry, given that mutation rates will be low (but this remains to be tested). The only index studied for this kind of marker locus is individual heterozygosity.
Infinitealleles and infinitesites models: Under the infinitealleles model (IAM; Kimura and Crow 1964), each mutation event produces a new allele at rate u per generation. The infinitesites model (ISM; Kimura 1969) shares the same property but is more specifically designed for DNA sequences. In this model the number of nucleotide sites in the sequence is assumed to be so large that each new substitution occurs at a site that has not mutated before. The total mutation rate of the sequence is u = lμ, where l is the sequence length in base pairs and μ is the substitution rate per nucleotide site per generation (see above). The first index considered is sequence heterozygosity H, which takes the value 0 when the two alleles of the marker are strictly identical in sequence and 1 when there is at least one different nucleotide site between the two allelic states. Note that analytical derivations of heterozygosity under the ISM and the IAM with mutation rate u are strictly equivalent. For neutral DNA sequences, the number of nucleotide differences between the two alleles of an individual, denoted by p, is also used as an alternative index.
The partial selfing model and singlelocus heterozygosityfitness correlation: We investigate a simple scenario in which HFC relates to the actual mating system of the population. This scenario examplifies shortterm inbreeding. Consider a large population of size n at inbreeding equilibrium, with freely recombining loci and one marker locus with a given mutation model. A proportion S of offspring is produced by selfing at each generation, whereas a proportion 1 − S comes from outcrossing events. Each individual is characterized by the number of generations of selfing in its pedigree, k, starting from the most recent outcrossing event (0 ≤ k ≤ ∞). The population is thus partitioned into inbreeding classes C_{k}, 0 ≤ k ≤ ∞, each consisting of individuals having the same inbreeding level (i.e., the same k). Note that, under this model, there is an infinite number of inbreeding classes, although in practice only classes with low k are largely represented. In the absence of selection, inbreeding class k has frequency Pr(C_{k}) = (1 − S)S^{k}. Selection against homozygous genomes reduces the frequency of classes with high k. The effect of selection thus resembles a reduction in the selfing rate S to a lower value S_{sel} that can be computed as explained in David 1999 (appendix a). In practice, David (1999) shows that frequencies of inbreeding classes computed using S_{sel} instead of S provide a good approximation of the frequencies with selection. In what follows, we therefore take into account selection against homozygous genomes simply by replacing the “raw” S value by the value of S_{sel}.
Let us define the probabilities of IIS within a population: Q_{0} for a pair of genes drawn from the same individual and Q_{1} for genes from different individuals. In what follows, moments are expressed as functions of Q_{0} and Q_{1} that are derived in appendix a under each mutation model. The inbreeding coefficient in class k, f_{k}, is defined by
Moments of heterozygosity H are simply expressed as functions of probabilities of IIS:
Note that cov(H, W) is null for S = 0 and S = 1, since there is then no variance in the inbreeding level among individuals. Equations 4b, 5b, and 6b can now be used to derive the correlation coefficient between heterozygosity and the fitness trait from (1). The maximum correlation coefficient (ρ_{max}), assuming no withininbreeding class variance for the fitness trait (
Extension of the partial selfing model to other genotypic indices: We first focus on the squared difference in repeat units for a microsatellite marker under the stepwise mutation model of mutation. Cov(d^{2}, W) is derived as in (6a), replacing H by d^{2}. Moments of d^{2} are obtained as detailed in appendix a (Equation A7). This yields
We then focus on another genotypic index: the number of nucleotide differences for a neutral sequence under the infinitesites model of mutation. The correlation coefficient between p and W under the ISM is computed similarly to the correlation coefficient between d^{2} and W. Moments of p are derived in appendix a (Equation A5) and Cov(p, W) is as in (8). We obtain
The population admixture model (hybridization): We next consider a simple situation where the HFC is related to longterm processes. We assume a large, random mating population of size N at mutationselection equilibrium, which splits into two randomly mating finite subpopulations of size n ⪡ N. Let the two subpopulations diverge for a long time τ (assumed to be equal to N generations for the sake of simplicity) without any gene flow between them and then merge into a single, infinite, panmictic population. Given enough divergence time, different deleterious mutations inherited from the ancestral population will reach different frequencies in each subpopulation. We investigate the genotypephenotype relationship in the resulting mixed population. Each individual is characterized by the probability x that the two alleles it has at a given locus originate from the same subpopulation, depending on the number g of panmictic generations following the admixture. The quantity 1 − x can be interpreted as the individual “degree of outbreeding,” i.e., the “disparity between the genome of the two parents” (Coulsonet al. 1999).
Just after the admixture (g = 1) the correlation is due to the coexistence of two inbreeding classes, C_{w} and C_{b}, with equal probabilities 1/2. C_{w} (w stands for “within”) is composed of individuals with both parents originating from the same subpopulation (x = 1) while C_{b} individuals (b stands for “between”) have one parent in each of the two subpopulations (x = 0). C_{w} individuals are expected to have lower fitness than C_{b} individuals as they are more likely to be homozygous for their deleterious mutations. They also have more similar alleles at a neutral marker locus than C_{b} individuals.
Each generation after the contact, new inbreeding classes C_{x} appear as a consequence of recombination. At generation g after the contact, assuming free recombination for the sake of simplicity,
Let Q_{w} and Q_{b} be the probability of identityinstate of genes from C_{w} and from C_{b}, respectively. f_{x}, the inbreeding coefficient of class x, is defined as
Here W_{0} is the value of the fitness trait for C_{b} individuals (i.e., with x = 0) and the inbreeding load β, the difference in fitness for C_{b} and C_{w} individuals, measures heterosis between the two subpopulations before the mixture. The random effect ε with mean 0 and variance
Moments of heterozygosity H are simply expressed as functions of probabilities of IIS:
The covariance term is calculated using the partition into inbreeding classes,
The maximum correlation coefficient is
The model can be generalized to other cases, as for the partial selfing model. A summary of the analytical expressions for the squared correlation coefficients obtained under both population models is provided in Table 1.
Numerical parameters for the models: Analytical results derived in the previous sections were explored numerically using Mathematica 3.0 programs (Wolfram 1996). The ranges of values explored for the various parameters of the model are presented below.
Mutation rates: To analyze the impact of mutation on heterozygosityfitness correlations with respect to different mutation models, a large range of mutation rates is investigated, from 10^{−9} to 10^{−2} per locus and per generation. Mutation rates ranging from 10^{−6} to 10^{−2} are investigated to compare d^{2} and H under the SMM, as this is the accepted range for a microsatellite locus. Assuming a substitution rate equal to 10^{−9} per site per generation, total mutation rates ranging from 10^{−9} to 10^{−6} correspond to sequences of realistic length (few to several hundreds of base pairs). These mutation rates were thus considered when comparing p and H under the ISM.
Population parameters: Under the partial selfing model the standard situation modeled is a large population (n = 10^{3} individuals) with intermediate selfing rate (S = 0.4). The effect of a change in the population size (from 10^{3} to 10^{6}) and the selfing rate (from 0 to 1) is analyzed. Under the admixture model the standard assumptions are small subpopulations (n = 10^{2}) descended from a large ancestor population (N = 10^{4}) after a long divergence time (τ = 10^{4} generations). The effect of the subpopulation size (from 10^{2} to 10^{4}) and the joint influence of the ancestor population size and the divergence time in generations (assumed to be equal for the sake of simplicity and ranging from 10^{4} to 10^{6}) are investigated. The inbreeding load does not need to be estimated here as we focus on maximum correlation coefficients, which do not depend on this parameter (see Table 1).
RESULTS
We first focus on the effects of the mutation models and mutation rates on the heterozygosityfitness correlation for both population scenarios, evaluating the effect of population parameters. We then analyze the influence of genotypic indices on the genotypefitness correlation obtained. Finally the form of the relation between the genotype and the fitness (the expected value of the fitness trait as a function of d^{2}) is discussed under the SMM.
Impact of marker mutation on the heterozygosityfitness correlation: Under the partial selfing model, the maximal correlation coefficient increases with the marker mutation rate for all mutation models studied (Figure 1A). This is not surprising since this correlation is an increasing function of marker gene diversity [(1 − Q_{1}) in Equation 7a], which itself increases with the mutation rate irrespective of the mutation model (see appendix a). For low mutation rates (<10^{−5}), the correlation coefficient is never higher than a few percent, in agreement with most experimental data (David 1998). Mutation models have a low impact except for high mutation rates (>10^{−5}), which are realistic only for microsatellite markers. For high mutation rates the correlation coefficient increases with the value of K in the Kalleles model (KAM) and is maximal for the IAM. The SMM behaves like a KAM with large K (K > 10, data not shown).
A different picture is obtained under the admixture model, as the maximum correlation coefficient exhibits a maximum for an intermediate mutation rate (Figure 1B). This “optimal mutation rate” is of the order of 1/τ (the divergence time; Figure 1B and data not shown). The correlation coefficient obtained under the admixture model may be much higher than that obtained under the partial selfing model. Under the partial selfing model the correlation is limited by the fact that all genotypes, be they homozygotes or heterozygotes, are represented in all inbreeding classes although the frequency of heterozygotes is halved in each generation of selfing. In contrast, under the admixture model, parameter values can be found such that all homozygotes are in class C_{w} and all heterozygotes are in class C_{b}, which allows the maximum correlation (in the absence of environmental variance in fitness) to approach unity. However, assuming freely recombining loci, this maximum correlation is halved each generation following the admixture (Equation 18). Furthermore, in natural populations, withininbreeding class variance for the fitness trait may weaken the real correlation [Equation 7b, replacing σ^{2}(f_{k}) by σ^{2}(x)]. Mutation models rank in the same order as those under the partial selfing model with regard to the correlation coefficient, except with very high mutation rates (>10^{−3}) for which the SMM provides a stronger correlation than the IAM. Ultimately, for both population models studied, marker mutation influences HFC primarily through its effect on marker diversities: (1 − Q_{1}) for the partial selfing model (Equation 7a) and (1 − Q_{w}) and (1 − Q_{b}) for the admixture model (Equations 18).
Impact of population parameters on the heterozygosityfitness correlation: Under the partial selfing model, the correlation coefficient increases with the population size and the selfing rate, until S is very close to 1 (Figure 2), whatever the mutation rate and model (data not shown). Ohta and Cockerham (1974) obtained qualitatively similar results for the effect of the selfing rate on associative overdominance in an infinite partially selfing population. However, Charlesworth (1991) found that associative overdominance is maximal for intermediate selfing rates. This discrepancy probably results from the measure of associative overdominance chosen. Charlesworth (1991) studied the apparent selection coefficient for heterozygotes at the neutral marker, which is expressed in phenotypic units and thus depends on the inbreeding load. The maximal heterozygosityfitness correlation coefficient used here is independent of the inbreeding load β, unlike the covariance term, which scales with β (Equation 6b). Using the latter as a measure of associative overdominance, we obtain a maximum associative overdominance for intermediate selfing rates (data not shown), just as found by Charlesworth (1991). Taking into account withinpedigree variance for the fitness trait (σ_{ε} > 0 in Equation 4b) the correlation is weakened and becomes an increasing function of the inbreeding load (Equation 7b).
In the admixture model, the correlation coefficient decreases with the subpopulation size, particularly when the mutation rate is high (Table 2). The correlation increases when the ancestral population size and the divergence time increase simultaneously, particularly for low mutation rates. The effect of the inbreeding load β is similar to that in the partial selfing scenario (Equation 7b).
Impact of the genotypic index on the correlation with the fitness trait: The squared difference in repeat units vs. heterozygosity under the stepwise mutation model: Under the partial selfing model, for low mutation rates, d^{2} and heterozygosity H are equally correlated with fitness (Figure 3A). Increasing the mutation rate strongly increases the correlation with heterozygosity but not with d^{2} (Figure 3A), irrespective of population size or selfing rate (data not shown). This can be explained by the distribution of the expected value of the fitness trait conditioned on the value of d^{2}. Using (C3) and (C4) in appendix c we found that all heterozygotes (nonnull d^{2} values) correspond numerically to the same expected fitness (e.g., W = 0.89 with W_{0} = 1, β = 1, u = 10^{−4}, and other parameters as in Figure 3A), which is higher than that of homozygotes (d^{2} = 0) (e.g., W = 0.72 with W_{0} = 1, β = 1, u = 10^{−4}, and other parameters as in Figure 3A). Therefore, heterozygosity is as informative as d^{2}, and it provides a stronger correlation due to its lower variance.
Under the admixture model, however, the correlation between heterozygosity and fitness is weakened for high mutation rates (Figures 1B and 3B) whereas the correlation with d^{2} is almost insensitive to the mutation rate (Figure 3B). With a small subpopulation size (n = 100) and mutation rate (u = 10^{−4}), the correlation between the fitness trait and d^{2} is weaker than that with heterozygosity, as in the partial selfing model (data not shown). In contrast, when the subpopulation size and the mutation rate are large enough (so that nu ⪢ 1) d^{2} provides a better correlation than H (Figure 3B). Again, this can be explained by the distribution of the expected value of the fitness trait conditional on d^{2}, derived according to appendix c (Figure 4). Indeed in the first case (nu ⪡ 1), d^{2} has little value for predicting fitness once it is >1. However, when nu ⪢ 1, the expected value of fitness increases progressively with the value of d^{2}, so that in this case d^{2} is far more informative than H.
The number of nucleotide differences vs. heterozygosity under the infinitesites model: Under the partial selfing model, given the low mutation rate assumed for a neutral sequence, the same correlation with fitness is obtained whether heterozygosity H or the number of nucleotide differences p is considered (Figure 5A). Only very high mutation rates (>10^{−4}), corresponding to unreasonably long sequences, would introduce a difference, with a higher correlation with H than with p (data not shown).
Under the population mixture model, identical correlation coefficients are obtained with p and H except with large subpopulation size (>10^{3}) and large mutation rates (>10^{−6}; data not shown). In this case p is slightly less correlated to fitness than H (Figure 5B). Again, only very high mutation rates (>10^{−4}) would result in a higher correlation with p than with H (data not shown). Unsurprisingly, these results are qualitatively similar to that obtained in the comparison of d^{2} and H under the SMM with equivalent mutation rates.
The expected value of the fitness trait as a function of d^{2}: The expected shape of the genotypephenotype relationship is not linear (see Figure 4), with fitness plateauing for high d^{2}values. This indicates that linear regression models previously used in empirical studies are not appropriate. We therefore suggest using a nonlinear regression based on the logistic equation
For the partial selfing model, δ(d^{2}) can be interpreted as a microsatellitebased indirect measure of inbreeding depression. For the admixture model it is interpreted as an indirect measure of heterosis between the two subpopulations before the mixture. Numerical data obtained using appendix c were fitted to this nonlinear model, for a large range of parameters defined in Tables 3 and 4 (data not shown). This model explains 95–100% of the variance in fitness obtained for both population models studied (neglecting withininbreeding class variance for the fitness trait). It thus provides a satisfactory description of the genotypephenotype relationship on the basis of three parameters. Under the partial selfing model, δ(d^{2}) increases with the marker mutation rate, the selfing rate, and the population size (Table 3). Under the mixture model, δ(d^{2}) increases with the marker mutation rate and when the divergence time and the population size increase simultaneously. It decreases with the subpopulation size (Table 4).
DISCUSSION
A unified framework for studying genotypefitness relationships: To our knowledge, mutational processes at marker loci have not previously been incorporated into analytical models of associative overdominance, although they are recognized as relevant for empirical issues (Pembertonet al. 1999). Evaluating heterozygosityfitness correlations requires the evaluation of probabilities of identityinstate, which depend both on the demographic scenario assumed and on mutational processes and have been extensively used to address population structure issues (reviewed in Rousset 2001). This could be done in a unified framework using the same formalism as previously used (Rousset 1996). It is equivalent to a formalism based on the computation of distributions of coalescence times but does not require an explicit computation of such distributions. Under this formalism, neglecting withininbreeding class variance in fitness for a given fitness model, the maximal correlation coefficient between a fitness trait and a neutral genotype does not depend on the inbreeding load (the reduction in fitness associated with complete inbreeding), as already reported for other demographic scenarios (Bierneet al. 2000b). There is therefore no need for an explicit expression of the inbreeding load, which depends on the genetic architecture of the fitness trait and on the population scenario. Rather than exploring the continuum of possible historical scenarios exhaustively, we chose to exemplify a shortterm (with only very recent coalescence events) and a longterm (with coalescence events deeper in the pedigree of individuals) inbreeding process with, respectively, the partial selfing and the admixture model.
The fitness model we used (Equations 3 and 14) was first proposed by Morton et al. (1956). Neglecting purging selection and disequilibrium between selected loci and assuming additive effects among fitness loci (or multiplicative effects, considering logfitness rather than fitness), we obtained the Morton model for fitness under various genetic architectures and various population scenarios (Charlesworth and Charlesworth 1987, 1999; Bierneet al. 2000b; Whitlocket al. 2000).
Variance for the fitness trait may exist within an inbreeding class, particularly when recombination is limited, due to the segregation of chunks of chromosomes rather than independent loci. In natural populations, environmental effects may also increase this variance. If this is large enough compared to genetical effects, HFC would be weakened and also would depend on the inbreeding load, which should then be estimated. Empirically estimating within and betweeninbreeding class variance for the fitness trait might be difficult in natural populations, but would be particularly useful to assess the validity of our model.
Choosing appropriate marker genes: Microsatellites have recently been used as a tool to infer fitness differences due to variation in the level of inbreeding between individuals (Pembertonet al. 1999). Are microsatellites good markers to address such issues?
We have shown that HFC due to partial selfing increases with the marker mutation rate (and thus with marker variability), making microsatellites good markers in such a situation. Partial selfing is an example of shortterm inbreeding generated by the current mating system. One generation of outcrossing reduces individual inbreeding to zero, so that only selfing since the most recent outcrossing event (i.e., a few generations) needs to be taken into account. However, shortterm inbreeding can also be due to recent demographic episodes such as bottlenecks. A sudden reduction in population size produces random inbreeding among individuals, i.e., mating between kin due to chance, as opposed to systematic inbreeding due to the mating system (Malécot 1969). Bierne et al. (2000b) analyzed theoretically the situation of a population experiencing a recent and drastic bottleneck sustained over a few generations. They found strong HFC, but it is transient due to the rapid homogenization of the inbreeding level among individuals after a few generations of a sustained bottleneck. Bierne et al. (2000b) did not use an explicit mutation model for the marker genes. However, they found that highly variable markers produce larger HFC than less variable ones, consistent with our result under partial selfing. Heterozygosityfitness correlations increase with marker diversity and thus with marker mutation rates, whenever they are related to shortterm inbreeding. Microsatellite loci seem therefore ideal to investigate fitness consequences of shortterm inbreeding.
However, when the origin of HFC is longterm inbreeding (e.g., admixture of anciently diverged small subpopulations), the correlation increases up to a certain mutation rate but then decreases (Figure 1B). Using highly mutable markers can thus lead to confusing results under a longterm inbreeding scenario. The explanation is that even highly “inbred” genotypes can have heterozygous marker loci if marker mutation rates are sufficiently high. Marker heterozygosity is therefore not a good index of individual fitness. There is thus an “optimal mutation rate,” which decreases as the timescale of the population inbreeding history increases. In conclusion, one should use molecular markers whose mutational process scales roughly with the assumed cause of inbreeding. However, empirical studies provide only imprecise estimates of mutation rates and a temporal scale of inbreeding. In practice we suggest only that empirical studies should avoid using markers that mutate at an obviously inappropriate rate.
Mutation models of marker genes play little role in HFC, except when the mutation rate is high. In the infinitealleles model there is no homoplasy, whereas the Kalleles model is characterized by increasing homoplasy with decreasing Kvalues. The comparison between IAM and KAM, and among various Kvalues, suggests that homoplasy is a limiting factor for HFC, irrespective of the population scenario. Homoplasy represents a loss of information about identitybydescent, making marker genotypes less closely related to the inbreeding level of individuals. As the homoplasy generated by the stepwise mutation model is different from that produced by a KAM, the results obtained for these two mutation models cannot be directly compared.
Choosing an appropriate genotypic index for microsatellite data: It has been suggested recently that d^{2} may be a more appropriate index than heterozygosity when analyzing fitness consequences of inbreeding (Coulsonet al. 1999; Pembertonet al. 1999). Heterozygosity should be adequate to distinguish inbred and outbred individuals, whereas d^{2} may give the additional opportunity to detect “highly outbred” individuals having high d^{2} values, i.e., large coalescence time for their two microsatellite alleles, vs. “moderately outbred” individuals (Coulsonet al. 1999).
In a red deer population from the Isle of Rum, Coulson et al. (1998) found a positive correlation between birth weight and mean d^{2}, although not between birth weight and heterozygosity. They argued that this could be due to the demographic history of the population, namely population divergence followed by mixing (our admixture model), which is supported by historical data. However, we have shown that, in our admixture model, at least, the conditions under which d^{2} would lead to a better correlation than H in the admixture model are restrictive (a very high marker mutation rate and large subpopulation sizes and divergence time) and furthermore they do not seem to apply to what is known of the red deer population of the Isle of Rum (Pembertonet al. 1999). Divergence times of ~100 generations, with large subpopulation size, parameters that appear realistic for such a population, always lead to better predictive power of H as compared to d^{2} (not shown). Note that the fact that there are only two classes of inbreeding and fitness in our admixture scenario just after the contact does not necessarily imply that a binary variable such as H will perform better than a quantitative variable such as d^{2}. If the proportion of high fitness individuals (class C_{b}) within a class of d^{2} increases progressively with d^{2} (Figure 4, with u = 10^{−2}) rather than in a stepwise way (in most cases), d^{2} performs better than H.
One can suggest another explanation for the stronger correlation with fitness obtained using d^{2} as compared to H. First, we cannot exclude the possibility that a more sophisticated demographic scenario may generate this result under less restrictive conditions than the admixture scenario formulated here, although we do not have a precise idea of what such a scenario should be. Alternatively, the results of Coulson et al. (1998, 1999) may be due to low homozygosity in their population, making heterozygosity an inadequate genotypic index when compared to d^{2}. Finally these results could simply be due to chance, i.e., random variation. Recently, several studies have found fitnessrelated traits, such as adult breeding success in the red deer (Slateet al. 2000) or survival and parasite resistance in the Soay sheep (Coltmanet al. 1999), to be more strongly correlated with heterozygosity than with d^{2}. In a captive wolf population with known inbreeding levels, Hedrick et al. (2001) found that d^{2} was less predictive of the known inbreeding coefficient than microsatellite heterozygosity. This supports the hypothesis that d^{2} is generally not more powerful than heterozygosity to detect fitness consequences of inbreeding.
It seems generally unlikely that d^{2} will provide a higher correlation with fitness than H. For low mutation rates, heterozygosity and d^{2} are equivalent under both our population models, and for high mutation rates, heterozygosity always provides a higher correlation than d^{2} with fitness under the partial selfing model. Under the admixture model only, d^{2} is more correlated with a fitness trait than heterozygosity under restrictive conditions: large mutation rate, divergence time, and subpopulation size. These conditions are expected to be associated with a relatively low inbreeding load (β) in the mixed population, because there is a small probability of fixation, or of sufficient change in frequency, of a deleterious mutation in a large subpopulation. The inbreeding load does not influence the maximum correlation (Equation 18a) but affects the actual correlation, taking into account withininbreeding class variance for the fitness trait. Given that sources of variation in fitness other than inbreeding must exist in natural populations (generating σ_{ε} in Equations 4b and 15b), this small inbreeding load weakens the actual correlation (see Equation 7b), thus reducing the chance of detecting it.
Finally, a conclusion of our analysis is that there seems to be little theoretical reason to use d^{2} instead of H when analyzing correlations between microsatellite genotype and fitness, even for longterm inbreeding scenarios. When doing so, however, one should choose an appropriate regression model to analyze genotypefitness correlations (Equations 19 and 20) rather than the linear model used in previous studies.
Acknowledgments
We thank Nicolas Bierne, Sylvain Glémin, Philippe Jarne, Sally Otto, Josephine Pemberton, John Slate, and John Thompson for stimulating discussions and helpful comments on a previous draft of the manuscript. This work was supported by funds from the Centre National de la Recherche Scientifique (CNRS; Jeune Equipe program) to P. Jarne.
APPENDIX A: IDENTITYINSTATE IN THE PARTIAL SELFING MODEL
We derive probabilities of identityinstate within and between individuals under the infinite sites, the stepwise mutation, and the Kalleles mutational models (see Rousset 1996). Generating functions of allele size difference (under the SMM) or of the number of nucleotide differences (under the ISM) are also derived and used to obtain the moments of d^{2} and of p.
IAM and ISM: Recursions for the IIS probabilities under the IAM or the ISM are
The recursions on ψ are then mathematically similar to those for IBD (Rousset 1996):
SMM: To derive Q_{0} and Q_{1} under the SMM it is useful to introduce ψ_{0} (resp. ψ_{1}), the generating function of the probability p_{j}_{,0} (resp. p_{j,1}), that two randomly chosen alleles within an individual (resp. between individuals) differ by j steps. Because steps can take either nonnegative or negative values,
By differentiating the generating functions ψ one obtains the various moments of d^{2}:
KAM: The recursions for IIS probabilities are
The equilibrium values of the system are
APPENDIX B: IDENTITYINSTATE IN THE MIXTURE MODEL
We derive probabilities of identityinstate under the ISM, SMM, and KAM (see Rousset 1996) focusing on C_{w}, C_{b}, and C_{x} individuals. Generating functions of allele size difference (under the SMM) or of the number of nucleotide differences (under the ISM) are also derived, as they are useful to derive moments of d^{2} and of p. In contrast to the partial selfing model there is no need to distinguish identity within (subscript 0) and between (subscript 1) individuals, as random mating is assumed.
IAM and ISM: The probability of identity (IIS and IBD are confounded here) for two alleles randomly drawn from the ancestral random mating population of size N at inbreeding equilibrium is Q_{a} = γ/(2N(1 − γ) + γ) with γ = (1 − u)^{2}. Let two subpopulations from the ancestral population diverge without gene flow for τ generations. Just before the admixture, the probabilities of identity of two randomly chosen alleles are Q_{w}(τ) = Q_{eq} + (γ(1 − 1/(2n)))^{τ}(Q_{a} − Q_{eq}) within a subpopulation and Q_{b}(τ) = γ^{τ}Q_{a} between subpopulations, with Q_{eq} = γ/(2n(1 − γ) + γ). After the admixture, neglecting mutations since the admixture event and considering an infinite mixed population for the sake of simplicity,
Under the ISM we obtain analogous expressions for generating functions, as we did in appendix a. Let ψ_{w} (resp. ψ_{b}) be the generating function of the probability p_{j,w} (resp. p_{j,b}) that two randomly chosen alleles that originated from the same subpopulation (respectively from the two subpopulations) differ by j nucleotide substitutions [
Neglecting mutations since the admixture event and considering an infinite mixed population for the sake of simplicity, we finally obtain
SMM: Again, we obtain analogous expressions for generating functions, as we did in appendix a. The generating functions ψ_{w}, ψ_{b}, and ψ_{0x} are defined as in the previous section, but considering “steps” (varying from −∞ to + ∞) instead of “nucleotide substitutions” (varying from 0 to + ∞). Recursions on ψ(z) are identical to that in the previous section (Equations B2) except r(z) = (1 − u + 1/2(uz + u/z))^{2}. By definition, we have p_{0,w} = Q_{w} and p_{0,b} = Q_{b}. Furthermore, for all j ≠ 0, p_{−j,w} = p_{+j,w} = ½ Pr(d^{2} = j^{2}C_{w}) and p_{−j,b} = p_{+j,b} = ½ Pr(d^{2} = j^{2}C_{b}). Therefore
Using inverse Fourier transforms,
KAM: In the ancestral random mating population of size N at inbreeding equilibrium we have
APPENDIX C: THE EXPECTED FITNESS CONDITIONAL ON d^{2}
Under the partial selfing scenario, Equation 3 yields
Pr(d^{2} = j^{2}) is derived from the L_{j} transform of the generating function ψ_{0} (appendix a, SMM). According to (A6), Q_{0} = Pr(d^{2} = 0) = L_{0}(ψ_{0}(e^{ix})) and similarly, for all j ≠ 0,
Finally,
Using the same method in the population admixture model we obtain
At the first generation after the admixture (g = 1) we have
Footnotes

Communicating editor: D. Charlesworth
 Received May 23, 2001.
 Accepted October 1, 2001.
 Copyright © 2001 by the Genetics Society of America