Viability, behavior, and color expression in the offspring of matings between common wall lizard Podarcis muralis color morphs.

Color polymorphisms are widely studied to identify the mechanisms responsible for the origin and maintenance of phenotypic variability in nature. Two of the mechanisms of balancing selection currently thought to explain the long-term persistence of polymorphisms are the evolution of alternative phenotypic optima through correlational selection on suites of traits including color and heterosis. Both of these mechanisms can generate differences in offspring viability and fitness arising from different morph combinations. Here, we examined the effect of parental morph combination on fertilization success, embryonic viability, newborn quality, antipredator, and foraging behavior, as well as inter-annual survival by conducting controlled matings in a polymorphic lacertid Podarcis muralis, where color morphs are frequently assumed to reflect alternative phenotypic optima (e.g., alternative reproductive strategies). Juveniles were kept in outdoor tubs for a year in order to study inter-annual growth, survival, and morph inheritance. In agreement with a previous genome-wide association analysis, morph frequencies in the year-old juveniles matched the frequencies expected if orange and yellow expressions depended on recessive homozygosity at 2 separate loci. Our findings also agree with previous literature reporting higher reproductive output of heavy females and the higher overall viability of heavy newborn lizards, but we found no evidence for the existence of alternative breeding investment strategies in female morphs, or morph-combination effects on offspring viability and behavior. We conclude that inter-morph breeding remains entirely viable and genetic incompatibilities are of little significance for the maintenance of discrete color morphs in P. muralis from the Pyrenees.

Population polymorphism refers to a situation where 2 or more discrete phenotypic forms coexist in sympatry among individuals of the same sex and age (i.e., morphs), with the rarer morph being too frequent to be solely the result of recurrent mutation (Ford 1945; Huxley 1955). Often, this phenotypic variation takes the form of distinct color morphs that are genetically inherited, fixed in adults, and not condition-dependent (i.e., color polymorphism; Galeotti et al. 2003; Roulin 2004; Mckinnon and Pierotti 2010). Explaining the maintenance of color polymorphism in natural populations remains one of the major challenges in evolutionary biology as it requires some form of selective balance to counteract the eroding effects of differential morph fitness and drift (Galeotti et al. 2003; Roulin 2004; Gray and McKinnon 2007).

Color morphs often differ in features other than color (Forsman et al. 2008; Mckinnon and Pierotti 2010). Such co-variations are thought to have originated via disruptive, correlational selection favoring optimal trait combinations in the different morphs, each representing alternative adaptive peaks (i.e., alternative strategies; Sinervo and Svensson 2002; Shuster and Wade 2003; Gray and McKinnon 2007). Consider the classical example of the peppered moth Biston betularia, with 2 alternative morphs conferring a cryptic advantage in different habitats. Genetic correlations between color morph and behavior (i.e., choice of resting sites) would be expected to arise, as moths that choose backgrounds matching their own color survive, whereas those with the wrong preference perish (Kettlewell 1955; Majerus 1998; Sinervo and Svensson 2002; Cook et al. 2012). Color polymorphic strategies can persist over long timescales only if they have the same average fitness over time (determined by background availability in the Peppered moth example), which may occur under balancing selection, including frequency-dependent or spatio-temporally varying selection (Zamudio and Sinervo 2003; Roulin 2004; Gray and McKinnon 2007; Mckinnon and Pierotti 2010; Wellenreuther et al. 2014). However, the observation of alternative color-coded strategies within a population does not imply that the polymorphism is balanced and stable in time. In fact, the evolution of alternative strategies may have (albeit not always) consequences for reproductive compatibility between the morphs (Roulin 2004; Gray and McKinnon 2007; Mckinnon and Pierotti 2010). This is because correlational selection is expected to favor the evolution of assortative mating patterns and genetic mechanisms ensuring the co-segregation of optimal trait combinations. Such correlations are often mediated through linkage disequilibria (which can sometimes create genetic divergence and hence, post-zygotic incompatibilities between the morphs), or pleiotropy (which is not expected to affect morph compatibility; Denoël et al. 2001; Sinervo and Svensson 2002; Gray and McKinnon 2007; Pryke and Griffith 2009; Svensson et al. 2009; Mckinnon and Pierotti 2010; Bastiaans et al. 2014; Küpper et al. 2015; Svensson 2017). For instance, in the polymorphic Gouldian finch Erythrura gouldiae, where 2 Z-linked color morphs correlate with 2 alternative competitive strategies in males (Pryke and Griffith 2006), offspring produced by genotypically mixed mating pairs shows consistent inviability effects from egg to adulthood (Pryke and Griffith 2009).

Heterosis, that is, the situation in which heterozygosity at the genes involved in the expression of color polymorphism confers a fitness advantage, represents an alternative mechanism of balancing selection by which sympatric color morphs may be maintained in natural populations (Roulin 2004; Roulin and Bize 2007; Gratten et al. 2008, 2010; Johnston et al. 2013; Timberlake 2013; Wellenreuther et al. 2014). Increased fitness in heterozygotes may result from heterozygosity at a single locus with pleiotropic effects, or from epistatic interactions among multiple heterozygous loci (i.e., when polymorphisms are governed by several tightly linked genes with major phenotypic effects; Gray and McKinnon 2007; Jamie and Meier 2020). Heterosis thus differs from disruptive selection in that genetically intermediate (i.e., heterozygous) individuals are favored, and inter-morph reproductive barriers are never expected to arise. In fact, if heterozygosity at the genes coding for the color polymorphism (or at genes linked to them) provides a fitness advantage, selection is likely to favor the evolution of disassortative mating preferences and offspring of same-morph parents may suffer a fitness cost (Krüger et al. 2001; Tuttle 2003; Roulin 2004). In the color polymorphic white-throated sparrow Zonotrichia albicollis, genetic incompatibilities caused by a strongly deleterious recessive allele (i.e., a large chromosomal inversion) are thought to explain the pronounced pattern of disassortative pairing observed in the wild, which plays a key role in maintaining the polymorphism (Thorneycroft 1975; Houtman and Falls 1994; Tuttle 2003). However, it is worth noting that if the heterozygous advantage concerns only viability selection (i.e., survival to adulthood), the color polymorphism may be maintained even if the morphs show random or even maladaptative mating preferences (Krüger et al. 2001; Roulin and Bize 2007).

