## Abstract

Under several hypotheses for the evolutionary origin of imprinting, genes with maternal and reproductive effects are more likely to be imprinted. We thus investigate the effect of genomic imprinting in single-locus diallelic models of maternal and fertility selection. First, the model proposed by Gavrilets for maternal selection is expanded to include the effects of genomic imprinting. This augmented model exhibits novel behavior for a single-locus model: long-period cycling between a pair of Hopf bifurcations, as well as two-cycling between conjoined pitchfork bifurcations. We also examine several special cases: complete inactivation of one allele and when the maternal and viability selection parameters are independent. Second, we extend the standard model of fertility selection to include the effects of imprinting. Imprinting destroys the “sex-symmetry” property of the standard model, dramatically increasing the number of degrees of freedom of the selection parameter set. Cycling in all these models is rare in parameter space.

BIOLOGISTS have long recognized the fundamental importance of the maternal environment to the developing organism (see Wade 1998 for an overview). Although we have as yet no comparable data from mammals, it is clear from DNA microarray studies on *Drosophila melanogaster* that a significant proportion of the maternal genome is expressed in offspring. For example, Arbeitman *et al*. (2002) found that 1212 different RNA transcripts present in the first hours of development were maternally deposited during oogenesis. Crucially, this maternal genetic effect was compounded by standard genetic expression since all but 27 of these same genes were subsequently transcribed from the embryo's own copies.

The list of mammalian maternal-effect genes is small, although growing rapidly. Perhaps the best-known example is the mouse chromosome 7 gene *Mater*, which has no known phenotypic effect on its bearers, except that females homozygous for a null mutant are sterile because the maternally derived protein is necessary for normal embryonic development (Tong *et al.* 2000). Indeed, most of our examples of mammalian maternal-effect genes affect the early development of the offspring of homozygous females: the offspring of *Zar1* null females, for example, usually die at the two-cell stage (Wu *et al.* 2003). The mouse gene *stella* (also called *PGC7*) exhibits a maternal effect, but the paternal contribution is also relevant: when mated to *stella*-deficient males, *stella*-deficient females have no live pups, but when mated to wild-type males they produce about a third of the usual number of offspring (Payer *et al.* 2003).

Genes that are active early in mammalian development are also likely targets for genomic imprinting, according to a number of explanations for the evolutionary origin of imprinting (Spencer 2000). For instance, the genetic-conflict hypothesis argues that growth-affecting genes active during fetal development may be agents of genetic conflict and thus become imprinted (Haig 1992). The ovarian time-bomb hypothesis (Varmuza and Mann 1994; see also Iwasa 1998) posits that genes essential for the initiation of embryogenesis will be imprinted. Thus mammalian genes with a strong maternal effect may also be those more likely to be subject to imprinting.

Indeed, both the murine paternally expressed gene-1 (*Peg1*) and its human ortholog (*PEG1*) are imprinted (Kaneko-Ishino *et al.* 1995; Kobayashi *et al.* 1997) and the former is known to have strong direct and maternal effects (Lefebvre *et al.* 1998). Mice paternally inheriting a mutation of *Peg1* (also known as mesoderm-specific transcript, *Mest*), exhibited growth retardation and reduced survival, and adult females were unable to successfully raise pups, irrespective of the pups' own genotype (Lefebvre *et al.* 1998). Similarly, mutations in the mouse paternally expressed gene-3 (*Peg3*) affect several features of fetal and postnatal development in offspring, as well as various aspects of maternal care in mothers (Curley *et al.* 2004).

Moreover, mutations at loci involving imprinting—such as those involved in the maintenance of DNA methylation, crucial to the monoallelic expression that defines imprinting—will almost always exhibit a maternal effect if the usual imprinting pattern is disrupted. For instance, heterozygous offspring of female mice homozygous for a deletion in the DNA methyl-transferase-1 (*Dnmt1*) locus showed biallelic expression at the *H19*, which is normally expressed maternally, and died during gestation (Howell *et al*. 2001).

