Negative assortative mating and maintenance of shell colour polymorphism in Littorina (Neritrema) species

Colour polymorphism is a widespread phenomenon in natural populations of several species. In particular, it is especially common on marine gastropod species from the genus Littorina. Recently, it has been argued that intrapopulation shell colour polymorphism in Littorina fabalis could be caused by negative frequency-dependent sexual selection via a mechanism of mate choice (indirectly estimated through negative assortative mating). Here we try to determine the existence of negative assortative mating in three species of the subgenus Neritrema (L. fabalis, L. obtusata, L. saxatilis) that share a similar shell colour polymorphism, to ascertain if this mechanism could represent an ancestral character in this subgenus that could be contributing to the maintenance of the colour polymorphism observed in each species. We collected or reanalysed from previous studies a sample of mating pairs of these three species from seven locations from NW Spain and NE Russia and estimated assortative mating using the IPSI index. Our results suggest that all species and populations show a systematic tendency towards negative assortative mating when shell colour is grouped in the broad categories ‘light’ and ‘dark’. Moreover, a more detailed analysis of each colour separately suggests that shell colour may not be the main target of assortative mating, but perhaps a physically-linked trait to the real target of selection. This hypothesis opens interesting new lines of research in Littorina snails.


Introduction
Colour is a key trait in evolution, usually contributing to individual survival (Olsson et al. 2013;Rolán-Alvarez et al. 2015a) and reproduction (Cordero 1989;Hugall and Stuart-Fox 2012), and therefore potentially influencing fitness. In addition, colour is a relatively easy trait to detect and measure, and so it has been thoroughly studied in several taxa (McKinnon and Pierotti 2010;Hugall and Stuart-Fox 2012;Svensson 2017;Cuthill et al. 2017). In particular, shell colour has been extensively studied in gastropod snails (reviewed in Williams 2017).
Of special interest in studies of colour is colour polymorphism, in which there exist various discrete and genetically determined colour phenotypes within a population (White and Kemp 2016). Colour polymorphisms are expected to be infrequent, given that directional natural selection and genetic drift tend to erode genetic variability over time, therefore reducing observed polymorphisms (Lewontin 1974;Nielsen 2005). However, under certain circumstances selection can maintain a polymorphism in a population, and even despite the tendency of genetic drift to fixation, a polymorphism can be transitorily maintained from a few to several generations depending on effective population size and gene flow (Wright 1977;Clarke 1979;Kimura 1983;Hartl and Clark 2006). This makes it rather difficult to identify the true causal factor responsible for the evolutionary persistence of colour polymorphism within a population. One of the most effective types of selection able to maintain a polymorphism in wild populations is negative frequencydependent selection (Ayala and Campbell 1974;Fitzpatrick et al. 2007;Svensson and Connallon 2019). Negative frequency-dependent selection occurs when the fitness of a particular genotype correlates with its frequency, in a manner that a decrease in frequency increases its fitness relative to alternative genotypes (Ayala and Campbell 1974). Fitness will peak as frequency falls below a certain threshold, hence protecting the polymorphism from being lost in the population.
The genus Littorina has been shown to be a particularly useful model organism in the field of evolutionary ecology (Rolán-Alvarez et al. 2015a). Besides speciation being the main focus of the evolutionary research carried out in these species (e.g. Morales et al. 2019;Costa et al. 2020;Galindo et al. 2020), the maintenance of intrapopulation colour polymorphism has also been studied (Rolán-Alvarez and Ekendahl 1996;Rolán-Alvarez et al. 2012, 2015bJohannesson and Butlin 2017;reviewed in Rolán-Alvarez et al. 2015a). Recently, it has been shown that shell colour polymorphism in populations of the intertidal gastropod Littorina fabalis (Turton 1825) can be maintained by negative frequency-dependent sexual selection via mate choice (Estévez et al. 2020). Another two species of the subgenus Neritrema (genus Littorina), L. saxatilis (Olivi 1792) and L. obtusata (Linnaeus 1758) have also shown colour polymorphism within certain populations (Reid 1996;Rolán-Alvarez et al. 2015a). These three species are direct developers (Reid 1996), lacking a planktotrophic larval stage, and therefore showing a potential for fast adaptation to local conditions (reviewed in Johannesson 2003;Galindo and Grahame 2014;Rolán-Alvarez et al. 2015a). Another important characteristic in these (dioecious) Littorina species is the possibility to observe mating pairs, which allows to readily estimate assortative mating in the wild (reviewed in Rolán-Alvarez et al. 2015a). In this context, assortative mating is understood as the phenotypic correlation between mates (Jiang et al. 2013;Janicke et al. 2019), which can be measured by isolation indexes in the case of qualitative traits, like shell colour (reviewed Gilbert and Starmer 1985;Rolán-Alvarez and Caballero 2000). Assortative mating for colour and its evolutionary consequences have been widely studied in several species over the last few years (Calderon et al. 2010;Pryke 2010;De Lanuza et al. 2013;Gade et al. 2016;Hedrick et al. 2016;Yang et al. 2016;Fargevieille et al. 2017). Besides, frequent biases in the estimation of assortative mating as well as its potential solutions have been recently identified (Rolán-Alvarez et al. 2015b).
Although there are alternative explanations (see Taborsky et al. 2009;Kopp et al. 2018), assortative mating in simple organisms such as marine gastropods is presumably caused by mate choice (Ng et al. 2019;Estévez et al. 2020). On the other hand, negative assortative mating can be directly related to negative-frequency dependent selection, as it has been shown both theoretically (Pusey and Wolf 1996;Hedrick et al. 2016Hedrick et al. , 2018 and empirically in several species (Takashi and Hori 2008;Field and Barrett 2012;Holman et al. 2013;Hedrick et al. 2016Hedrick et al. , 2018Estévez et al. 2020). This, in turn, makes it plausible to assume that negative assortative mating is a relatively good proxy for both mate choice and negative-frequency dependent selection. Here we study assortative mating for shell colour in new populations of L. fabalis and two other related species with a similar colour polymorphism, L. obtusata and L. saxatilis. The aim of this study is to determine the existence of a negative assortative mating pattern in these two new species, as well as to corroborate its presence in the new populations of L. fabalis. Given that previous results strongly suggest a causal link between negative assortative mating and both mate choice and negative frequencydependent selection, we suggest that if a negative assortative mating pattern is found in these two new species, an ancestral plesiomorphic behavioural trait for mate choice may exist within the subgenus Neritrema, which encompasses L. fabalis, L. obtusata and L. saxatilis. This could pave the way to design future studies to directly test the existence of mate choice and negative frequency-dependent selection in the Littorina clade.

Study species and trait
We selected three colour polymorphic species from the genus Littorina (L. obtusata, L. fabalis and L. saxatilis), subgenus Neritrema. This subgenus is the only one with non-planktotrophic development within the genus, and it is formed by ten species with a Northern distribution (Reid 1996;Reid et al. 2012). Within this subgenus, we find the group of flat periwinkles composed by Littorina fabalis and L. obtusata, which are two sibling species (split 0.85 Mya, Sotelo et al. 2020) that live on intertidal brown macroalgae (e. g. Fucus sp.). L. fabalis feeds on diatoms and other microalgae that grow on the seaweed surface, while L. obtusata feeds directly on the macroalgae (Reid 1996). They also differ in habitat, L. fabalis commonly dwelling lower on the shore while L. obtusata lives in the upper part of the intertidal range. In certain localities their distribution partially overlaps, whereas in others only one species is found depending on the environmental conditions of the shore. In NW Spain they can be observed in true sympatry living on the same Fucus vesiculosus plants ( Rolán-Alvarez and Ekendahl 1996). On the other hand, they are both oviparous, laying egg masses on the algae from which juveniles hatch as miniature snails (Reid 1996).
The rough periwinkle, L. saxatilis (split with L. fabalis 2.5 Mya, Panova et al. 2011), is the only ovoviviparous species in this subgenus (Reid 1996). The embryos are developed in a brood pouch inside the female until eventually hatching as miniature snails. L. saxatilis lives on the rock surfaces across the whole intertidal range, feeding on diatoms, microalgae and detritus (Otero-Schmitt et al. 1997). The three species have very similar distribution across NE Atlantic coasts (Reid 1996) and in many shores they can be found overlapping, although their microhabitats might differ.
Shell colour in these species has been typically studied by visual (qualitative) codification (reviewed in Reid 1996; Rolán-Alvarez et al. 2015a) but at least in one case (L. fabalis from Abelleira) the qualitative classes have been validated by reflectance spectrometry (Rolán-Alvarez et al 2015b). Shell colour variation in the three species has been systematically attributed to one or a few genes (reviewed in Rolán-Alvarez et al. 2015a). Shell colours have not been commonly associated to selective pressures, but there are a few exceptions (reviewed in Reid 1996;Rolán-Alvarez et al. 2015a). For example, unknown predators cause colour-based selective mortality in Littorina sp. from Iceland (Johannesson and Ekendahl 2002), and yellow-painted L. obtusata specimens were cooler and survived better than black-painted ones (Phifer-Rixey et al. 2008). On the other hand, survival differences between colours in L. saxatilis has been associated with physiological differences to resist salinity gradients (Sokolova and Berger 2000).

Experimental methods
To study assortative mating for shell colour in L. fabalis, L. obtusata and L. saxatilis, the experimental design consisted in the capture of mating pairs and the surrounding unmated individuals in natural populations following previous studies (Rolán-Alvarez and Ekendahl 1996;Rolán-Alvarez et al. 2012, 2015bEstévez et al. 2020). These three species show a similar mating strategy, where males usually perform an active search of females, following them by their mucus trails, and typically mate by size similarity (Saltin et al. 2013;Ng et al. 2019). In the field, a mating pair can be detected when the male is located over the female in a certain position on the right-hand side of the shell (with the axis of coiling of the two snails forming a 'V'-like pattern) and the penis is inserted in the mantle cavity of the female (Fig. 1A, B). Associated to every mating pair we collected the four closest unmated individuals of the same species (typically within 5 cm of the mating pair), and these six snails were defined as an individual microarea for further analyses (see Fig. 2 in Estévez et al. 2020). These samplings were done for several days during the summer for each year since snail density and matings increase during this season (Rolán-Alvarez and Ekendahl 1996).
All the snails were taken alive to the laboratory and stored at − 20° C. The mating pairs were dissected to discard malemale matings from the analyses, as they are known to be quite frequent in L. saxatilis (Rolán-Alvarez et al. 1999). Following these dissections, we discarded around 10% for L. saxatilis and 12% for flat periwinkles. Then, the background shell colour of each individual was scored using a colour printed model following Rolán-Alvarez et al. (2015b) and classified into the categories described in Fig. 1. The banding pattern present in some shells was ignored in further analyses, as its frequency differed between species and it  Estévez et al. 2020) these colours were classified into two categories: light (white, yellow, orange) and dark (olive, brown). In addition, in this study we have also performed the analyses based on the individual colours, comparing each colour against the rest, on a binary manner (see Table S1).