Lizards are an excellent model group in which to study the evolutionary dynamics associated with color polymorphism. Many phylogenetically distant taxa have convergently evolved discrete color morphs which are heritable (often involving few autosomal genes with major effects in the few species examined to date; Olsson et al. 2007; Sinervo et al. 2010; Haisten et al. 2015; Rankin et al. 2016; McLean et al. 2017; Stuart-Fox et al. 2020). Morphs may be functionally relevant to crypsis or thermoregulation, and/or correlated to other physiological or behavioral traits under selection (Zamudio and Sinervo 2003; Olsson et al. 2009; Baird 2013; Lattanzio and Miles 2014; Ortega et al. 2015; Yewers et al. 2016; Yewers 2017; Fernández et al. 2018; Stuart-Fox et al. 2020). Furthermore, polymorphic lizards sometimes represent different points in the balance between color polymorphism maintenance and incipient isolation. Whereas there is compelling evidence for the long-term maintenance of color polymorphism via frequency-dependent selection in some cases (i.e., Uta stansburiana;Sinervo and Lively 1996; Sinervo et al. 2000), color-assortative pairing, limited gene flow, and weak post-zygotic barriers to inter-morph breeding have also been described (Rosenblum 2006; Gray and McKinnon 2007; McLean and Stuart-Fox 2014; Mclean et al. 2014a, 2014b; Dong et al. 2019; McLean et al. 2020). For instance, in the side-blotched lizard U. stansburiana, where 3 morph-specific sexual behaviors (in males) and 2 alternative breeding strategies (in females) persist cycling in frequency over time (Sinervo and Lively 1996; Sinervo et al. 2000), females have been found to show morph-assortative mate choice preferences and decreased egg viability in disassortative pairings (Bleay and Sinervo 2007; Lancaster et al. 2014).

In wall lizards (genus Podarcis), a widely distributed group present in most of Europe, northwestern Africa, and the Mediterranean islands, several species exhibit a striking ventral color polymorphism (Arnold et al. 2007; Huyghe et al. 2007; Sacchi et al. 2007; Runemark et al. 2010; Marshall and Stevens 2014; Brock et al. 2020a, 2020b). One such species (which is also the most widely distributed) is the European common wall lizard Podarcis muralis. In this species, adults of both sexes may show up to 5 alternative ventral color morphs: 3 uniform (pure) morphs, that is, orange (o), white (w), and yellow (y), and 2 mixed-morph mosaics combining orange and white (ow) or yellow and orange (yo) (Fig. 1; Sacchi et al. 2007, 2013; Pérez i de Lanuza et al. 2013, 2019). Published spectra show that chromatic differences between pure color morphs are due mainly to variation in reflectance between 400 and 600 nm (i.e., discriminable by humans), and further empirical evidence confirms that the morphs are perceived as chromatically distinct by conspecifics (Pérez i de Lanuza et al. 2014, 2018a). No study has established morph inheritance, but recent research suggests that orange and yellow color expressions are caused by recessive homozygosity at 2 separate loci in the regulatory regions of 2 genes associated with pterin (sepiapterin reductase; SPR) and carotenoid (beta-carotene oxygenase; BCO2) metabolism, respectively (Andrade et al. 2019). Because of their ventral position and conspicuousness to the species’ visual system (Pérez i de Lanuza and Font 2015; Pérez i de Lanuza et al. 2018a), at least the pure color morphs of P. muralis have often been thought to represent the visible mark of an underlying set of alternative phenotypic optima (Sinervo et al. 2007; Calsbeek et al. 2010; Scali et al. 2013, 2016; Sacchi et al. 2017a; Mangiacotti et al. 2019). Although recent evidence from a mesocosm experiment does not support the existence of alternative strategies concerning socio-sexual and spatial behavior (Abalos et al. 2020), research on P. muralis color morphs has detected inter-morph differences in several (mostly physiological) traits relevant to fitness such as immune response, hormonal profile, and prevalence of infection by parasites, hence suggesting the existence of alternative trade-offs in resource allocation and life-history strategies in the different morphs (Galeotti et al. 2007, 2010; Calsbeek et al. 2010; Zajitschek et al. 2012; Sacchi et al. 2017a, 2017b). Female color morphs from Italian polymorphic populations have also been found to show differential breeding investment, and clutch viability and newborn mass were also found to vary depending on parental morph combination (Galeotti et al. 2013). Lastly, field observations suggest an assortative morph bias in pair formation (68% of male–female pairs observed in the wild are homomorphic; Pérez i de Lanuza et al. 2013, 2016).

Taken together, these findings suggest that disruptive, correlational selection coupled with assortative mating may have favored different co-adapted trait complexes in each morph while perhaps introducing some degree of genetic divergence between the color morphs. These differences, though minor (Andrade et al. 2019), could still have important implications on offspring viability (i.e., the probability that an organism may develop normally and survive long enough to reproduce) if they play an organizational role in early development; Sinervo and Svensson 2002; Gray and McKinnon 2007; Mckinnon and Pierotti 2010; Svensson 2017; Jamie and Meier 2020). If this is the case, we would expect to find decreased embryonic viability, newborn quality, and survival to adulthood in the offspring of disassortative matings. In contrast, higher viability and overall performance in the offspring of disassortative pairings may suggest the existence of heterosis in the genes coding for the expression of P. muralis polymorphic coloration. The genetic basis of P. muralis color polymorphism prevents genotype to be directly inferred from color morph. However, we know that the offspring of disassortative matings will have higher probability of being heterozygous for either loci, and that the offspring produced in crosses involving white lizards will be the most variable genotypically (Tables 1 and 2). The magnitude of genetic incompatibilities is often larger in complex traits involving the expression of many genes and closely related with fitness (Mousseau and Roff 1987; DeRose and Roff 1999; Ala-Honkola et al. 2009). Beyond fertilization success or egg viability, behavior represents the end product of many genes interacting with other genes and with the environment, and is also believed to play a crucial role in determining survival into adulthood in juvenile lizards (Civantos and Forsman 2000; Clobert et al. 2000; Miles 2004). Thus, searching for morph combination effects on behaviors directly related with survival into adulthood such as foraging ability and predator avoidance seems key to evaluate the existence of any potential incompatibilities affecting offspring fitness.

Table 1.

