The effects of Pleistocene glaciations on the phylogeography of Melitaea cinxia ( Lepidoptera : Nymphalidae )

Partial (600 bp) sequences of mitochondrial cytochrome oxidase I (COI) gene were used to infer the phylogeography of Melitaea cinxia (Lepidoptera: Nymphalidae) across the entire distributional range of the species, encompassing north Africa and Eurasia. Cladistic analysis of 49 distinct haplotypes (haplotype and nucleotide diversity were 0.95 and 0.027, respectively) revealed strong phylogeographic structure in M. cinxia, characterised by four major clades: Morocco; Western (Iberia, France, Italy); Central (central and northern Western Europe, Balkans, Greece, Anatolia, Levant); and Eastern (eastern Baltic, Urals, Iran, Siberia, China); separated by average pairwise distances of beween 2 and 6 percent. This pattern is consistent with the location of southern glacial refugia in the Iberian, Italian and Balkan peninsulas, as well as multiple eastern refugia. The Western clade is further structured into south-central Iberian, northern Iberian (and French) and southern Italian sub-clades; and the Eastern clade into Near Eastern and Far Eastern sub-clades; with weaker phylogeographical concordance within the Central clade, except for a large area in central and northern Western Europe which is monomorphic for COI haplotype. The Baltic and eastern Europe have been primarily colonized by the Far Eastern sub-clade, rather than the Central (Balkan) clade, highlighting the importance of including Near and Far Eastern populations in phylogeographic studies of Palearctic species. Maps showing the extent of clades and sub-clades suggest several regions of secondary contact and possible hybridization. Interspecific comparison of representative M. cinxia haplotypes supports a monophyletic origin of all M. cinxia. 675 * Corresponding and current address: Laboratory of Genetics, Department of Biology, University of Turku, 20014 Turku, Finland; e-mail: niklas.wahlberg@utu.fi We take a phylogeographic approach (Avise, 2000) to investigate the biogeographic history of the ecologically well-studied butterfly Melitaea cinxia (L.), the Glanville fritillary. M. cinxia is found on dry meadows containing its host plants Plantago and Veronica (both in Plantaginaceae) from the Atlantic coast to the Baikal Lake region in Siberia (Hanski, 1999) and to northwestern China (I. Hanski, pers. comm.). It is distributed between approximately 35°N and 60°N, though in the southern parts of its range it is restricted to higher elevations. There is a population in the Atlas Mountains of Morocco in northern Africa, which is often considered to be subspecifically distinct (Tolman & Lewington, 1997). Otherwise, M. cinxia is known for its morphological invariability. Although some populations have been described as subspecifically distinct (Tuzov, 2000; Churkin & Kolesnichenko, 2002), many authors consider all Eurasian populations of M. cinxia to belong to the nominotypical subspecies (e.g. Higgins, 1941; Tolman & Lewington, 1997). MATERIAL AND METHODS We sampled Melitaea cinxia from as wide a range of localities as possible (Table 1). Individuals were mainly sampled as adults, but in some cases larvae were collected and either reared to adults or used directly in DNA extraction. We extracted DNA mainly from two legs of dried specimens according to protocols described previously (Zimmermann et al., 2000; Wahlberg et al., 2003). We amplified the first half (600 bp) of the cytochrome oxidase subunit I (COI) gene in the mitochondrial genome from all the sampled individuals using the primer pair LCO and HCO (Folmer et al., 1994). The cycling profile for the PCR reactions was 95°C for 2 min, 35 cycles of 94°C for 30 s, 47°C for 30 s, 72°C for 1 min and a final extension period of 72°C for 10 min. PCR fragments were then sequenced using either an ABI 377 Automated Sequencer or a Beckman-Coulter CEQ8000 capillary sequencer (Bromma, Sweden) and a dye terminator cycle sequencing kit. Sequences were aligned by eye as the length of this protein-coding gene is highly conserved and insertions or deletions have not been observed in melitaeines (Wahlberg & Zimmermann, 2000; Zimmermann et al., 2000). Each unique haplotype sequence is available from GenBank


