Frequency-dependent Batesian mimicry maintains colour polymorphism in a sea snake population

Evolutionary theory suggests that polymorphic traits can be maintained within a single population only under specific conditions, such as negative frequency-dependent selection or heterozygote advantage. Non-venomous turtle-headed sea snakes (Emydocephalus annulatus) living in shallow bays near Noumea in New Caledonia exhibit three colour morphs: black, black-and-white banded, and an intermediate (grey-banded) morph that darkens with age. We recorded morph frequencies during 18 consecutive years of surveys, and found that the numbers of recruits (neonates plus immigrants) belonging to each morph increased in years when that morph was unusually rare in the population, and decreased when that morph was unusually common. Thus, morph frequencies are maintained by negative frequency-dependent selection. We interpret the situation as Batesian mimicry of highly venomous sea snakes (Aipysurus, Hydrophis, Laticauda) that occur in the same bays, and range in colour from black-and-white banded to grey-banded. Consistent with the idea that mimicry may protect snakes from attack by large fish and sea eagles, behavioural studies have shown that smaller fish species in these bays flee from banded snakes but attack black individuals. As predicted by theory, mimetic (banded) morphs are less common than the cryptically-coloured melanic morph.

Polymorphic traits present a challenge to simplistic models of natural selection. Two morphs are unlikely to confer exactly the same fitness, and hence we would expect natural selection against the less-fit morph to eliminate the polymorphism 1 . An extensive literature examines why this prediction is not always met and shows that multiple morphs can be maintained by processes such as heterozygote advantage, negative frequency-dependent selection, differential habitat use based on genotype or phenotype, correlational selection, or differential selection on males and females [1][2][3][4] . Some of the most dramatic cases of polymorphism involve multiple colour-morphs that function as Batesian mimics of sympatric toxic species. In particular, some butterfly populations contain individuals with distinctive colour patterns that strongly resemble those of multiple co-occurring "models" that are avoided by visually-hunting predators 5,6 .
Although Batesian mimicry likely is widespread, most examples of this phenomenon in reptiles involve species-level resemblances between harmless mimics and toxic models, as in the coral-snake system in South America 7,8 and the viper-colubrid system in Europe 9,10 . In contrast, studies of intraspecific colour polymorphism in reptiles mostly have focused on sex-specific morphs in lizards, maintained by asymmetric advantages in intrasexual combat among males (e.g., "rock-paper-scissors" dynamics 11 ) or by sex-based differences in microhabitat use that affect the link between colour patterns and fitness 12 . The best empirical evidence for balancing selection on colour morphs in natural populations comes from longterm studies on invertebrates 13 .
Although sea snakes are diverse and abundant in tropical oceans, logistical obstacles have discouraged detailed research on these animals 14 . In the present study, we examine a species that is highly variable in colour both among and within populations: the turtle-headed sea snake Emydocephalus annulatus 15 . In the populations that we study in the IndoPacific archipelago of New Caledonia, most individuals are jet-black but some are brightly banded in black-and-white whereas others are intermediate [15][16][17] . To clarify the mechanisms that maintain this polymorphism, we analyse data from an 18-year mark-recapture program to explore the distribution of colour morphs with respect to sex and body size, and the degree to which morph frequencies remain relatively constant through time and space. If recruitment is negatively frequency-dependent, we expect that an increased proportion www.nature.com/scientificreports/ of snakes of a particular morph will be followed by a reduction in relative rate of recruitment of that morph, thus keeping morph frequencies stable.