Theoretical allelic frequencies in the gametes of P. muralis lizards according to their color morph

Color morph	Possible genotypes	Segregation				Gametes
		Locus O		Locus Y		oy	oY	Oy	OY
		o	O	Y	Y
O	1/3 ooYY	1	0	1/3	2/3	1/3	2/3	0	0
	2/3 ooYy
Y	1/3 OOyy	1/3	2/3	1	0	1/3	0	2/3	0
	2/3 Ooyy
W	1/9 OOYY	1/3	2/3	1/3	2/3	1/9	2/9	2/9	4/9
	2/9 OOYy
	2/9 OoYY
	4/9 OoYy
YO	1/1 ooyy	1	0	1	0	1	0	0	0

Genotypic frequencies in the parental lizards are calculated assuming panmixia and allelic frequencies in the gametes are calculated assuming Mendelian segregation in the 2 unlinked and bi-allelic loci described in Andrade et al. (2019).

Table 2.

Theoretical offspring genotypic frequencies per parental morph combination

Morph cross	Possible offspring genotypes									Heterozygous for locus		Double homozygous
	ooYY	ooYy	ooyy	OoYY	OoYy	Ooyy	OOYY	OOYy	OOyy	O	Y
O × O	0.44	0.44	0.11	—	—	—	—	—	—	0	0.44	0.55
O × W	0.15	0.11	0.07	0.30	0.30	0.07	—	—	—	0.67	0.41	0.22
O × Y	—	0.22	0.11	—	0.44	0.22	—	—	—	0.67	0.67	0.11
W × O	0.15	0.11	0.07	0.30	0.30	0.07	—	—	—	0.67	0.41	0.22
W × W	0.05	0.05	0.01	0.20	0.20	0.05	0.20	0.20	0.05	0.44	0.44	0.31
W × Y	—	0.07	0.07	—	0.30	0.11	—	0.30	0.15	0.41	0.67	0.22
Y × O	—	0.22	0.11	—	0.44	0.22	—	—	—	0.67	0.67	0.11
Y × W	—	0.07	0.07	—	0.30	0.11	—	0.30	0.15	0.41	0.67	0.22
Y × Y	—	—	0.11	—	—	0.44	—	—	0.44	0.44	0	0.55
Offspring morph	O	O	YO	W	W	Y	W	W	Y	—	—	—

Genotypic frequencies in the parental lizards are calculated assuming panmixia (see Table 1) and offspring genotypic frequencies are calculated by multiplying the gametic frequencies in Table 1. The comparison of assortative crosses between orange and yellow lizards (0.55 probability of being homozygous for both loci) with the offspring of yellow–orange crosses (0.11 probability of being a double homozygote) can be particularly informative about the existence of any potential homozygote advantage or heterosis. In addition, notice that doubly heterozygous (OoYy) offspring (in which heterosis should be stronger) cannot be produced in orange-orange and yellow-yellow crosses.

Here, we tested for morph combination effects on fitness using controlled matings among all-pair combinations of P. muralis pure color morphs, an ecologically relevant question considering the high prevalence of color-assortative mating observed in wild populations. Our objective was 4-fold: (1) to compare primary reproductive parameters (i.e., clutch size and juvenile mass) among female morphs in order to evaluate the existence of alternative breeding strategies; (2) to evaluate the potential existence of advantages to mating with respect to morph at both the pre- and post-zygotic levels by comparing fertilization success, embryonic viability, and newborn quality in matings between different parental morph combinations; (3) to test for morph combination effects on viability selection (i.e., survival to adulthood) by assessing the antipredator response and foraging ability of newborn lizards, as well as raising them in a common garden and keeping track of their inter-annual growth and survival; and (4) to describe morph inheritance in the year-old juveniles.

Materials and Methods

Lizard capture and female hibernation