There is, consequently, a natural link between loci subject to imprinting and those exhibiting maternal effects. The population genetic effects of such an association, however, have yet to be fully explored. In this article we make a start by examining the mathematical properties of a simple model of maternal selection at a locus subject to imprinting. We do so by incorporating the effect of imprinting into the two-allele single-locus model propounded by Gavrilets (1998; see also Spencer 2003), which describes the population-genetic consequences of fitness differences among both the maternal and the zygote's own gene products.

Genetic conflict may also be manifested in fertility selection, where the relative number of offspring is a property of the maternal and paternal phenotypes. It thus makes sense to investigate the consequences of incorporating imprinting into the standard model of fertility selection at a single diallelic locus (Bodmer 1965).

Although imprinting is often depicted as the inactivation of either the paternal or the maternal copy of a gene, there is usually significant variation among different tissues. Murine insulin-like growth factor 2 (*Igf2*), for example, is maternally inactive in most tissues during embryogenesis, but has standard biallelic (Mendelian) expression in two structures associated with the central nervous system, the choroid plexus and leptomeninges (DeChiara *et al.* 1991; Pedone *et al.* 1994). Isoform 2 of human *PEG1/MEST* is not imprinted in most organs, but exhibits substantial variation among individuals in the level of imprinting in the placenta (McMinn *et al.* 2006). Thus, in its most general form, imprinting is the differential expression of paternally and maternally inherited genes. One important consequence for population genetics is that reciprocal heterozygotes at an imprinted locus will have different mean phenotypes that are not necessarily the same as either homozygote (Spencer 2002).

## MODELS OF MATERNAL SELECTION

### General maternal selection model:

Consider a single locus with two autosomal alleles, *A*_{1} and *A*_{2}, in a randomly mating, dioecious population, in which the effects of mutation and genetic drift are negligible. Imprinting means that we must distinguish between maternally and paternally derived alleles: by *A _{i}A_{j}* we mean a genotype with a maternally derived

*A*allele and a paternal

_{i}*A*. Suppose

_{j}*w*is the fitness of individuals of genotype

_{ijkl}*A*with genotype

_{k}A_{l}*A*mothers. Since

_{i}A_{j}*k*=

*i*or

*j*and

*i*,

*j*,

*k*,

*l*∈ {1, 2}, there are 12 different fitness parameters, as shown in Table 1. If

*g*is the postselection frequencies of adults with genotypes

_{ij}*A*(with ), then the recursion equations for these frequencies in the following generation are(1)in which(2)are the respective frequencies of

_{i}A_{j}*A*

_{1}and

*A*

_{2}, and , the population's mean fitness, is the sum of the right-hand sides of (1) so that the iterated frequencies () also add to one. To our knowledge, Equations 1 are not formally equivalent to those previously used to describe any other population genetic system. Note that the normalization of Equations 1 means that the relative rather than absolute sizes of the fitness parameters (the

*w*'s) determine the dynamical behavior and so just 11 parameters are needed to describe the system.

_{ij}The model exhibits a number of interesting properties, including those possessed by Gavrilets' (1998) model (since the latter is a special case of the imprinting model with *w _{ijkl}* =

*w*=

_{ijlk}*w*=

_{jikl}*w*). Thus, two distinct polymorphic equilibria (

_{jilk}*i.e*., values of such that for all

*i*and

*j*and with each inequality strict for at least one value) may be locally stable, or an internal equilibrium and one or both fixations may simultaneously be so. This result implies that at least five distinct equilibria (not all of which are stable, of course) may be feasible (

*i.e*., real solutions with ) for certain parameter values, and, indeed, we give an example in Table 2.

We have not been able to analytically solve this model to find all possible equilibria, but not one of 10^{5} randomly generated fitness sets (*i.e.*, sets of 12 pseudorandom *w _{ijkl}* values drawn independently from the uniform distribution between 0 and 1) possessed >5 distinct, feasible equilibria. Moreover, we did not find any sets of fitnesses that afforded ≥3 stable polymorphic equilibria in 10