Statistical analyses
Assortative mating for qualitative traits can be estimated by several indices (reviewed in Rolán-Alvarez and Caballero 2000), from which the I PSI index is the most robust under different sampling scenarios (Pérez-Figueroa et al. 2005). This index can be estimated in a hypothetical population with two colours (brown, B; yellow, Y) as, Being PSI BB the ratio between the observed number of brown female-brown male mating pairs and the expectation  Table S1), and so only the L. obtusata pie chart is presented. Russia: 1 (WS-1) and 2 (WS-2). Spain: 3 (Portela), 4 (Soutoxuste), 5 (Cangas), 6 (Castelete) and 7 (Abelleira) of this type of mating pair under random mating (using exclusively the mated individuals to estimate random mating; see further details in Carvajal-Rodríguez and Rolán-Alvarez 2006 and references therein). The remaining PSI values are calculated accordingly. The I PSI index ranges from -1 (maximum negative assortative mating) to 1 (maximum positive assortative mating), being 0 complete random mating. This index has been used successfully in previous studies of assortative mating in Littorina (Rolán-Alvarez et al. 2015b; Estévez et al. 2020). The I PSI index was calculated for the trait shell colour, categorized as a binary trait in two different analyses. The first one (pooling analysis), grouped colours into light (white, yellow, orange) and dark (olive, brown) as mentioned in the experimental design and following previous studies. The second analysis (individual analysis) considered every colour separately against the rest of the colours grouped. These two analyses were carried out for each species and population datasets. The I PSI index gives an estimate of assortative mating for qualitative traits similar to Pearson's correlation coefficient for quantitative traits (Rolán-Alvarez et al. 2015b). However, these authors also showed that the existence of trait heterogeneity at a smaller scale than the scale of mate choice would bias assortative mating towards positive values (i. e. the scale-of-choice effect sensu Rolán-Alvarez et al. 2015b). Therefore, to reduce the potential bias towards positive values of assortative mating, we also performed an analysis subdividing each dataset (species/population) into three different homogeneous groups of microareas (mating pair + four unmated individuals) based on the colour frequencies of the microareas (pooling analysis: light vs dark; individual analysis: one colour vs the rest) as follows: low frequency (0-33.3% light or other colour), intermediate frequency (