Methods
Study site and species. We work in two small shallow (< 3 m deep) bays beside the city of Noumea, in New Caledonia (22°16′ S, 166°26′ E). The substrate is dominated by live coral, coral rubble, boulders and sand, with sea snakes primarily using coral substrates 16,18 . Water temperatures average around 27 °C in midsummer (February) and 23 °C in midwinter (August: https:// www. seate mpera ture. org/ austr alia-pacifi c/ new-caled onia/ noumea. htm). Turtle-headed sea snakes (Emydocephalus annulatus) are medium-sized (to 800 mm snout-vent length [SVL]) hydrophiine elapids ("true" sea snakes 19 ) that are entirely aquatic. Adult females produce a litter of about two large offspring (300 mm SVL) on an approximately biennial cycle 20 , usually in May-June. Courtship and mating occur primarily in the cooler months (June-August), with ovulation in spring (September-October) followed by a long gestation. Both sexes mature at about 500 mm SVL, and at two to three years of age 20 . The unusual diet of this species (eggs of demersal-spawning fishes 21 ) has driven major changes in the trophic apparatus. Venom plays no role in subduing prey, resulting in a considerable reduction in size of the fangs and venom glands, potency of the venom, and ability to gape the mouth open widely. Thus, unlike almost all other sea snakes, this species is essentially non-venomous 19 .
Snakes in our study populations occur in three colour morphs: jet-black (melanic), brightly black-and-white banded, and intermediate (grey-banded: see Fig. 1). All snakes were classified to colour morph by the same observer (RS), except that CG scored snakes in 2021. To quantify inter-rater reliability, we asked 10 volunteers (colleagues and students) to classify 30 photographs of snakes into morphs. Participants were told that snakes occur in three morphs (black, black-and-white-banded, and intermediate) but were not told how many snakes fell into each category. We provided 10 randomly-selected photographs of each morph. www.nature.com/scientificreports/ In previous work, we have shown that darker colours increase algal fouling (because algae preferentially settle onto dark backgrounds), increasing hydrodynamic resistance and thus, reducing swimming speeds 23 . However, black skin also confers an advantage: it facilitates elimination of trace-element pollutants from the body via sloughing, because melanin binds to those trace elements, and because of the higher sloughing frequency induced by higher rates of algal fouling 15 . Detailed surveys have not revealed any other ecological traits associated with colour morph: for example, rates of growth and survival are similar between banded and melanic individuals 23 , as is microhabitat use 16 . The mechanistic basis of colour polymorphism in E. annulatus is not known, but studies on other snake species have revealed a genetic basis to intraspecific colour variation, with simple Mendelian inheritance 22 . Surveys. Every January from 2004 to 2021, we snorkelled (in groups of 3 to 14 people) through three coralsubstrate sites within two adjacent bays: the Baie des Citrons and Anse Vata (see 19 for site details). All snakes sighted were captured by hand, and taken to a nearby laboratory for measuring, weighing and microchipping prior to release at the site of capture. Analyses using the program MARK suggest that we captured 24 to 80% of resident snakes at each site each year (mean = 54%). These slow-moving snakes are relatively easy to see, especially when they are active in shallow water, so a snake's colour is unlikely to have affected its probability of capture.
During processing, we scored the colour of each snake into three categories (see above and Fig. 1). Pregnancy and litter sizes were assessed by palpation of the large developing embryos. We handled and scored each individual only once in each year, but many snakes were captured in multiple years. The current analysis is based on 2,673 captures of 1,427 snakes.

Statistical analysis.
We used generalized linear models implemented in JMP Pro 15.0 and SAS 9.4 to analyse data on colour morph frequencies and their correlates. We used nominal logistic regression (with a glogit link function) to assess the effects of sex body size, site, and year on the frequencies of each colour morph, using only the first capture per individual to avoid pseudoreplication. To examine variation through space and time in the proportion of snakes of each colour morph, we used models with number of individuals of each colour morph as the dependent variable, and independent variables of morph and site (to examine spatial variation) or morph and year (to examine temporal variation). Interaction terms were included and a negative binomial distribution was used to fit the model. We were unable to include both Site and Year in the same models, because not all sites were studied over the same period.
To examine survival rates for snakes of each morph, we pooled the data from Baie des Citrons south and Anse Vata (the two sites for which 17-year mark-recapture data were available) for Cormack-Jolly-Seber analysis implemented in MARK software 24 . We used colour morph as a grouping variable with three levels and ran 16 models where survival and recapture parameters were constant, varying over time, varying among colour morphs or varying interactively with time and morphs. We then compared the fit of the 16 models by ranking their corrected Akaike (AICc) values. We considered models whose AICc values differed by < 2.0 to have equivalent support.
To evaluate ontogenetic changes in colour pattern, we examined mark-recapture records to quantify the numbers of cases in which a snake remained the same colour as in earlier years, versus was classified as a different morph.
To assess frequency-dependence in recruitment, we compared the proportion of snakes of each colour morph within a given population to the proportion of snakes of that colour in the following year. We conducted this analysis using a model with study site as a factor, % of snakes in a specific colour morph as a covariate (plus the interaction between these factors), and change in % of snakes of that colour between that year and the next as the dependent variable. We ran separate analyses for each of the three colour morphs, using a Gaussian distribution.
Although proportions facilitate interpretation of patterns in the data, ratio measures can introduce statistical problems 25 . We thus conducted an alternative analysis for density-dependence that used data on actual counts: the changes in numbers of snakes of each colour morph were included as three independent variables, and the number of snakes of each colour morph in the preceding year as the dependent variable. We used a generalized mixed model for these analyses, incorporating Site as a random effect and using a Gaussian distribution.

Results
How discrete are the three colour morphs? Analysis of polymorphic traits requires unambiguous definitions of each state. Our survey of inter-rater reliability gave scores from 10 observers for each of 30 snakes. Observers agreed with each other, and with our a priori classification, for 296 of those 300 evaluations (98.7%). One person scored an intermediate-phase snake as black-and-white banded, and three people scored an intermediate-phase snake as black.
In the case of E. annulatus, then, the melanic and black-and-white banded states are very different (never confused with each other in the above survey), but the intermediate (grey-banded) phase is more problematic (mis-classified in 4 of 100 cases in the survey above). The three phases are very different in juvenile (youngof-the-year) individuals, but increasing melanism of intermediate-phase snakes renders them increasingly less distinguishable from melanic snakes in later age classes. That transition is evident from SVL distributions of the three morphs, with the intermediate morph most common in smaller snakes (Fig. 2)   Negative frequency-dependence. Years with a higher-than-usual proportion of snakes of a particular colour morph were generally followed by lower-than-usual rates of recruitment of snakes of that morph. For all three morphs, the change in proportions was negatively related to the prior proportion of snakes exhibiting that   Fig. 4). We note, however, that these three tests are not statistically independent because the frequency of one morph influences the relative frequency of the others.
In an alternative analysis using data on counts rather than proportions and using Site as a random effect, the number of snakes of a given colour morph was negatively related to the change in numbers of that morph during the following year, but not to changes in numbers of other morphs. Thus, the number of banded snakes was negatively associated with the change in number of banded individuals (  www.nature.com/scientificreports/ but not with changes in the numbers of banded or melanic snakes (F 1,35.08 = 0.18, P = 0.67 and F 1,35.21 = 0.77, P = 0.39 respectively).

Discussion
Independent of their sex, about 80% of the turtle-headed sea snakes in our study populations are melanic (black), but around 10% are brightly banded in black-and white. The remaining 10% exhibit medium (grey) bands when young, but darken to resemble the melanic morph as they grow older. Those proportions are similar among populations (despite very little gene flow between sites 26 ) and were consistent through time across the entire study area. However, relative proportions of individuals of the three colour morphs shifted through time within each site, with negative frequency-dependence: as soon as one morph became more common than usual, recruitment of that morph declined and hence, maintained the status quo (Fig. 4). What biological processes maintain numerical stability of the three colour morphs? Many of the processes identified in mathematical models, and documented in previous studies, do not apply to our study system. For example, there is no suggestion of differential selection on males and females (sex is not linked to colour), nor of ecological differences between morphs (detailed surveys have revealed no such divergences in habitat or diet 15 ). Asymmetric advantages in interference competition or sexual rivalry between individuals of different morphs (as in some lizards 11 ) cannot apply, because males do not engage in agonistic interactions, and do not discriminate among female colour morphs when selecting targets for courtship 27 . Heterozygote advantage also seems unlikely, because survival rates of intermediate-phase snakes (which may or may not be heterozygotes) were no higher than were those of melanic or black-and-white banded individuals.
We suggest that the stable polymorphism seen in our harmless sea snakes is due to Batesian mimicry of sympatric deadly snakes, reducing rates of predation from visually-hunting fishes and raptorial birds. Specifically, the melanic morph may be cryptic, as well as conferring advantages in eliminating pollutants from the body 15 . In contrast, resemblance in colour pattern between venomous snakes and the banded colour phases of E. annulatus may protect those putative mimics. In support of that hypothesis, we note: (1) annual recruitment to our study populations depends primarily on rates of survival of neonates (changes in total numbers from one year to the next are largely due to changes in the numbers of yearling snakes that survive 28 ); (2) annual recruitment is higher after years when windy conditions during winter curtail predation 28 ; (3) our study sites contain seven species of highly venomous sea snakes 13 Fig. 1 for photographs of some of these presumed mimic-model similarities). Because juvenile sea snakes tend to be more brightly coloured than conspecific adults, the strongly-banded models are similar in size to the putative mimic species, E. annulatus; (4) predation likely is a significant source of mortality, based on reports of sea snakes in the diets of large fishes and birds 29-34 ; (5) many E. annulatus in our study populations have injuries to their tails, presumably the result of foraging strikes by predators (Shine et al., unpubl. data); (6) small nest-defending fishes in our study sites attack melanic snakes but ignore or avoid banded and striped snakes 35 .
Although the evidence for the Batesian mimicry hypothesis remains circumstantial, it fits well with our experience during fieldwork. Not infrequently, we have mistaken one of the banded morphs of E. annulatus for another snake species, and vice versa. The same mistake is evident in published work: for example, Fig. 1 of Norman et al. 's 36 paper on mimicry ironically includes a photograph of a banded E. annulatus misidentified as a sea krait (Laticauda), the taxon both they and we have identified as a likely model for mimetic morphs. The range of colour patterns in putative models is increased by the occurrence of snakes covered with algae, rendering bands less distinct; and limited penetration of light into deeper water reduces chromatic cues 18 . Especially in deeper water, or when cloudy water and strong currents limit visibility, a harmless snake may be avoided even if it partially mimics a dangerous species.
On a broader spatial scale, the shallow inshore bays where we have studied E. annulatus are among a minority of places where melanism is common in this species, perhaps reflecting advantages of darker colouration in eliminating trace-element pollutants 15 . In most other areas, most individuals are banded or blotched in colour 15 , consistent with the likelihood of higher predation rates on snakes, as well as lower pollutant loads, in these less disturbed habitats.
Sea snakes may be important components of model-mimic complexes in other contexts as well. For example, some octopuses dynamically adopt a banded appearance that resembles a venomous sea snake 36 . Famously, snake-eels (e.g., Leiuranus semicinctus, Myrichthys colubrina) have evolved a similar banded-snake-mimicking pattern 37,38 . Dudgeon & White 39 likewise interpret the morphology and brightly banded colours of juvenile zebra sharks (Stegostoma fasciatum) as mimicry of banded sea snakes. These examples suggest multiple independent evolutionary origins of banded sea-snake colouration in marine organisms, presumably because of the fitness benefits available through mimicry. In at least some cases, predators strongly avoid aposematically coloured seasnake models 40 . Additionally, sea snakes may mimic other sea snakes. Most species of sea snakes possess potent venom and a high proportion of them are banded or blotched in colour 19 . That consistency in colour might be due to multiple processes including phylogenetic conservatism, convergent adaptation to marine habitats, and/ www.nature.com/scientificreports/ or Mullerian mimicry. That is, deadly sea snakes may benefit from resembling other deadly species, facilitating avoidance by predators (as suggested for other types of venomous snakes 41,42 ). If so, the colour polymorphism we have documented in E. annulatus may be part of a much broader and more complex model-mimic system. Future work could test the interpretation of Batesian mimicry by looking into predator-prey interactions. For example, we predict that larger fishes (like the smaller ones we have already studied 35 ) would direct fewer attacks towards banded morphs than towards melanic morphs, at least in a situation where banded morphs are relatively infrequent. Artificial models of snakes, differing in colour pattern, might allow an empirical test of that prediction. More generally, we need to overcome the logistical constraints of studying predator-prey interactions underwater, if we are to more robustly test the idea that colour patterns in sea snakes are driven by Batesian mimicry.

Data availability
Data has been deposited in the Dryad digital repository and is available at https:// doi. org/ 10. 5061/ dryad. x3ffb g7mq.