^{6}randomly generated sets, each with 100 random initial genotype frequencies (drawn using the broken-stick method). These simulation results indicate that 5 may well be the maximum number of feasible equilibria. (Up to 10 equilibria arise from solving for all

*i*and

*j*, but many of these are complex roots or unfeasible solutions.)

The two fixation equilibria (*A*_{1} fixed, *g*_{11} = 1, *g*_{12} = *g*_{21} = *g*_{22} = 0; *A*_{2} fixed, *g*_{11} = *g*_{12} = *g*_{21} = 0, *g*_{22} = 1) always exist. The former is locally stable if(3)which, in the absence of imprinting, collapses to the condition found by Gavrilets (1998). The analogous condition for local stability of the *A*_{2} fixation is(4)These inequalities reveal the importance of the asymmetry of imprinting in even simple matters, in this case the stability of monomorphisms. The greater a certain imprinting effect—the difference between fitnesses of heterozygotes inheriting the rare allele from the two sorts of heterozygous mothers (*e.g.*, for the fixation of *A*_{1}, those heterozygotes inheriting *A*_{2} from their mothers and having respective fitnesses *w*_{1221} and *w*_{2121})—the less likely (other things being equal) that monomorphism is locally stable. It is interesting that it is the maternal difference between otherwise identical heterozygous offspring that is crucial, rather than the difference between the offspring themselves.