In total, we captured 135 female and 44 male adult lizards by noosing from 14 polymorphic localities (>1 km apart) in the eastern Pyrenees (Supplementary Table S1). In each of these localities, we captured lizards showing each of the 3 pure color morphs. We only captured lizards with a snout-to-vent length (SVL) larger than 56 mm, a reliable cut-off point for adult morph expression in P. muralis polymorphic populations from eastern Pyrenees (Pérez i de Lanuza et al. 2013). We captured 60 females at the end of the 2018 breeding season (September 5–8) and transferred them to the Station d’Ecologie Théorique et Expérimentale (Moulis, France), where they were kept under natural conditions in 45 outdoor plastic tubs (170 cm diameter, 60 cm high; Supplementary Figure S1). Each tub housed a group of 3–5 females coming from the same locality. Tubs were covered with anti-bird netting to prevent predation while allowing local invertebrates to get inside the tubs and serving as a primary food source for the lizards. To increase habitat complexity, in addition to the naturally occurring vegetation inside the tubs, we provided a perforated construction brick to serve as shelter, 2 wood logs, and a water dish to collect rainwater. We visited the tubs monthly to trim the weeds and provide additional food (Tenebrio molitor larvae and Acheta domesticus adults). On March 2019, we captured 75 additional females from the same localities, re-captured the females within the outdoor tubs, and installed all of them in individual hibernation pods consisting of a plastic bowl (12 cm diameter, 8 cm high) covered with a perforated lid and filled with moist coco husk (1:2 coco: water by weight. Because environmental conditions at Moulis are milder than at the Cerdanya plateau, we hibernated both batches of females to (1) stimulate reproductive physiology (in order for vitellogenesis to be completed, females need to hibernate ˂5 °C for a period determined by the duration of winter at the species natural range; Gavaud 1983, 1991; Carretero 2007), (2) promote synchronicity between both batches of females, and (3) delay the onset of their cycle until early May, which typically marks the start of the P. muralis breeding season at the Cerdanya plateau. Hibernation pods were thus stored for 9 weeks in a dark climatic chamber (1st week: 10 °C/2nd to 8th weeks: 5 °C/9th week: 10 °C). On 8 May 2019, we took the females out of hibernation and installed them into the outdoor tubs for acclimation, whereas we captured the males at the same localities as the females. Two days before the onset of the experiment, we measured SVL (0.1 mm) and mass (±0.01 g) of each lizard with a ruler and a spring balance (Pesola, Schindellegi, Switzerland). We also removed ∼5 mm from the tail tip of each individual and preserved the tissue in 90% ethanol for genetic analyses.

Staged mating trials and juvenile husbandry

To study the effect of parental morph combination on reproductive output and juvenile viability, we carried out staged mating trials by releasing 131 females (45o, 44w, and 42y) and 44 males (14o, 15w, and 15y) into 44 outdoor plastic tubs (Supplementary Figure S1). We designed mating trials to guarantee fertilization of triads of females (one for each pure morph) by a single male. In each tub, we released 3 females (1o: 1w: 1y) followed 5 days later (14 May 2019) with 1 male showing either orange, white, or yellow ventral coloration (N = 14 tubs for each male ventral color). Due to female post-hibernation mortality, the yellow female morph was absent in the 2 remaining tubs (a: 2o♀, 1w♀, 1w♂; b: 1o♀, 1w♀, 1y♂). To reduce any possible effects of female size asymmetry, time of capture, and prior social interactions, we allowed a maximum SVL difference of 2 mm among females housed in the same tub, minimized co-habitation of females captured at different times, and only put lizards together in the same experimental tub if they had been captured at least 300 m apart. Prior to release, we marked females with a nontoxic permanent marker (Edding®751) to facilitate individual recognition without the need to recapture them (Supplementary Figure S2; Ferner and Plummer 2016; Ngo et al. 2016). In total, the males and females spent 30 and 35 days, respectively, inside the experimental tubs. During this period, we monitored male–female interactions by visiting the outdoor tubs twice per day and recording the identity of the male–female pairs we sighted either (1) co-perching (i.e., lying in physical contact) or (2) engaged in copulatory behavior (i.e., intertwined in the doughnut-shaped mating posture typical of lacertid lizards; Weber 1957; Carpenter and Ferguson 1977; in Den Bosch and Zandee 2001).

On 15 June 2019, we released the males at their capture location (previously determined using a GPS device), and housed females in the laboratory until oviposition (see Supplementary Appendix S1). We lost 9 clutches due to failure to retrieve the females from the outdoor tubs before oviposition (or locate the clutch). In addition, 10 females did not produce a clutch (due to a lack of fecundation or suboptimal conditions). For the remaining 112 females, we counted the number of fertile (452) and infertile (105) eggs within each clutch by noting the presence of a calcified shell and vascularization 48 h after oviposition (Köhler 2006). We incubated the 452 fertile eggs in plastic cups filled with moist coco husk (1:2 coco: water by weight) and covered it with a perforated lid at a constant temperature of 28 °C (Van Damme et al. 1992; Braña and Ji 2000) in an incubator (Memmert GmbH + Co.KG, Schwabach, Germany). On hatching, each of the 424 born juveniles was measured (SVL), weighted, sampled for DNA, temporarily marked, and housed in 55 × 38 × 28 cm high glass terraria according to their mother’s tub assignment. For 28 embryos that died before hatching, we obtained DNA samples via dissection of the eggs. Juvenile lizards were fed small live prey (pinhead A. domesticus and flightless Drosophila hydei) and dusted with vitamins. To confirm exclusive paternity from the experimental male and to determine parentage in those clutches found in the outdoor tubs, we conducted genetic analyses of paternity based on 6 known microsatellite loci described in P. muralis (Richard et al. 2012; Heathcote et al. 2014). Thus, we could reliably assign paternity to every offspring examined (strict: 271 juveniles, relaxed: 41 juveniles; see Supplementary Appendix S1). On 18 August 12019, after conducting lab-staged behavioral tests on a subset of the juveniles (see below), we released 372 juveniles (52 juveniles died shortly after hatching) in groups of 10–12 unrelated individuals into 33 outdoor plastic tubs for their long-term maintenance. Prior to release, we permanently marked the juveniles on the ventral scales using a disposable medical cautery unit (Ekner et al. 2011).

Juvenile behavioral tests

To evaluate prospective fitness, 3–6 days after hatching, we conducted lab-staged behavioral tests on 2 juveniles from each clutch. Offspring mass being often considered a proxy of newborn quality (Stamps 1988; Civantos and Forsman 2000; Galeotti et al. 2013), we chose the heaviest and lightest from each clutch to encompass the offspring mass range of each female. These tests were repeated 1 year later for the 50 surviving juveniles found in the outdoor tubs in July 2020. Each test comprised 3 stages each lasting 5 min: acclimation period, antipredator trial, and foraging trial (see Supplementary Video S1). We conducted tests during the lizards’ natural daily period of activity (10:00 to 13:30 and 16:00 to 19:00 h local time). In total, we conducted 247 behavioral tests. We started each test by transferring a juvenile lizard inside a plastic cup into an experimental arena consisting of a 50 × 25 × 40 cm high glass terrarium with all except the front wall covered with cardboard, and a piece of filter paper as substrate. After the 5-min acclimation period, in which the lizard could move freely across the arena, we simultaneously introduced a shelter in the arena and simulated the presence of an avian predator using an artificial plastic model with the size, shape, and coloration of an Eurasian magpie (Pica pica), a generalist native predator that feeds on wall lizards (Martin and López 1990). The model magpie was placed in front of the arena, at a distance of 5 cm from the arena’s front wall, and at the end of a 75 cm pole. The shelter (consisting of 2 superimposed tiles separated along the edges by thin pieces of felt on 3 sides so as to leave a single entrance to a 90 × 90 × 12 mm hollow) was positioned in the middle and in contact with the back wall, with its entrance facing the front wall (Supplementary Figure S3). After 5 min, we removed the magpie model and introduced a novel prey item (one <8 mm T. molitor larvae) into the arena. In every trial, we positioned the larvae in the juvenile’s field of vision and at a distance of 4–8 cm from its head (or the shelter’s entrance, whenever the juvenile was sheltered). We stopped the trial 5 min after the prey item was introduced. To minimize the impact of human presence on the lizard’s behavior, all trials were filmed remotely with a photographic camera (Canon®, EOS 7D, Tokyo, Japan). During trials, the arena was illuminated by 2 40 W lamps (providing a uniform 25 °C temperature in the arena), and a full spectrum lamp (Reptistar 5.0) equipped with a high-frequency ballast (Quicktronic, Osram, Munich, Germany). The experimental arena and shelter were thoroughly cleaned with alcohol and rinsed with water after every test. We played back each filmed test and a single researcher recorded: (1) the latency to enter the shelter once presented with the magpie model; (2) the latency to bite the novel prey item; and (3) the time elapsed between prey capture and consumption (prey handling time).

Inter-annual survival and morph inheritance

In early July 2020, we re-captured the juvenile lizards installed in our outdoor tubs. We also searched for marked lizards in adjacent patches of lizard habitat (∼1 km from the study site) to account for any potential escape, but none were found. We measured SVL and mass in each of the 50 surviving juveniles. We also assigned them visually to a color morph according to their throat coloration (Pérez i de Lanuza et al. 2013). We then used the spreadsheet available in Rankin et al. (2016) to calculate the morph frequencies expected in the offspring of each morph combination and compare them to the observed frequencies. Following Andrade et al. (2019), we considered 2 different loci each with 2 alleles (O/o, Y/y), and a recessive homozygous base for orange and yellow color expression. Genotypic frequencies in the parental lizards were obtained assuming panmixia. To increase sample size, we added data from 17 juvenile lizards of known parental morph combinations which were conceived during a previous experiment (see Abalos et al. 2020) and also maintained in the same outdoor tubs for a similar period (September 2018 to July 2019). We compared observed and expected offspring morph frequencies for each parental morph combination using G likelihood-ratio tests.

Statistical analyses

To account for the random factors in our experimental design, we ran generalized linear mixed models (GLMMs) using the lme4 package (Bates et al. 2014) in R (R Core Team 2019). Model selection was conducted using backward single-term deletions (P < 0.05) of the saturated model followed by model comparisons via likelihood ratio tests (at α = 0.05). To examine the effect of parental morph combination on clutch size and juvenile mass, we used linear mixed models (LMMs) including female capture locality, capture year, and experimental tank, as random factors. To investigate the effect of parental morph combination on inter-annual survival, growth, and performance during the behavioral tests, we used logistic mixed models including the same random factors (in addition to mother identity for the models on behavior). All numerical variables were centered and scaled before running the models (Schielzeth 2010). For each fitted model, we calculated marginal pseudo-r (following Nakagawa et al. 2017), checked the absence of multicollinearity among predictors, influential data points, and graphically explored that residuals conformed to homoscedasticity and normality assumptions (when pertinent) using the performance package in R (Lüdecke et al. 2020). For logistic models, we checked that there was a linear relationship between the logit of the outcome and each of the predictor variables. We also checked for overdispersion using the PsychHelperFunctions package (Huff 2020). For further details on the statistical analyses, see Supplementary Appendix S1.

Results

Effect of parental morph combination on reproductive output

There was no morph bias in body size or mass in the adult lizards (Supplementary Table S2). In the experimental tubs, males and females were frequently observed perching in close proximity (N = 1002) or engaged in copulatory behavior (N = 102), but color morph played no role on the patterns of male–female association (co-perching: χ = 8.33, P = 0.40; copulatory behavior: χ = 1.64, P = 1.0; Supplementary Figure S4). The levels of fertilization success (i.e., fertilized eggs/clutch size) were 81.4% and embryonic viability (i.e., viable juveniles/fertilized eggs) was 93.7%. Egg incubation lasted 37.7 ± 0.4 (x± 95% CI henceforth) days with a mean clutch size of 4.6 ± 0.2 eggs. Clutch size correlated significantly with female weight and juvenile mass (LMM: female weight: χ = 15.88, P < 0.001; juvenile mass: χ = 20.23, P < 0.001; pseudo-r = 0.32), but not with female morph, male morph, or the interaction between the two (P > 0.2; Supplementary Table S3). Similarly, we found no difference in clutch size between homomorphic and heteromorphic pairings (P > 0.2; Supplementary Figure S5 and Table S4). The absolute number of viable juveniles produced was also independent of parental morph and morph combination (Supplementary Tables S3 and S4). However, we found a significant morph bias in fertilization success (GLMM binomial: female morph × male morph: χ = 19.50, P < 0.001; Supplementary Table S3), which was mainly driven by the large number of infertile clutches laid by yellow females paired with orange males (Figure 2). Eggs produced by white homomorphic pairs (which showed the highest fertilization success) were 1.61 (odds ratio [OR] 95% CI = 0.99, 2.79) times more likely to be fertile than eggs produced by Y–O (♀–♂) pairs. We also found significantly lower fertilization success in heteromorphic pairs with respect to homomorphic matings (GLMM binomial: χ = 9.46, P = 0.002; Supplementary Table S4), although the OR for this difference was small (OR [95% CI] = 1.11 [0.86, 1.44]). Embryonic viability was unaffected by parental morph combination (Figure 2), and showed only a weak negative relationship with female mass (GLMM binomial: female morph × male morph: P > 0.5; female mass: P = 0.012; Supplementary Tables S3 and S4). Juvenile mass was negatively correlated with clutch size (GLMM: χ = 21.75, P < 0.001, pseudo-r = 0.21; Supplementary Table S4) but unrelated to parental morph or morph combination (P > 0.4; Supplementary Tables S3 and S4; Figure 3). According to our model, each additional egg laid by a female entailed a decrease of 0.02 ± 0.01 g in the average juvenile mass of its clutch, ∼5% of the average mass of juveniles in our sample (0.38 ± 0.01 g). We did not find evidence of differential solutions to the trade-off between offspring quantity and quality among female pure morphs.

Figure 2.

(A) predicted probability that an egg is fertilized (i.e., fertilization success) according to morph combination. (B) predicted probability that a fertilized egg turns into a viable juvenile (i.e., embryonic viability). Circles mark the position of the mean and bars extend to the 95% CI. Numbers and numbers within brackets correspond to the number of females and eggs within each category, respectively.

Figure 3.

(A) boxplots showing average juvenile mass by female (g) according to morph combination. Boxes indicate the interquartile range (IQR, 50% of data). Horizontal lines represent the median and bars extend to 1.5 times the IQR. A jittered dot cloud represents each value in our dataset. (B) scatter plot showing the negative relationship between average juvenile mass and clutch size. Colored lines and shaded areas correspond to the linear regression slope and 95% CI for each female morph.

Color variation in the ventral surface of adult P. muralis lizards. Both sexes may show up to 3 pure color morphs (white, orange, and yellow), or 2 intermediatemosaics combining ow or yo scales (i.e., mixed-morphs). In males, the color covers the entire ventral surface whereas in females it is sometimes restricted to the throat. Photos by Javier Abalos.

Effect of parental morph combination on juvenile behavior

In 2019, we run 197 behavioral tests on the offspring of 99 different females representing 42 of the outdoor experimental tubs (see Figure 4 for information on sample size per morph combination). Nearly half of the juveniles entered the shelter on presenting them with the magpie model (98 juveniles, 49.7%, mean latency = 133 ± 19 s), whereas two-thirds of the juveniles captured and consumed the novel prey item (133 juveniles, 67.5%, mean latency = 57 ± 11 s). We observed a positive relationship between juvenile mass and performance in the viability tests: the probability that a juvenile entered the shelter doubled with each 0.2 g increase in juvenile mass (GLMM binomial: χ = 4.12, P = 0.042, pseudo-r = 0.03; Supplementary Table S5; Figure 4), whereas latency until the novel prey item was captured (GLMM: χ = 22.67, P < 0.001, pseudo-r = 0.16; Supplementary Table S5; Figure 4) and consumed (GLMM: χ = 22.10, P < 0.001, pseudo-r = 0.16; Supplementary Table S6) decreased exponentially with growing juvenile mass. However, we did not find any significant effect of parental color morph (or morph combination) on the juveniles’ antipredator response or foraging ability (P > 0.2, Supplementary Tables S5 and S6). In July 2020, we repeated the behavioral tests with the 50 surviving juveniles, corresponding to the offspring of 34 different females (23 of the outdoor experimental tubs). Yearling lizards entered the shelter and ate the prey more frequently than the animals tested few days after hatching (GLMM binomial: enter: 73.9%, χ = 9.43, P = 0.004; eats: 87.8%, χ = 7.23, P <0.007; Supplementary Figure S6). They also showed longer latencies to enter the shelter (GLMM: χ = 16.16, P < 0.0001), but no difference in the latency to consume the prey (Supplementary Figure S6). Again, we did not find any significant effect of parental morph combination on the juvenile’s antipredator response or foraging ability (P > 0.05; Supplementary Table S7).

Figure 4.

(A) Barplot showing the proportion of juveniles that entered the shelter during antipredator trials grouped by parental morph combination. Numbers inside the bars correspond to sample size (e.g., 10 out of the 22 juveniles fathered by orange homomorphic pairs entered inside the shelter). (B) Positive relationship between the probability that a juvenile entered the shelter and juvenile mass. Black line represents individual probabilities predicted by the logistic mixed model. Blue line and shaded area correspond to the smoothed tendency line and its 95% CI. (C) Barplot showing the proportion of juveniles that captured and ate the novel prey item (i.e., a small T. molitor larva), grouped by parental morph combination. Numbers inside the bars correspond to sample size. (D) Scatterplot showing the negative relationship between juvenile mass (g) and the time elapsed from the introduction of the novel prey item until the juvenile bit the larva (i.e., latency in seconds).

Survival, inter-annual growth, and color expression

Short-term survival (3–4 weeks after hatching) averaged 87.7% (372 lizards out of 424) and was positively related with juvenile mass (GLMM binomial: χ = 9.39, P = 0.002, pseudo-r = 0.09), but not with parental morph (P > 0.1; Supplementary Table S8; Figure 5). Inter-annual survival (1 year after hatching) averaged 13.4% (50 lizards out of 372) and was significantly related with juvenile mass and male color morph (GLMM binomial: juvenile mass: χ = 10.43, P = 0.001; male morph: χ = 8.57, P = 0.014; Supplementary Table S8; Figure 5). Namely, offspring sired by yellow males had 1.44 (OR [95% CI] = 0.49, 3.17) higher odds of staying alive until July 2020 than juveniles fathered by orange morph males (Figure 5). Lastly, juvenile lizards gained an average of 2.19 ± 0.24 g in 1 year, but growth was unrelated to any of the examined predictors (Supplementary Table S8). Observed morph frequencies in subadult lizards were largely congruent with the frequencies expected for each morph combination (assuming 2 loci each with 2 alleles and a recessive, homozygous base for orange and yellow expression; P = 0.297; Supplementary Table S9; Figure 6).

Figure 5.

Predicted probability that a juvenile survived according to parental morph combination (A, C) and mass at eclosion (B, D), measured at 2 different points in time: 3–4 weeks after hatching (September 2019: A, B), and 1 year later (July 2020: C, D). Circles mark the position of the mean and bars extend to the 95% CI. Numbers on top of the error bars correspond to the surviving number of juveniles in each morph combination category. Black line represents individual probabilities predicted by the logistic mixed model. Blue line and shaded area correspond to the smoothed tendency line and it is 95% CI. There is a positive relationship between survival and growing juvenile mass. In addition, juveniles sired by orange male lizards showed significantly lower levels of inter-annual survival than juveniles sired by yellow morph lizards (although mortality and dispersion could have been conflated in our experimental tubs).

Figure 6.

(A) Morph assignment for 67 subadult P. muralis lizards 1 year after hatching. Numbers correspond to sample size per morph combination. (B) Barplot comparing the observed morph frequencies (O) per parental morph combination with the frequencies expected (E) in a panmictic population if the expression of orange and yellow (instead of white) ventral coloration required recessive homozygosity at 2 separate loci, each with 2 alleles.

Discussion

In this study, we examined the effect of parental morph combination on fertilization success, embryonic viability, newborn mass, and juvenile behavior (i.e., antipredator response and foraging ability) in common wall lizards (P. muralis), a species where ventral color morphs are often thought to reflect alternative life-history trade-offs (Calsbeek et al. 2010; Sacchi et al. 2017a, 2017b; Mangiacotti et al. 2019). By keeping the juveniles in outdoor tubs for ∼1 year, we could also evaluate the effect of parental morph combination on inter-annual growth and survival. Overall, our results do not support the existence of strong advantages to either morph-assortative or disassortative mating. In fact, although we did find a higher fertilization success of females when paired assortatively, parental morph combination was mostly irrelevant to the different stages of development we examined (i.e., pre-hatching, post-hatching, survival to adulthood).

Both pre-copulatory (e.g., inability to copulate) and post-copulatory (e.g., sperm rejection) barriers to reproduction may result in biased fertilization. However, no sex/morph combination was underrepresented in our recording of inter-sexual interactions, which argues against pre-copulatory mechanisms as an explanation for our results. Adaptive sperm-sorting strategies by which females (typically nonmate-limited in the wild) would eschew the sperm of less-preferred male morphs have been invoked to explain morph-biased patterns of fertilization in U. stansburiana, both in the wild and in a mating experiment where pairing females with a single male resulted in infertility (Calsbeek and Sinervo 2002; Lancaster et al. 2014). For 2 reasons we think that our results in P. muralis are unlikely to be explained by sperm-sorting strategies favoring assortative mating. The difference we found in fertilization success was mainly driven by the high prevalence of infertile clutches in yellow females mated to orange males (5 out of a total of 9 infertile clutches), in contrast with the almost perfect fertilization success of white homomorphic pairs. However, the converse crossing (i.e., orange female with yellow male) showed no abnormal levels of fertilization success, as expected if limited gene flow between these 2 morphs hampered hybridization (Gray and McKinnon 2007; Pryke and Griffith 2009; Lancaster et al. 2014). In addition, because the offspring of white homomorphic pairings ought to be the most genotypically diverse (Tables 1 and 2), finding the highest levels of fertilization success in this crossing (instead of in orange or yellow homomorphic crossings) cannot be interpreted as evidence of a strong homozygote advantage.

In agreement with Andrade et al. (2019), our results argue against the hypothesis that P. muralis morphs may reflect complex suites of co-adapted genes under linkage disequilibria, although do not preclude the possibility that correlational selection mediated by pleiotropy (which also results in the co-segregation of traits but does not create genetic divergence and incompatibilities) may underlie color polymorphism in P. muralis. This is because linkage disequilibria are predicted to erect mainly post-zygotic barriers between the morphs, and hence significantly decrease embryonic viability (which was not observed in our experiment) rather than fertilization success (Sinervo and Svensson 2002; Gray and McKinnon 2007; Mckinnon and Pierotti 2010). For instance, conducting controlled crosses in U. stansburiana (in which color morphs reflect multi-trait reproductive strategies in both males and females), Lancaster et al. (2014) reported that disassortative pairing entails lower probabilities of hatching in fertilized eggs, which they interpret as a form of outbreeding depression causing increased embryonic mortality. Likewise, significantly lower hatchability was documented in the fertilized eggs of disassortative pairings between 2 of the color morphs exhibited by the Gouldian finch (Pryke and Griffith 2009). In both of these species, the combined effect of assortative mating preferences and linkage disequilibria between co-adapted gene complexes within each morph (i.e., correlational selection) seems to have introduced significant genetic divergence and reproductive isolation between morphs (Bleay and Sinervo 2007; Pryke and Griffith 2007). In contrast, our results suggest that despite the high prevalence of morph-assortative pairings observed in wild populations of P. muralis (Pérez i de Lanuza et al. 2013, 2016), inter-morph breeding remains entirely viable. This agrees with the recent finding of only minor genetic divergence between color morphs in polymorphic P. muralis populations from eastern Pyrennees, which suggests inter-morph gene flow must be sufficiently high to prevent the build-up of strong genetic differentiation between the morphs (Andrade et al. 2019). In short, although we found lower fertilization success in Y♀–O♂ heteromorphic crossings (which remains to be confirmed using a larger sample), overall our results suggest that color polymorphism is not promoting post-zygotic barriers in P. muralis from eastern Pyrennees, whereas arguing against a fitness advantage of homomorphic matings as an explanation for the high prevalence of morph-assortative pairings in natural populations (Pérez i de Lanuza et al. 2013, 2016). In fact, the genetic basis of P. muralis color polymorphism (2 unlinked and bi-allelic loci)—which prevents assortative crosses between pure morphs to yield genotypically homogeneous offspring—may partially explain why assortative mating did not generate reproductive barriers, and why morph-biased mating preferences are unlikely to evolve in this species.

We found that juvenile mass at birth (but not parental morph combination) influences juvenile behavior, with heavier juveniles showing the most adaptive response more often and with shorter latencies in lab-staged antipredator and foraging ability tests. Offspring mass was also found to correlate positively with inter-annual survival. This is in agreement with evidence from several studies in lizards, suggesting that large body size increases survival early in life (Sinervo et al. 1992; Clobert et al. 2000; Svensson and Sinervo 2000; Warner and Andrews 2002; Uller and Olsson 2010). The selective advantage of larger juvenile lizards (i.e., bigger-is-better hypothesis) is often assumed to result from their better performance at foraging and predator avoidance (Ferguson and Fox 1984; Stamps and Krishnan 1997; Civantos and Forsman 2000; Miles 2004; Massot and Aragón 2013). Here we present experimental evidence in support of this assumption. The observable consequences of genetic incompatibilities are often larger in traits involving the expression of many genes and closely related to fitness (DeRose and Roff 1999). This is because dominance variance is expected to be greater in such traits, as additive genetic variance has been decreased by directional selection (Mousseau and Roff 1987; Roff 1997; Ala-Honkola et al. 2009). For instance, research on poeciliid fish shows that genetic stress negatively affects male reproductive behavior and the area of a carotenoid ornament, but not other morphological characters more indirectly related to fitness (Sheridan and Pomiankowski 1997; Ala-Honkola et al. 2009). The expression of adaptive antipredator and foraging behavior depends on the interaction between many genes with early environmental conditions and is reasonably believed to correlate directly with survival (and hence, fitness) in juvenile lizards (Civantos and Forsman 2000; Clobert et al. 2000; Miles 2004). For this reason, we deem our finding of more adaptive behavior in heavier juveniles (irrespective of their parental morph combination) a compelling argument against the existence of morph incompatibilities in P. muralis.

These results are in stark contrast with those of a previous mating experiment in P. muralis from northern Italy where parental morph combination was found to significantly impact reproductive output, and female color morphs were found to show alternative breeding investment strategies (Galeotti et al. 2013). In that study, orange and yellow females showed an overall higher fertilization success and embryonic viability when paired assortatively, whereas white morph females performed better when paired with orange males. In addition, yellow morph females tended to produce many light offspring, and white females produced few heavy juveniles. Similar morph-specific strategies have also been found in females of the common lizard (Zootoca vivipara), a close relative of the genus Podarcis also showing a complex ventral color polymorphism (Vercken et al. 2007; Vercken and Clobert 2008; but see Cote et al. 2008). Whereas our results agree with these studies on the negative relationship between juvenile mass and clutch size (a general life-history trade-off previously demonstrated in lizards; Galán 1997; Ji and Braña 2000; Galeotti et al. 2013; Mitchell et al. 2018), we failed to replicate their results concerning color polymorphism (Figure 3 andSupplementary Figure S7). It is possible that these contrasting results are biologically meaningful. Most of the evidence suggesting physiological or behavioral morph differences in P. muralis comes from studies on the Southern Alps lineage, which is only distantly related to the Western European lineage of the Pyrenees (Giovannotti et al. 2010; Schulte et al. 2012; Gassert et al. 2013; Salvi et al. 2013; Yang et al. 2018). Recent studies in polymorphic organisms have reported correlations between local environmental factors and morph population frequencies (McLean and Stuart-Fox 2014; McLean et al. 2014a, 2014b, 2015), incorporating climate and other abiotic components of natural selection into the study of the adaptive value of color polymorphism (Roulin 2004; Svensson 2017; Svensson et al. 2020). In P. muralis, both climatic and demographic factors (i.e., sex ratio) may affect morph composition, suggesting the need to investigate the environmental dependence of morph fitness in future research (Pérez I de Lanuza et al. 2017, 2018b; Pérez i de Lanuza and Carretero 2018). In fact, we deem the finding that P. muralis color morphs do not show strong genetic incompatibilities an important motivation for new studies tackling the environmental dependence of morph fitness in natural populations (see also Abalos et al. 2020).

Although we did not find evidence of a strong morph combination effect on reproductive output, we found lower inter-annual survival in the offspring of orange morph males. Calsbeek et al. 2010 reported lower levels of inter-annual survival for free-roaming orange morph lizards inhabiting a stone church and a cemetery at <1 km from where our experimental tubs were positioned. Having also found larger body sizes and higher levels of parasite infection in orange morph lizards, Calsbeek et al. (2010) argued that selection may have favored a higher investment in dominance and reproduction in the orange morph, at the cost of a weaker immune function and shorter lifespan than the other morphs. Alternative life-history strategies subject to environmentally dependent fitness have been described in other polymorphic systems such as the tawny owl (Strix aluco), where gray-colored individuals live longer and produce more offspring than the brown morph under cold environmental conditions, but the situation is reversed under warmer conditions (Brommer et al. 2005). Unfortunately, both ours and Calsbeek et al. (2010) results may conflate mortality with dispersal, as no experimental constrain to dispersal was established when monitoring the free-roaming population, and we cannot completely discard that some of our juveniles may have escaped the outdoor tubs through the overgrown vegetation during theCOVID-19 lockdown. Although in both cases the study site was surrounded by unsuitable lizard habitat and researchers failed to find marked lizards in the adjacent area, these results do not allow us to disentangle the relative effects of escape and mortality on the observed pattern of inter-annual re-sighting. Still, the fact that 2 independent studies found a similar bias in inter-annual re-sighting concerning the orange morph suggests the potential value of investigating differences in life-history traits such as life expectancy, longevity, or dispersal among P. muralis color morphs in future research.

Although based on a limited sample size, the observed morph frequencies among the subadult lizards allow us to reject certain mechanisms of inheritance that control morph expression in other polymorphic lizards. For instance, a mechanism based on a single locus with 3 co-dominant alleles (e.g., U. stanburiana; Sinervo et al. 2010; Olsson et al. 2007) entails that yellow–orange offspring can only result from Y–O matings and that W–W crossings should always produce white morph lizards, 2 conditions that are not met in our sample (Figure 6). Likewise, a mechanism based on 2 separate loci (each with 2 alleles) where presence of the dominant O and Y allele results in color expression (Ctenophorus decresii; Rankin et al. 2016) would restrict the outcome of certain inter-morph matings in a way incompatible with our observations (e.g., W–Y crossings could never result in orange offspring). An alternative model of inheritance for C. decresii was proposed in Stuart-Fox et al. (2020) where 1 autosomal, bi-allelic locus determines whether individuals express color (orange and/or yellow) instead of grey coloration. According to our results, this model also seems implausible for P. muralis because it prevents the production of offspring expressing the dominant phenotype in matings between phenotypically recessive individuals (e.g., white morph offspring in matings between colored adults, if the putative white allele was dominant over the color allele). In contrast, observed morph frequencies were congruent with those expected under the assumption that yo colorations are controlled by 2 separate loci (each with 2 alleles) where recessive homozygosity of the y and o alleles results in color expression (Andrade et al. 2019). In fact, orange and yellow subadults were primarily produced in O–O and Y–Y homomorphic crossings (respectively), whereas heteromophic and W–W homomorphic crossings resulted in a wider range of differently colored offspring (Figure 6). Several corollaries follow these results. White coloration showing genetic dominance over orange and yellow offers a simple explanation for the overall higher prevalence of this morph in natural P. muralis populations (Pérez i de Lanuza et al. 2017, 2018b; Andrade et al. 2019). Dominance relationships may also have important implications on polymorphism maintenance, as recessivity lowers exposure to selection hence allowing for the expression of alternative (i.e., specialist) phenotypes only under particular conditions. In turn, this highlights the need to genotypically characterize P. muralis lizards when studying morph frequencies and selection in natural populations. Lastly, these results evidence taxonomic diversity in the genetic architecture of lizard color polymorphism. Assessing the evolutionary causes and consequences of such diversity may be a fruitful avenue of research for future comparative studies (Stuart-Fox et al. 2020).

Altogether, our results suggest that, at least in the examined populations (i.e., eastern Pyrenees), P. muralis color morphs do not reflect alternative breeding investment strategies in females and parental morph combination does not influence offspring viability, behavior, and survival. Although we found lower fertilization success in Y–O crossings, we list several reasons why this result should be confirmed in new studies before inferring its possible consequences for the long-term maintenance of color polymorphism. Additional mating experiments crossing adults of known genotype could more precisely evaluate the potential existence of heterosis and/or pleiotropic effects in the genes governing P. muralis color polymorphism, as well as explore the possible environmental dependence of morph combination effects on offspring fitness. Overall, our results suggest that inter-morph breeding remains entirely viable in P. muralis while maintaining the uncertainty about which mechanism is preventing any of the color morphs from going extinct.

Click here for additional data file.