Results and discussion
Even though directional natural selection and genetic drift would preclude the persistence of genetic diversity (Lewontin 1974;Nielsen 2005), colour polymorphisms exist in many invertebrate and vertebrate species. The same colour polymorphism can even sometimes be found in closely related species, suggesting the persistence of a shared evolutionary mechanism maintaining polymorphism across species radiation (Jamie and Meier 2020). We have studied a group of marine gastropod species within the subgenus Neritrema (genus Littorina) that share a similar shell colour polymorphism ( Fig. 2; Table S1). The results of the pooling analyses of assortative mating (light vs dark) using mating pairs from natural populations of the three studied species (Littorina fabalis, L. obtusata and L. saxatilis) show a consistent negative assortative mating across all populations and Table 1 Estimates of assortative mating (I PSI ) for shell colour in natural populations of L. fabalis, L. obtusata and L. saxatilis The studied localities for each species and the number (N) of microareas (mating pair + unmated individuals) are indicated. The first estimate of I PSI takes into account the colours grouped into dark (olive, brown) or light (white, yellow, orange). The second estimate is carried out for each separate colour vs the rest of the colours (e.g. white vs yellow/orange/olive/brown). For each I PSI the standard deviation is shown ? P < 0.10; *P < 0.05; **P < 0.01; ***P < 0.001 a From Estévez et al. (2020) 1996;Rolán-Alvarez et al. 2012, 2015bEstévez et al. 2020) and are also similar to those found in other species (Takashi and Hori 2008;Field and Barrett 2012;Holman et al. 2013;Hedrick et al. 2016Hedrick et al. , 2018. Even though our results cannot unravel the evolutionary mechanisms responsible for the observed patterns, they strongly suggest that negative assortative mating for shell colour in Neritrema species is caused by mate choice and therefore it could be at least partially responsible for the maintenance of the colour polymorphism within population via negative frequency-dependent sexual selection. There are at least two different pieces of evidence that support this claim, as we have already suggested: first, theoretical analysis have shown that any mate-choice based negative assortative mating can render a negative frequencydependent sexual selection that could be able to maintain the polymorphism (Pusey and Wolf 1996;Hedrick et al. 2016Hedrick et al. , 2018. Second, recently Estévez et al. (2020) have shown that this relationship between negative assortative mating, mate choice and negative frequency-dependent selection holds for L. fabalis (see similar evidences for other model organisms: Takashi and Hori 2008;Field and Barrett 2012;Holman et al. 2013;Hedrick et al. 2016Hedrick et al. , 2018. These two lines of evidence, in addition to the results presented in this study, suggest that an ancestral plesiomorphic behavioural mechanism could exist for mate choice within the Neritrema subgenus. We find the ancestral behavioural mechanism hypothesis plausible, since, from an evolutionary point of view, it is reasonable to expect similar behavioural mechanisms operating under similar ecological conditions in closely related species. This is akin to what in evolutionary studies of cognition has been dubbed the "evolutionary parsimony" principle (de Waal 1999). Thus, if L. fabalis, L. obtusata and L. saxatilis show a similar behavioural pattern, live in the same ecosystem, and they all share a (relatively) recent common ancestor, it is reasonable to expect a common behavioural mechanism behind this behavioural pattern. Nevertheless, this hypothesis, although plausible, needs to be experimentally tested in the future before concluding its validity. In any case, there are a few potential caveats in the present study that need to be addressed to dispel doubts about our results. First, the I PSI and similar indexes allow us to estimate assortative mating, and are frequently used as a proxy of mate choice, at least in laboratory conditions (reviewed in Gilbert and Starmer 1985;Rolán-Alvarez and Caballero 2000;Pérez-Figueroa et al. 2008). However, their capability to estimate mate choice has been discussed and several sources of bias have been identified (Rolán-Alvarez and Caballero 2000;Pérez-Figueroa et al. 2008). In any case, when they have been exclusively used to estimate assortative mating in the laboratory, these indices performed appropriately Rolán-Alvarez et al. 2015b). On the other hand, it has been shown that correlation indexes could be potentially inflated towards positive values when estimated directly in the wild, especially in species with low mobility, as well as when the study species shows great variability (e. g. colour) at a scale smaller than the scale at which the samples were collected (Rolán-Alvarez et al. 2015b). This bias could be corrected by estimating the I PSI in homogeneous groups (e. g. by similar colour), as it has been shown by simulations (Rolán-Alvarez et al. 2015b) and has been carried out here. In fact, this correction has been used previously in cases of negative and positive assortative mating (Rolán-Alvarez et al. 2015b;Ng et al. 2016;Estévez et al. 2018Estévez et al. , 2020. Thus, the correction used here does not correct, a priori, towards either negative or positive values of assortative mating. For example, the same correction was used to estimate the real scale at which mate choice occurs in two species of marine gastropods (Estévez et al. 2018), one characterised by positive assortative mating (Echinolittorina malaccana) and the other by negative assortative mating (L. fabalis), highlighting that this method is suitable to estimate both types of assortative mating. Therefore, the correction that we applied in this study assures a more reliable estimate of assortative mating in natural populations.
Second, another potential bias of our results could be caused by the actual observation of the mating pairs during the collection of the samples. For example, mating pairs in flat periwinkles (L. fabalis, L. obtusata) are commonly found among a dense canopy of Fucus vesiculosus, typically dark brown, allowing for a possible sampling bias towards the most visible mating pairs (light-light, dark-light). Although, to avoid this potential bias, we performed the search of mating pairs very close to the algae (~ 30 cm) and for each sampled microarea we checked all the individual algae, so it was unlikely to miss dark-dark pairs more often than the rest. Moreover, this bias in the sampling procedure has been rejected in previous studies as light-light mating pairs, which are highly noticeable on the Fucus canopy, appear at significantly lower frequencies than expected under random mating (Rolán-Alvarez et al. 2015b;Estévez et al. 2020). This potential sampling bias is even more unlikely in L. saxatilis, given that in this case all mating pairs are collected on the rock surface, as this species grazes on microalgae that grows on the bare rock. Therefore, in the light of the current results it seems plausible a hypothesis where a mating strategy of negative assortative mating for shell colour could represent an ancestral character within the subgenus Neritrema.
Recently, Estévez et al. (2020) have investigated the behavioural mechanism responsible for the negative assortative mating observed in L. fabalis from Abelleira (NW Spain) and the White Sea (NE Russia) populations. They found indirect evidence that mate choice is involved in assortative mating, with males showing preference for females with different shell colours. This, in turn, points to negative frequency-dependent selection as the evolutionary mechanism maintaining colour polymorphism in this species. However, their results did not provide any evidence regarding the potential contribution of colour to fitness, so that it would be favoured by mate choice. Their results neither clarify if colour, and not any other genetically linked trait, is the real target of selection. Here, we present a complementary result that could help to suggest the relationship between shell colours and mating in these species. In the analysis of each colour separately (Table 1 and Supplementary Table S2), we found an interesting pattern given that, although in general all colours showed the exact same trend as in those cases in which they are grouped in dark and light colours, there was an exception. The olive shell colour in L. fabalis changes from a pattern of random (positive but not significant estimates) mating in Spain to a pattern of negative assortative mating in the Russian populations (Table 1 and Supplementary Table S2; Fig. 3). This result is less probable if the shell colour itself is the target of the putative behavioural mechanism causing the assortative mating, but it could be expected if the colour is linked to another adaptive trait, physically or through pleiotropy. A relevant question is, therefore, what circumstances make a particular mating behaviour to be present in one location but not in others. One hypothesis is that the colour genes are associated with several other genes in a chromosomal inversion. In such scenario, different alleles of the colour gene could be associated with different chromosome variants (inversed or not). In support of this view is the fact that colour vision has not been detected and olfactory cues seems more relevant than visual cues in determining movement in aquatic gastropods (Seyer 1992;Wyeth 2019). Actually, Cephalopods, the animals with the most complex visual system within the phylum Mollusca, have only one visual pigment and therefore cannot see in colour (Nilsson 2013). If cephalopods do not have colour vision despite having a complex visual system, it seems implausible that gastropods could have it. Because of this, previous hypotheses about mate choice for shell colour in gastropods were based on the possibility that they can smell different colours when males follow mucus trails of females in the context of mate searching (see Estévez et al. 2020). Actually, a relationship between assortative mating for colour forms and sex chromosomes have been claimed in one finch species (Pryke 2010). Chromosomal inversions have been recently detected and associated with presumably adaptative traits in several organisms (Wellenreuther and Bernatchez 2018), including littorinids (Westram et al. 2018;Faria et al. 2019), and in at least one case has been related to negative assortative mating too (Hedrick et al. 2018). This could mean that the trait responsible for the pattern of negative assortative mating would be linked to the colour gene within the inversion, but could not be the colour itself. Although speculative, this hypothesis can be tested in the future and it opens a potentially fruitful new line of research.