One novel behavior is the potential for long-period cycling of genotype frequencies, an example of which is shown in Figures 1 and 2. Further analysis (see appendix) revealed that this dynamical behavior was due to a pair of supercritical Hopf bifurcations. In a supercritical Hopf bifurcation, a stable equilibrium point becomes unstable and is encircled by an attracting, invariant closed curve. A numerical investigation of this example revealed that the first Hopf bifurcation occurs at *w*_{1111} ≈ 0.0254 and the second at *w*_{1111} ≈ 0.3395. The asymptotic dynamics on this closed curve can be either periodic or aperiodic and both types of behaviors are exhibited in this case. With standard biallelic expression (*i.e.*, Gavrilets' 1998 model), the only cycles known are of period 2 (Spencer 2003). Consequently, the mean fitness, , need not be maximized (Figure 1a).

The bifurcation diagrams of some of the examples of 2-cycles are also worthy of note: in no cases did we find further bifurcations, giving 4-cycles. Indeed, as shown in Figure 3, 2-cycling appeared and then disappeared as one parameter was varied. This behavior has the appearance of two conjoined pitchfork bifurcations, one reversed, with their prongs aligned. So far as we know, this is also a novel finding in any population genetics model. A number of special cases deserve further analysis.

### Multiplicative maternal selection model:

This special case further assumes that the selective pressures of the maternal effects and ordinary viability selection are independent, as, for example, when selection occurs at two separate stages in the life cycle of each individual, the first as the result of its mother's phenotype and the second as the result of its own phenotype. These effects thus act multiplicatively, and so(5)for *i* = 1, 2, 3. Equations 1 thus become(6)Normalization means that just six parameters—three *m*'s and three *v*'s—now specify the dynamical and equilibrial behavior of the system. Again, up to two distinct polymorphic equilibria may be locally stable, but no cases of cycling were found.

Moreover, since(7)only two of Equations 6 are truly independent. In fact, they can be rewritten in terms of the frequencies, *p*_{f} and *q*_{f} (*p*_{m} and *q*_{m}) of maternally (paternally) derived *A*_{1} and *A*_{2} alleles in zygotes after the maternal-effect selection has acted, in the same way as in the absence of imprinting (Spencer 2003). Since(8)and *p*_{m} = *p*, we have(9)and(10)which are the recursions for different selection pressures on males and females at a single diallelic locus (Pearce and Spencer 1992). The fitness of *A _{i}A_{j}* males is

*v*and that of females,

_{ij}*m*, so (as in the absence of imprinting; Spencer 2003) the maternal component of selection effectively acts on females only.

_{ij}v_{ij}In contrast to this result, the formal equivalence between the multiplicative case of maternal selection and fertility selection in the absence of imprinting (Gavrilets 1998) does not extend to our model with imprinting (see below).

### Complete paternal inactivation:

For many imprinted loci (*e.g.*, *Igf2-r* in rodents), the paternal copy of a gene is effectively silenced. If this silencing occurs to genes in both the mother and offspring, then there are just four distinct fitnesses: *w _{ijkl}* = α

_{ik}(

*i*,

*k*∈ {1, 2}). As before, this allows Equations 1 to be simplified and, indeed, they may be reduced to just two independent recursions:(11)where(12)The natural interpretation of

*p*

_{f}and

*p*

_{m}is a little different from the multiplicative case, however. This time,

*p*

_{f}is the postselection frequency of

*A*

_{1};

*p*

_{m}is the preselection frequency.

This system affords just three equilibria: two trivial fixations and one potential polymorphism given by the pseudo-Hardy–Weinberg form , , and , in which(13)and . The fixation of *A*_{1} is locally stable provided that and that of *A*_{2} if . These fixations can, of course, be simultaneously locally stable. If both these inequalities are satisfied, the polymorphic equilibrium is feasible but unstable; if both are reversed the polymorphic equilibrium is feasible and may be locally stable. If the polymorphic equilibrium remains unstable, genotype frequencies oscillate between two values (see Figure 4, for an example). No cases of cycles of length >2 were found in 10^{7} numerical examples with the four α_{ij} values independently sampled from the uniform distribution between zero and one.

### Complete maternal inactivation:

We now investigate the case in which the maternal allele is completely silenced in both the mother and her offspring. Maternal effects mean that this counterpart to the previous model is not its formal equivalent, in contrast to most pairs of population genetic models of maternal and paternal inactivation (see, *e.g.*, Pearce and Spencer 1992). Supposing that *w _{ijkl}* = β

_{jl}, Equations 1 reduce to(14)Note that the different β

_{jl}values appear a different number of times in Equations 14 than the corresponding α

_{ik}'s in this form of the equations for complete paternal inactivation and, consequently, the system cannot be rewritten in terms of two variables. Nevertheless, as in the complete paternal inactivation case, for various fitnesses the two fixations and a polymorphism are possible equilibria. In contrast, however, oscillations in genotype frequency do not appear to occur: none were found in 10

^{7}numerical cases, each with four β

_{jl}values independently drawn from

*U*[0, 1]. Nor, in spite of this apparent simplicity of the system, have we been able to find an analytical expression for the polymorphic equilibrium, and we cannot show that it is unique. Nevertheless, some 10

^{7}random numerical examples, each with 100 random initial genotype frequencies, failed to reveal any cases with more than one locally stable polymorphism.

### Maternal selection only:

This case assumes that selective differences are due only to the mother's genes. This assumption means that *w _{ijkl}* = γ

_{ij}and Equations 1 reduce to Equations 6 with

*v*= 1 and

_{ij}*m*= γ

_{ij}_{ij}for all

*i*and

*j*. Thus, this system behaves as if no selection acts on males and viability selection with fitnesses γ

_{ij}acts on females. Like the case of complete paternal inactivation, this system affords just three equilibria: two trivial fixations and one potential pseudo-Hardy–Weinberg polymorphism, in which(15)The stability of the equilibria is simple: fixation of

*A*

_{1}is locally stable provided that , fixation of

*A*

_{2}is stable if , and the polymorphism is stable if both these inequalities are reversed—amounting to mean heterozygote advantage—and this condition simultaneously guarantees feasibility.

## MODELS OF FERTILITY SELECTION

### General fertility selection model:

Using the same conventions as for the models of maternal selection, let us suppose that *f _{ijkl}* is the fertility of the cross of

*A*females with

_{i}A_{j}*A*males and that

_{k}A_{l}*A*individuals have viability

_{m}A_{n}*v*. This parameterization of 16 fertility and 4 viability parameters is simply the imprinting version of the model of fertility selection proposed by Bodmer (1965) and leads to the recursion equations for the genotype frequencies in the following generation of(16)

_{mn}At first glance, Bodmer's (1965) model apparently requires nine fertility parameters, corresponding to the 3 × 3 = 9 possible matings between the three different genotypes. But as Feldman *et al.* (1983) pointed out, reciprocal matings between unlike parental genotypes produce the same proportions of offspring genotypes. Consequently, the fertility parameters corresponding to these three crosses always appear together in the recursions, and their average, rather than their individual values, is what matters. This “sex-symmetry” property means there are just six independent parameters in Bodmer's model, since the iterations are unchanged if we set each of these three pairs of fertilities to the same value as their arithmetic means.

Fertility selection with imprinting—Equations 16—does not display this property, however, because imprinting destroys the symmetry of reciprocal matings. For instance, an *A*_{1}*A*_{1} × *A*_{1}*A*_{2} cross, with fertility *f*_{1112}, produces one-half *A*_{1}*A*_{2} offspring, but an *A*_{1}*A*_{2} × *A*_{1}*A*_{1} cross, with fertility *f*_{1211}, gives none. Hence, the parameter *f*_{1112} appears in the recursion for the frequency of *A*_{1}*A*_{2}, , but *f*_{1211} does not, and so both values matter.

The local stability condition for the fixations again reveals the importance of imprinting. For instance, the condition for local stability of the fixation of *A*_{2} is(17)In the absence of imprinting, this inequality collapses to simple heterozygote disadvantage.

Not surprisingly, the large number of parameters of the above model allows some unusual dynamical behavior. Cycling is possible (since it is in the nonimprinting case; Doebeli and de Jong 1998) but it must be quite rare: we found no cycling in 10^{7} simulations.

## DISCUSSION

The models above demonstrate that when selection acts on loci that engender maternal genetic effects and that are subject to genomic imprinting, novel genotype-frequency dynamics may arise. These behaviors include oscillations between two distinct polymorphic values, as well as longer-period cycling lasting many generations due to Hopf bifurcations (also known as Andronov–Hopf bifurcations). This last phenomenon is particularly interesting because, as one fitness parameter is continuously altered (as in a standard bifurcation diagram), these cycles appear at once (at the Hopf bifurcation), rather than as the culmination of a sequence of bifurcations. Moreover, these cycles disappear at a second Hopf bifurcation, as the fitness parameter is further changed. We know of no other single-locus population-genetic models exhibiting such behavior, although single Hopf bifurcations do arise in models of (i) constant viability selection and recombination for two loci each with two alleles (Hastings 1981; Akin 1982, 1983) and (ii) constant selection and multilocus mutation with selfing (Yang and Kondrashov 2003).

The way in which the cycling between two polymorphic values occurs (for certain fitness values) in some of the models is also noteworthy. In all cases examined, a single locally stable polymorphic equilibrium becomes unstable at a pitchfork bifurcation, as one parameter is varied (*e.g.*, Figures 3 and 4). The population then oscillates between two polymorphic values. In the general model of maternal effects with viability selection, however, further changes in the parameter sometimes led to these two values and the unstable equilibrium coalescing into a locally stable equilibrium again (Figure 3).

Nevertheless, although these mathematical properties are interesting, examination of the various models with randomly assigned fitnesses showed that cycling of all types was rare in parameter space. Moreover, the examples in Figures 1 and 3 exhibit strong interactions between maternal and offspring genotypes that generate large differences in fitnesses between, say, the same offspring genotypes with different maternal genotypes. Whether or not evolution disproportionately favors such parameters is a separate question, of course (Spencer and Marks 1988; Marks and Spencer 1991), but the biological importance of cycling in these models is certainly open to question. Nevertheless, any random allele frequency changes due to genetic drift are unlikely to prevent cycling (unless they lead to the fixation of one allele) because in all cases these cycles are attracting.

Two of the special cases—complete paternal inactivation and maternal selection only—possess another interesting feature: the “pseudo-Hardy–Weinberg” form of the sole polymorphic equilibrium even though the recursions cannot be reduced to those in allele frequencies. Previous models concerned with the evolution of imprinting from standard expression (Spencer *et al.* 1998, 2004) have found a similar result, but these latter models all involved both imprintable and unimprintable alleles.

Pearce and Spencer (1992) investigated the effect of imprinting on standard models of viability selection. They found that in all cases—except for that of different selection pressures on males and females—the models were formally equivalent to models without imprinting (but with suitably adjusted viabilities). In contrast, none of the models derived above have any formal equivalence to known models without imprinting. Moreover, the properties of some models without imprinting are destroyed by imprinting: the sex-symmetry property of fertility selection (Feldman *et al.* 1983), for example, does not hold in the presence of imprinting. This work thus adds to the growing literature that shows how standard biallelic Mendelian expression permits several simplifications in population-genetic and quantitative-genetic models (Spencer 2002). Thus, although the number of imprinted genes is small (Morison *et al*. 2005), their very existence illuminates our understanding of population-genetic and other processes.

## APPENDIX: HOPF BIFURCATIONS IN THE GENERAL MATERNAL SELECTION MODEL

The structure shown in the bifurcation diagram of Figure 1b suggests a supercritical Hopf bifurcation since the equilibrium point bifurcates not via period doubling cascade, but immediately into a region where the dynamics are more complicated. The oscillatory behavior illustrated in Figure 1a reinforces the suggestion that a Hopf bifurcation has occurred. In a supercritical Hopf bifurcation, a stable equilibrium becomes unstable and is surrounded by an attracting invariant circle as the parameter is varied. The dynamics on this circle may be periodic or aperiodic, but in either case are oscillatory in nature.

The Hopf bifurcation theorem for maps describes the necessary conditions for the existence of such a bifurcation. Let with , , and *F* three times differentiable. For a Hopf bifurcation to occur an equilibrium point *x*_{0} and a parameter value must exist such that the linearization of with respect to *x* at has a pair of complex conjugate eigenvalues with *a* and *b* nonzero and . (There is one additional restriction that if *a =* then .) Moreover, one must check that as the parameter is varied monotonically about , the quantity goes from a value <1 to a value >1 or vice versa.

The algebraic complexity of this model has so far precluded an analytical verification of the hypotheses of the Hopf bifurcation theorem. However, we have numerically verified these conditions using the software Maple 10. The program (Dorn 2005) works as follows:

Increment the parameter

*w*_{1111}from 0 to 0.4.For each parameter value perform the following:

Numerically compute the critical point.

Compute the four eigenvalues of the linearization at this critical point.

Compute the modulus for each of these eigenvalues.

Plot each of the moduli against the corresponding parameter value.

Figure A1 shows the result of this process. The top sequence plots the moduli of the complex conjugate pair of eigenvalues as a function of *w*_{1111}. The horizontal line at *y* = 1 represents the threshold that must be crossed for a Hopf bifurcation to occur. Note that this sequence crosses this line at ∼*w*_{1111} = 0.02 and again at *w*_{1111} = 0.33. Further refinements of this program yielded more accurate bifurcation values of *w*_{1111} = 0.0254 and *w*_{1111} = 0.3395. For these two approximate bifurcation values we manually verified that these eigenvalues satisfy the conditions on *a* and *b* given above. Finally, the crossing of the sequence through the line *y* = 1 at both *w*_{1111} = 0.0254 and *w*_{1111} = 0.3395 verifies the final condition stated above.

We note for completeness that there is one other sequence plotted in Figure A1. This represents the real, nonzero eigenvalue of the linearization. The fourth eigenvalue equals zero for all parameter values due to the normalization of Equations 1.

## Acknowledgments

We thank Mike Paulin for discussions about the models and Ken Miller for assistance with the figures. Two anonymous reviewers also provided helpful suggestions. Additional thanks are due to Gustavus Adolphus College and the Allan Wilson Centre at Massey University for funding and hosting the sabbatical for T.L. that facilitated this collaboration. Financial support for this work was provided by the Marsden Fund of the Royal Society of New Zealand contract U00-315 (H.G.S.).

## Footnotes

↵2

*Present address:*Department of Mathematics, Oregon State University, Corvallis, OR 97331-4605.Communicating editor: M. K. Uyenoyama

- Received February 22, 2006.
- Accepted June 4, 2006.

- Copyright © 2006 by the Genetics Society of America