INTRODUCTION
Quaternary climatic oscillations, punctuated by the Pleistocene glaciations, caused massive changes to the distributions of species in the Palaearctic (Hewitt, 2000(Hewitt, , 2004;;Schmitt, 2007).Repeated cycles of demographic contraction and expansion, into and out of pockets of climatically protected regional refugia, combined with individual ecologies, has led to varied patterns of present day phylogeographic concordance among species (Hewitt, 1999).These patterns indicate the existence of multiple glacial refugia, large differences in the ages of intraspecific divergences, and differential rates of post-glacial expansion among refugial populations.Isolation among refugial populations promotes genetic and phenotypic differentiation as a result of independent adaptation to local environments and genetic drift, with consequences for reproductive isolation between discrete refugial lineages and the creation of hybrid zones where diverged lineages come into secondary contact (Hewitt, 1999).Thus, the Pleistocene glaciations are likely to have increased the effective size of many populations at the species level, but the same process typically reduces genetic diversity at the regional level, usually in proportion to the distance from the parent refugium (e.g.Bernatchez & Wilson, 1998;Schmitt & Seitz, 2002;Michaux et al., 2003) and depending on the mode of contraction and expansion (Ibrahim et al., 1996).
Phylogeographic studies of Palaearctic species are burgeoning but it remains the case that current phylogeographic paradigms for Europe are limited by the relatively small number of comprehensive case studies (Hewitt, 1999(Hewitt, , 2004)).Moreover, the overwhelming majority of intraspecific phylogeographic studies of Palaearctic species have focused on European populations of species that often have much wider distributions (Taberlet et al., 1998;Hewitt, 1999), meaning that our understanding of the history of Palaearctic species is based upon patterns found in only one-fifth of the landmass of this biogeographical region.This is especially a problem for insects, for which the fossil and subfossil record is very poor (with the exception of beetles), and for which inference of evolutionary history is currently based entirely on molecular data, particularly DNA sequences from the mitochondrial genome.Butterflies are no exception and the few published studies have been entirely restricted to western Europe (e.g., Joyce & Pullin, 2001;Schmitt & Seitz, 2001;Aagaard et al., 2002;Schmitt & Seitz, 2002;Habel et al., 2005;Schmitt et al., 2005;Weingartner et al., 2006), with the exception of one study on Aglais urticae (Vandewoestijne et al., 2004), and another on Hesperia comma (Forister et al., 2004) focused on North America but with some sampling of eastern and western Eurasia.
We take a phylogeographic approach (Avise, 2000) to investigate the biogeographic history of the ecologically well-studied butterfly Melitaea cinxia (L.), the Glanville fritillary.M. cinxia is found on dry meadows containing its host plants Plantago and Veronica (both in Plantaginaceae) from the Atlantic coast to the Baikal Lake region in Siberia (Hanski, 1999) and to northwestern China (I.Hanski, pers. comm.).It is distributed between approximately 35°N and 60°N, though in the southern parts of its range it is restricted to higher elevations.There is a population in the Atlas Mountains of Morocco in northern Africa, which is often considered to be subspecifically distinct (Tolman & Lewington, 1997).Otherwise, M. cinxia is known for its morphological invariability.Although some populations have been described as subspecifically distinct (Tuzov, 2000;Churkin & Kolesnichenko, 2002), many authors consider all Eurasian populations of M. cinxia to belong to the nominotypical subspecies (e.g.Higgins, 1941;Tolman & Lewington, 1997).

MATERIAL AND METHODS
We sampled Melitaea cinxia from as wide a range of localities as possible (Table 1).Individuals were mainly sampled as adults, but in some cases larvae were collected and either reared to adults or used directly in DNA extraction.We extracted DNA mainly from two legs of dried specimens according to protocols described previously (Zimmermann et al., 2000;Wahlberg et al., 2003).We amplified the first half (600 bp) of the cytochrome oxidase subunit I (COI) gene in the mitochondrial genome from all the sampled individuals using the primer pair LCO and HCO (Folmer et al., 1994).The cycling profile for the PCR reactions was 95°C for 2 min, 35 cycles of 94°C for 30 s, 47°C for 30 s, 72°C for 1 min and a final extension period of 72°C for 10 min.PCR fragments were then sequenced using either an ABI 377 Automated Sequencer or a Beckman-Coulter CEQ8000 capillary sequencer (Bromma, Sweden) and a dye terminator cycle sequencing kit.Sequences were aligned by eye as the length of this protein-coding gene is highly conserved and insertions or deletions have not been observed in melitaeines (Wahlberg & Zimmermann, 2000;Zimmermann et al., 2000).Each unique haplotype sequence is available from GenBank (Table 1).
Initially all the M. cinxia sequences were analyzed together to identify unique haplotypes in the program Arlequin 3 (Excoffier et al., 2005).The phylogenetic relationships of the unique haplotypes were then analyzed using parsimony and rooting with the outgroup Melitaea aurelia (Wahlberg & Zimmermann, 2000;Wahlberg et al., 2005b).Most parsimonious cladograms were searched for from the equally weighted and unordered data matrix using a heuristic search algorithm of NONA 2.0 (Goloboff, 1998).The heuristic searches were conducted with 10-1000 random addition replicates using tree bisectionreconnection (TBR) branch swapping with up to twenty trees held during each step.Trees were drawn with the aid of Winclada (Nixon, 2002).Clade support was evaluated using bootstrap with 1000 pseudoreplicates, with 5 random addition replicates within each pseudoreplicate.In addition, a minimumspanning network was constructed for all individuals using the program Arlequin 3 (Excoffier et al., 2005).
To verify the monophyly of M. cinxia haplotypes, a subset of ten haplotypes, representing the main lineages, was analysed with all the homologous Melitaea sequences from the study by Wahlberg & Zimmermann (2000).The resulting Melitaea tree was rooted using the related Chlosyne acastus from the same study.

RESULTS
The 120 individuals of M. cinxia from 58 localities yielded 49 haplotypes, of which 34 were found only in one individual (Table 1).The 600 bp showed 77 polymorphic sites with 71 observed transitions and 11 observed transversions.Haplotype diversity (h) was 0.95 (±0.01) and nucleotide diversity ( ) was 0.027 (±0.013).The uncorrected pairwise differences of the haplotypes varied widely from 0.1% to more than 6%, with the three haplotypes from Morocco being most diverged from the others (Table 2).One haplotype (Haplotype 25, Table 1) was particularly widespread, being found in Ukraine, northern Italy, Denmark, Sweden, Poland, Belgium and northern France.A few other haplotypes were found in more than one collection locality, although generally from localities relatively close to each other.Haplotype 4 was found in southern, central and northern Spain; haplotype 19 in northern Spain and eastern France; and haplotype 43 in Saratov and northern Iran.The greatest geographical separation between two identical haplotypes was for haplotype 31, found in the Lake Baikal region and in Åland, Finland.
Of the 77 variable sites, 44 were parsimony informative within M. cinxia.Cladistic analysis of the 49 M. cinxia haplotypes resulted in 408 equally parsimonious trees (length 178 steps, CI = 0.68, RI = 0.89).The low level of homoplasy (as shown by the CI and RI values) indicates that there is strong structure in the data set.Indeed the resulting trees are impervious to assumptions of transition:transversion weightings.Of the twelve inferred transversions ten are synapomorphic (or autapomorphic) with no reversals.Two are homoplastic and occur independently in different clades.Weighting transversions up to 10 times transitions has no effect on the topology or number of trees found.
There are four major clades in the haplotype tree (Fig. 2), which correspond to clearly definable geographical areas (Fig. 3).The Western clade contains haplotypes found mainly on the Iberian peninsula but also in    Finland to the Baikal Lake region and down into northwestern China and northern Iran.The haplotypes found in Morocco form their own lineage.
The Eastern clade appears to be the sister group of the Central clade, with the Western clade being sister to these two in the most parsimonious reconstruction of the data set.While this arrangement receives only moderate (79%) bootstrap support, the sister relationship of the Eastern and Central clades is supported by a unique nonsynonymous change at site 505 which leads to an amino acid change of an isoleucine to a valine.There are 8 other inferred amino acid changes, but these are autapomorphic or define small clades of two or three haplotypes.It can be noted that a neighbour joining algorithm using distances estimated with a variety of models (uncorrected p-values, the HKY85 model and the LogDet model) places the Western and Central clades as sister groups.However, such a scenario would entail independent acquisition of the amino acid change described above, or alternatively a single acquisition and a subsequent reversal to isoleucine in the Western clade.We prefer the more parsimonious scenario of a single acquisition of the amino acid change with no subsequent reversal.
A minimum-spanning network of the M. cinxia haplotypes (Fig. 4) shows the same patterns as the phylogeny, with the four major groups of haplotypes (Morocco, Western, Central and Eastern) clearly separated from each other by 9 or more inferred changes.Because of unequal sampling effort among localities the size of haplotype nodes in Fig. 4 does not accurately reflect the abundance of each haplotype in M. cinxia at the species level.This caveat applies principally to haplotypes 20, 32 and 38, respectively from Isle of Wight, Åland and Klucze, the three largest single locality samples.On the other hand, the large nodes for haplotypes 4 and 25, which include samples from many geographically separated localities, are genuinely common haplotypes at the species level.
The minimum-spanning network (Fig. 4) and subclades in Fig. 2 reveal further phylogeographic structuring within the major clades, particularly the Western and Eastern clades.A colour-coded map (Fig. 5) shows the distribution of these sub-clades.The Western clade has three distinct groups which we designate "Iberia", "France" and "southern Italy", with France being the sister group to Iberia (Fig. 2).Samples in the Iberia subclade originate almost exclusively from Spain, whereas the France sub-clade is predominantly composed of samples from southern and eastern France but also includes samples from the northeastern corner of Spain (and one central Spanish sample), the coastal Italy-France border, and the Isle of Wight.The southern half of Italy, including Sicily, appears to be monophyletic at the subclade level, though two samples falling into the southern Italian sub-clade originate from northern Italy and the northeast Italy-France border.
The Eastern clade consists of two geographically structured sub-clades, one occurring predominantly in the southwestern distribution of the clade ("Ural" sub-clade), and the other covering a very extensive area from the far east to the northwest of the clade ("Siberia" sub-clade).The Ural sub-clade contains haplotypes from the southern Ural Mountains, Saratov (Russia), northern Iran, Poland, Lithuania and the island of Hiiumaa (Estonia); the Siberia sub-clade contains haplotypes from Buryatia (Lake Baikal region), northwestern China and, surprisingly, the northern Estonian coast, Åland (Finland) and the Stockholm archipelago (Sweden).The Central clade does not exhibit any strong internal structure (Fig. 4) and shows the least amount of variation (

DISCUSSION
The COI haplotypes sequenced for this study support a monophyletic origin of all Melitaea cinxia, and have revealed a high degree of phylogeographic structuring across the species range.Taking the evidence given in Fig. 2, the deepest and earliest divergence occurred between Eurasian populations and populations in the North African Atlas mountains.The next split is likely to have occurred between the ancestors of the Western clade and Central + Eastern clades, closely followed by the final major divergence between the ancestors of Central and Eastern clade populations.Interestingly, the deep split between European and North African populations has been found in other butterfly species, such as Pararge aegeria (Weingartner et al., 2006) and two species of Polyommatus (Wiemers & Fiedler, 2007).However, the COI distance between the two populations of P. aegeria   is less than half (about 2%, Weingartner et al., 2006) the distance we have found for M. cinxia (more than 5%, Table 2).Weingartner et al. (2006) suggest that the ancestor of the genus Pararge colonised North Africa over 5 million years ago and subsequently recolonised Europe about 1 million years ago.In the case of M. cinxia, we cannot say much about the direction of colonisation since it may be that there is no extant sister species of M. cinxia, but rather that M. cinxia is sister to a large clade of Eurasian species related to M. diamina, M. arcesia and M. athalia (Wahlberg & Zimmermann, 2000;Wahlberg et al., 2005b).In any case, it appears that the split between populations in North Africa and Eurasia occurred about twice as much earlier in M. cinxia than in P. aegeria (assuming a similar rate of molecular evolution in COI).Alternatively, the rate of molecular evolution could be very different in these two species of Nymphalidae, leading to different estimates of relative ages of divergence.Absolute times of divergence in M. cinxia can only be investigated in a study of the entire genus Melitaea using age estimates from a recent study on the subfamily Nymphalinae (Wahlberg, 2006).
The Western and Central clades correspond well to having originated in the previously hypothesized glacial refugia in the Iberian, Italian, and Balkan Peninsulae (Hewitt, 1999).In the Western clade, all the populations from southern Italy appear to be sister to the populations from the Iberian peninsula and France, but in northern Italy haplotypes from France sub-clade, Central clade, and southern Italy sub-clade are all present.This pattern of haplotypes in Italy suggests that the characteristic southern Italian sub-clade diverged from the Iberian-French lineage in a southern Italian refugium during one or more glacials, and has subsequently expanded northwards, while French and Balkan haplotypes have colonized northern Italy from the west and east respectively.The other division of the Western clade into Iberia and France sub-clades supports the existence of two separate refugia for M. cinxia in southwestern Europe, one in southern Spain and another along coastal regions of southern France or northeastern Spain, consistent with mounting evidence of multiple glacial refugia in Iberia (Bella et al., 2006;Gómez & Lunt, 2007;Schmitt, 2007).
The geographic position of the Central clade is concordant with the idea of a Balkan refuge, although additional more southerly refuges in Greece, Turkey, and the Levant cannot be ruled out.Although there is much less divergence among COI lineages in the Central clade compared to the other two major clades, there is a degree of phylogenetic separation between haplotype 25 from central and northwestern Europe (northern France, Germany, Belgium, northern Italy, Denmark, southern Sweden, Poland, Hungary, Ukraine) and the remaining Central clade hap-682 lotypes in southeastern Europe, indicating that central and northwestern Europe was colonized by a single, probably Balkan lineage, while southeastern populations are derived from a diversity of refugial lineages.
The role of temperate Asia in the Pleistocene biogeographical history of Palaearctic insect species has not been studied in detail.Our results indicate that the area has been important in the evolutionary history of M. cinxia.Our sampling is rather sparse compared to European populations (Fig. 3), yet we find internal phylogeographic structure in the Eastern clade (Fig. 2), corresponding to a Caspian Sea refuge and a Far Eastern refuge.Given the relative geographical proximity of the Baltic region to the European refugia, it was a surprise to find that current populations of M. cinxia in Åland, Finland, and eastern Sweden are closely related to Far Eastern populations, in contrast to several other animal and plant species in this region that are typically derived from a Balkan refugium (Hewitt, 1999).Possible explanations include a faster rate of post-glacial expansion from a Far Eastern refugium, associated with a more continuously favourable habitat, perhaps combined with pre-adaptation to colder northern climates.Another interesting aspect of the Far Eastern expansion is the retention of genetic diversity at what we presume to be the northwestern colonizing front of this clade.The sparse sampling across Near and Far Eastern Asia does not allow critical assessment of haplotype diversity along hypothetical colonization routes, but the presence of eight distinct Siberia clade haplotypes in the Baltic region implies a relatively even demographic expansion.This pattern of mitotype diversity is in marked contrast to the significantly reduced COI variation in north European populations of the saproxylic beetle, Phyto kolwensis (Painter et al., 2007), also derived from eastern Eurasia.The far greater extent, during the last glacial maximum, of steppe habitat suitable for M. cinxia, compared to the much more restricted distribution of old spruce forest, on which P. kolwensis is dependent, is likely to have been a key factor underlying this difference (Grichuk, 1984;Tarasov et al., 2000).
The haplotype maps (Fig. 3 and 5) give some idea of the position of secondary contact zones between the major clades, and also between sub-clades.The Western (both the Iberia and France sub-clades) and Central clades meet in northeastern France and northern Italy (southern Italy and France sub-clades), where contact between France and southern Italy sub-clades is also likely.The main regions of secondary contact between Iberia and France sub-clades are in Catalonia and central Spain.In eastern Europe we detected Central and Eastern (Siberia and Urals sub-clade) haplotypes in one Polish population, and haplotypes from both Eastern sub-clades in Lithuania, Estonia, and eastern Sweden.Eastern Ukraine, eastern Turkey, and Kazakhstan, where we have not sampled, represent further regions of potential secondary contact.This study of mitochondrial haplotypes gives no information on the extent of hybridization and nuclear introgression between the clades and sub-clades, but pro-vides a basis for designing hybrid zone studies, including fitness consequences of hybridization in this species.
M. cinxia is a widely distributed species but has suffered severe declines in western Europe since the 1960s, largely as a result of agricultural intensification and urbanization (Hanski & Kuussaari, 1995).In this context, this phylogeographic study has some bearing on the conservation of this species, particularly as regards phylogenetic uniqueness of regional populations and reintroduction or augmentation programmes.For example, we note that the remnant UK population on the Isle of Wight is monomorphic for COI and is derived from the France sub-clade, but that surviving Belgian and German populations belong to the Central clade.These observations suggest that the now extinct populations in mainland UK and Holland were derived from France sub-clade and Central clade, respectively.
The COI haplotype variation in M. cinxia surpasses any other Eurasian butterfly species studied to date (Vandewoestijne et al., 2004;Wahlberg et al., 2005a;Weingartner et al., 2006).In contrast to the highly mobile Aglais urticae (Vandewoestijne et al., 2004), this variation is shown to have a strong phylogeographic signal, consistent with short-range dispersal (Hanski et al., 1994;Hanski, 1999) and host plant specificity (Kuussaari et al., 2000) in M. cinxia.Our results suggest that M. cinxia has had continuously large, widespread populations over the past 2 million years, and that local population dynamics (Hanski, 1999;Ehrlich & Hanski, 2004), coupled with multiple isolated glacial refugia, have allowed many haplotypes to persist through repeated cycles of demographic expansion and contraction.
France, Italy and Isle of Wight.The Central clade is distributed throughout central Europe from southern Sweden down to northern Italy and across to Turkey and Lebanon.The Eastern clade stretches out from central Sweden across

Fig. 5
Fig. 5 additionally gives an impression of where secondary contact between sub-clades is likely to occur.

Fig. 2 .
Fig. 2. Strict consensus of 408 most parsimonious trees for Melitaea cinxia haplotypes inferred for 600 bp of COI gene.Major lineages coded with symbols: square = Morocco, triangle = Western clade, star = Central clade, and circle = Eastern clade.Numbers above branches are bootstrap proportions.

Fig. 3 .
Fig. 3. Map showing localities of origin of specimens.Symbol codes according to Fig. 2: square = Morocco, triangle = Western clade, star = Central clade, and circle = Eastern clade.Overlapping symbols indicate where more than one major clade is represented in the same locality.

Fig. 5 .
Fig. 5. Geographic distribution of COI sub-clades in Europe, as indicated by haplotype tree (Fig. 2) and minimum-spanning network (Fig. 4).Colour codes are the same as in Fig. 4: dark blue = Iberia sub-clade, light blue = France sub-clade, purple = southern Italy sub-clade, green = Central clade, orange = Urals sub-clade, yellow = Siberia sub-clade.Multi-coloured circles indicate the presence of haplotypes from more than one clade or sub-clade at the same locality.The sub-clade designation of eastern samples not included in this map is Urals for Urals and Iran, and Siberia for Baikal and China.

TABLE 1 .
Geographical origin and collection year of Melitaea cinxia individuals used in this study, listed in order of COI haplotype number.

Table 2 )
. Fig 3 indicates regions of secondary contact between the major clades;