High altitude phylogeograph of selected Moroccan herptofauna

High Altitude Phylogeography of Selected Moroccan Herpetofauna Mafalda Cotrim Roberto Barata DEPARTAMENTO DE BIOLOGIA FACULDADE DE CIÊNCIAS DA UNIVERSIDADE DO PORTO Porto, 2013 High altitude phylogeography of selected Moroccan herpetofauna Mafalda Cotrim Roberto Barata Tese submetida à Faculdade de Ciências da Universidade do Porto para obtenção do grau de Doutor em Biologia DEPARTAMENTO DE BIOLOGIA FACULDADE DE CIÊNCIAS DA UNIVERSIDADE DO PORTO Porto, 2013 Nota Prévia Nos termos do nº 2 do artigo 8º do Decreto-Lei nº 388/70, foram incluídos em alguns capítulos desta dissertação os resultados de trabalhos já publicados ou em publicação. Em todos estes trabalhos, a candidata participou na obtenção, análise e discussão dos resultados, bem como na elaboração da publicação, embora sejam resultado de colaborações. A instituição de origem da candidata foi a Faculdade de Ciências da Universidade do Porto (FCUP), tendo o trabalho sido realizado sob orientação do Professor Doutor David James Harris, Professor Convidado da FCUP e Investigador do Centro de Investigação em Biodiversidade e Recursos Genéticos (CIBIO). A instituição de acolhimento foi inicialmente a Universitat de Biologia de Barcelona (UB) e posteriormente o Instituto de Biología Evolutiva do Consejo Superior de Investigaciones Científicas – Universitat Pompeu Fabra (CSIC-UPF), sempre sob a co-orientação do Professor Doutor Salvador Carranza, Investigador das referidas instituições. O trabalho laboratorial foi realizado no CIBIO, na UB e no CSIC-UPF. Este trabalho foi apoiado pela Fundação para a Ciência e a Tecnologia através da atribuição de uma bolsa de doutoramento de referência SFRH/BD/41488/2007, co-financiada pelo POPH/FSE. A todos os que nas Montanhas nos receberam e que do pouco que tinham partilharam, E a todos os outros que transformaram as 9 viagens a Marrocos numa aventura Resumo atenção no planeamento da conservação. Além disso, o uso de uma taxonomia integrativa é avidamente aconselhada, tendo em conta contudo, que os processos de especiação são complexos e que diferentes tipos de variação não são adquiridos necessariamente ao mesmo tempo nem na mesma ordem. Os resultados inesperados obtidos na filogenia dos Chalcides realçam as limitações e potenciais problemas do uso de um único caractere em taxonomia. Estes resultados salientam a importância do uso de múltiplas ferramentas no estudo da biodiversidade e dos processos evolucionários que lhe dão origem. Além disto, o estudo de taxa com um habitat limitado, isolado e de difícil acesso é importante para o conhecimento da biodiversidade real. A surpreendente biodiversidade observada nos répteis de Marrocos é, provavelmente, um exemplo da realidade em muitos outros habitats semelhantes, alertando para a necessidade do estudo dos mesmos. Abstract ABSTRACT Quantifying biodiversity is very important to understand the evolutionary history of life on Earth but specially to try to slow down or even reverse the loss of diversity that we are facing. Delimiting species, despite being a controversial issue, is of major importance since species are the basic units in areas such as ecology, biogeography and evolution and has serious implications for conservation biology. Montane systems, with high levels of endemisms and reduced gene flow between isolated habitats are, generally, poorly known and their species are particularly sensitive to climatic flutuations. This sensitivity issue is due to the fact that their small window of tolerances to temperature and elevation ranges can restrict their ability to persist. The specific goals of this thesis aimed to increase the knowledge about high altitude reptiles from Morocco. This included to study their distribution, levels of genetic, morphological and ecological diversity and the evolutionary relations between related taxa. With all this we intended to contribute to the knowledge about evolution history and diversity of high altitude endemic reptiles, not just from Morocco but also to see how patterns observed here reflect the patterns observed in other geographic regions. First, the known distributions of Quedenfeldtia and Atlantolacerta species were extended and some questions about the distributions and ecological requirements of the two Quedenfeldtia species were clarified. Secondly, high levels of genetic diversity were identified in the two high altitude specialists’ species, Q. trachyblepharus and Atlantolacerta andreanskyi accompanied by low levels of morphological variation. This cryptic diversity is concordant with previous results obtained for other North African reptiles, and high altitude reptiles from Europe. However the high diversity patterns found in these species are more likely to the “paleoendemics” found in southern African Mountains. Specifically in Atlatolacerta genus almost all sampled populations were demonstrated to be different isolated lineages, being proposed to be classified as six different species. The analyses of the three Chalcides species recovered a complex pattern of relationship that questions the current taxonomy. Furthermore, the present results question the taxonomy of various Chalcides species that is based mostly on colour patterns. Finally, this study reinforces the importance of assessing cryptic diversity as this provides information on diversity and speciation processes deserving special consideration in conservation planning. Additionally, the use of an integrative taxonomy is highly recommended, however, it is always necessary to bear in mind that, speciation is a complex process and different kinds of variation are not achieved necessarily at the same time or order. Abstract The unexpected results obtained in Chalcides species highlight the limitation and potential problems of using only one character in taxonomy. These results highlighted the importance of the use of multiple tools in the study of biodiversity and the evolutionary processes that gives rise. Additionally, the study of taxa with a limited and isolated habitat of difficult access is important to be aware of the real biodiversity. The amazing biodiversity observed in Moroccan reptiles, probably, is an example of what happen in several other similar habitats. Resumé RESUMÉ La quantification de la diversité est très importante pour comprendre l'histoire évolutive de notre planète mais surtout pour ralentir, voire inverser, l'actuelle perte de la diversité. La délimitation des espèces, en dépit d'être une question controversée, est d'une importance majeure car les espèces sont les unités de base dans des domaines tels que l'écologie, biogéographie et l'évolution et à de graves conséquences pour la biologie de la conservation. Les systèmes de montagne, avec des niveaux élevés d' endémismes et le flux génétique réduit entre habitats isolés, sont généralement mal connus; de plus les espèces sont particulièrement sensibles aux fluctuations climatiques à cause de leur étroite fenêtre de tolérances aux fluctuations de température (comme une élevation) qui peut limiter leur capacité à persister. Les objectifs spécifiques de cette thèse étaient d'augmenter notre connaissance sur les reptiles de haute altitude du Maroc. Cela comprenait la distribution, la diversité génétique, morphologique et écologique et les relations entre taxons apparentés. Nous avions l'intention par ce travail de contribuer à la connaissance de l'histoire de l'évolution et de la diversité de reptiles endémiques de haute altitude, non seulement du Maroc mais aussi voir comment les tendances observées ici sont reflétées dans d'autres régions géographiques. Tout d'abord, les distributions connues des espèces Quedenfeldtia et Atlantolacerta ont été élargies et les quelques doutes sur les distributions et les exigences écologiques des deux espèces de Quedenfeldtia ont été clarifiées. Deuxièmement, la diversité génétique élevée a été estimée sur les espèces spécialisées de haute altitude Q. trachyblepharus et Atlantolacerta andreanskyi accompagnée d'un faible niveau de variation morphologique.Cette diversité cryptique est concordante avec les résultats antérieurs obtenus pour les reptiles d'Afrique du Nord et de haute altitude en Europe. Cependant, ces patrons de diversité élevée sont plus susceptibles d'être "paleoendemics", comme celui (ceux) trouvés dans les montagnes de l'Afrique du Sud. Plus précisément, dans les populations d' Atlatolacerta, la quasi-totalité de l'échantillon appartenaient à différentes lignées isolées, qui ont été classés comme six espèces différentes. Les analyses des trois espèces de Chalcides montrent une relation complexe et floue qui remettent en question la taxonomie actuelle des Chalcides qui est basée la plupart du temps sur les modèles de couleur. Finalement, les résultats de cette étude renforcent l'importance d'évaluer la diversité cryptique car ils fournissent des informations sur la diversité et les processus de spéciation et méritent une attention particulière pour la planification de la conservation. En outre, l'utilisation d'une taxonomie integrée est extrêmement recommandé, cependant, il est toujours Resumé nécessaire de rappeler que la spéciation est un processus complexe et différents types de variation ne sont pas atteints en même temps ou ordre. Les résultats inattendus obtenus chez les espèces Chalcides mettent en évidence le problème de la seule utilisation de seul caractère en matière de taxonomie. Ces résultats mettent en évidence l'importance de l'utilisation de plusieurs outils dans l'étude de la biodiversité et les processus évolutifs qui lui donne lieu. En outre, l'étude des taxons dans un habitat limité et isolé, à l'accès difficile, est important pour connaître la biodiversité réelle. L'incroyable biodiversité observée chez les reptiles du Maroc, sans doute, est un exemple de ce qui se passe dans plusieurs autres habitats similaires. 15 TABLE OF CONTENTS Resumo * Summary * Résumé Page CHAPTER 1. GENERAL INTRODUCTION 23 Section 1.1. Quantifying biodiversity – Species definition and delimitation 27 1.1.1. Quantifying biodiversity 27 1.1.2. The problematic of species concept 27 1.1.3. Species delimitation 29 Section 1.2. Mountains as centers of speciation 31 Section 1.3. Study area: the Atlas Mountains, Morocco 33 1.3.1. Morocco 33 1.3.2. The Atlas Mountains 33 Section 1.4. Study group: high altitude reptiles from Morocco 37 1.4.1. Quedenfeldtia spp.(Boettger, 1883) 38 1.4.2. Atlantolacerta andreanskyi (Werner, 1929) 42 1.4.3. Chalcides spp. (Laurenti, 1768) 46 Section 1.5. Tools and methods 49 1.5.1. Molecular methods 50 1.5.2. Morphology 51 1.5.3. Ecological niche modelling 53 Section 1.6. Objectives and organization of the thesis 53 1.6.1. Objectives 53 1.6.2. Organization and thematic of the thesis 54 References 57 CHAPTER 2. CRYPTIC DIVERSITY IN QUEDENFELDTIA SPP. (BOETTGER, 1883)65 ARTICLE 1. 67 Barata M. Perera A. Martínez-Freiría F. and Harris D.J. 2012. Cryptic diversity within the Moroccan endemic day geckos Quedenfeldtia (Squamata: Gekkonidae): a multidisciplinary approach using genetic, morphological and ecological data. Biological Journal of Linnean Society, 106(4): 828-850. Table of contents 16 CHAPTER 3. CRYPTIC DIVERSITY IN ATLANTOLACERTA ANDREANSKYI 101 ARTICLE 2. 105 Barata M. Carranza S. and Harris D.J. 2012. Extreme genetic diversity in Atlantolacerta andreanskyi (Werner, 1929): A mountain cryptic species complex. BMC Evolutionary Biology, 12: 167. ARTICLE 3. 137 Barata M. Perera A. and Harris D.J. (submitted). Cryptic diversity in the Moroccan high altitude lacertid Atlantolacerta andreanskyi (Werner, 1928): a taxonomical assessment. CHAPTER 4. PHYLOGENETIC RELATIONSHIPS OF THREE CHALCIDES SPP.185 ARTICLE 4. 187 Barata M. Geniez P. Carranza S. and Harris D.J. (in preparation). Complex estimates of phylogenetic relationships between three species of Chalcides skinks from Morocco. CHAPTER 5. NEW OBSERVATIONS OF AMPHIBIANS AND REPTILES IN MOROCCO 201 ARTICLE 5. 203 Barata M. Perera A. Harris D.J. Van Der Meijden A. Carranza S. Ceacero F. García-Muñoz E. Gonçalves D. Henriques S. Jorge F. Marshall J.C. Pedrajas L. and Sousa P. 2011. New observations of amphibians and reptiles in Morocco, with a special emphasis on the Eastern Region. Herpetological Bulletin, 116: 4-14. CHAPTER 6. GENERAL DISCUSSION 219 Section 6.1. Sampling high altitude reptiles in Morocco 223 Section 6.2. Identifying levels of genetic variation and detecting cryptic diversity 225 Section 6.3. Clarifying the relation between Chalcides montanus,C. polylepis and C. manueli 230 Section 6.4. Final Remarks 231 Section 6.5. Future Perspectives 233 References 235 Table of contents 17 INDEX OF FIGURES Page Chapter 1 Figure 1. 28 The 34 hotspots identified by Conservation International in 2005 (Mittermeier et al. 2004). Figure 2. 32 Mediterranean Basin dissecation representation (6 Ma) and mountain chains from Hsu et al. (1973). Figure 3. 34 Map of Morocco with the identification of some important features, to frame the study area. Figure 4. 35 Origin of the Appalachian Orogen, a result of three separate continental collisions involving the North American continent and the collision of the African and North American continents during the Alleghenian Orogeny at the end of the Paleozoic (USGS – United States Geological Survey). Figure 5. 36 Representation of the movements of European and African continents, around 35 and 15 Mya. Figure 6. 39 Distribution map of Quedenfeldtia species based in Bons and Geniez (1996). Figure 7. 40 Photographs of males from the Quedenfeldtia species. Q. trachyblepharus and Q. moerens. Figure 8. 40 Previous and current higher order classification of extant Gekkota (Gamble et al. 2008). Figure 9. 41 Gecko phylogeny with time-calibrated using a Bayesian uncorrelated relaxed clock from Gamble et al. (2010). Figure 10. 43 Distribution map of Atlantolacerta andreanskyi based on Bons and Geniez (1996). Figure 11. 44 ML tree of a reanalysis of the mtDNA data set of Fu (2000), (cytochrome b, cytochrome oxidase I and 12S rRNA + 16S rRNA) (adapted from (Arnold et al. 2007). Figure 12. 45 Bayesian phylogenetic tree of the Lacertini based on mitochondrial DNA sequence (cytochrome band 12S rRNA) adapted from Arnold et al. (2007). Figure 13. 46 Adapted tree from Carranza et al. (2008) showing the detailed phylogenetic relationships in the Western clade. Figure 14. 47 Distribution map of three Chalcides species based on Bons and Geniez (1996). Figure 15. 48 Photo from a Chalcides montanus specimen (photo from Gabriel Martínez). Figure 16. 48 Photo from a Chalcides polylepis specimen from near Guelmin (photo from Gabriel Martínez). Figure 17. 49 Photo from Chalcides manueli specimen, from Essaouira (photo from Philippe Geniez). Figure 18. 52 Representation of standard measurements and scales counting used in lizards (Kaliontzopoulou et al. 2007). Figure 19. 52 Some examples of differentiation in colour patterns in Quedenfeldtia species. Table of contents 18 Chapter 2 Article 1 Figure 1. 71 Study area location, toponomies used in text, and distribution of the Quedenfeldtia species. Figure 2. 81 Trees derived from a Bayesian partitioned analysis for the mitochondrial DNA (12S rRNA, ND4 and tRNA’s) and nuclear sequences (Rag1, ACM4, MC1R and PDC). Figure 3. 84 Scatterplots of the first two axes of the multivariate Principal Component Analysis (PCA) and Multiple Correspondence Analysis (MCA). Figure 4. 87 Average and standard deviation of probability of occurrence and areas of probable sympatry for Quedenfeldtia moerens and Q. trachyblepharus in North-West Africa. Figure 5. 88 Response curves for the most related ecogeographical variables to the distribution of Quedenfeldtia moerens and Q. trachyblepharus in North-West Africa. Chapter 3 Article 2 Figure 1. 109 Atlantolacerta andreanskyi distribution map. Figure 2. 112 Trees resulting from partitioned Bayesian analysis, mtDNA tree (12S, ND4 and flanking tRNA-His), nuclear concatenated tree (RAG1, ACM4, MC1R, PDC and C-MOS), concatenated tree from the combined mitochondrial and nuclear DNA data and Species tree from mitochondrial and nuclear DNA data from the Bayesian Inference of Species Trees (STARBEAST). Figure 3. 114 Parsimony networks corresponding to MC1R (A), RAG1 (B), C-MOS (C), ACM4 (D) and PDC (E) nDNA sequence variation from all the populations. Figure 4. 115 Population structure estimation. Each individual is represented by a thin vertical line, which is partitioned into K coloured segments that represent the individual’s estimated membership fractions in K clusters. The bigger vertical divisions separate individuals from different populations. Article 3 Figure 1. 141 Atlantolacerta andreanskyi phylogenetic trees, mtDNA tree (12S and ND4) and nuclear (MC1R, PDC, ACM4, CMOS and RAG1), adapted from Barata et al. (2012a). Figure 2. 143 Distribution map of Atlantolacerta andreanskyi, populations used in the current study and the known distribution of the species as available in Bons and Geniez (1996). Figure 3. 148 Variation in multivariate size (mSIZE) and iso-corrected linear measurements in males and females of the A. andreanskyi lineages included in this study Figure 4. 149 Canonical discriminate function analysis (CDFA) of morphometric variables, for males and females. Figure 5. 152 Variation in pholidotic variables (log-transformed values), of males and females of the A. andreanskyi lineages, included in this study. Figure 6. 153 Canonical discriminate function analysis (CDFA) of pholidotic variables, for males and females. Figure 7. 154 Multiple Correspondence analysis (MCA) of male and female colour pattern. Table of contents 19 Figure 8. 177 Representative specimen from each lineage, J. Sirwa, Oukaimeden, Tizin Tichka, Outabati, J. Azourki and J. Ayache. Chapter4 Figure 1. 191 Distribution map of Chalcides samples, locations of the samples used in this study and the distribution from Bons and Geniez (1996). Figure 2. 194 Tree resulting from partitioned Bayesian analysis from mitochondrial DNA (12S rRNA and Cytb). Figure 3. 195 Tree resulting from partitioned Bayesian analysis from nuclear DNA (MC1R). Figure 3. 196 Parsimony network corresponding to MC1R gene fragment. Chapter 5 Figure 1. 206 Map of Morocco with the distribution of the sampling localities presented in this study. Figure 2. 211 Distribution map and photographs of Tarentola deserti,Stenodactylus sthenodactylus, Chalcides ocellatus,Trogonophis wiegmanni. Figure 3. 212 Distribution map and photographs of Leptotyphlops macrohynchus,Scutophis moilensis, Spalerosophis dolichospilus, Telescopus tripolinanus. INDEX OF PHOTOGRAPHS Cover – Morocco 2008, by Mafalda Barata 1 Cover – Moroccan man, J. Awlime base, 2009,by Fátima Jorge 5 Chapter 1 –High Atlas, Morocco 2008, by Mafalda Barata 23 Chapter 2 –Quedenfeldtia moerens, Morocco,by Mafalda Barata 65 Chapter 3 –Atlantolacerta andreanskyi, Oukaimeden, 2009,by Salvador Carranza 101 Article 2 – A. andreanskyi, 2009 by Salvador Carranza 103 –Landscape, J. Awlime, 2010, by Mafalda Barata 103 Article 3 –A. andreanskyi, 2010,by Dianna Steiner 135 Chapter 4 –Chalcides polylepis, Sidi Yahia, 2007,by Ana Perera 185 Chapter 5 –Landscape, Jebel Awlime, 2011,by Mafalda Barata 201 Chapter 6 –Moroccan children, 2009,by Mafalda Barata 219 Back Cover –Landscape; J. Sirwa, 2008 by Sónia Ferreira 26 General introduction 27 1.1. Quantifying Biodiversity - Species Definition and Delimitation 1.1.1. Quantifying Biodiversity Human beings always felt the need of understanding the processes that rules everything around them, and this is the main reason for the evolution of any kind of Science. The classification of biodiversity dates back at least to the 1700s when Linnaeus (1707-1778) developed a system of naming, ranking and classifying organisms, which is still in use today, his Systema Naturae. After Aristoteles classification of all known organisms in two groups, Kingdoms Plantae and Animalia, Linnaeus created the basis of the modern scientific systematics and his ideas on classification still influence all biologists even the ones that disagree with the roots of his ideas. On the other hand, the establishment of the fundamental ideas of evolutionary biology took place in 1859, with the publication of Darwin’s book “On the Origin of Species”, even though some of the ideas were older. Nowadays, the loss of biodiversity, most of it driven by human activities (habitat destruction, pollution and introduction of exotic species), is increasing (Begon et al. 2006), and this fact weighs more than the simple curiosity in the need of knowledge. The attempt of slowing down the loss of biodiversity and reverse the processes that lead to extinctions, highlight the necessity of quantifying diversity accurately. According to Myers (2000), biodiversity hotspots are defined as areas containing exceptional concentrations of endemic species and facing exceptional loss of habitats. Following this definition, the identification of those areas is the first and most important step to prevent biodiversity loss (Myers 2003). Recently, Mittermeier et al. (2004) updated to 34 (Fig. 1), the initial hotspot list of Conservation International (CI) that identified 25 terrestrial areas of the world for priority conservation. Nowadays, the knowledge of biodiversity remains unsatisfactory due to Linnean and Wallacean shortfalls or, in other words, there are many species that have not been formally described, and geographical distributions of most species are still poorly understood and usually contain many gaps (Whittaker et al. 2005). 1.1.2. The Problem of Species Concept Species are the fundamental basic units for studies of ecology, evolution, systematic and conservation biology (Wiens 1999). However the seed of the discussion about their definition was growing long before Darwin (Britton 1908; Wilkins 2009), and even now it has not reached a general consensus (reviewed by(Mayden 1997; de Queiroz 1998; Harrison 1998). Darwin’s book “On the origin of species” (1859) further heated the discussion about the species concept, bringing uncertainty for some biologists due to Darwin’s explanation of the General introduction 28 evolutionary process as a gradual process, but without giving details about how one species gives rise to two (Bailey 1896). Figure 1. The 34 hotspots identified by Mittermeier et al. (2004) (Conservation International, 2005). Therefore all species definitions are, in some way incomplete, since they are static concepts that attempt to define a continuous process that occurs in very different organisms in completely diverse habitats. Some of the most used are: Biological species concept (Dobzhansky 1935; Mayr 1942; 1963): has historically been the most widely used and accepted definition. Ernest Mayr defined species as "groups of actually or potentially, interbreeding natural populations that are reproductively isolated from other such groups." This definition allows the existence of subspecies. Ecological species concept (Vanvalen 1976; Ridley 1993): Ridley defined species as a “group of organisms that is exploiting or adapted to a set of resources - niche”. This group evolves separately from all other groups outside its range (Andersson 1990). Evolutionary species concept (Simpson 1961; Wiley 1978; Templeton 1989): Wiley defined species as a “single lineage of ancestor-descendant populations, which maintain its identity from other such lineages and which has its own evolutionary tendencies and historical fate” Cohesion species concept (Templeton 1989): “A species is the most inclusive group of organisms having the potential for genetic and/or demographic exchangeability”. General introduction 29 Phylogenetic species concept (Cracraft 1983; 1989): “species is an irreducible (basal) cluster of organisms, diagnosable distinct from other such clusters, and within which there is a parental pattern of ancestry and descent”. It implies monophyly (commonly inferred from possession of shared derived character states), exclusive coalescence of alleles (all alleles of a given gene descend from a common ancestral allele not shared with those of other species) and diagnosability of qualitative, fixed differences according to different authors (de Queiroz 1998). Cracraft (2002) argue that the question “what is a species?” remains the most important of the “seven great questions of systematic biology”, which is true especially because after clarifying this question, biologists can more easily focus on the problem of species delimitation. According to de Queiroz (2007), one of the main problems is that the issue of species delimitation has long been confused with the definition of species itself, leading to disagreement in the methods to define boundaries and numbers of species. The same author (de Queiroz (2007) defends that this problem is not as “serious as it appears”, because despite the differences between the several species concepts, they have an essential conceptual agreement which provides the basis for a “unified species concept”. This common element associates species to separately evolving metapopulation lineages, or more specifically, an ancestor-descendent series (Simpson 1961; Hull 1980). de Queiroz (2007) finalizes this unified species concept saying that a species is not an entire metapopulation lineage but only a segment of that lineage in a way that species derive from other species. The only property needed to delimit species would be detecting a segment of a metapopulation lineage evolving separately (de Queiroz 1998). Other evidences mentioned in the several different species concepts, like reproductive isolation, reciprocal monophyly, phenetic distinguishability or occupation of a distinct niche or adaptive zone will be useful lines of evidence important to determine species delimitation but fail to be determinants of a species existence (de Queiroz 2007). 1.1.3. Species Delimitation Defining the species boundaries and describing new ones is important to various biological fields including biogeography, ecology, evolutionary biology, and conservation (Sites and Marshall 2003; Agapow 2005). Assuming the general metapopulation lineage species concept explained above (section 1.1.2), speciation involves lineage separation and divergence that can lead to reproductive isolation, ecological divergence, morphological distinctness and reciprocal monophyly (de Queiroz 1998; Harrison 1998). These criteria were and still are the ones used by systematists as evidence to delimit species, and these criteria can arise at different times and in different order during lineage formation processes (de Queiroz 2007). General introduction 30 Species delimitation becomes extremely difficult, especially when any of the criteria are not achieved, but this is expected in recent or adaptive radiations, as speciation is a continuous process (Wake 2006). The principal advantage of the separation of species concept from species delimitations criteria is to release the species delimitation studies from the controversy that is surrounding the species concept. Consequently, such studies can concentrate on investigating all evidence relevant to the recognition of evolutionary independent metapopulations lineages (de Queiroz 2007). According to Knowles and Carstens (2007) species lineages can be delimited long before reciprocal monophyly is achieved as several recent works have shown that gene genealogies give information about the history of species split (gene trees), despite the presence of incomplete lineage sorting (Degnan and Salter 2005; Maddison and Knowles 2006; Carstens and Knowles 2007). Various authors agree that a substantial amount of time is needed for observing reciprocal monophyly after the initial divergence of species (Hudson and Coyne 2002; Hudson and Turelli 2003). The use of one or two mitochondrial genes to reconstruct the “species tree” has been largely criticized due to the differences found between gene trees and species trees (Maddison 1997; Degnan and Rosenberg 2009; Edwards 2009). The absence of recombination that makes this molecule so attractive to perform phylogenetic research is also one of its most important limitations, differing from the species tree not only due to stochastic processes of lineage sorting but also due to events of hybridization and introgression. Furthermore, this molecule only reflects the evolution of female lineage, that can differ from the species history (Rosenberg and Nordborg 2002)see examples in(Shaw 2002). Back in 1988, Pamilo and Nei suggested that it was better to combine information from different loci than to add more samples, but the Bayesian methods based on the coalescent theory makes possible the use of the gene trees to construct a species tree instead of simply concatenating all the information (Heled and Drummond 2010). New emerging methods to delimit species based on the coalescent theory are changing the species delimitation field. These methods use recent theoretical models that combine gene trees and species phylogenies from multilocus sequence data. The methodology starts from a sample of genes and trace backwards in time to infer the demographic history from the population since the most recent common ancestor of the sampled genes. Several recent works used this approach to determine species boundaries and identify new lineages (Heled and Drummond 2010; Leaché and Fujita 2010; Yang and Rannala 2010). The recent conceptual advances in integrative taxonomy (Padial et al. 2010) support that the approach to deal with the problem of delimiting species is the use of multiple and complementary disciplines for a consistent identification of species. The use of multiple lines General introduction 31 of evidence, using different data types and diagnostic methods to get the most relevant information is achieving a consensus between biologists (Sites and Marshall 2004; de Queiroz 2007; Knowles and Carstens 2007; Leache et al. 2009). However, some uncertainty in the species boundaries of recent lineages may remain due to incomplete separation, secondary introgression, sampling deficiencies, or disagreement of criteria. Furthermore, the degree of congruence that different characters must show to consider a population or a group of populations as a separate species, splits integrative taxonomists (Padial et al. 2010). Following these points of reasoning, to better understand the processes involved in species delimitation requires a model with particular characteristics. Firstly, distinct sources of data (for example, morphological and molecular data) must be obtainable. Reptiles are ideal for this kind of studies, since the alpha taxonomy, although complex is not particularly challenging, they can be found in large numbers, and genetic markers are widely available (see the review in(Camargo et al. 2010). Additionally, preferably, the geographic locality should be discreet so that variation can be assessed across the range of the model organisms. Islands are often employed for this purpose, as “natural laboratories” to study evolution. Analogously, in this thesis, selected taxa of reptiles from Mountains have been used as model organisms to study the evolution history and diversity of the reptiles living on those conditions. Mountains share many of the advantages that make islands attractive settings for evolutionary studies – as Mayr (1967) states, “islands may demonstrate certain biological phenomena almost with the clarity of test-tube experiments; indeed every island is an experiment on its own” and mountains can be considered, for their lizard fauna, as islands. 1.2. Mountains as Centres of Speciation Mountains, like most other isolated areas, are generally species deprived but rich in endemisms, as is typical of islands (Darwin 1859). Isolation is the key feature that makes these habitats so attractive models for evolutionary and ecological studies (Emerson 2002). A variety of mechanisms could have produced montane restrictions in non-volant species such as reptiles. Mountains in northern latitudes, which were especially affected by the climatic fluctuations of the Pleistocene, are likely to have relatively recent fauna, possibly having been colonized by cold-tolerant fauna from lower regions that moved into the mountains during warmer periods. Moreover, the formation of temporary corridors of suitable habitat may have allowed pre-adapted mountain forms to invade new habitats. In other areas, lineages may have simply adapted into montane forms during the orogenesis of the mountain system. In this case, the geological age of the mountains will be correlated with the age of the species. Alternatively, mountain taxa may have been restricted or even displaced by other taxa in nearby lowlands and due to this developed into mountain specialists (Carranza et al. 2004). These different mechanisms will imprint identifiable characteristics in the genetic General introduction 32 make-up of the populations. However, it is possible that during colder climates, such as those that existed during glacial periods of the Pleistocene, these populations extended to lower elevations and more continuous ranges permitted gene flow between populations that were otherwise isolated by regions of unsuitable habitat (Wiens 2004). Nevertheless, high mountain populations may have persisted in separate refugia, resulting in high levels of genetic, and perhaps phenotypic, divergence (Garcia-Paris et al. 2000; Bowie et al. 2006; Cadena et al. 2007). High altitude specialist species are particularly sensitive to climatic flutuations due to the fact that their small window of tolerances to temperature and elevation ranges can restrict their ability to persist in, or disperse across, different habitats (Janzen 1967; Ghalambor et al. 2006; Deutsch et al. 2008; McCain 2009). As many of the species that inhabit mountain biodiversity hotspots are endemic to a single mountain, or limited number of adjacent mountains they face higher levels of extinction risk (Gifford and Kozak 2011). The mountains around the Mediterranean basin (Fig. 2), such as the Pyrenees,the Alps, the Balkans and Rhodope Mountains and the Atlas Mountains, have been glaciated on several occasions through the Quaternary (Hughes et al. 2006). Due to this fact and to its location at the interface between the North Atlantic Ocean, the western Mediterranean Sea and the Sahara Desert (Hughes 2008), the Atlas Mountains are an especially interesting setting for evolutionary studies. Figure 2. Mediterranean Basin desiccation representation (6 Ma) and mountain chains from Hsu et al. (1973). Pyrenees Alps Balkans Rhodope Mountains Atlas Mountains General introduction 33 1.3. Study Area: The Atlas Mountains, Morocco 1.3.1. Morocco North Africa is a model region for phylogeographic studies of Montane fauna. Morocco (Fig. 3) is in the northwest extreme of the African continent and includes an area of 710 850 Km2, with the Sahara Desert occupying half of it. Morocco together with Algeria and Tunisia form an area known as the Maghreb. The geology, weather, relief, fauna and flora make this area very different from the remaining continent and, in various aspects, more similar to Southern Europe. Morocco suffers the influence of the Mediterranean Sea, Atlantic Ocean and the Sahara climate. All these characteristics and diversity are responsible for the enormous richness and high level of endemism found in the herpetofauna of Morocco. Of the 104 species that constitute the herpetofauna, 22 are endemic. This makes Morocco the richest country in the Occidental Paleartic area in terms of herpetofauna (Bons and Geniez 1996). These combine to mean that the area contains considerable diversity but has definable boundaries and includes various features that can be used as calibration points to try to date phylogeographic breaks. During the mid-Tertiary, the collision between the Eurasia and Africa plates resulted in a diverse, complex and unusual geographic and topographical scenario. A great diversity of confined climates arose and during the mid-Pliocene to Pleistocene this region was repeatedly affected by alternated humid and arid phases (Le Houerou 1997). The principal vegetation in Morocco are the forests, distributed in the Tingitana Peninsula, Rif, Middle Atlas, Beni Snassen, the plateaux of Debdou and Jerada, some parts of the North face of the High Atlas and Jebel Sirwa, the Souss plain and the occidental extreme of the Anti-Atlas. The principal tree species are Quercus ilex, Tretraclinis articulata, Argania spinosa, Quercus suber, Juniperus sp., Cedrus libanotica atlântica, Pinus pinaster and Abies maroccana. Other areas are, generally, covered by shrubs (Bons and Geniez 1996). 1.3.2. The Atlas Mountains The Atlas Mountains are a mountain range that reaches 2500 Km through Morocco, Algeria and Tunisia. In Morocco, the Atlas Mountains are divided in three main mountains belts, Middle Atlas (max. 3340 m a.s.l.), High Atlas (max. 4167 m a.s.l.) and Anti-Atlas (max. 2531 m a.s.l.), the last extending in a NE/SW direction. The highest peak is the Toubkal, in the High Atlas, with 4167 m; this and other higher mountains are covered with snow during winter and spring and these mountains constitute a very important reservoir of water, where the most important rivers come from, including the Moulouya River (520 Km), Sebou river (458 Km), Oum-er-Rbia river (555 Km) and Tennsitt river (270 Km) (Schleich et al. 1996). General introduction 34 Figure 3. Map of Morocco with the identification of some important features, to frame the study area. This Mountain range separates Morocco into two bioclimatic regions, acting as a barrier to the dispersal of several Mediterranean species coming from the north, and also to species coming from the Sahara (Bons and Geniez 1996). The Atlas Mountains were formed in three subsequent events of Earth’s history. The first tectonic deformation was in the Paleozoic (around 300 Mya) and involved only the Anti-Atlas formation. The Anti-Atlas chain is believed to have initially been formed as part of Alleghanian orogeny (Hatcher 2008) and as a result of the collision between America and Africa (Fig. 4). At that time, North America was part of the super-continent Euroamerica, while Africa was part of Gondwana. The collision between the two super-continents formed the super-continent Pangaea, which comprised all major continental landmasses. In the Mesozoic Era (before 65 Mya) a second event took place. This event consisted in an extension of the Earth’s crust that rifted and separated the American and African continents. Most of the rocks that form the High Atlas today were deposited under the ocean during this period. In the Tertiary (between 68 to 1.8 Mya), approximately 35 Mya the landmasses of Europe and Africa collided in the area where is currently the Strait of Gibraltar in the Mesinian, and the mountain chains that today comprises the Atlas uplifted. This tectonic convergence was the responsible for the closure of the Strait of Gibraltar with the consequent elimination of much of the original General introduction 35 Figure 4. Origin of the Appalachian Mountains, a result of three separate continental collisions involving the North American continent and the collision of the African and North American continents during the Alleghenian Orogeny at the end of the Paleozoic (USGS –United States Geological Survey.) Tethys and the formation of some of the more than twenty mountains ranges that today encircles the Mediterranean Basin (Fig. 5A), as the High Atlas, Alps and Pyrenees. From roughly 15 Mya ago (Fig. 5B), Africa continued moving northwards but, for the first time in many millions of years, it also started to move westwards clockwise, producing the opening of the Red Sea, and the formation of more mountains in Turkey, Southern Europe and North Africa and slowly the connection between the Mediterranean and the Atlantic Ocean at its western end was closed again (Teixell et al. 2005; Missenard et al. 2006). General introduction 42 The analyses of Rato and Harris (2008) indicated the paraphyly of Saurodactylus with S. mauritanicus, and S. brosseti closer to Teratoscincus przewalskii than to S. fasciatus. Based on this fact, they proposed to divide the Saurodactylus in two genera, with mauritanicus and brosseti remaining in Saurodactylus, and creating a new monotypic genus to fasciatus. Such a finding shows the importance of extensive sampling within other groups, such as Quedenfeldtia, to fully assess their phylogenetic relationships. In general, geckos present high genetic diversity with very conservative morphology (Gamble et al. 2008; Rato and Harris 2008; Perera and Harris 2010), recent studies show that two North African genera (Quedenfeldtia and Saurodactylus, both endemic to Morocco) are basal to the American Sphaerodactylidae family that diverged by vicariance and dispersal events after fragmentation of Gondwana (Gamble et al. 2008; Gamble et al. 2010). In consequence, such North African endemics as Quedenfeldtia are expected to retain high levels of genetic diversity because of their old evolutionary histories (Busack 1986; Rato and Harris 2008). 1.4.2. Atlantolacerta andreanskyi (Werner, 1929) The Atlas Dwarf Lizard, Atlantolacerta andreanskyi, is a lacertid lizard endemic from the highest peaks (2400 to 3800 m) of western and central parts of the High Atlas Mountains (Fig. 10) (Bons and Geniez 1996; Schleich et al. 1996). They can be found in alpine meadows, scree, amongst boulders, and in areas of thorn cushion vegetation and thickets near small watercourses or plateaux in the top of the mountains that retain some water from rain or snowmelt (Geniez 2005); authors personal observation). It is a small lacertid lizard similar to a half-grown Zootoca vivipara in size and pattern, with a light brown middorsal region and dark brown flanks. This species does not have accentuated sexual dimorphism: the male’ head is relatively larger and the body shorter, its extended foreleg reaches the anterior border of the eye while in females only its posterior border is reached (Schleich et al. 1996). The systematics of Atlantolacerta andreanskyi has been historically complex. Initially it was placed in several different genera and subgenera within the Lacertini subtribe, including Zootoca (Pasteur and Bons 1960), Lacerta part II (Arnold 1973), Lacerta incertae sedis (Guillaume 1987), and Lacerta sensu lato (Arnold 1989). Later, (Volobouev et al. 1990), based on cytochemical methods, suggested that A. andreanskyi and L. vivipera (within the subgenus Zootoca), considered sister species, presented similar patterns. The taxonomic situation was revised when Arnold et al. (2007) assigned andreanskyi as a basal member of Eremiadini. A reanalysis of the mtDNA dataset of Fu (2000) obtained the same results (Fig 11). General introduction 43 Figure 10. Distribution map of Atlantolacerta andreanskyi based on Bons and Geniez (1996). This position was unexpected given that this species lacks the synapomorphies that characterize most other Eremiadini, namely a derived condition of the ulnar nerve and the presence of a fully developed armature in the hemipenis, which has folded lobes when retracted. Because of its probably basal position, without close relationship to any other genus of Eremiadini and its distinctive morphology, the High Atlas endemic was described as a new monotypic genus: Atlantolacerta by (Arnold et al. 2007). Besides the basal position, they suggest that A. andreanskyi is between 16 and 12 Mya (Fig. 12), an old origin inside the Lacetidae. Pavlicev and Mayer (2009) latter confirmed the phylogenetic relationships and generic status of A. andreanskyi, using a combined analysis of mitochondrial and nuclear sequences (C-MOS and RAG1). The different populations of A. andreanskyi present an apparently disjunct distribution (Bons and Geniez 1996; Schleich et al. 1996), and this situation is similar to an archipelago, with the different “islands” being represented by mountaintops unconnected due to areas of unsuitable habitat below 2400 m. As a result of this scenario, minimal gene flow is currently expected between the different populations; even though it is not known how the different climatic events occurred during the Miocene and Pleistocene have affected this species. Despite this, some aspects of the biology of A. andreanskyi are already well known (Busack 1987; Carretero et al. 2006), although all available information comes from only one population, at Oukaimeden. The genetic structure of the different populations, as well as the relationships between them have never been assessed before. General introduction 44 Quedenfeldtia Figure 11. ML tree of a reanalysis of the mtDNA data set of Fu (2000), based on 4522 bp (1026 bp of cytochrome b, 1048 bp of cytochrome oxidase I and 2448 bp of the ribosomal genes 12S rRNA + 16S rRNA) (adapted from(Arnold et al. 2007). General introduction 45 Quedenfeldtia Figure 12. Bayesian phylogenetic tree of the Lacertini based on mitochondrial DNA sequence (cytb and 12S rRNA) adapted from Arnold et al. (2007). These analyses support the origin of the Lacertini in the Mid-Miocene (between 12 and 16 Mya). Atlantolacerta andreanskyi, a member of the Eremiadini, is sister to the Lacertini. General introduction 46 1.4.3. Chalcides spp. (Laurenti 1768) Chalcides is a genus of the Scincidae family, comprising approximately 24 species with a wide distribution, from Southern Europe, North Africa to Somalia and Kenya, Turkey, Iraq, Arabia, coastal Iran and Pakistan (Carranza et al. 2008). Pasteur (1981) suggested that Morocco was an important evolutionary centre for this genus, probably because most of the species occur there and in surrounding areas, including several endemisms. This genus has a typical elongated body, round in cross sections and limbs short or reduced (Schleich et al. 1996). Many of the species are morphologically similar and difficult to identify. The taxonomy of Chalcides has been revised, all or in part, by several authors (Boulenger 1887; 1890; 1896; Boulenger 1920; Lanza 1957; Pasteur 1981; Caputo 1993; Mateo et al. 1995; Greenbaum 2005; Greenbaum et al. 2006) However it was difficult to estimate a phylogeny based on morphological features and there are still uncertainty in species boundaries and in the relationships between them (Fig. 13;(Carranza et al. 2008). Figure 13. Adapted tree from Carranza et al. (2008) showing the detailed phylogenetic relationships in the Western clade.Note that C. montanus is nested with a paraphyletic C. polylepis. General introduction 47 In this study we focused on the relationships between three endemic Chalcides species from Morocco, Chalcides montanus,Chalcides polylepis and Chalcides manueli (distribution shown in Fig. 14). Figure 14. Distribution map of three Chalcides species based on Bons and Geniez (1996). Black dots represent C. montanus, white dots represent C. polylepis and red dots C. manueli. Chalcides montanus Werner, 1931 is the Chalcides found at higher altitude (until 2500 m) in North Africa, normally found in cold and humid mountain regions, with low bushes. It is a small Chalcides (Fig. 15) with white parallel lines on the neck and a yellow venter. Morphologically it is similar to C. polylepis (but smaller), with which it is sympatric (Schleich et al. 1996). Before Caputo and Mellado (1992) attributed specie rank to this form, it was considered a subspecies of Chalcides ocellatus.Chalcides lanzai was long considered a subspecies of C. montanus (Schleich et al. 1996), but mtDNA confirmed that they are unrelated (Carranza et al. 2008). General introduction 48 Figure 15. Photo from a Chalcides montanus specimen (photo from Gabriel Martínez). Chalcides polylepis Boulenger 1990, also similar to C. ocellatus is one of the largest of the genus (Fig. 16) with a larger body and head, and ocelli forming parallel lines on the upper side. Unlike C. montanus, this species is usually found in hillsides or flat country up to 2000 m (Schleich et al. 1996). In a recent phylogeny, three samples of C. montanus from one locality appear within the C. polylepis clade, making it paraphyletic (Carranza et al. 2008). The authors suggested two possible explanations: C. montanus might be just a highland form of C. polylepis, with a different altitudinal habitat, or the second possibility is that C. montanus mitochondrial DNA was received from C. polylepis through introgression. However this issue remained unresolved as Carranza et al (2008) used only mtDNA evidence. Figure 16. Photo from a Chalcides polylepis specimen from near Guelmin (photo from Gabriel Martínez). Chalcides manueli Hediger 1935, is a small Chalcides also similar to C. ocellatus but without any conspicuous dorsal pattern, see Fig. 17 (Schleich et al. 1996). It was known from localities in the coast between Essaouira and Sidi Ifni (Bons and Geniez 1996; Carranza et al. 2008). Recently, specimens morphologically identified as C. montanus from Jebel Sirwa and General introduction 49 Tizin Tichka where genetically (mtDNA, 12S rRNA) closer to C. manueli (Harris et al. 2010; Barata et al. 2011) from Sidi Ifni and sister taxa to C. manueli from Essaouira, confirming that taxonomy and relationships within Moroccan Chalcides are complex and in need of revision. This species is listed as vulnerable in the IUCN red list, due to its small range distribution being, in fact, only known for 11 localities (Joger et al. 2006). Figure 17. Photo from Chalcides manueli specimen, from Essaouira (photo from Philippe Geniez). Although this thesis was initially focussed only on high altitude reptile species from Morocco, given the clearly conflicting evidence between morphological characters and mtDNA for these three species of Chalcides, sampling was extended so that greater numbers of all three could be included in an improved estimate of phylogenetic relationships, and to assess fully if introgression was indeed taking place. 1.5. Tools and methods Before molecular techniques where developed, morphological characters were the primary source of information used in the study of taxonomy. Different species were, for a long time, classified based only on different morphological characters (e.g. see(Schleich et al. 1996). However, use of morphological characters can underestimate species diversity due to cryptic species (Baker and Bradley 2006), species can look very similar or even be identical but are reproductively isolated. Furthermore, unrelated taxa can acquire similar appearance as a consequence of convergent evolution or mimicry. On the other hand, sampling a few genetic markers does not reveal patterns that are identifiable using morphological, behavioural or ecological information, that could be crucial to separate two species (Knowles and Carstens General introduction 50 2007). A classical example is the case of the cichlid fish species from some African Lakes that show large functional diversity but have limited genetic variation. The hybrids are viable and fertile (Seehausen et al. 2003). In such situations, morphology is still a powerful tool to study ecology and behaviour and can be used to support the characterization of new species. Recently, another method that is earning the confidence of biologists, especially in ecology but even in species delimitation area, is ecological niche modelling. 1.5.1. Molecular Methods The methods for delimiting species changed dramatically in the seventies, when molecular techniques become more widespread (Avise et al. 1979). Currently, the amount of DNA sequences published is still increasing exponentially, and the software available for analysing this data is also improving (Knowles and Carstens 2007; Kubatko and Degnan 2007). DNA has become the most popular source of data for reconstructing phylogenies (Harris 1999), and even for understanding the mechanisms acting in populations or lineages. However, despite the initial idea that molecular data would readily resolve all phylogenies (Diaz 2007), these new methods, besides keeping old methodological problems, brought new ones, maintaining a constant search for increasing objectivity (Suárez-Díaz and Anaya-Muñoz 2008). Most early studies relied only mtDNA and although this genome is a very useful tool, it has particular problems. The non-recombinant characteristic can have important limitations that are now widely recognized. The first one is that the analyses of mtDNA correspond to the study of a single locus that reflects the history of that molecule rather that the species history (Zhang and Hewitt 2003). This can be due to the effects of natural selection and introgression. Furthermore, mtDNA only reflects the evolutionary history of the female lineage, which in some situations can be completely different from the population or species history (see(Magri et al. 2006). Currently, despite the relative ease of obtaining sequence data from multiple loci, the most appropriate methods used to analyse this data are still debatable. Although some authors argue that standard methods that concatenate multigene data are enough for accurate phylogenetic estimates (Chen and Li 2001; Rokas et al. 2003), many studies revealed that gene histories could be different from the species histories (Kolaczkowski and Thornton 2004; Mossel and Vigoda 2005; Kubatko and Degnan 2007). This lack of concordance between gene trees and species trees can result from diverse processes such as coalescence, gene flow, selection, hybridisation, and gene duplication (Maddison 1997; Kubatko and Degnan 2007; Edwards 2009). As a result of this, information from different unlinked genetic markers (mitochondrial and nuclear) is thus necessary for delimiting evolutionary lineages, as well as for establishing phylogenetic relationships. General introduction 51 Apart from that, a widespread problem in the use of mtDNA alone is introgression between different taxa. This event has being observed in several different groups as reptiles (Pinho et al. 2008), insects (Zakharov et al. 2009), and mammals (Alves et al. 2006; Boratynski et al. 2011) and can be confounded with ancestral polymorphism, showing how evolutionary inferences can be misrepresented when they are based on single locus. Bayesian analysis was successively introduced to calculate trees and to test hypothesis (Huelsenbeck and Ronquist 2001; Drummond and Rambaut 2007). Recently, several new methods were developed in order to calculate species tree from multiple loci in a coalescent framework (see(Blair and Murphy 2011). The well-known benefits of considering the stochasticity of genetic processes promoted the development of coalescent based approaches (Kubatko and Degnan 2007; Heled and Drummond 2008). Coalescent methods (mentioned before in section 1.1.3) do not uses gene genealogies directly to infer demographic history. Instead, this methods use those genealogies as a nuisance parameter to get estimations of biogeographically informative parameters (Hey and Machado 2003). Instead of focusing on gene trees, this method uses data from multiple loci and multiple specimens from population to estimate the “species tree”, an estimate of the history of divergence (Belfiore et al. 2008; Brumfield et al. 2008). On the other hand, these coalescence-based methods are still inadequately applied in non-model species phylogenies, especially at the interface between populations and species, where the threat of incomplete lineage sorting is greatest (Degnan and Rosenberg 2009), and where lack of prior information about the organisms and surrounding environment is often lacking. 1.5.2. Morphology Morphology is a branch of science dealing with the form and structure of the organisms and their specific structural descriptions as shape, colour, pattern and structure of internal parts as bones and organs. Most taxa differ morphologically from other taxa but there are exceptions like cryptic species and convergent evolution as mentioned before (section 1.5). For a long time and before molecular techniques were developed, taxonomy was based only on morphological information. There are several fields of morphology, including the study of the patterns of the locus of structures within the body plan of an organism (comparative morphology), the relationship between form and the function (functional morphology), effects of external factors upon the morphology of organisms under experimental condition (experimental morphology), the study of the internal structures of an organism (anatomy), and the external appearance of the organisms (eidonomy). There are several studies regarding the morphology of reptiles, such as indicators in ecotoxicology (Amaral et al. 2012), detection of cryptic species (Kaliontzopoulou et al. 2012), parasitology (Carretero et al. 2011; Maia et al. 2011), geometric morphometrics to General introduction 58 Boulenger G.A. 1896. A list of the reptiles and batrachians collected by Dr Ragazzi in Shoa and Eritrea. Annali del Museo Civico di Storia Naturale 16(2nd Series): 546-554. Bowie R.C.K. Fjeldså J. Hackett S.J. Bates J.M. and Crowe T.M. 2006. Coalescent models reveal the relative roles of ancestral polymorphism, vicariance, and dispersal in shaping phylogeographical structure of an African montane forest robin. Molecular Phylogenetics and Evolution 38(1): 171-188. Britton N.L. 1908. The taxonomic aspect of the species question. American Naturalist 42: 225-242. Brumfield R.T. Liu L. Lum D.E. and Edwards S.V. 2008. Comparison of Species Tree Methods for Reconstructing the Phylogeny of Bearded Manakins (Aves: Pipridae, Manacus) from Multilocus Sequence Data. Systematic Biology 57(5): 719-731. Busack S.D. 1986. Biogeographic Analysis of the Herpetofauna Separated by the Formation of the Strait of Gibraltar. National Geographic Research 2: 17-36. Busack S.D. 1987. Notes on the biology of Lacerta andreanszkyi (Reptilia: Lacertidae). AmphibiaReptilia 8: 231-236. Cadena C.D. Klicka J. and Ricklefs R.E. 2007. Evolutionary differentiation in the Neotropical montane region: Molecular phylogenetics and phylogeography of Buarremon brush-finches (Aves, Emberizidae). Molecular Phylogenetics and Evolution 44(3): 993-1016. Camargo A. Sinervo B. and Sites J.W. 2010. Lizards as model organisms for linking phylogeographic and speciation studies. Molecular ecology 19(16): 3250-3270. Caputo V. 1993. Taxonomy and evolution of the Chalcides chalcides complex (Reptilia, Scincidae) with description of two new species. Museo Regionale di Scienze Naturali Bollettino 11(1): 47-120. Caputo V. and Mellado J. 1992. A new species of Chalcides (Reptilia: Scincidae) from northeastern Morocco. Boll. Zool. 59. Carranza S. Arnold E.N. and Amat F. 2004. DNA phylogeny of Lacerta (Iberolacerta) and other lacertine lizards (Reptilia: Lacertidae): did competition cause long-term mountain restriction? Systematics and Biodiversity 2(01): 57-77. Carranza S. Arnold E.N. Geniez P. Roca J. and Mateo J.A. 2008. Radiation, multiple dispersal and parallelism in the skinks, Chalcides and Sphenops (Squamata: Scincidae), with comments on Scincus and Scincopus and the age of the Sahara Desert. Molecular Phylogenetics and Evolution 46(3): 1071-1094. Carranza S. Arnold E.N. Mateo J.A. and Geniez P. 2002. Relationships and evolution of the North African geckos, Geckonia and Tarentola (Reptilia: Gekkonidae), based on mitochondrial and nuclear DNA sequences. Molecular Phylogenetics and Evolution 23(2): 244-256. Carretero M.A. Perera A. Harris D.J. Batista V. and Pinho C. 2006. Spring diet and trophic partitioning in an alpine lizard community from Morocco. African Zoology 41(1): 113-122. Carretero M.A. Roca V. Larbes S. Ferrero A. and Jorge F. 2011. Intestinal Helminth Parasites of Wall Lizards, Podarcis vaucheri Complex (Sauna: Lacertidae) from Algeria. Journal of Herpetology 45(3): 385-388. Carstens B.C. and Knowles L.L. 2007. Estimating Species Phylogeny from Gene-Tree Probabilities Despite Incomplete Lineage Sorting: An Example from Melanoplus Grasshoppers.Systematic Biology 56(3): 400-411. Cheddadi R. and Bar-Hen A. 2009. Spatial gradient of temperature and potential vegetation feedback across Europe during the late Quaternary. Climate Dynamics 32(2-3): 371-379. Cheddadi R. Lamb H.F. Guiot J. and van der Kaars S. 1998. Holocene climatic change in Morocco: a quantitative reconstruction from pollen data. Climate Dynamics 14(12): 883-890. Chen F.C. and Li W.H. 2001. Genomic divergences between humans and other hominoids and the effective population size of the common ancestor of humans and chimpanzees. American Journal of Human Genetics 68(2): 444-456. Cracraft J. 1983. Species concepts and speciation analysis. Current Ornithology 1: 159-187. Cracraft J. 1989. Speciation and its ontology: the empirical consequences of alternative species concepts for understanding patterns and processes of differentiation in Speciation and its consequences by D. Otte and J. A. Endler. Sunderland, Mass, Sinauer Associates. Cracraft J. 2002. The Seven Great Questions of Systematic Biology: An Essential Foundation for Conservation and the Sustainable Use of Biodiversity. Annals of the Missouri Botanical Garden 89: 127-144. Darwin C. 1859. On the origin of species by means of natural selection or the Preservation of Favoured Races in the Struggle for Life.1St ed. London, John Murray. de Queiroz K. 1998. The General Lineage Concept of Species, Species Chteria, and the Process of Speciation: A Conceptual Unification and Terminological Recommendations in Endless General introduction 59 forms: Species and speciation by D. J. Howard and S. H. Berlocher. New York, Oxford University Press. de Queiroz K. 2007. Species Concepts and Species Delimitation. Systematic Biology 56(6): 879-886. Degnan J.H. and Rosenberg N.A. 2009. Gene tree discordance, phylogenetic inference and the multispecies coalescent. Trends in ecology & evolution 24(6): 332-340. Degnan J.H. and Salter L.A. 2005. Gene tree distribution under coalescent process. Evolution 59(1): 24-37. Deutsch C. Tewksbury J. Huey R. Sheldon K. Ghalambor C. Haak D. and Martin P. 2008. Impacts of climate warming on terrestrial ectotherms across latitude. Proc Natl Acad Sci U S A105(18): 6668-6672. Diaz E.S. 2007. The rhetoric of informational molecules: Authority and promises in the early study of molecular evolution. Science in Context 20(4): 649-677. Dobzhansky T. 1935. A critique of the species concept in biology. Philosophy of Science 2: 344-355. Donnellan S.C. Hutchinson M.N. and Saint K.M. 1999. Molecular evidence for the phylogeny of Australian gekkonoid lizards. Biological Journal of the Linnean Society 67(1): 97-118. Drummond A.J. and Rambaut A. 2007. BEAST: Bayesian evolutionary analysis by sampling trees. BMC Evolutionary Biology 7(1): 214. Edwards S.V. 2009. Is a New and General Theory of Molecular Systematics Emerging? Evolution 63(1): 1-19. Elith J. and Leathwick J.R. 2009. Species Distribution Models: Ecological Explanation and Prediction Across Space and Time. Annual Review of Ecology Evolution and Systematics 40: 677-697. Emerson B.C. 2002. Evolution on oceanic islands: molecular phylogenetic approaches to understanding pattern and process. Molecular Ecology 11: 951-966. Fahd S. Barata M. Benítez M. Brito J.C. Caro J. Carvalho S. Chirosa M. Feriche M. Herrera T. Márquez-Ferrando R. Nesbitt D. Pleguezuelos J.M. Reques R. Rodríguez M.P. Santos X. Sicilia M. and Vasconcelos R. 2007. Presencia de la víbora hocicuda Vipera latastei en el Atlas Medio (Marruecos) y otras citas herpetológicas para la región. Boletín de la Asociación Herpetológica Española 18: 26–34. Fonseca M.M. Brito J.C. Paulo O.S. Carretero M.A. and Harris D.J. 2009. Systematic and phylogeographical assessment of the Acanthodactylus erythrurus group (Reptilia: Lacertidae) based on phylogenetic analyses of mitochondrial and nuclear DNA. Molecular Phylogenetics and Evolution 51(2): 131-142. Fonseca M.M. Brito J.C. Rebelo H. Kalboussi M. Larbes S. Carretero M.A. and Harris D.J. 2008. Genetic variation among spiny-footed lizards in the Acanthodactylus pardalis group from North Africa. African Zoology 43(1): 8-15. Fu J.Z. 2000. Toward the phylogeny of the family Lacertidae - Why 4708 base pairs of mtDNA sequences cannot draw the picture. Biological Journal of the Linnean Society 71(2): 203-217. Gamble T. Bauer A.M. Colli G.R. Greenbaum E. Jackman T.R. Vitt L.J. and Simons A.M. 2010. Coming to America: multiple origins of New World geckos. Journal of Evolutionary Biology 24(2): 231-244. Gamble T. Bauer A.M. Greenbaum E. and Jackman T.R. 2008. Evidence for Gondwanan vicariance in an ancient clade of gecko lizards. Journal of Biogeography 35(1): 88-104. Garcia-Barros E. Gurrea P. Lucianez M.J. Cano J.M. Munguira M.L. Moreno J.C. Sainz H. Sanz M.J. and Simon J.C. 2002. Parsimony analysis of endemicity and its application to animal and plant geographical distributions in the Ibero-Balearic region (western Mediterranean). Journal of Biogeography 29(1): 109-124. Garcia-Paris M. Good D.A. Parra-Olea G. and Wake D.B. 2000. Biodiversity of Costa Rican salamanders: Implications of high levels of genetic differentiation and phylogeographic structure for species formation. Proceedings of the National Academy of Sciences of the United States of America 97(4): 1640-1647. Geniez P. 2005. Lacerta andreanskyi in 2006 IUCN Red List of Threatened Species by www.iucnredlist.org. Ghalambor C.K. Huey R.B. Martin P.R. Tewksbury J.J. and Wang G. 2006. Are mountain passes higher in the tropics? Janzen's hypothesis revisited. Integrative and Comparative Biology 46(1): 5-17. Gifford M.E. and Kozak K.H. 2011. Islands in the sky or squeezed at the top? Ecological causes of elevational range limits in montane salamanders. Ecography 35(3): 193-203. General introduction 60 Greenbaum E. 2005. Systematics of West African skinks in the Chalcides thierryi group: Composition, distribution, and redescription of types. African Journal of Herpetology 54(1): 17-29. Greenbaum E. Campbell A.C. and Raxworthy C.J. 2006. A Revision of Sub-Saharan Chalcides (Squamata: Scincidae), with Redescriptions of Two East African Species. Herpetologica 62(1): 71-89. Guillaume C.P. 1987. Les petites Lacertidés du bassim méditerranéen occidental (genera Podarcis et Archaeolacerta essentiellement). Sur quelques problèmes d'ordre systématique et viogéographique. Montpellier, 474. Guisan A. and Zimmermann N.E. 2000. Predictive habitat distribution models in ecology. Ecological Modelling 135(2-3): 147-186. Han D. Zhou K. and Bauer A.M. 2004. Phylogenetic relationships among gekkotan lizards inferred from C-mos nuclear DNA sequences and a new classification of the Gekkota. Biological Journal of the Linnean Society 83(3): 353-368. Harris D.J. 1999. Molecular systematics and evolution of lacertid lizards. Natura Croatica 83(3): 161180. Harris D.J. Perera A. Barata M. Tarroso P. and Salvi D. 2010. New distribution notes for terrestrial herpetofauna from Morocco. North-Western Journal of Zoology 6(2): xx-xx. Harrison R.G. 1998. Linking evolutionary pattern and process, the relevange of species concept for the study of speciation in Endless forms, species and speciation by D. J. Howard and S. H. Berlocher, Oxford University Press. Hatcher R.D.J. 2008. Tracking Lower-to-Mid-to-Upper Crustal Deformation Processes through Time and Space through Three Paleozoic Orogenies in the Southern Appalachians Using Dated Metamorphic Assemblages and Faults. Geological Society of America 40(6). Heled J. and Drummond A.J. 2008. Bayesian inference of population size history from multiple loci. BMC Evolutionary Biology 8. Heled J. and Drummond A.J. 2010. Bayesian inference of species trees from multilocus data Molecular Biology and Evolution 27(3): 570-580. Hewitt G. 2000. The genetic legacy of the quaternary ice ages. Nature 405: 907-913. Hewitt G.M. 1999. Post-glacial re-colonization of European biota. Biological journal of the Linnean Society 68: 87-112. Hewitt G.M. 2001. Speciation, hybrid zones and phylogeographyor seeing gemes in space and time. Molecular ecology 10: 537-549. Hey J. and Machado C.A. 2003. The study of structured populations [mdash] new hope for a difficult and divided science. Nat Rev Genet 4(7): 535-543. Hoogmoed M.S. 1974. Echsen aus Nordafrika. Herpetologishe Impressionen aus Südmarokko Aquarien Magazine 8(7): 304-310. Hsu K.J. 1978. Messinian Salinity Crisis - Evidence of Late Miocene Eustatic Changes in World Ocean. Naturwissenschaften 65(3): 151-151. Hsu K.J. Ryan W.B.F. and Cita M.B. 1973. Late Miocene Desiccation of Mediterranean. Nature 242(5395): 240-244. Hudson R.R. and Coyne J.A. 2002. Mathematical consequences of the genealogical species concept. Evolution 56(8): 1557-1565. Hudson R.R. and Turelli M. 2003. Stochasticity overrules the "three-times rule": Genetic drift, genetic draft, and coalescence times for nuclear loci versus mitochondrial DNA. Evolution 57(1): 182-190. Huelsenbeck J.P. and Ronquist F. 2001. MrBayes: Bayesian inference of phylogeny. Bioinformatics 17: 754-755. Hughes P.D. 2008. The glacial history of the Atlas Mountains, Morocco. Quaternary Newsletter 116: 48-50. Hughes P.D. Woodward J.C. and Gibbard P.L. 2006. Quaternary glacial history of the Mediterranean mountains. Progress in Physical Geography 30(3): 334-364. Hull D.L. 1980. Individuality and Selection. Annual Review of Ecology and Systematics 11: 311-332. Huntley B. and Birks J.D.S. 1983. An Atlas of Past and Present Pollen Maps for Europe: 0-13,000 years ago. Cambridge, Cambridge University Press. Janzen D.H. 1967. Why Mountain Passes Are Higher in Tropics. American Naturalist 101(919): 233- &. Joger U. Slimani T. El Mouden H. and Geniez P. 2006. Chalcides manueli. In: IUCN 2011. IUCN Red List of Threatened Species Version 2011.2. General introduction 61 Johnston T.H. 1924. The relation of climate to the spread of prickly pear. Transactions of the Royal Society of South Australia 48: 269-295. Kaliontzopoulou A. Carretero M.A. and Liorentel G.A. 2007. Multivariate and geometric morphometrics in the analysis of sexual dimorphism variation in Podarcis lizards. Journal of Morphology 268(2): 152-165. Kaliontzopoulou A. Carretero M.A. and Llorente G.A. 2010a. Intraspecific ecomorphological variation: linear and geometric morphometrics reveal habitat-related patterns within Podarcis bocagei wall lizards. Journal of Evolutionary Biology 23(6): 1234-1244. Kaliontzopoulou A. Carretero M.A. and Llorente G.A. 2010b. Sexual dimorphism in traits related to locomotion: ontogenetic patterns of variation in Podarcis wall lizards. Biological Journal of the Linnean Society 99(3): 530-543. Kaliontzopoulou A. Carretero M.A. and Llorente G.A. 2012. Morphology of the Podarcis wall lizards (Squamata: Lacertidae) from the Iberian Peninsula and North Africa: patterns of variation in a putative cryptic species complex. Zoological Journal of the Linnean Society 164(1): 173-193. Knowles L.L. and Carstens B.C. 2007. Delimiting species without monophyletic gene trees. Systematic Biology 56(6): 887-895. Kolaczkowski B. and Thornton J.W. 2004. Performance of maximum parsimony and likelihood phylogenetics when evolution is heterogeneous. Nature 431(7011): 980-984. Kubatko L.S. and Degnan J.H. 2007. Inconsistency of phylogenetic estimates from concatenated data under coalescence. Systematic Biology 56(1): 17-24. Lanza B. 1957. Su alcuni "Chalcides" del Marocco (Reptilia: Scincidae). Monit. Zool. Italiano 65: 8598. Le Houerou H.N. 1997. Climate, flora and fauna changes in the Sahara over the past 500 million years. Journal of Arid Environments 37(4): 619-647. Leaché A.D. and Fujita M.K. 2010. Bayesian species delimitation in West African forest geckos (Hemidactylus fasciatus). Procedings of Royal Society B xx: xx-xx. Leache A.D. Koo M.S. Spencer C.L. Papenfuss T.J. Fisher R.N. and McGuire J.A. 2009. Quantifying ecological, morphological, and genetic variation to delimit species in the coast horned lizard species complex (Phrynosoma). Proceedings of the National Academy of Sciences of the United States of America 106(30): 12418-12423. Lee M.S.Y. 2004. The molecularisation of taxonomy. Invertebrate Systematics 18(1): 1-6. Maddison W.P. 1997. Gene trees in species trees. Systematic Biology 46(3): 523-536. Maddison W.P. and Knowles L.L. 2006. Inferring phylogeny despite incomplete lineage sorting. Systematic Biology 55(1): 21-30. Magri D. Vendramin G.G. Comps B. Dupanloup I. Geburek T. Gomory D. Latalowa M. Litt T. Paule L. Roure J.M. Tantau I. van der Knaap W.O. Petit R.J. and de Beaulieu J.L. 2006. A new scenario for the Quaternary history of European beech populations: palaeobotanical evidence and genetic consequences. New Phytologist 171(1): 199-221. Maia J.P.M.C. Harris D.J. and Perera A. 2011. Molecular Survey of Hepatozoon Species in Lizards from North Africa. Journal of Parasitology 97(3): 513-517. Martinez-Freiria F. Santos X. Pleguezuelos J.M. Lizana M. and Brito J.C. 2009. Geographical patterns of morphological variation and environmental correlates in contact zones: a multiscale approach using two Mediterranean vipers (Serpentes). Journal of Zoological Systematics and Evolutionary Research 47(4): 357-367. Mateo J.A. Geniez P. and Bons J. 1995. Saurians of the genus Chalcides Laurenti 1768 (Reptilia, Scincidae) in Morocco, I. review and distribution. Revista Española de Herpetología 9: 7-36. Mayden R.L. 1997. A hierarchy of species concepts: the denouement in the saga of the species problem in Species: The units of diversity by M. F. Claridge, H. A. Dawah and M. R. Wilson. London, Chapman and Hall, 381-423. Mayr E. 1942. Systematics and the Origin of Species from the Viewpoint of a Zoologist. . New York, Columbia University Press. Mayr E. 1963. Animal Species and Evolution. Cambridge, Belknap Press of Harvard University Press. Mayr E. 1967. The challenge of island faunas. Australian Natural History 15: 359-374. McCain C.M. 2009. Global analysis of bird elevational diversity. Global Ecology and Biogeography 18(3): 346-360. Missenard Y. Zeyen H. de Lamotte D.F. Leturmy P. Petit C. Sebrier M. and Saddiqi O. 2006. Crustal versus asthenospheric origin of relief of the Atlas Mountains of Morocco. Journal of Geophysical Research-Solid Earth 111(B3). General introduction 62 Mittermeier R.A. Gil P.R. Hoffmann M. Pilgrim J. Brooks T. Mittermeier C.G. Lamoreux J. and G.A.B. d.F. 2004. Hotspots Revisited: Earth’s Biologically Richest and Most Endangered Ecoregions. Mexico City, CEMEX. Mossel E. and Vigoda E. 2005. Phylogenetic MCMC algorithms are misleading on mixtures of trees. Science 309(5744): 2207-2209. Mouret V. Guillaumet A. Cheylan M. Pottier G. Ferchaud A.L. and Crochet P.A. 2011. The legacy of ice ages in mountain species: post-glacial colonization of mountain tops rather than current range fragmentation determines mitochondrial genetic diversity in an endemic Pyrenean rock lizard. Journal of Biogeography 38(9): 1717-1731. Myers N. 2003. Biodiversity hotspots revisited. Bioscience 53(10): 916-917. Myers N. Mittermeier R.A. Mittermeier C.G. da Fonseca G.A.B. and Kent J. 2000. Biodiversity hotspots for conservation priorities. Nature 403(6772): 853-858. Padial J.M. Miralles A. De la Riva I. and Vences M. 2010. The integrative future of taxonomy. Front Zool 7: 16. Pamilo P. and Nei M. 1988. Relationships between Gene Trees and Species Trees. Molecular Biology and Evolution 5(5): 568-583. Pasteur G. 1981. A Survey of the Species Groups of the Old World Scincid Genus Chalcides.Journal of Herpetology 15(1): 1-16. Pasteur G. and Bons J. 1960. Catalogue des reptiles actuels du Maroc, Vol. 18: Série Zoologie, Travaux de l´Institut scientifique chérifien. Pavlicev M. and Mayer W. 2009. Fast radiation of the subfamily Lacertinae (Reptilia: Lacertidae): History or methodical artefact? Molecular Phylogenetics and Evolution 52(3): 727-734. Perera A. and Harris D.J. 2010. Genetic variability within the Oudri’s fan-footed gecko Ptyodactylus oudrii in North Africa assessed using mitochondrial and nuclear DNA sequences. Molecular Phylogenetics and Evolution 54: 634-639. Perera A. Vasconcelos R. Harris D.J. Brown R.P. Carretero M.A. and Perez-Mellado V. 2007. Complex patterns of morphological and mtDNA variation in Lacerta perspicillata (Reptilia; Lacertidae). Biological Journal of the Linnean Society 90(3): 479-490. Pinho C. Harris D.J. and Ferrand N. 2008. Non-equilibrium estimates of gene flow inferred from nuclear genealogies suggest that Iberian and North African wall lizards (Podarcis spp.) are an assemblage of incipient species. BMC Evolutionary Biology 8. Rato C. and Harris D.J. 2008. Genetic variation within Saurodactylus and its phylogenetic relationships within the Gekkonoidea estimated from mitochondrial and nuclear DNA sequences. Amphibia-Reptilia 29(1): 25-34. Rhoujjati A. Cheddadi R. Taieb M. Baali A. and Ortu E. 2010. Environmental changes over the past c. 29,000 years in the Middle Atlas (Morocco): A record from Lake Ifrah. Journal of Arid Environments 74(7): 737-745. Ridley M. 1993. Evolution. Cambrigde, Blackwell Science. Rokas A. Williams B.L. King N. and Carroll S.B. 2003. Genome-scale approaches to resolving incongruence in molecular phylogenies. Nature 425(6960): 798-804. Rosenberg N.A. and Nordborg M. 2002. Genealogical trees, coalescent theory and analysis of genetic polymorphisms. Nature 3: 380-391. Schleich H.H. Kastle W. and Kabisch K. 1996. Amphibians and Reptiles from North Africa. Königstein, Germany., Koeltz Scientific Publications. Seehausen O. Koetsier E. Schneider M.V. Chapman L.J. Chapman C.A. Knight M.E. Turner G.F. van Alphen J.J.M. and Bills R. 2003. Nuclear markers reveal unexpected genetic variation and a Congolese-Nilotic origin of the Lake Victoria cichlid species flock (vol 270, pg 129, 2003). Proceedings of the Royal Society of London Series B-Biological Sciences 270(1533): 2637-2638. Shaw K.L. 2002. Conflict between nuclear and mitochondrial DNA phylogenies of a recent species radiation: What mtDNA reveals and conceals about modes of speciation in Hawaiian crickets. Proceedings of the National Academy of Sciences of the United States of America 99(25): 16122-16127. Sillero N. 2011. What does ecological modelling model? A proposed classification of ecological niche models based on their underlying methods. Ecological Modelling 222(8): 1343-1346. Simpson G.G. 1961. Principles of Animal Taxonomy. New York, Columbia University Press. Sites J.W. and Marshall J.C. 2003. Delimiting species: a Renaissance issue in systematic biology. Trends in Ecology & Evolution 18(9): 462-470. Sites J.W. and Marshall J.C. 2004. Operational criteria for delimiting species. Annual Review of Ecology Evolution and Systematics 35: 199-227. General introduction 63 Suárez-Díaz E. and Anaya-Muñoz V.H. 2008. History, objectivity, and the construction of molecular phylogenies. Studies in History and Philosophy of Science Part C: Studies in History and Philosophy of Biological and Biomedical Sciences 39(4): 451-468. Teixell A. Ayarza P. Zeyen H. Fernàndez M. and Arboleya M.L. 2005. Effects of mantle upwelling in a compressional setting: the Atlas Mountains of Morocco. Terra Nova 17: 456-461. Templeton A.R. 1989. The meaning of species and speciation: A genetic prespective in Speciation and its Consequences by D. Otte and J. A. Endler. Sunderland, Sinauer, 3-27. Tolley K.A. Tilbury C.R. Measey G.J. Menegon M. Branch W.R. and Matthee C.A. 2011. Ancient forest fragmentation or recent radiation? Testing refugial speciation models in chameleons within an African biodiversity hotspot. Journal of Biogeography 38(9): 1748-1760. Townsend T.M. Larson A. Louis E. and Macey J.R. 2004. Molecular phylogenetics of Squamata: The position of snakes, Amphisbaenians, and Dibamids, and the root of the Squamate tree. Systematic Biology 53(5): 735-757. Underwood G. 1954. On the classification and evolution of geckos. Proceedings of the Zoological Society of London 124(3): 469-492. Vanvalen L. 1976. Ecological Species, Multispecies, and Oaks. Taxon 25(2-3): 233-239. Volobouev V. Pasteur G. Bons J. Guillaume C.P. and Dutrillaux B. 1990. Sex chromosome evolution in reptiles: divergence between two lizards long regarded as sister species Lacerta vivipara and Lacerta andreanskyi. Genetica 83(1): 85-91. Wake D.B. 2006. Problems with species: Patterns and processes of species formation in salamanders. Annals of the Missouri Botanical Garden 93(1): 8-23. Whittaker R.J. Araujo M.B. Paul J. Ladle R.J. Watson J.E.M. and Willis K.J. 2005. Conservation Biogeography: assessment and prospect. Diversity and Distributions 11(1): 3-23. Wiens J.J. 1999. Polymorphism in systematics and comparative biology. Annual Review of Ecology and Systematics 30: 327-362. Wiens J.J. 2004. Speciation and ecology revisited: Phylogenetic niche conservatism and the origin of species. Evolution 58(1): 193-197. Wiley E.O. 1978. Evolutionary Species Concept Reconsidered. Systematic Zoology 27(1): 17-26. Wilkins J.S. 2009. The first biological species concept : Evolving Thoughts. Retrived (2009-09-24). Yang Z. and Rannala B. 2010. Bayesian species delimitation using multilocus sequence data. PNAS 107(20): 9264-9269. Zakharov E.V. Lobo N.F. Nowak C. and Hellmann J.J. 2009. Introgression as a likely cause of mtDNA paraphyly in two allopatric skippers (Lepidoptera: Hesperiidae). Heredity 102(6): 590-599. Zhang D.X. and Hewitt G. 2003. Nuclear DNA analysis in genetic studies of population: practice, problems and prospects. Molecular ecology 12: 563-584. General introduction 64 65 CHAPTER 2. CRYPTIC DIVERSITY IN QUEDENFELDTIA SPP.(BOETTGER, 1883) Mafalda Barata, Morocco, 2008 ARTICLE 1. Barata M. Perera A. Martínez-Freiría F. and Harris D.J. 2012. Cryptic diversity within the Moroccan endemic day geckos Quedenfeldtia (Squamata: Gekkonidae): a multidisciplinary approach using genetic, morphological and ecological data.Biological Journal of Linnean Society. 106(4): 828-850. 66 67 Cryptic diversity within the Moroccan endemic day geckos Quedenfeldtia (Squamata: Gekkonidae): a multidisciplinary approach using genetic, morphological and ecological data MAFALDA BARATA1,2,3* ANA PERERA1, FERNANDO MARTÍNEZ-FREIRÍA1and D. JAMES HARRIS1,2 1CIBIO-UP, Centro de Investigação em Biodiversidade e Recursos Genéticos, Campus Agrário de Vairão, 4485-661 Vairão, Portugal 2Departamento de Biologia, Faculdade de Ciências da Universidade do Porto, 4099-002 Porto, Portugal 3CSIC-UPF, Institute of Evolutionary Biology, Passeig Marítim de la Barceloneta, 37-49, E - 08003 Barcelona, Spain Abstract Quedenfeldtia (Boettger, 1883) is a genus of diurnal geckos, endemic to the Atlas Mountains in Morocco, with two species being recognized: Quedenfeldtia moerens and Quedenfeldtia trachyblepharus.Quedenfeldtia moerens is found across a wide variety of habitats, from sea level to 3000 m a.s.l., whereas Q. trachyblepharus occupies exclusively high mountain regions reaching up to 4000 m a.s.l. This differentiation, offers an interesting model for study biogeographical patterns and evolutionary scenarios in a North African endemic. Analysis of two mitochondrial (12S rRNA and ND4) and four nuclear (ACM4, MC1R, PDC, and Rag1) DNA markers revealed high genetic variation, consistent with other recent phylogeographical studies, and with the two currently described species. However, within each species, a subdivision into two groups with geographical consistence was found. Multivariate morphological analyses confirmed the existence of two main phenotypes, whereas ecological niche modelling identified various environmental variables associated with the distribution of each species, and helped to predict occurrences outside the confirmed ranges. The results obtained in the present study indicate the possible existence of additional ‘cryptic’ species within this genus, a condition found in many North African reptiles, and particularly common in geckos. In general, North African montane fauna appears to reflect the occurrence of diverse palaeoendemics, as seen in Central Africa Mountain systems, rather than the pattern of recent postglacial recolonization observed in Europe. © 2012 The Linnean Society of London, Biological Journal of the Linnean Society, 2012, 106(4), 828–850. Keywords: Atlas Mountains – cryptic species – ecological modelling – mitochondrial DNA – morphology – North Africa – nuclear DNA – Quedenfeldtia moerens –Quedenfeldtia trachyblepharus. Cryptic diversity in Quedenfeldtia spp. 74 2.2.2. Phylogenetic analysis Mitochondrial genes (12S rRNA and ND4 with flanking tRNAs) and nuclear gene fragments (Rag1, ACM4, MC1R and PDC) were concatenated in two independent datasets and the most appropriate evolutionary model for each gene was calculated using jModeltest (Posada 2008) under the Akaike information criteria. Two different phylogenetic approaches, Bayesian inference (BI) and Maximum Likelihood (ML) were implemented to analyse the phylogenetic relationships within the genus. The BI analyses for the mitochondrial and nuclear datasets were performed using MrBayes v.3.1.2 (Huelsenbeck and Ronquist 2001) with partitions defined by gene, codon position for ND4, and without partitions. Partitioned Bayesian analysis was performed following Brandley et al. (2005), and Bayesian Factors (BF) were used to identify the best partitioning strategy for the concatenated dataset. All analyses started with randomly generated trees and ran for 2 million generations, saving one tree every 1000 generations. Two independent replicates were conducted and inspected for consistency to check for local optima (Huelsenbeck and Bollback 2001). Tracer v.1.4 (Drummond and Rambaut 2007) was used to verify if stationarity had been reached, both in terms of likelihood scores and parameters (ESS – Effective sample size - higher than 200 were considered acceptable), and the first 500 trees were discarded. A majority-rule consensus tree was generated from the remaining 19500 trees. ML analysis was executed using RaxML version 7.0.4 (Stamatakis 2006) for concatenated data under the GTR+G+I model. Bootstrapping (1000 pseudo-replicates) was used to evaluate the stability of nodes of the phylogenetic trees (Felsenstein 1985) for ML analysis. Uncorrected distances (p-distance), between and within clades were calculated in Mega v.3.0 (Kumar et al. 2004). Any heterozygotes present in the nuclear data were treated as ambiguous (N) in all analyses, since we sequenced only a few samples from each lineage and the analysis of nuclear genes was based on concatenated sequences. 2.2.3. Species delimitation Bayesian species delimitation was conducted using the program Bayesian Phylogenetics and Phylogeography, BPP v.2.0b (Rannala and Yang 2003; Yang and Rannala 2010) using the four nuclear loci. This method accommodates both species phylogeny and lineage sorting due to ancestral polymorphism. A gamma prior G(1, 10) was used on the population size parameters (qs). The age of the root in the species tree (t0) was assigned the gamma prior G (1, 10), while the other divergence time parameters were assigned the Dirichlet prior (Yang and Rannala 2010): equation 2). We used algorithm 0 with the finetuning parameter = 15.0 in order to ensure adequate rjMCMC mixing. This involves specifying a reversible jump algorithm to achieve dimension matching between species delimitation models with different Cryptic diversity in Quedenfeldtia spp. 75 number of parameters. Each species delimitation model was assigned equal prior probability and each analysis was run at least twice to confirm consistency between runs. The guide tree plays an important role in the result of the species delimitation model (Leaché and Fujita 2010); therefore we used the guide tree: ((Qt_Oukaimeden; Qt_JSirwa); (Qm_North; Qm_South)) based on the estimate of relationships from both mtDNA and nDNA trees (see results). 2.3. Morphological analysis 2.3.1. Variables quantification In total, 111 adult males (M) and 95 adult females (F) from five different localities were analyzed: Oukaimeden (23M, 23F), Jebel Sirwa (22M, 21F), Agoudal (17M, 12F), Ida-ouBuzia (28M, 21F) and Tafraoute (21M, 18F). Seven linear measurements and eight pholidotic characters including the ones identified by Arnold (1990) for their taxonomic value, and five colour pattern variables were taken from each individual. Details of the variables collected can be found in Appendix 1. All linear measurements were recorded in the field by the same person (AP), to the nearest 0.01 mm, in order to minimize measurement error. Pholidotic and colour pattern variables were retrieved from digital pictures by the same person (MB) at least twice and the mean value was recorded. Pholidotic variables were treated as continuous when they had more than five categories (UPLAB, SUBLAB, SNEY, PRECL and LAM, Appendix 1), while variables with a lower number (DS, SLC_1 and SLC_2) were considered as categorical. All bilateral variables were taken from the right side of the animal when possible. In cases where data could not be collected due to member amputations or poor quality of pictures (3.4%, n = 85) they were replaced by the group mean. 2.3.2. Statistical analysis Body measurements were log-transformed and checked for homoscedasticity (Levene's test) and normality (Shapiro-Wilks test) assumptions. Since linear measurements were highly correlated to the body size measured as snout-vent length (SVL, Pearson correlation in all cases p < 0.01), we used an isometric correction (Somers 1986) to separate the relative contribution of size and shape to the total variation analysed. Thus, the isometric body size (mSIZE) was used as a multivariate representation of size, whereas the combination of the remaining isometric-size corrected variables were used as representation of shape (Kaliontzopoulou et al. 2010). In order to obtain a general approximation regarding the level of differentiation between sexes, species and populations, we performed multivariate analysis of the variance Cryptic diversity in Quedenfeldtia spp. 76 (MANOVA) on shape (represented by the iso-corrected linear measurements) and pholidosis (only continuous variables) using SPECIES, SEX, POPULATION, nested in SPECIES, and their interaction, as factors. Since normality and homoscedasticity assumptions were not always met, non-parametric permutational MANOVAs were used (Anderson 2001) using the function adonis, implemented in the R package Vegan (Oksanen et al. 2012), and based on 1,000 permutations of the Euclidean distance matrix between group means. The same nonparametric procedure was used in univariate analysis of the variance (ANOVA) to estimate differences in SIZE and all continuous variables using the same factor design. Regarding categorical variables (colour pattern and non-continuous pholidotic variables), we plotted them in order to visualize their relative frequencies, and differences between SPECIES and POPULATION were compared using chi-squared tests. In order to investigate the source of morphological variation in Quedenfeldtia at a multivariate level, the morphological relationships among populations and the relative contribution of each dataset to their variability, we performed a Principal Component Analysis (PCA) on linear measurements (iso-corrected variables) and pholidosis (only continuous ones). For colour pattern, we used a Multiple Correspondence Analysis (MCA) to assess, which characters contributed most to the differentiation between populations. Finally, since we were interested in investigating the existence of further morphological differences between the two current Quedenfeldtia species besides the general colour pattern and pholidotic characteristics described by (Arnold 1990), we performed Canonical Discriminant Function Analysis (CDFA) on linear measurements (using both SIZE and SHAPE datasets) and continuous pholidotic variables. This allowed assessment at a multivariate level of which variables were the major contributors to the differentiation between species, and to create canonical discriminant functions to calculate the probability of classifying correctly the individuals based on them. We used the leave-one-out option to cross-validate the classification results. Since this procedure (Jacknife prediction) generates individual classifications using discriminate functions based on all observations except the given case, it provides a more accurate estimate of the classification values. All statistical analyses were performed using SPSS v.17.0 and R (R Development Core Team,(2009). Significance level was considered at p < 0.05. 2.4. Ecological niche modelling analyses 2.4.1. Data sampling A total of 92 and 35 localities for Q. moerens and Q. trachyblepharus, respectively, were gathered from published (Bons and Geniez 1996; Harris et al. 2010; Barata et al. 2011) and unpublished data (Table 1) and were initially considered for ecological modelling purposes. Cryptic diversity in Quedenfeldtia spp. 77 However, in order to avoid spatial bias, and thus decrease the level of spatial autocorrelation in species presence, a randomly process of removing localities from clusters of species occurrence was performed (for details see(Brito et al. 2011). The Nearest Neighbour Index was used to assess the degree of data clustering, along of this process, finally obtaining values of 0.97 and 0.89 for Q. moerens and Q. trachyblepharus, respectively, which indicate a random distribution for both species. Spatial analyses were done with the “Spatial Analyst” extension of ArcGIS (ESRI 2006).Following this, a total of 44 and 20 no clustered localities for Q. moerens and Q. trachyblepharus, respectively, were used for running ecological models (model samples), whereas 48 and 15 localities for Q. moerens and Q. trachyblepharus, respectively, were used for a secondary test of ecological models (validation samples). 2.4.2. Statistical analysis Choosing a correct extension for the study area is a complex but necessary requirement for ecological niche modelling (VanDerWal et al. 2009; Anderson and Raza 2010). In presencebackground-modelling techniques, such as MaxEnt (Phillips et al. 2006), different extensions may vary the nature of pseudoabsences randomly chosen by the algorithm, affecting the model output (VanDerWal et al. 2009). Therefore, calibration exercises with different sizes of the study area were performed (VanDerWal et al. 2009; Anderson and Raza 2010), and finally an area of 150 km of buffer from all localities was chosen to develop ecological models (Fig. 1). Subsequently, predictions were projected to a larger area of 837,034 km2 covering all of Morocco and adjacent areas of Algeria, Western Sahara and Mauritania (Fig. 1). Three sets of environmental and ecogeographical variables (hereafter EGV) were selected for the ecological models according to their meaning to the ecology and distribution of other reptiles (Brito et al. 2008; Kaliontzopoulou et al. 2008; Martinez-Freiria et al. 2008; Brito et al. 2011). These sets included one topographical grid (USGS, 2006) that was used to derive the variable Slope, with the “Slope” function of ArcGIS; 19 climate grids (Hijmans et al. 2005); and a land cover grid from the years 2004-2006 (Bicheron et al. 2008). In order to convert the categorical land cover EGV into a continuous variable, one binary grid was created for each habitat type that covered more than 2% of the study area. The Euclidean distance of each grid cell to the closest habitat-type was calculated for each individual habitat grid (13 habitat types) using the “Euclidean Distance” tool of ArcGIS (Brito et al. 2011). Finally, spatial correlation among all EGVs was tested using the “Band Collection Statistics” tool of ArcGIS and only 11 EGVs (r < 0.750) were chosen for modelling purposes (Appendix 2). Cryptic diversity in Quedenfeldtia spp. 78 Ecological niche models were developed with the MaxEnt v.3.3 software (Phillips et al. 2006). The Maximum Entropy modelling approach requires only presence data as input, performs well in comparison to other methods, and has been used successfully in ecological niche-based modelling of many species distributions (Elith et al. 2006; Elith et al. 2011). Model samples and EGVs were imported into MaxEnt and a total of 40 and 20 model replicates for Q. moerens and Q. trachyblepharus, respectively, were run in the buffer area with random seed and 80% training / 20% testing data partition in each run. Samples for each replicate were chosen by bootstrap allowing sampling with replacement. Models were run with auto-features (Phillips et al. 2006), and the Area under the Curve (AUC) of the receiveroperating characteristics (ROC) plots was taken as a measure of individual model fit (Fielding and Bell 1997). Finally, each model replicate was projected to the projection area and the individual model replicates in the projection area were added to generate a mean forecast of probability of species occurrence (Marmion et al. 2009). Standard deviation between individual model probabilities of presence was used as an indication of prediction uncertainty (Buisson et al. 2010; Carvalho et al. 2010). Mean models were reclassified according to the minimum training presence thresholds given by MaxEnt, which includes all presence data in suitable cells, to display areas of probable absence and presence for each species. To evaluate model quality, the model and validation samples were intersected with the threshold models to calculate the percentage of correct classification of presences in each probability category for each species. This technique was also used to assign a classification to the “indeterminate” Quedenfeldtia samples derived from Bons and Geniez (1996) (n = 26). Identification of areas of probable sympatry between both species was determined by the overlap of threshold models using the “Raster Calculator” function of ArcGIS. The importance of each EGV for explaining the distribution of the species was determined by its average percent contribution to the models (Brito et al. 2008; Martinez-Freiria et al. 2008; Brito et al. 2011). The relation between occurrence of the species and EGVs was determined by the visual examination of response curves profiles from univariate models (Phillips et al. 2006). Similar profiles for a given EGV were taken as an indication of parallel relationships between the occurrence of these species and the range of variation of the EGV (Brito et al. 2008; Martinez-Freiria et al. 2008; Brito et al. 2011). This indicates also the possible occurrence of sympatry and eventual competition within the range of values of the EGV equally selected by both species. Conversely, a distinct profile among species was taken as an indication of divergent relationships and possible allopatry. Finally, niche overlap and similarity between both species was tested using ENMtools v.1.3 (http://enmtools.blogspot.com). First, niche overlap was quantified from models generated for the two species using Schoener’s D (Schoener 1968), I statistic (Warren et al. 2008) and Cryptic diversity in Quedenfeldtia spp. 79 relative rank metrics (RR,(Warren and Seifert 2011), which range from 0 (no overlap) to 1 (total overlap). The overlaps for I and D are calculated by taking the difference between species in suitability scores at each grid cell, after suitabilities are standardized so that they sum to 1 over the geographic space being measured (Warren et al. 2008; 2010; Warren and Seifert 2011). The relative rank statistic is an estimate of the probability that the relative ranking of any two patches of habitat is the same for the two models, irrespective of the quantitative difference in suitability estimates (Warren and Seifert 2011). Then, an identity test, using 50 randomized pseudoreplicates, was carried out to test the hypothesis of whether models generated from the two species are more different than expected if they were drawn from the same underlying distribution (Warren et al. 2008; 2010). By comparing the observed values of D, I and RR to the null distribution obtained using the identity test; it is possible to ask whether models produced by the two species are statistically significantly different (Warren et al. 2010). 3. Results 3.1. Genetic analysis 3.1.1. Mitochondrial DNA The mitochondrial data consisted of 12S rRNA (407 bp), ND4 and tRNAs (692 bp) fragments totalling 1099 bp. A total of 42 haplotypes were retrieved from the 42 individuals of Quedenfeldtia, 16 from Q. trachyblepharus and 26 from Q. moerens. Best-fit models of nucleotide substitution for each gene and the concatenated fragment were: 12S rRNA - GTR+G; ND4 1st position - HKY+G; ND4 2nd position - HKY; ND4 3rd position - HKY+G, tRNAs - K80+G, and concatenated data - GTR+I+G. The division of the data into five partitions had an improvement effect in the -lnLk (BF = 10.596; BF > 10 in favour of the hypothesis being tested (Kass and Raftery 1995; Brandley et al. 2005) and thus the partitioned analysis was preferred to estimate the phylogeny. ML and BI analyses were mostly congruent and most estimates of relationships were well supported (Fig. 2A). The analyses divided the samples in two main groups concordant with the two current species, Q. moerens and Q. trachyblepharus. Quedenfeldtia moerens was further divided in two geographically congruent groups. One well-supported lineage (BI – 1 and ML - 92) grouped the samples from the South, and the other included the populations from the North, although the latter was paraphyletic with respect to Q. moerens from the South. The lineage that included the samples of Q. moerens from the South was further divided into two sub-lineages corresponding to the central and southern populations sampled. Quedenfeldtia trachyblepharus was also divided in two lineages, one that included all Cryptic diversity in Quedenfeldtia spp. 80 samples from Oukaimeden and Toubkal and a second that included the remaining samples, all from the Jebel Sirwa area. All lineages exhibited high differentiation between them, with 13% (p-distance) between species and 8% and 9.5% between sub-lineages of Q. moerens and Q. trachyblepharus, respectively. Interestingly, there was low diversity within groups, up to a maximum of circa 3%, although sample size was limited (8 to 16 individuals) 3.1.2. Nuclear DNA For the combined four nuclear genes from 23 individuals, 13 from Q. moerens and 10 from Q. trachyblepharus were analysed (Fig. 2B). The four partial genes Rag1 (471 bp), ACM4 (399 bp), MC1R (690 bp) and PDC (372 bp) were concatenated, giving a total length of 1932 base pairs. Best-fit models of nucleotide substitution for each gene fragment were: Rag1 - HKY+G; ACM4 - HKY+I; MC1R - HKI+I+G and PDC - K80+I, and for the concatenated data - GTR+I+G. ML and BI analysis using nuclear markers were mostly congruent (Fig. 2B). The analyses of individual nuclear fragments (not shown) were not discordant but some did not resolve the two groups of Q. moerens (ACM4) or Q. trachyblepharus (Rag1) and one of the genes (PDC) showed the same polyphyly in the Q. moerens groups, as the estimate of relationships based on mtDNA tree. Quedenfeldtia moerens samples were divided in two lineages that coincide with the mtDNA results but the two groups were monophyletic. However, for the nDNA there was no differentiation within the Southern populations. Further, there was no support for the monophyly of Q. trachyblepharus. Moreover, results differ from the mitochondrial analysis, since individuals from Toubkal were well differentiated from Oukaimeden and almost identical to the remaining Q. trachyblepharus localities. The nuclear divergence (p-distance) between the two species where 1.5%, and between the lineages, 0.5% between Quedenfeldtia trachyblepharus from Oukaimeden and J. Sirwa and 0.8% between Quedenfeldtia moerens from North and South. 3.1.3. Bayesian species delimitation When assuming four lineages, Bayesian species delimitation analysis (BPP) strongly supports the guide tree, as found in previous analyses (four clades), with speciation probability ≥ 0.99 on all nodes (Guide tree with posterior probability for presence of nodes: ((Qt_Oukaimeden; Qt_JSirwa)#0.991; (Qm_North; Qm_South)#1.0)#1.0) We used only one guide tree, the one obtained in both mtDNA and nDNA phylogenetic analysis, since this resulting tree was very simple with no ambiguous relationships. Following Leaché and Fujita (2010) the use of random trees, with an artificial increase of sister species, has a negative impact on the result. Cryptic diversity in Quedenfeldtia spp. 81 Figure 2. Trees derived from a Bayesian partitioned analysis for the mitochondrial DNA (12S rRNA, ND4 and tRNA’s) on the left (A) and nuclear sequences (Rag1, ACM4, MC1R and PDC) on the right (B), using the models described in the text. Bayesian posterior probabilities (0-1) and bootstrap values (> 50%) for ML (1-100) are indicated near the branches. Between the trees are grouped the two species (Q. trachyblepharus and Q. moerens) and the clades found inside each species (Oukaimeden, J. Sirwa, South, North and in mtDNA the South clade is divided in two smaller groups, represented by the sampled population, Tafraoute and Ida-ou-Bouzia). Cryptic diversity in Quedenfeldtia spp. 82 3.2. Morphological analysis Details on the basic statistics of the two species by locality and sex can be found in Appendix 2. 3.2.1. Continuous variables Our results from the MANOVAs showed general differences between SPECIES, SEX and POPULATION in size and shape, although interaction between factors was not always significant (Table 2). SEX was the main factor explaining the variation in size (R2= 0.27, Table 2), while SPECIES and POPULATION were the most explicative factors on shape (Table 2). Regarding differentiation among body measurements (using size-corrected variables), males had wider heads (HW), while females exhibited a larger trunk length (ILL) (Table 2). The degree of sexual dimorphism (interaction SEX*SPECIES) was, in general, similar in the two species (except for HW and FLL, Table 2) although at the intraspecific level (interaction SEX*POPULATION, the latter nested to SPECIES) only differences in multivariate size were observed (Table 2). Regarding the morphological differentiation between the two species, Q. moerens showed significantly larger multivariate size and longer limbs (HLL and FLL) and head (HL), while Q. trachyblepharus exhibited wider heads (HW) and longer trunks (ILL; Table 2). Moreover, there were also differences between populations at the intraspecific level in all the characters analysed (Table 2 and Appendix 1). For pholidosis, the results indicated differences between SPECIES, SEX and POPULATION, as well as their interactions (Table 2), although the most contributing factor to the variation in pholidosis was POPULATION (R2= 0.18, Table 2). In fact, populations differed in all the pholidotic characters analyzed, except for the number of supralabial scales (UPLAB; Table 2). Such differences were also evident at the species level (Table 2). Differences in morphological and pholidotic traits were reflected in a relative grouping of the individuals in populations when analysing variation at the multivariate level using a PCA (Table 3 and Fig. 3A). Individuals from the two species appeared differentiated across the first two axes, although a high overlap was observed, especially between females of Q. trachyblepharus and Q. moerens from the South-Central locality of Ida-ou-Bouzia (Fig. 3A). Individuals from the other two Q. moerens populations from North (Agoudal) and South (Tafraoute) were better differentiated (Fig. 3A). Variation across the first axes was mostly explained by shape-related variables in males and females, including limbs length (both FLL and HLL), trunk length (ILL) and head width (HW, only in females) the most contributing ones (Table 3). However, number of lamellae scales (LAM) was also an important factor for the variation observed across PC1 in both sexes. Regarding PC2 and PC3, size contributed more across PC2 in males and PC3 in females, while SNEY was also important to explain the Cryptic diversity in Quedenfeldtia spp. 83 Table 2. Results of the non-parametric (M) ANOVAs on the effect of species, sex and population (as factor nested in species) and their interactions on size (isometric size), shape (remaining iso-corrected linear measurements) and pholidosis (only continuous). For each factor and dataset analysed, degrees of freedom (df) and results of the R2, F-value and level of significance (p) are provided. Analyses were based on 999 permutations. See materials and methods for more details. Numbers in bold indicates statistically significant values. species sex pop sex*species sex*pop Residuals Total df 1 1 3 1 3 197 206 MANOVA R20.16 0.27 0.17 0.00 0.03 0.37 1.00 F88.06 143.26 30.28 2.61 4.81 p 0.001 0.001 0.001 0.114 0.006 R20.18 0.02 0.21 0.01 0.01 0.58 1.00 F62.30 5.14 23.68 1.74 1.30 p 0.001 0.001 0.001 0.121 0.198 R20.09 0.02 0.18 0.03 0.03 0.66 1.00 F25.46 4.77 18.12 7.57 2.87 p 0.001 0.001 0.001 0.001 0.002 ANOVA R20.38 0.01 0.09 0.00 0.01 0.50 1.00 Body F 152.30 3.96 12.38 0.54 1.95 measurements p 0.001 0.054 0.001 0.464 0.113 R20.32 0.02 0.10 0.01 0.00 0.55 1.00 F112.34 6.96 11.65 4.58 0.18 p 0.001 0.022 0.001 0.033 0.914 R20.02 0.01 0.07 0.00 0.03 0.87 1.00 F4.61 2.74 4.92 0.64 2.08 p 0.034 0.092 0.001 0.423 0.094 R20.01 0.00 0.32 0.00 0.01 0.66 1.00 F2.87 0.37 32.25 0.83 0.75 p0.080 0.526 0.001 0.336 0.508 R20.12 0.06 0.19 0.01 0.01 0.61 1.00 F39.72 18.89 20.83 1.74 1.36 p 0.001 0.001 0.001 0.195 0.244 R20.22 0.00 0.29 0.01 0.01 0.47 1.00 F90.53 0.07 40.93 4.39 1.50 p 0.001 0.790 0.001 0.036 0.235 R20.21 0.00 0.39 0.00 0.01 0.39 1.00 F107.35 2.14 66.62 0.01 0.93 p 0.001 0.153 0.001 0.928 0.426 ANOVA UPLAB R20.02 0.03 0.01 0.03 0.04 0.88 1.00 Pholidosis F 4.06 6.22 0.50 6.23 2.61 p 0.041 0.013 0.683 0.020 0.058 SNEY R20.01 0.01 0.16 0.06 0.00 0.76 1.00 F3.69 2.68 13.45 15.41 0.06 p0.062 0.098 0.001 0.002 0.978 PRECL R20.10 0.04 0.14 0.05 0.09 0.57 1.00 F35.06 14.51 16.69 18.87 10.18 p 0.001 0.001 0.001 0.001 0.001 SUBLAB R20.02 0.00 0.23 0.01 0.01 0.73 1.00 F5.05 0.38 20.49 3.11 0.61 p 0.03 0.56 0.00 0.07 0.60 LAM R20.27 0.01 0.28 0.00 0.01 0.43 1.00 F123.46 3.66 42.15 0.16 1.77 p 0.00 0.06 0.00 0.70 0.15 HH ILL FLL HLL size shape pholidosis SVL HW HL Cryptic diversity in Quedenfeldtia spp. 90 species complex. One lineage (Oukaimeden) comprises the samples from only two localities, Oukaimeden and Toubkal, both in the High Atlas and separated by 10 km from each other. The other clade (J. Sirwa), grouped the remaining specimens of Q. trachyblepharus from the eastern High Atlas (Aguelmous), and an isolated mountain in the southern High Atlas (Jebel Sirwa and El Jazib n-Triri).However, nuclear markers only separated Oukaimeden (only two sequences available) from the remaining localities, with samples from Toubkal being very similar to the samples from the J. Sirwa clade. This discrepancy between the markers again highlights the limitations of using only mtDNA for phylogeographic studies (Godinho et al. 2008), and also seems to confirm that the complexity of these species may reflect their “melting pot” characteristics (Canestrelli et al. 2010; Canestrelli et al. 2012), although incomplete lineage sorting in the nuclear markers is also a possibility. Morphologically the two populations were considered Q. trachyblepharus (Arnold 1990), suggesting a conservative morphology. The analysis of colour patterns (MCA) did not reveal differences between them, although some biometric traits such as head dimensions or trunk length discriminate partially both populations, but only in males. Such discrimination, as also seen in Q. moerens, may be indicative of different selective pressures on sexual dimorphism. We also found in most of the individuals from Oukaimeden the presence of a double scale under the fourth hind limb toe, which was only found occasionally in J. Sirwa. It is necessary to assess this character in other populations to determine how well it corresponds with the two genetically distinct forms identified. The main four lineages, two within each species, are genetically very different from each other which can be a result of a long-term isolation while the variation found between them suggests that they diverge approximately at the same time. Following the time frame proposed by Gamble et al. (2010), this could be as long ago as 15-17 Myr, and would therefore precede the formation of the Atlas mountains 9.0 Myr (Gomez et al. 2000; Babault et al. 2008). Thus, these species can reasonably be considered as “palaeoendemics”, mirroring the situation in central African mountains much more closely than the patterns observed in Iberian montane species (Mouret et al. 2011; Tolley et al. 2011). The Bayesian species delimitation using the nuclear markers supports the existence of the four clades with high speciation probability (≥ 0.99%), indicating the possible presence of cryptic species. The magnitude of divergence between the subclades (8% in Q. moerens and 9.5% in Q. trachyblepharus, ND4), was similar to that found in other groups that were suggested to be cryptic species (Brown et al. 2002; Perera and Harris 2010). Quedenfeldtia trachyblepharus showed high levels of divergence in a geographically restricted area, yet very low variation within each clade (1%). This seems to suggest that the current distribution is the result of a rapid expansion from very small populations or a Cryptic diversity in Quedenfeldtia spp. 91 bottleneck effect. Assessment of additional populations would be needed to confirm this, but again it seems to indicate that rather than a pattern of expansion from a single primary refugium, as is often observed in European herpetofauna (eg(Rowe et al. 2006), in these species many small but distinct refugia existed during the last glaciations leading to maintenance of high genetic diversity between lineages but with limited variation within them. Previous studies (Bons and Geniez 1996; Schleich et al. 1996; Sindaco and Jeremcenko 2008) showed that while Q. moerens had a more widespread distribution from the south, almost from the desert to the mountains in central East Morocco, Q. trachyblepharus was restricted to high mountain altitudes of the High Atlas. However, many records for these species were considered undetermined, particularly in the High Atlas (Bons and Geniez 1996); Fig. 1). Results from ecological niche models identified areas for the occurrence of both species mostly in concordance with previous studies (Bons and Geniez 1996) and successfully assigned a probable classification for most of the indeterminate localities (classification > 76%). Moreover, ecological models showed suitable cells for the occurrence of Q. moerens in areas where the species was not reported. Historical reasons (e.g expansion from the south and/or barriers to dispersal) could have prevented Q. moerens from reaching these areas, but fieldwork should be conducted to investigate its possible presence. This is particularly important given the case of Podarcis lizards in North Africa, where models indicated possible presence in extra-range regions that were later confirmed after prospection (Kaliontzopoulou et al. 2008). Probable areas for the occurrence of both species were almost allopatric and environmental variables related to the distribution of both species were different or similar but with mostly different profiles in the response curves, which also suggests allopatry. For instance, Q. moerens occur in flat and dry areas with temperatures below 25ºC during summer and close to sparse vegetation whereas Q. trachyblepharus occurs in sloping and humid areas with high levels of precipitation during winter. Nevertheless, there are specific environmental conditions where both species could coexist, with a large area of habitat in the High Atlas that could be suitable for both species. However, no regions of sympatry have been found in the field to date. Possibly Q. trachyblepharus is isolated as a result of competition with Q. moerens, since Q. moerens is limiting the south and north limits of Q. trachyblepharus, and according to Bons and Geniez (1996) is present between Oukaimeden and J. Sirwa, where there seems to be some kind of geographic barrier between these two clades. Probably these two species split and then adapted to different environmental conditions, temperate and drier climate for Q. moerens and high and humid mountains in Q. trachyblepharus. Following this hypothesis, when the temperature increased after the last glacial maxima Q. moerens Cryptic diversity in Quedenfeldtia spp. 92 expanded its distribution upwards, forcing Q. trachyblepharus to move to higher altitudes, where long-term isolation coupled with a reduction in population sizes (bottlenecks) might have promoted the high inter-population and low intra-population differentiation observed today. The presence of the two species in the central area of the distribution, between Q. moerens from the South and Q. trachyblepharus from Oukaimeden, indicates this could be where the species originated. 5. Conclusions The extreme genetic variability found within the genus Quedenfeldtia, including within the two species, is more than the level used to distinguish species in various other groups. However, geckos generally show high genetic variability (Harris et al. 2004a; Rato and Harris 2008), and species definition is still a controversial issue (de Queiroz 2007). Indications are that both currently accepted species might mask additional units that could deserve recognition as full species. Within Q. moerens subgroups show more phenotypic variation, probably related to the different selective pressures resulting from the wider range of habitats they occupy. On the other hand Q. trachyblepharus, despite the limited distribution, show higher genetic diversity but low morphological differentiation between populations. Our results indicate that, for this genus, phylogeographic patterns are much more similar to those recovered from montane species in Central Africa – highly divergent palaeoendemics - rather than the predominantly genetically uniform pattern observed in European montane reptiles. Some evidence of a “melting pot” scenario in which mtDNA patterns do not fully coincide with nuclear markers highlights the need for further assessments prior to any alterations to taxonomy. However, assessment of morphological data indicates some characters that might be useful to discriminate genetically divergent lineages in the field. In the Pyrenean rock lizard, Iberolacerta bonnali, phylogeographic patterns tended to reflect recolonization history rather than current habitat (Mouret et al. 2011). Our modelling approach indicates that, at least at a two-taxon level, patterns in Quedenfeldtia can be associated with habitat. It may be that within these two major lineages, genetic patterns better reflect refugial areas, but this can only be tested by sampling extensively across the region. It will be important in the future therefore to increase the sampling of Q. trachyblepharus to other high mountain localities, as well as investigate the possible existence of contact areas between the two species, in order to assess in more detail the evolutionary history of this group. Cryptic diversity in Quedenfeldtia spp. 93 Acknowledgements MB is supported by the FCT grant SFRH/BD/41488/2007, AP is supported by the FCT grant SFRH/BPD/26546/2006 and FMF is supported by the FCT grant SFRH/BPD/69857/2010. This work was funded by FCT grant PTDC/BIA-BEC/105327/2008. Thanks to all the colleagues from CIBIO that helped with fieldwork in Morocco, and to H. El Mouden, T. Slimani and M. Znari for their help. All individuals were captured under the permit of the Haut Commisariat aux Eaux and Forets of Morocco (HCEFLCD/DLCDPN/DPRN/DFF Nº 14/2010). We thank four anonymous reviewers for their helpful comments. References Agapow M. 2005. Species: demarcation and diversity in Phylogeny and Conservation by A. Purvis, J. L. Gittleman and T. Brooks. Cambridge, UK, Cambridge University Press, 57–75. Anderson M.J. 2001. A new method for non-parametric multivariate analysis of variance. Austral Ecology 26(1): 32-46. Anderson R.P. and Raza A. 2010. The effect of the extent of the study region on GIS models of species geographic distributions and estimates of niche evolution: preliminary tests with montane rodents (genus Nephelomys) in Venezuela. Journal of Biogeography 37(7): 13781393. Arévalo E. Davis S.K. and Sites J.W. 1994. Mitochondrial-DNA Sequence Divergence and Phylogenetic-Relationships among 8 Chromosome Races of the Sceloporus-Grammicus Complex (Phrynosomatidae) in Central Mexico. Systematic Biology 43(3): 387-418. Arnold E.N. 1990. The two species of Moroccan day-geckoes, Quedenfeldtia (Reptilia: Gekkonidae). Journal of Natural History 24(3): 757-762. Babault J. Teixell A. Arboleya M.L. and Charroud M. 2008. A late cenozoic age for longwavelength surface uplift of the atlas mountains of Morocco. Terra Nova 20(2): 102-107. Barata M. Harris D.J. and Castilho R. 2008. Comparative phylogeography of northwest African Natrix maura (Serpentes: Colubridae) inferred from mtDNA sequences. African Zoology 43(1): 1-7. Barata M. Perera A. Harris D.J. Van Der Meijden A. Carranza S. Ceacero F. García-Muñoz E. Gonçalves D. Henriques S. Jorge F. Marschall J.C. Pedrajas L. and Sousa P. 2011. New observations of amphibians and reptiles in Morocco, with a special emphasis on the eastern region. Herpetological Bulletin 116: 4-14. Bauer A.M. de Silva A. Greenbaum E. and Jackman T. 2007. A new species of day gecko from high elevation in Sri Lanka, with a preliminary phylogeny of Sri Lankan Cnemaspis (Reptilia, Squamata, Gekkonidae). Zoosystematics and Evolution 83(S1): 22-32. Bicheron P. Defourny P. Brockmann C. Schouten L. Vancutsem C. Huc M. Bontemps S. Leroy M. Achard F. Herold M. Ranera F. and Arino O. 2008. Globcover: Products Description and Validation Report., Medias-France and Postel. Blackburn D.C. and Measey G.J. 2009. Dispersal to or from an African biodiversity hotspot? Molecular ecology 18(9): 1904-1915. Boettger O. 1883. Die reptilien und amphibien von Marocco II. Travaux de lnstitut Scientifique Chdrifien, Sbrie Zoologique 18: 1-130. Bons J. 1959. Les lacertiliens du Sud-Ouest Marocain, travaux de institut. Scientifique Cherifien, Série Zoologique 18: 1-130. Bons J. and Geniez P. 1996. Amphibiens et reptiles du Maroc (Sahara Occidental compris) Atlas Biogéographique. Barcelone, Asociación Herpetologica Espanola. Brandley M.C. Schmitz A. and Reeder T.W. 2005. Partitioned Bayesian Analyses, Partition Choice, and the Phylogenetic Relationships of Scincid Lizards. Systematic Biology 54(3): 373-390. Brito J.C. Fahd S. Geniez P. Martinez-Freiria F. Pleguezuelos J.M. and Trape J.F. 2011. Biogeography and conservation of viperids from North-West Africa: An application of ecological niche-based models and GIS. Journal of Arid Environments 75(11): 1029-1037. Brito J.C. Santos X. Pleguezuelos J.M. and Sillero N. 2008. Inferring evolutionary scenarios with geostatistics and geographical information systems for the viperid snakes Vipera latastei and Vipera monticola.Biological Journal of the Linnean Society 95(4): 790-806. Cryptic diversity in Quedenfeldtia spp. 94 Brown R.P. Suárez N.M. and Pestano J. 2002. The Atlas mountains as a Biogeografical divide in North-West Africa: evidence from mtDNA evolution in the Agamid lizard Agama impalearis. Molecular Phylogenetics and Evolution 24: 9 (324-332). Buisson L. Thuiller W. Casajus N. Lek S. and Grenouillet G. 2010. Uncertainty in ensemble forecasting of species distribution. Global Change Biology 16(4): 1145-1157. Busack S.D. 1986. Biogeographic Analysis of the Herpetofauna Separated by the Formation of the Strait of Gibraltar. National Geographic Research 2: 17-36. Canestrelli D. Aloise G. Cecchetti S. and Nascetti G. 2010. Birth of a hotspot of intraspecific genetic diversity: notes from the underground. Molecular ecology 19(24): 5432-5451. Canestrelli D. Sacco F. and Nascetti G. 2012. On glacial refugia, genetic diversity, and microevolutionary processes: deep phylogeographical structure in the endemic newt Lissotriton italicus. Biological Journal of the Linnean Society 105(1): 42-55. Carranza S. Romano A. Arnold E., N. and Sotgiu G. 2008. Biogeography and evolution of European cave salamanders, Hydromantes (Urodela: Plethodontidae), inferred from mtDNA sequences. Journal of Biogeography 35: 724-738. Carvalho S.B. Brito J.C. Crespo E.J. and Possingham H.P. 2010. From climate change predictions to actions - conserving vulnerable animal groups in hotspots at a regional scale. Global Change Biology 16(12): 3257-3270. de Queiroz K. 2007. Species Concepts and Species Delimitation. Systematic Biology 56(6): 879-886. Drummond A.J. and Rambaut A. 2007. BEAST: Bayesian evolutionary analysis by sampling trees. BMC Evolutionary Biology 7(1): 214. Elith J. Graham C.H. Anderson R.P. Dudik M. Ferrier S. Guisan A. Hijmans R.J. Huettmann F. Leathwick J.R. Lehmann A. Li J. Lohmann L.G. Loiselle B.A. Manion G. Moritz C. Nakamura M. Nakazawa Y. Overton J.M. Peterson A.T. Phillips S.J. Richardson K. Scachetti-Pereira R. Schapire R.E. Soberon J. Williams S. Wisz M.S. and Zimmermann N.E. 2006. Novel methods improve prediction of species' distributions from occurrence data. Ecography 29(2): 129-151. Elith J. Phillips S.J. Hastie T. Dudik M. Chee Y.E. and Yates C.J. 2011. A statistical explanation of MaxEnt for ecologists. Diversity and Distributions 17(1): 43-57. ESRI 2006. ArcMap 9.2. Redlands, California, USA, Environmental Systems Research Institute, Inc. Felsenstein J. 1985. Confidence-Limits on Phylogenies - an Approach Using the Bootstrap. Evolution 39(4): 783-791. Fielding A.H. and Bell J.F. 1997. A review of methods for the assessment of prediction errors in conservation presence/absence models. Environmental Conservation 24(1): 38-49. Fjeldsa J. and Lovett J.C. 1997. Geographical patterns of old and young species in African forest biota: the significance of specific montane areas as evolutionary centres. Biodiversity and Conservation 6(3): 325-346. Fonseca M.M. Brito J.C. Paulo O.S. Carretero M.A. and Harris D.J. 2009. Systematic and phylogeographical assessment of the Acanthodactylus erythrurus group (Reptilia: Lacertidae) based on phylogenetic analyses of mitochondrial and nuclear DNA. Molecular Phylogenetics and Evolution 51(2): 131-142. Fonseca M.M. Brito J.C. Rebelo H. Kalboussi M. Larbes S. Carretero M.A. and Harris D.J. 2008. Genetic variation among spiny-footed lizards in the Acanthodactylus pardalis group from North Africa. African Zoology 43(1): 8-15. Gamble T. Bauer A.M. Colli G.R. Greenbaum E. Jackman T.R. Vitt L.J. and Simons A.M. 2010. Coming to America: multiple origins of New World geckos. Journal of Evolutionary Biology 24(2): 231-244. Gamble T. Bauer A.M. Greenbaum E. and Jackman T.R. 2008a. Evidence for Gondwanan vicariance in an ancient clade of gecko lizards. Journal of Biogeography 35(1): 88-104. Gamble T. Bauer A.M. Greenbaum W. and Jackman T.R. 2008b. Out of the blue: a novel, transAtlantic clade of geckos (Gekkota, Squamata). Zoologica Scripta 37(4): 355-366. Godinho R. Crespo E.G. and Ferrand N. 2008. The limits of mtDNA phylogeography: complex patterns of population history in a highly structured Iberian lizard are only revealed by the use of nuclear markers. Mol Ecol 17(21): 4670-4683. Gomez F. Beauchamp W. and Barazangi M. 2000. Role of the Atlas Mountains (northwest Africa) within the African-Eurasian plate-boundary zone. Geology 28(9): 775-778. Harris D.J. Batista V. and Carretero M.A. 2004a. Assessment of genetic diversity within Acanthodactylus erythrurus (Reptilia : Lacertidae) in Morocco and the Iberian Peninsula using mitochondrial DNA sequence data. Amphibia-Reptilia 25(2): 227-232. Cryptic diversity in Quedenfeldtia spp. 95 Harris D.J. Batista V. Carretero M.A. and Ferrand N. 2004b. Genetic variation in Tarentola mauritanica (Reptilia : Gekkonidae) across the Strait of Gibraltar derived from mitochondrial and nuclear DNA sequences. Amphibia-Reptilia 25(4): 451-459. Harris D.J. Carretero M. Perera A. Perez-Mellado V. and Ferrand N. 2003. Complex patterns of genetic diversity within Lacerta (Teira)perspicillata: Preliminary evidence from 12S rRNA sequence data. Amphibia-Reptilia 24(3): 386-390. Harris D.J. Perera A. Barata M. Tarroso P. and Salvi D. 2010. New distribution notes for terrestrial herpetofauna from Morocco. North-Western Journal of Zoology 6(2): 309-315. Hijmans R.J. Cameron S.E. Parra J.L. Jones P.G. and Jarvis A. 2005. Very high resolution interpolated climate surfaces for global land areas. International Journal of Climatology 25(15): 1965-1978. Huelsenbeck J.P. and Bollback J.P. 2001. Empirical and hierarchical Bayesian estimation of ancestral states. Systematic Biology 50(3): 351-366. Huelsenbeck J.P. and Ronquist F. 2001. MrBayes: Bayesian inference of phylogeny. Bioinformatics 17: 754-755. Kaliontzopoulou A. Brito J.C. Carretero M.A. Larbes S. and Harris D.J. 2008. Modelling the partially unknown distribution of wall lizards (Podarcis) in North Africa: ecological affinities, potential areas of occurrence, and methodological constraints. Canadian Journal of ZoologyRevue Canadienne De Zoologie 86(9): 992-1001. Kaliontzopoulou A. Carretero M.A. and Llorente G.A. 2010. Intraspecific ecomorphological variation: linear and geometric morphometrics reveal habitat-related patterns within Podarcis bocagei wall lizards. Journal of Evolutionary Biology 23(6): 1234-1244. Kass R.E. and Raftery A.E. 1995. Bayes Factors. Journal of the American Statistical Association 90(430): 773-795. Katoh K. Misawa K. Kuma K. and Miyata T. 2002. MAFFT: a novel method for rapid multiple sequence alignment based on fast Fourier transform. Nucleic Acids Research 30(14): 30593066. Kluge A.G. 1987. Cladistic relationships in the Gekkonoidea (Squamata, Sauria): Zoology. University of Michigan, Museum of UMMZ. Kocher T.D. Thomas W.K. Meyer A. Edwards S.V. Paabo S. Villablanca F.X. and Wilson A.C. 1989. Dynamics of mitochondrial-DNA evolution in animals - Amplification and sequencing with conserved primers. Proceedings of the National Academy of Sciences of the United States of America 86(16): 6196-6200. Kumar S. Tamura K. and Nei M. 2004. MEGA3: Integrated software for molecular evolutionary genetics analysis and sequence alignment. Briefings in Bioinformatics 5(2): 150-163. Leaché A.D. and Fujita M.K. 2010. Bayesian species delimitation in West African forest geckos (Hemidactylus fasciatus). Procedings of Royal Society B 277(1697): 3071-3077. Loveridge A. 1947. Revision of the African lizards of the family Gekkonidae.Bulletin of the Museum of Comparative Zoology 98: 1-469. Marmion M. Parviainen M. Luoto M. Heikkinen R.K. and Thuiller W. 2009. Evaluation of consensus methods in predictive species distribution modelling. Diversity and Distributions 15(1): 59-69. Martinez-Freiria F. Sillero N. Lizana M. and Brito J.C. 2008. GIS-based niche models identify environmental correlates sustaining a contact zone between three species of European vipers. Diversity and Distributions 14(3): 452-461. Medail F. and Quezel P. 1999. The phytogeographical significance of SW Morocco compared to the Canary Islands. Plant Ecology 140(2): 221-244. Mouret V. Guillaumet A. Cheylan M. Pottier G. Ferchaud A.L. and Crochet P.A. 2011. The legacy of ice ages in mountain species: post-glacial colonization of mountain tops rather than current range fragmentation determines mitochondrial genetic diversity in an endemic Pyrenean rock lizard. Journal of Biogeography 38(9): 1717-1731. Myers N. Mittermeier R.A. Mittermeier C.G. da Fonseca G.A.B. and Kent J. 2000. Biodiversity hotspots for conservation priorities. Nature 403(6772): 853-858. Oksanen J. Blanchet F.G. Kindt R. Legendre P. Minchin P.R. O’Hara R.B. Simpson G.L. Solymos P. M. H. H. Stevens and Wagner H. 2012. Vegan: Community Ecology Package. R package version 1.17-7. http://CRAN.R-project.org/package=vegan. Perera A. and Harris D.J. 2010. Genetic variability within the Oudri’s fan-footed gecko Ptyodactylus oudrii in North Africa assessed using mitochondrial and nuclear DNA sequences. Molecular Phylogenetics and Evolution 54: 634-639. Cryptic diversity in Quedenfeldtia spp. 96 Phillips S.J. Anderson R.P. and Schapire R.E. 2006. Maximum entropy modeling of species geographic distributions. Ecological Modelling 190(3-4): 231-259. Pinho C. Harris D.J. and Ferrand N. 2007. Contrasting patterns of population subdivision and historical demography in three western Mediterranean lizard species inferred from mitochondrial DNA variation. Molecular ecology 16(6): 1191-1205. Pinho C. Rocha S. Carvalho B.M. Lopes S. Mourao S. Vallinoto M. Brunes T.O. Haddad C.F.B. Goncalves H. Sequeira F. and Ferrand N. 2010. New primers for the amplification and sequencing of nuclear loci in a taxonomically wide set of reptiles and amphibians. Conservation Genetics Resources 2: 181-185. Posada D. 2008. jModelTest: Phylogenetic model averaging. Molecular Biology and Evolution 25(7): 1253-1256. R D.C.T. 2009. R: A language and environment for statistical computing. R Foundation for Statistical Computing. Vienna, Austria, available at: http://www.R-project.org. Rannala B. and Yang Z. 2003. Bayes estimation of species divergence times and ancestral population sizes using DNA sequences from multiple loci Genetics 164: 1645-1656. Rato C. Carranza S. Perera A. Carretero M.A. and Harris D.J. 2010. Conflicting patterns of nucleotide diversity between mtDNA and nDNA in the Moorish gecko, Tarentola mauritanica.Molecular Phylogenetics and Evolution 56(3): 962-971. Rato C. and Harris D.J. 2008. Genetic variation within Saurodactylus and its phylogenetic relationships within the Gekkonoidea estimated from mitochondrial and nuclear DNA sequences. Amphibia-Reptilia 29(1): 25-34. Rowe G. Harris D.J. and Beebee T.J.C. 2006. Lusitania revisited: a phylogeographic analysis of the natterjack toad Bufo calamita across its entire biogeographical range. Molecular Phylogenetics and Evolution 39: 335–346. Sambrook J. Fritsch E.F. and Maniatis T. 1989. Molecular cloning: a laboratory manual.3nd edt. New York, Cold Sring Harbor Laboratory Press. Sattler T. Bontadina F. Hirzel A.H. and Arlettaz R. 2007. Ecological niche modelling of two cryptic bat species calls for a reassessment of their conservation status. Journal of Applied Ecology 44(6): 1188-1199. Schleich H.H. Kastle W. and Kabisch K. 1996. Amphibians and Reptiles from North Africa. Königstein, Germany., Koeltz Scientific Publications. Schoener T.W. 1968. Anolis lizards of Bimini: resource partitioning in a complex fauna. Ecology Letters 49: 704-726. Sindaco R. and Jeremcenko V.K. 2008. The reptiles of the Western Palearctic. Annotated checklist and distributional Atlas of the turtles, crocodiles, amphisbaenians and lizards of Europe, North Africa, Middle East and Central Asia, 580 S, Edition Belvedere, Monografie della Societas Herpetologica Italica. Somers K.M. 1986. Multivariate Allometry and Removal of Size with Principal ComponentsAnalysis. Systematic Zoology 35(3): 359-368. Stamatakis A. 2006. RAxML-VI-HPC: Maximum likelihood-based phylogenetic analyses with thousands of taxa and mixed models. Bioinformatics 22(21): 2688-2690. Tolley K.A. Tilbury C.R. Measey G.J. Menegon M. Branch W.R. and Matthee C.A. 2011. Ancient forest fragmentation or recent radiation? Testing refugial speciation models in chameleons within an African biodiversity hotspot. Journal of Biogeography 38(9): 1748-1760. VanDerWal J. Shoo L.P. Graham C. and William S.E. 2009. Selecting pseudo-absence data for presence-only distribution modeling: How far should you stray from what you know? Ecological Modelling 220(4): 589-594. Vidal N. and Hedges S.B. 2004. Molecular evidence for a terrestrial origin of snakes. Proceedings of the Royal Society B-Biological Sciences 271: S226-S229. Warren D.L. Glor R.E. and Turelli M. 2008. Environmental Niche Equivalency Versus Conservatism: Quantitative Approaches to Niche Evolution. Evolution 62(11): 2868-2883. Warren D.L. Glor R.E. and Turelli M. 2010. ENMTools: a toolbox for comparative studies of environmental niche models. Ecography 33(3): 607-611. Warren D.L. and Seifert S.N. 2011. Ecological niche modeling in Maxent: the importance of model complexity and the performance of model selection criteria. Ecological Applications 21(2): 335-342. Weiss S. and Ferrand N. 2007. Current perspectives in phylogeography and the significance of South European refugia in the creation and maintenance of European Biodiversity in Phylogeography of Southern European Refugia by S. Weiss and N. Ferrand. Dordrecht, Netherlands, Springer, 341-357. Cryptic diversity in Quedenfeldtia spp. 97 Yang Z. and Rannala B. 2010. Bayesian species delimitation using multilocus sequence data. PNAS 107(20): 9264-9269. APPENDIX Table A1. Morphological, pholidotic and colour patterns variables included in the study. 1-Body measurements: SVLsnout-vent length, from the tip of the snout until the cloaca. HLhead length, from the tip of the snout to the posterior ear cavity. HWtotal head width at its widest part at the level of the temporal region. HHhead height from occipital to jaws ILLinter-limb length from the posterior edge of forelimb insertion to the anterior edge of hindlimb insertion. HLLhind-limb length from the longest toe to the base of the hindlimb. FLLforelimb length from the longest toe to the base of the forelimb. 2-Scales: Continuous variables UPLABnumber of upper labial scales on the right side until the limit of the mouth opening counted on the right side of the head. SUBLABnumber of sub-labial scales until the limit of the mouth opening counted on the right side of the head. SNEYnumber of linear scales between the eye and the nostril counted on the right side of the head. PRECLnumber of "pre-cloacal" scales, counted from leg axis to cloaca LAMnumber of non-divided enlarged side to side lamellae on the fourth right hind toe Discrete variables DSpresence of a divided subdigital scale on the fourth hind toe: no (0), yes (1) SLC_1contact between nostril and first supralabial scale: no (0), yes (1) SLC_2supralabial scale fused with postnasal or rostral scale: no (0), yes (1) 3-Colouration: Ventral VHBGular colourwhite (0), yellow (1) VHP-Head pigmentationwhite (0), Black (1) VBB-Trunk ventral pigmentationwhite (0), yellow (1) Dorsal DB-Colourationbrown uniform (0), two colours (1), uniform with dots (2) DO-Ocelliwithout (0), 1 to 3 (1), more than three (2) Cryptic diversity in Quedenfeldtia spp. 98 Males (N =28) Females (N=21) Males (N =17) Females (N=12) Males (N =21) Females (N=19) Males (N=22) Females (N=21) Males (N=23) Females (N=23) SVL 47.00 ± 1.71 44.45 ± 2.09 44.26 ± 1.06 43.46 ± 1.31 47.45 ± 2.58 43.13 ± 2.64 45.29 ± 2.14 43.07 ± 1.89 47.54 ± 1.61 45.8 ± 1.04 43.00 – 50.00 40.00 – 49.00 42.00 – 46.00 41.00 – 46.00 42.00 – 52.50 39.00 – 47.50 41.00 – 48.50 40.00 – 46.50 44.00 – 50.0 44.00 – 48.00 HL 9.38 ± 0.34 8.74 ± 0.35 8.51 ± 0.36 8.14 ± 0.34 9.33 ± 0.50 8.51 ± 0.31 8.97 ± 0.35 8.64 ± 0.62 9.26 ± 0.34 8.88 ± 0.31 8.77 - 10.13 7.95 - 9.28 7.91 - 9.04 7.83 - 8.99 8.26 - 10.12 7.96 - 9.25 8.38 - 9.61 7.62 - 9.56 8.41 - 9.82 8.27 - 9.46 HW 11.08 ± 0.59 10.32 ± 0.61 10.32 ± 0.35 9.90 ± 0.39 11.14 ± 0.59 10.06 ± 0.59 10.25 ± 0.47 9.40 ± 0.39 10.28 ± 0.40 9.82 ± 0.36 9.92 - 12.05 9.07 - 11.12 9.65 - 11.23 9.19 - 10.49 9.86 - 12.06 8.95 - 11.06 9.30 - 10.94 8.76 - 10.08 9.33 - 10.85 9.12 - 10.44 HH 5.97 ± 0.20 5.65 ± 0.28 5.30 ± 0.27 5.08 ± 0.20 5.58 ± 0.36 5.14 ± 0.35 5.20 ± 0.33 4.87 ± 0.39 5.47 ± 0.23 5.22 ± 0.26 5.52 - 6.34 5.08 - 6.15 4.95 - 5.95 4.79 - 5.44 4.95 - 6.41 4.24 - 5.53 4.64 - 5.76 4.17 - 5.56 4.97 - 5.84 4.87 – 5.90 ILL 21.06 ± 1.46 20.75 ± 1.88 18.84 ± 1.01 19.20 ± 0.75 21.05 ± 1.17 19.50 ± 1.44 19.34 ± 1.29 19.41 ± 0.89 21.43 ± 1.08 21.46 ± 1.40 18.64 – 24.40 16.77 - 23.22 17.26 - 20.99 18.03 - 20.42 19.10 - 23.79 17.01 – 22.10 16.98 - 21.65 18.23 – 21.50 19.00 – 23.60 19.11 - 24.62 FLL 19.10 ± 0.83 17.95 ± 1.02 19.09 ± 0.52 18.36 ± 0.45 21.20 ± 1.26 19.08 ± 1.05 17.64 ± 1.00 16.59 ± 0.71 18.14 ± 1.09 16.61 ± 0.92 17.69 - 20.68 16.00 - 19.86 17.99 – 19.80 17.62 - 18.99 18.50 - 23.23 17.20 - 20.97 15.94 - 19.58 14.96 - 18.05 16.32 - 20.24 14.90 - 18.23 HLL 25.87 ± 1.24 23.47 ± 1.37 25.29 ± 0.85 23.69 ± 1.17 29.28 ± 1.86 25.95 ± 0.69 23.91 ± 1.08 22.05 ± 1.20 23.38 ± 1.24 21.72 ± 1.34 23.22 - 28.39 20.48 - 26.41 23.50 - 26.37 22.26 - 26.37 25.49 - 32.29 24.67 - 27.04 21.44 - 25.36 19.52 - 24.45 20.25 – 25.10 18.98 - 24.27 UPLAB 7.11 ± 0.63 6.57 ± 0.68 6.73 ± 0.66 6.89 ± 0.67 7.24 ± 0.7 6.58 ± 0.61 6.62 ± 0.65 6.67 ± 0.55 6.74 ± 0.62 6.73 ± 0.45 6 - 8 6 - 8 5 - 8 6 - 8 6 - 9 6 - 8 5 - 8 6 - 8 6 - 8 6 - 7 LOLAB 5.04 ± 0.43 4.48 ± 0.6 4.4 ± 0.47 4.78 ± 0.38 4.71 ± 0.64 4.37 ± 0.5 5.1 ± 0.61 5.06 ± 0.38 5.09 ± 0.42 5.27 ± 0.54 4 – 6 4 - 6 4 - 5 4 - 5 4 - 6 4 - 5 4 - 6 4 - 6 4 - 6 4 - 6 SNEY 11.88 ± 0.96 11.25 ± 0.94 10.82 ± 0.69 10.33 ± 0.55 11.17 ± 0.57 10.67 ± 0.51 10.53 ± 0.70 10.86 ± 0.43 11.00 ± 0.85 11.27 ± 1.05 10 – 14 10 - 13 9 - 12 9 - 11 10 - 12 10 - 12 9 - 12 10 - 12 10 - 13 10 - 14 PRECL 11.14 ± 1.38 11.00 ± 0.95 11.25 ± 1.44 10.56 ± 0.75 11.53 ± 1.07 12.11 ± 1.29 14.38 ± 1.46 11.55 ± 0.74 11.6 ± 0.63 11.37 ± 0.62 9 – 14 9 - 13 9 - 15 9 - 12 10 - 13 10 - 15 11 - 17 10 - 13 10 - 13 10 - 13 SUBLAB 7.79 ± 0.92 7.57 ± 0.75 6.94 ± 0.66 6.80 ± 0.83 6.40 ± 1.11 6.00 ± 0.75 6.89 ± 0.87 7.35 ± 0.95 7.45 ± 1.16 7.43 ± 1.07 6 – 9 6 - 9 5 - 8 5 - 8 5 - 9 5 - 7 5 - 9 6 - 10 5 - 9 5 - 9 LAM 9.32 ± 0.65 9.06 ± 0.67 10.31 ± 1.04 10.11 ± 0.67 9.95 ± 0.59 9.74 ± 0.56 9.53 ± 0.65 9.06 ± 0.86 7.91 ± 0.79 8.05 ± 0.47 8 – 11 8 - 10 9 - 13 9 - 11 9 - 11 9 - 11 8 - 11 8 - 11 7 - 10 7 - 9 DS 1.00 ± 0.00 1.00 ± 0.00 1.00 ± 0.00 1.00 ± 0.00 1.00 ± 0.00 1.06 ± 0.23 1.11 ± 0.29 1.12 ± 0.3 1.82 ± 0.39 1.9 ± 0.29 1 – 1 1 - 1 1 - 1 1 - 1 1 - 1 1 - 2 1 - 2 1 - 2 1 - 2 1 - 2 SLC 3.54 ± 0.88 3.43 ± 0.75 4.19 ± 0.73 4.4 ± 0.63 3.25 ± 0.62 3.06 ± 0.23 3.72 ± 0.81 4 ± 0.84 3.32 ± 0.55 3.24 ± 0.51 3 – 6 3 – 5 3 – 5 3 – 5 3 – 5 '3 - 4 3 - 5 3 - 5 3 - 5 3 - 5 Q. moerens Q. trachyblepharus Ida-ou-Bouzia Agoudal Tafraoute Jebel Sirwa Oukaimeden Table A2. Descriptive statistics for all the linear measurements and pholidotic variables of the different Quedenfeldtia localities analysed in this study. For each variable, mean ± standard deviation, minimum and maximum and sample size (N) per sex is shown. Cryptic diversity in Quedenfeldtia spp. 99 Table A3. Environmental factors used to model the distribution of Q. moerens and Q. trachyblepharus and their codes, units and range of variation. FINAL NOTE This article is formatted accordingly to the Journal where it was published, Biological Journal of Linnean Society. Code Name Units Range TAR Temperature Annual Range ºC 10.9 – 42.7 MTDQ Mean Temperature of Driest Quarter ºC 8.6 – 34.0 MTWQ Mean Temperature of Warmest Quarter ºC 10.4 – 36.2 PS Precipitation Seasonality (Coefficient of Variation) adimensional 19 - 122 PWQ Precipitation of Warmest Quarter mm 0 - 110 PCQ Precipitation of Coldest Quarter mm 2 - 573 Distance to closed to open (>15%) mixed broadleaved 0 – 15.561*105 and needleleaved forest (>5m) DVEG Distance to sparse (<15%) vegetation m 0 - 301024 DBAR Distance to bare areas m 0 - 69180 DWAT Distance to water bodies m 0 - 502573 SLOPE Slope (derived from altitude) degrees 0 – 31.42 DFOR m 106 Extreme genetic diversity in Atlantolacerta 107 Background An emerging pattern among European biotas is that the accentuated environmental instability that occurred during the Pleistocene did not lead to increased speciation rates, with many species and populations originating during the Miocene and proceeding through the Quaternary [1-2]. In many species, population fragmentation was triggered by the beginning of the Messinian Salinity Crisis, a short (600 000 years) but crucial period that occurred between 5.9 and 5.3 Mya during which the Mediterranean Sea desiccated almost completely, producing a general and drastic increase in aridity around the Mediterranean Basin [3-4]. As a result of this increased aridity, forests continued to be replaced by more open and arid landscapes forcing the mesic species to retreat to the moister Atlantic-influenced areas and to the mountainous regions, leading to high speciation in some groups [5-6]. Various studies have attempted to unravel the different roles that the global aridification at the end of the Miocene and the Pleistocene glacial cycles have played in the diversity and distribution of European faunas [7]. However, little is known about the effects that these climatic changes had on species living further South, in the African continent. Recent assessments of central African chameleons have uncovered evidence of long-isolated evolutionary histories, with the survival of palaeoendemics leading to considerable diversity [8]. In general, reptiles are excellent model organisms to assess the relative role that the PreQuaternary and Quaternary major climatic events have played in the origin, evolution and distribution of species [9]. Available data from some herpetofauna indicate that a similar pattern to the neighboring Iberian Peninsula exists in North Africa, with deep lineages originating at the end of the Miocene (Chalcides [10], Acanthodactylus [11-13], Podarcis [2, 14-15], Saurodactylus [16], Ptyodactylus [17], Salamandra [18], Pleurodeles [19]). However, the lack of informative nuclear markers in most of these studies may prevent the recovery of the true evolutionary history of the group [20-21], and makes it difficult to ascertain if these lineages correspond to species complexes or not. Since there is a strong likelihood of discordance between gene trees and species trees [22-24], information from different genetic markers (mitochondrial and nuclear) is thus necessary for delimiting evolutionary lineages, as well as for establishing phylogenetic relationships. Despite being key concepts in the fields of systematic and evolutionary biology, recognizing and delimiting species are highly controversial issues [e.g. 25-26]. Recognizing species is not only a taxonomic challenge, but is also essential for other biological disciplines such as biogeography, ecology and evolutionary biology [27], and has serious consequences for conservation biology and the design of effective conservation plans [28-29]. Delimiting species is also the first step towards discussing broader questions on evolution, biogeography, ecology or conservation. Recently, thanks to intellectual progress made in the field with the Extreme genetic diversity in Atlantolacerta 108 aim of identifying a common element among all the different species concepts, a single, more general, concept of species known as General Lineage Species Concept has been suggested [30]. This unified species concept emphasizes the common element found in many species concepts, which is that species are separately evolving lineages. Therefore, properties like reciprocal monophyly at one or multiple loci, phenotypic diagnosability, ecological distinctiveness, etc. are not part of the species concept but are used to assess the separation of lineages and to species delimitation [31]. This separation between species conceptualization and species delimitation and the proposal of a unified species concept has concentrated efforts in the development of new approaches for species delimitation, as for example with “integrative taxonomy” [32-33,]among others]. Under this new approach, species delineation is regarded as an objective scientific process that results in a taxonomic hypothesis. Therefore, the level of confidence in the taxonomic hypothesis supported by several independent character sets is much higher than for species supported by only one character [34]. Such an integrative view is especially useful in the case of taxonomic groups that are morphologically conservative, where cryptic species have probably been overlooked [17, 3536]. Normally, high altitude species carry signatures of the expansion and contraction cycles occurred during glacial and interglacial periods [37-39]. Because of this, they are of particular interest to study historical responses to climate change, since they are adapted to a small window of environmental changes, and usually present low tolerance to high temperatures [40]. In Europe, high altitude species often seem to have persisted through glacial periods by short movements to lower altitudes rather than to the classic "southern refugia" of lowland species. In this way current ranges may primarily reflect postglacial expansions [41]. However, it is not clear if the same phenomenon occurs in African montane taxa. Atlantolacerta andreanskyi (Werner, 1929) is a lacertid lizard endemic to the western and central parts of the High Atlas Mountains of Morocco. It is restricted to areas above 2400 m [42-43], where it is frequently found in the vicinity of small watercourses or plateaus in the top of the mountains that retain some water from rain or snowmelt. Habitat is normally screes and areas with boulders, meadows and, in particular, the base of cushion-like thorny plants in these places [42;]personal observation]. Although A. andreanskyi had initially been placed in several different genera within the subtribe Lacertini [44-48], recent phylogenetic analyses based on mitochondrial DNA and a combination of mitochondrial and nuclear markers [4950] suggest that A. andreanskyi is a member of the subtribe Eremiadini, and apparently sister to the remaining Eremiadini. This position would conform to this species lacking the synapomorphies that characterize most other Eremiadini, namely a derived condition of the ulnar nerve and the presence of a fully developed armature in the hemipenis, which has folded Extreme genetic diversity in Atlantolacerta 109 lobes when retracted. It is also distinctive within the Eremiadini regarding the presence of enlarged masseteric scale [49]. Because of its phylogenetic position, without close relationship to any other genus of Eremiadini and its distinctive morphology it was recently placed in a new monotypic genus, Atlantolacerta [49]. Atlantolacerta andreanskyi is distributed across 440 Km (straight line) of mountainous terrain, with the different populations presenting an apparently disjunct distribution [42-43;]see Fig. 1. As with many montane species, the situation observed in A. andreanskyi is similar to an archipelago, with the different “islands” being represented by mountaintops disconnected due to areas of unsuitable habitat below 2400 m. As a result of this scenario, minimal gene flow is currently expected between the different populations; however, it is not known how the different climatic events occurred during the Miocene and Pleistocene have affected this species. Even though some aspects of the biology of A. andreanskyi are already well known e.g.[51-52], the genetic structure of the different populations, as well as the relationships between the different populations have never been assessed before. Therefore, in order to shed some light on the previous questions and attempt to assess the evolutionary history of the species and identify the number of lineages, we sampled the distribution area of the species and performed several combined phylogenetic reconstructions and clustering analyses, using both mtDNA and nuclear markers. Figure 1.Atlantolacerta andreanskyi distribution map. The colour dots represent the localities of the populations sampled for this work. The white dots represent the distributions of the species by Bons and Geniez (1996).J. Awlime (1 yellow), Oukaimeden (2 red), Toupkal (3 orange), J. Sirwa (4 pink), Tizin Tichka (5 dark blue), Outabati (6 light blue), J. Azourki (7 light green) and J. Ayache (8 dark green). Extreme genetic diversity in Atlantolacerta 110 Results Mitochondrial genealogies A total of 1108 base pairs (bp) of concatenated mtDNA (12S rRNA 330 bp, ND4 592 bp and tRNA-His 186 bp) were obtained for 89 A. andreanskyi. The concatenated alignment of the ingroup sequences revealed 30 haplotypes (3 from Tizin Tichka, 7 from J. Ayache, 5 from J. Sirwa, 2 from Oukaimeden, 7 from J. Azourki, 2 from Outabati, 2 from Toubkal and 2 from J. Awlime) and contained 241 variable sites, of which 232 were parsimony informative. Analyses of the concatenated mtDNA data were mostly congruent (Fig. 2A). Seven wellsupported lineages were recovered from these analyses (pp > 0.95 and BP > 70%), corresponding to the populations from J. Awlime, J. Sirwa, Tizin Tichka, J. Azourki, Outabati, J. Ayache, and Oukaimeden and nearby Toubkal that were grouped together. Regarding the relationships among these clades, we could distinguish three main groups, Oukaimeden and Toubkal with J. Sirwa from the southern end of the distribution range; J. Ayache with Outabati from the northern distribution, and Tizin Tichka with J. Azourki from the central distribution range. The population from J. Awlime, from the extreme South of the range, is a genetically distinct lineage related to the northern group, although, both ML and BI analysis weakly support this topology (see Fig. 2A). All the populations present a low level of diversity in the mitochondrial DNA (uncorrected genetic distances 0 - 0.5% for the ND4+tRNA-His and 0 – 0.2% for the 12S; see Table 1), and a very high level of genetic divergence between populations (5.5 – 16.5% in the ND4+tRNAHis and 2.5 – 6.6% in the 12S; see Table 1). Nuclear genealogies A total of 77 specimens of A. andreanskyi were sequenced for five nuclear genes. The ACM4 was 447 bp long, presenting 47 haplotypes and 34 polymorphic sites, 33 of them parsimony informative; C-MOS was 534 bp long, with 32 haplotypes and 21 polymorphic sites, all of them parsimony informative; MC1R was 635 bp long, with 57 haplotypes and 36 variable sites, 35 of them parsimony informative; PDC was 441 bp long, with 60 haplotypes and 29 variable sites, 26 of them parsimony informative; RAG1 was 528 bp long, with 38 haplotypes and 19 variable sites, 18 of them parsimony informative. The differences in the genetic distances between the lineages are congruent with the geographic distance between them, supporting the grouping of the lineages in three main groups as seen in the analysis of mitochondrial sequences. The concatenated analyses of the 5 unphased nuclear markers are congruent with the results obtained in the mitochondrial DNA tree, although with some differences (Fig. 2B). Despite recovering the three main groups observed in the mtDNA analysis, according to the nuclear Extreme genetic diversity in Atlantolacerta 111 markers the J. Awlime population is not sister to the northernmost populations but branches off inside a polytomy with the westernmost lineages at the base of the tree. It is possible to distinguish some of the lineages, although in some cases they are not monophyletic. The J. Ayache population is monophyletic but makes Outabati paraphyletic. The same happens with Tizin Tichka, which makes the population from J. Azourki paraphyletic. The population from Oukaimeden is polyphyletic. Concatenated analysis (mtDNA and nDNA) The results of the ML and BI analyses of the mtDNA and nDNA (Fig. 2C) support the same seven lineages as recovered in the mitochondrial analysis, although in this case J. Awlime is sister to the central and northern lineages (Tizin Tichka, J. Azourki, Outabati, and J. Ayache) instead of being sister to only the northernmost lineages (Fig. 2A). As in the mtDNA analysis (Fig. 2A), the relationship of J. Awlime with the central and northern lineages is very poorly supported. This result was expected, given the higher resolving power of the mtDNA that contributed with 241 variable sites versus the 150 from the nDNA. Nuclear networks As show in Fig. 3 and Table 2, there is a moderate degree of haplotype sharing between populations, with most f them lacking private alleles for the nuclear genes analysed. Extreme genetic diversity in Atlantolacerta 112 1/99 1/100 1/100 1/73 1/100 1/100 1/100 0.86/73 0.59 1/99 1/100 0.99/92 1/100 1/92 1/100 1/100 1/100 0.99/73 1/100 0.59 1/100 1/66 1/99 1/99 0.88/92 1/99 1/92 0.72 0.81 1 0.8 0.55 0.82 0.85 1 0.99 D A B C 1 Oukaimeden J. Sirwa J. Azourki Tizin Tichka J. Ayache Outabati 1 0.95 1 0.84 *2.9 (1.0-5.6) *4.3 (1.4-7.8) *2.4 (0.8-4.4) *7.6 (4.3-11.9) *6.4 (3.1-10.2) Figure 2. Trees resulting from partitioned Bayesian analysis. (A) Mitochondrial DNA tree (12S, ND4 and flanking tRNA-His), (B) nuclear concatenated tree (RAG1, ACM4, MC1R, PDC and C-MOS), (C) Concatenated tree from the combined mitochondrial and nuclear DNA data. The partitions used the models described in the text. Bayesian posterior probabilities (0–1) and bootstrap values (> 50%) for ML (1–100) are indicated near the branches, (D) Species tree from mitochondrial and nuclear DNA data from the Bayesian Inference of Species Trees (STARBEAST). Clade posterior probabilities are shown to the left of the nodes, and divergence times and 95% intervals (calculated in BEAST using only ND4 +tRNA-His), to the right of the nodes. The trees were rooted using Podarcis bocagei,P. hispanica and P. carbonelli. The colours represent the different populations. Clustering analysis and individual assignment In our study, the obtained K differs with the combination between the ancestry model and the allele frequency model. When combined the No Admixture Model (ancestry model) with the Allele Frequencies Independent Model (allele frequency model) the best resulting K values were for K = 3: South (Oukaimeden, J. Sirwa, Toubkal and J. Awlime), center (Tizin Tichka and J. Azourki) and North (Outabati and J. Ayache) groups. With the other three combinations between the models the best result were for K = 6: J. Sirwa, Tizin Tichka, J. Azourki, Outabati, J. Ayache, and a group formed by Oukaimeden, Toubkal and J. Awlime (Fig. 4). Extreme genetic diversity in Atlantolacerta 113 Table 1. Genetic distances and divergence time estimate between populations. (A) Genetic distance (12S and ND4 +tRNA-His) between all the populations and (B) between main groups; and (C) divergence time estimates, calculated using BEAST with ND4 and tRNA-His. The diversity of each population is below the population's names. A Pop Tizin Tichka Oukaimeden J.Sirwa J.Ayache J.Azourki Outabati Toubkal J.Awlime p-distance (%) 0.1 0.4 0.3 0.5 0.2 0 0.4 0.1 12S / ND4 Tizin Tichka 13.1 12.7 14.5 10.5 15.3 12.9 13.6 0 Oukaimeden 47.7 15 13.2 16.1 1.7 13.2 0 J.Sirwa 4.2 2.8 16.1 12.7 16.5 7.5 11.6 0.2 J.Ayache 5.4 5.7 4.8 12.7 5.5 14.4 14.1 0.1 J.Azourki 4.3 4.3 3.8 6.6 14.2 13.2 13.1 0.1 Outabati 5.4 5.7 4.2 1.6 6 16 14 0 Toubkal 3.7 0.3 2.5 5.4 4 5.4 12.6 0 J.Awlime 4 4.7 4.5 5.1 6.4 5 4.3 0 B C Pop JAy+Out Tiz+JAz J.Awlime Ouk+JSi+Tou p-distance (%) 0.9 2.3 0 1.5 Beast Ma (95% HPD) ND4 12S / ND4 JAy+Out 13.7 13.4 14.6 North-South 7.6 (4.3-11.9) 2.9 JAw - Ouk 5.6 (2.5-9.7) Tiz+JAz 5.9 12.9 12.4 JAy+Out - Tiz+JAz 6.4 (3.1-10.2) 5.2 Ouk - JSi 2.9 (1.0-5.6) J.Awlime 5 5.2 12.2 Tiz - JAz 4.3 (1.4-7.8) 0.1 JAy - Out 2.4 (0.8-4.4) Ouk+JSi+Tou 5.1 4.1 4.6 Ouk - Tou 0.5 (0.1-1.2) 0.4 Species tree and divergence time estimates The results of the clustering analysis with K = 6 were used to define the species for the species tree analysis in STARBEAST. The tree inferred with information from mitochondrial and nuclear markers (phased) (figure 2D) recovered the same topology as in Fig. 2C, with all the relationships between the lineages supported by previous analyses. The divergence time estimates were calculated for the six populations (Table 2). High effective sample sizes were observed for all parameters in all BEAST analysis (posterior ESS Extreme genetic diversity in Atlantolacerta 114 Figure 3. Parsimony networks corresponding to MC1R (A), RAG1 (B), C-MOS (C), ACM4 (D) and PDC (E) nDNA sequence variation from all the populations. The colours used were the same as the used in the map (Figure 1) and trees (Figure 2), J. Awlime (yellow), Toubkal (orange), Oukaimeden (red), J. Sirwa (pink), Tizin Tichka (dark blue), J. Azourki (light blue), Outabati (light green), and J. Ayache (dark green). Lines represent a mutation step, circles represent haplotypes and dots missing haplotypes. The size of the circles is proportional to the number of alleles. values > 1000 for all four analyses) and assessment of convergence statistics in Tracer indicated that all analyses had converged. Maximum clade credibility tree for ND4+tRNA-His was identical in topology to those produced by Bayesian and ML analyses. According to the inferred dates resulted from BEAST (Fig. 2D), the two main mitochondrial lineages of A. andreanskyi (South versus central and North) split approximately 7.6 Ma (95% high posterior density (HPD) interval 4.3-11.9 Ma). The populations that are grouped in the three main clades (South, central and North) split approximately at the same time, being Tizin Tichka and J. Azourki the first to split at about 4.3 Ma (1.4-7.8), followed by Oukaimeden and J. Sirwa 2.9 Ma (1-5.6), and Outabati and J. Ayache 2.4 Ma (0.8-4.4). Tizin Tichka and J. Azourki diverged from Outabati and J. Ayache approximately 6.4 Ma (3.1-10.2). Extreme genetic diversity in Atlantolacerta 115 Table 2. Percentage of private alleles in all the populations and for each nuclear locus. PrivateAlleles (%) MC1R RAG1 C-MOS ACM4 PDC J. Awlime 33 50 0 67 0 J. Sirwa 96 12 0 42 92 Toubkal 50 0 50 100 50 Oukaimeden 41 33 70 29 57 Tizin Tichka 75 59 23 71 100 J. Azourki 60 100 60 84 90 Outabati 100 54 43 9 83 J. Ayache 92 85 57 80 20 Figure 4. Population structure estimation. Each individual is represented by a thin vertical line, which is partitioned into K coloured segments that represent the individual’s estimated membership fractions in K clusters. The bigger vertical divisions separate individuals from different populations. Populations are labeled below the figure. The colours used are the same used in Figure 1 and Figure 2. Discussion Extreme mtDNA diversity in A. andreanskyi Several recently published analyses of North African herpetofauna have revealed high levels of endemism and cryptic species [12, 14-15, 17]. In this analysis, the surprising result was the extreme diversity of mitochondrial DNA found between almost all the populations analysed. The genetic differentiation observed between populations (2.8% - 6.6% in 12S and 5.5% - 16.5% in ND4+tRNA-His) is similar and, in some cases, higher than the divergence found between Iberolacerta species (between 7.4% and 8.2% in the cytochrome b gene, [53]), a lacertid genus with most of its species occurring in the mountains of the Iberian Peninsula [41, 54]. Initially considered one species, there are now seven recognized species of Iberolacerta in the Iberian Peninsula. Genetic differentiation between these species is lower than between the different populations of A. andreanskyi. Although the mitochondrial phylogeny supports the existence of seven distinct groups, the clustering analysis only supports the existence of six lineages (J. Sirwa, Tizin Tichka, J. Azourki, Outabati, J. Ayache and a lineage formed by Oukaimeden, Toubkal and J. Awlime). Extreme genetic diversity in Atlantolacerta 122 GenBank acession codes 3/5 Specimen code Alleles Population Latitude Longitude Altitude 12S / ND4+tRNA-His / PDC / ACM4 / C-MOS / MC1R / RAG1 2612 2612a Oukaimeden (Ouk) 31.20426 -7.86705 2600 JX462077 / JX462170 / JX461613 / JX461967 / JX485224 / JX461781 / JX461435 2612b JX461614 / JX461968 / JX485225 / JX461782 / JX461436 2615 2615a Oukaimeden (Ouk) 31.20426 -7.86705 2600 JX462078 / JX462171 / JX461615 / JX461969 / JX485228 / JX461783 / JX461437 2615b JX461616 / JX461970 / JX485229 / JX461784 / JX461432 2616 2616a Oukaimeden (Ouk) 31.20426 -7.86705 2600 JX462079 / JX462172 / JX461617 / JX461971 / JX485226 / JX461785 / JX461439 2616b JX461618 / JX461972 / JX485227 / JX461786 / JX461440 1579 1579a Jebel Ayache (Jay) 32.53671 -4.79110 3043 JX462098 / JX462178 / JX461565 / JX461919 / JX485276 / JX461733 / JX461391 1579b JX461566 / JX461920 / JX485277 / JX461734 / JX461392 2552 2552a Jebel Ayache (Jay) 32.53671 -4.79110 3043 JX462099 / JX462177 / JX461589 / JX461943 / JX485278 / JX461757 / JX461413 2552b JX461590 / JX461944 / JX485279 / JX461758 / JX461414 2564 2564a Jebel Ayache (Jay) 32.53671 -4.79110 3043 JX462101 / JX462176 / JX461597 / JX461951 / JX485272 / JX461765 / JX461421 2564b JX461598 / JX461952 / JX485273 / JX461766 / JX461422 2608 2608a Jebel Ayache (Jay) 32.53671 -4.79110 3043 JX462102 / JX462175 / JX461611 / JX461965 / JX485270 / JX461779 / JX461433 2608b JX461612 / JX461966 / JX485271 / JX461780 / JX461434 2618 2618a Jebel Ayache (Jay) 32.53671 -4.79110 3043 JX462103 / JX462174 / JX461619 / JX461973 / JX485268 / JX461787 / JX461441 2618b JX461620 / JX461974 / JX485269 / JX461788 / JX461442 9189 9189a Jebel Ayache (Jay) 32.53671 -4.79110 3043 ... / ... / JX461675 / JX462043 / JX485282 / JX461859 / JX461509 9189b JX461676 / JX462044 / JX485283 / JX461860 / JX461510 9199 9199 Jebel Ayache (Jay) 32.53671 -4.79110 3043 JX462108 / JX462183 / ... / ... / ... / ... / ... 9255 9255a Jebel Ayache (Jay) 32.53671 -4.79110 3043 JX462110 / JX462179 / JX461681 / JX462051 / JX485290 / JX461867 / JX461515 9255b JX461682 / JX462052 / JX485291 / JX461868 / JX461516 9209 9209 Jebel Ayache (Jay) 32.53671 -4.79110 3043 JX462109 / JX462184 / ... / ... / ... / ... / ... 9191 9191a Jebel Ayache (Jay) 32.53671 -4.79110 3043 JX462106 / JX462181 / JX461677 / JX462045 / JX485292 / JX461861 / JX461511 9191b JX461678 / JX462046 / JX485293 / JX461862 / JX461512 9336 9336 Jebel Ayache (Jay) 32.53671 -4.79110 3043 JX462111 / JX462182 / ... / ... / ... / ... / ... 9193 9193a Jebel Ayache (Jay) 32.53671 -4.79110 3043 JX462107 / JX462188 / JX461679 / JX462047 / JX485288 / JX461863 / JX461513 9193b JX461680 / JX462048 / JX485289 / JX461864 / JX461514 2557 2557a Jebel Ayache (Jay) 32.53671 -4.79110 3043 ... / ... / JX461595 / JX461949 / JX485274 / JX461763 / JX461419 2557b JX461596 / JX461950 / JX485275 / JX461764 / JX461420 9145 9145a Jebel Ayache (Jay) 32.53671 -4.79110 3043 JX462104 / JX462185 / JX461673 / JX462041 / JX485286 / JX461857 / JX461507 9145b JX461674 / JX462042 / JX485287 / JX461858 / JX461508 5076 5076a Jebel Azourki (Jaz) 31.75847 -6.28826 2789 JX462120 / JX462206 / JX461647 / JX462005 / JX485296 / JX461821 / JX461475 5076b JX461648 / JX462006 / JX485297 / JX461822 / JX461476 5128 5128a Jebel Azourki (Jaz) 31.75847 -6.28826 2789 JX462126 / JX462207 / JX461665 / JX462029 / JX485298 / JX461845 / JX461495 5128b JX461666 / JX462030 / JX485299 / JX461846 / JX461496 5091 5091 Jebel Azourki (Jaz) 31.75847 -6.28826 2789 JX462122 / JX462208 / ... / ... / ... / ... / ... 5017 5017a Jebel Azourki (Jaz) 31.75847 -6.28826 2789 JX462113 / JX462209 / JX461635 / JX461989 / JX485308 / JX461805 / JX461459 5017b JX461636 / JX461990 / JX485309 / JX461806 / JX461460 5122 5122a Jebel Azourki (Jaz) 31.75847 -6.28826 2789 JX462125 / JX462210 / JX461661 / JX462023 / JX485300 / JX461839 / JX461491 5122b JX461662 / JX462024 / JX485301 / JX461840 / JX461492 Extreme genetic diversity in Atlantolacerta 123 GenBank acession codes 4/5 Specimen code Alleles Population Latitude Longitude Altitude 12S / ND4+tRNA-His / PDC / ACM4 / C-MOS / MC1R / RAG1 5105 5105a Jebel Azourki (Jaz) 31.75847 -6.28826 2789 JX462123 / JX462211 / JX461659 / JX462019 / JX485302 / JX461835 / JX461487 5105b JX461660 / JX462020 / JX485303 / JX461836 / JX461488 5072 5072a Jebel Azourki (Jaz) 31.75847 -6.28826 2789 JX462118 / JX462212 / JX461645 / JX462001 / JX485304 / JX461817 / JX461471 5072b JX461646 / JX462002 / JX485305 / JX461818 / JX461472 5037 5037a Jebel Azourki (Jaz) 31.75847 -6.28826 2789 JX462116 / JX462213 / JX461639 / JX461995 / JX485322 / JX461811 / JX461465 5037b JX461640 / JX461996 / JX485323 / JX461812 / JX461466 5011 5011a Jebel Azourki (Jaz) 31.75847 -6.28826 2789 JX462112 / JX462204 / JX461631 / JX461985 / JX485312 / JX461801 / JX461455 5011b JX461632 / JX461986 / JX485313 / JX461802 / JX461456 5034 5034 Jebel Azourki (Jaz) 31.75847 -6.28826 2789 JX462115 / JX462205 / ... / ... / ... / ... / ... 5080 5080 Jebel Azourki (Jaz) 31.75847 -6.28826 2789 JX462121 / JX462216 / ... / ... / ... / ... / ... 5025 5025a Jebel Azourki (Jaz) 31.75847 -6.28826 2789 JX462114 / JX462214 / JX461637 / JX461991 / JX485314 / JX461807 / JX461461 5025b JX461638 / JX461992 / JX485315 / JX461808 / JX461462 5043 5043a Jebel Azourki (Jaz) 31.75847 -6.28826 2789 JX462117 / JX462218 / JX461641 / JX461997 / JX485324 / JX461813 / JX461467 5043b JX461642 / JX461998 / JX485325 / JX461814 / JX461468 5073 5073 Jebel Azourki (Jaz) 31.75847 -6.28826 2789 JX462119 / JX462215 / ... / ... / ... / ... / ... 5111 5111 Jebel Azourki (Jaz) 31.75847 -6.28826 2789 JX462124 / JX462217 / ... / ... / ... / ... / ... 6016 6816a Outabati (Out) 32.17714 -5.33214 2441 JX462128 / JX462221 / JX461671 / JX462037 / JX485330 / JX461853 / JX461503 6816b JX461672 / JX462038 / JX485331 / JX461854 / JX461504 11754 11754a Outabati (Out) 32.17714 -5.33214 2441 JX462140 / JX462230 / JX461551 / JX461903 / JX485350 / JX461717 / JX461375 11754b JX461552 / JX461904 / JX485351 / JX461718 / JX461376 11746 11746a Outabati (Out) 32.17714 -5.33214 2441 JX462137 / JX462228 / JX461547 / JX461899 / JX485346 / JX461713 / JX461371 11746b JX461548 / JX461900 / JX485347 / JX461714 / JX461372 11743 11743a Outabati (Out) 32.17714 -5.33214 2441 JX462135 / JX462226 / JX461543 / JX461895 / JX485342 / JX461709 / JX461367 11743b JX461544 / JX461896 / JX485343 / JX461710 / JX461368 11717 11717a Outabati (Out) 32.17714 -5.33214 2441 JX462130 / JX462222 / JX461533 / JX461885 / JX485332 / JX461699 / JX461357 11717b JX461534 / JX461886 / JX485333 / JX461700 / JX461358 11755 11755a Outabati (Out) 32.17714 -5.33214 2441 JX462139 / JX462231 / JX461553 / JX461905 / JX485352 / JX461719 / JX461377 11755b JX461554 / JX461906 / JX485353 / JX461720 / JX461378 11727 11727a Outabati (Out) 32.17714 -5.33214 2441 JX462131 / JX462232 / JX461535 / JX461887 / JX485334 / JX461701 / JX461359 11727b JX461536 / JX461888 / JX485335 / JX461702 / JX461360 11752 11752a Outabati (Out) 32.17714 -5.33214 2441 JX462138 / JX462229 / JX461549 / JX461901 / JX485348 / JX461715 / JX461373 11752b JX461550 / JX461902 / JX485349 / JX461716 / JX461374 6643 6643 Outabati (Out) 32.17714 -5.33214 2441 JX462129 / JX462220 / ... / ... / ... / ... / ... 11741 11741a Outabati (Out) 32.17714 -5.33214 2441 JX462134 / JX462225 / JX461541 / JX461893 / JX485340 / JX461707 / JX461365 11741b JX461542 / JX461894 / JX485341 / JX461708 / JX461366 11734 11734a Outabati (Out) 32.17714 -5.33214 2441 JX462133 / JX462224 / JX461539 / JX461891 / JX485338 / JX461705 / JX461363 11734b JX461540 / JX461892 / JX485339 / JX461706 / JX461364 11745 11745a Outabati (Out) 32.17714 -5.33214 2441 JX462136 / JX462227 / JX461545 / JX461897 / JX485344 / JX461711 / JX461369 11745b JX461546 / JX461898 / JX485345 / JX461712 / JX461370 11733 11733a Outabati (Out) 32.17714 -5.33214 2441 JX462132 / JX462223 / JX461537 / JX461889 / JX485336 / JX461703 / JX461361 11733b JX461538 / JX461890 / JX485337 / JX461704 / JX461362 Extreme genetic diversity in Atlantolacerta 124 GenBank acession codes 5/5 Specimen code Alleles Population Latitude Longitude Altitude 12S / ND4+tRNA-His / PDC / ACM4 / C-MOS / MC1R / RAG1 6639 6639 Outabati (Out) 32.17714 -5.33214 2441 JX462127 / JX462219 / ... / ... / ... / ... / ... 3865 3865a Toubkal (Tou) 31.09415 -7.91367 2600 JX462142 / JX462236 / JX461627 / JX461981 / JX485360 / JX461797 / JX461451 3865b JX461628 / JX461982 / JX485361 / JX461798 / JX461452 13276 13276a Toubkal (Tou) 31.09415 -7.91367 2600 JX462143 / JX462237 / ... / JX461909 / JX485362 / JX461723 / JX461381 13276b ... / JX461910 / JX485363 / JX461724 / JX461382 5090 5090a Jebel Awlime (JAw) 30.81708 -8.86298 2967 JX462144 / JX462234 / JX461651 / JX462011 / JX485354 / JX461827 / JX461481 5090b JX461652 / JX462012 / JX485355 / JX461828 / JX461482 13179 13179a Jebel Awlime (JAw) 30.81708 -8.86298 2967 JX462146 / JX462235 / JX461555 / JX461907 / JX485358 / JX461721 / JX461379 13179b JX461556 / JX461908 / JX485359 / JX461722 / JX461380 5123 5123a Jebel Awlime (JAw) 30.81708 -8.86298 2967 JX462145 / JX462233 / JX461663 / JX462025 / JX485356 / JX461841 / JX461493 5123b JX461664 / JX462026 / JX485357 / JX461842 / JX461494 Extreme genetic diversity in Atlantolacerta 125 Sequences were aligned for each gene independently using the online version of MAFFT v.6 [91] with default parameters (gap opening penalty = 1.53, gap extension = 0.0) and FFT-NS-1 algorithm. Coding gene fragments (ND4,C-MOS,ACM4,RAG1,PDC and MC1R) were translated into amino acids and no stop codons were observed, suggesting that the sequences were all functional. Heterozygous individuals were identified based on the presence of two peaks of approximately equal height at a single nucleotide site. SEQPHASE [92] was used to convert the input files, and the software PHASE v2.1.1 to resolve phased haplotypes [93]. Default settings of PHASE were used except for phase probabilities that were set as ≥ 0.7 [94]. All polymorphic sites with a probability of < 0.7 were coded in both alleles with the appropriate IUPAC ambiguity code. Phased nuclear sequences were used for the structure analysis; networks and species tree analysis, and the unphased sequences for the phylogenetic analyses (see below). DnaSP [95] was used to calculate the number of haplotypes (h) and mutations (η). Mega v.3.0 [96] was used to estimate uncorrected p-distances and to obtain the number of variable and parsimony informative sites. Phylogenetic analyses Phylogenetic analyses were performed using maximum likelihood (ML) and Bayesian (BI) methods. JModelTest [97] was used to select the most appropriate model of sequence evolution under the Akaike Information Criterion [98]. ML analyses were performed with RAxML v.7.0.4 [99] with 100 random addition replicates. A GTR+I+G model was used and parameters were estimated independently for each partition (by gene). Reliability of the ML tree was assessed by bootstrap analysis [100] including 1000 replications. Bayesian analyses were performed with MrBayes v.3.1.2 [101] with best fitting models applied to each partition by gene and all parameters unlinked across partitions. The models selected for the different partitions were: 12S, GTR+I+G; ND4, GTR+G; tRNA-His, GTR+I+G; ACM4, HKY+I; CMOS, GTR+I+G; MC1R, HKY+I+G; PDC, GTR+I+G; and RAG1, GTR+I. Two independent runs of 5x106generations were carried out, sampling at intervals of 1000 generations producing 5000 trees. Convergence and appropriate sampling were confirmed examining the standard deviation of the split frequencies between the two simultaneous runs and the Potential Scale Reduction Factor (PSRF) diagnostic. Burn-in was performed discarding the first 1250 trees of each run (25%) and a majority-rule consensus tree was generated from the remaining trees. In both ML and BI alignment gaps were treated as missing data and the nuclear gene sequences were not phased. Extreme genetic diversity in Atlantolacerta 126 Nuclear Networks The genealogical relationships between the populations were assessed with haplotype networks for all the individual nuclear genes, constructed using statistical parsimony [102] implemented in the program TCS v 1.21 [103] with a connection limit of 95%. This analysis was made with the phased sequences. Haplotypes were coloured taking into account the population of origin. Population structure – Clustering analyses A model-based Bayesian clustering method was applied to all haplotypes using STRUCTURE v.2.3.2 [60, 104-105]. In this analysis, individuals are probabilistically assigned to either a single cluster (the population of origin), or more than one cluster (if there is admixture). STRUCTURE was run with haplotype information from the nuclear fragments independently. We ran our data with the all parameters combinations between the Ancestry Model and the Allele Frequency Model to compare the results. The genetic structure was forced to vary from K = 2 to K = 10 clusters, the latter corresponding to the number of geographic populations sampled plus two. STRUCTURE ran for 550 000 steps, of which the first 50 000 were discarded as burn-in. For each value of K ten independent replicates of the Markov Chain Monte Carlo (MCMC) were conducted. To detect the true number of clusters (K) we followed the graphical methods and algorithms outlined in Evanno et al. [61], with the comparison of the average posterior probability values for K (log likelihood; ln L) using the online version, STRUCTURE HARVESTER v0.6.5 (available at: http://taylor0.biology.ucla.edu/struct_ harvest/, April 2011). Species tree, and divergence time estimates Here we applied the coalescent-based species-tree approach implemented in STARBEAST [106] an extension of BEAST v1.6.1 [107] to test the origin and diversification patterns in Atlantolacerta, and to compare these results to those obtained from the ML and BI analyses of the concatenated dataset. This analysis needs a priori information regarding the species/populations delimitation and the species/populations assignation of the individuals in order to reconstruct the topology of the species tree. For this approach, we used the results obtained from previous clustering analyses to define the groups of individuals to be used as “species” (populations) in STARBEAST [106]. The clustering analysis supported the existence of six lineages, as Oukaimeden, Toubkal and J. Awlime were included in the same lineage. All five nuclear gene fragments, 12S and the fragment consistent of the ND4 and flanking tRNA-His were included in the analyses as 7 independent partitions. The phased dataset was Extreme genetic diversity in Atlantolacerta 127 used for the nuclear loci. The input file was formatted with the BEAUti utility included in the software package. We performed two independent runs of 1.5 x 108generations, sampling every 15 000 generations, from which 10% were discarded as burn-in. Models and prior specifications applied were as follows (otherwise by default): 12S - GTR+G; ND4 and tRNA-His - HKY+G; MC1R - HKY+I; ACM4 - HKY+I; C-MOS - GTR+I+G; RAG1 - HKY+I; PDC - GTR+I; Relaxed Uncorrelated Lognormal Clock (estimate); Yule process of speciation; random starting tree; alpha Uniform (0, 10). For all analyses implemented in BEAST, convergence for all model parameters was assessed by examining trace plots and histograms in Tracer v1.5 [108] after obtaining an effective sample size (ESS) > 200. The initial 10% of samples were discarded as burn-in. Runs were combined using LogCombiner, and maximum credibility trees with divergence time means and 95% highest probability densities (HPDs) were produced using Tree Annotator (both part of the BEAST package). Trees were visualized using the software FigTree v1.3.1 [109]. Several studies have already calculated divergence rates for reptiles, and particularly for lacertids [2, 15, 49]. Pinho et al. [15] used well-known and dated independent geological events in the Aegean [110] to estimate a maximum and minimum mutation rate for the ND4 mitochondrial fragment (and flanking tRNA-His) for the lacertid lizards of the genus Podarcis (0.0278 and 0.0174 mutation/site/million years, respectively). However, this was the only information available for our data, since we did not have any fossils or calibrations for nuclear markers. It is important to bear in mind that, in the absence of accurate calibration points in the phylogeny from external and independent data (fossil records, known biogeographic events, or paleoclimatic reconstructions) or as a result of the heterogeneity in the evolutionary rate between the calibrated and uncalibrated taxa, temporal estimates by means of molecular data could be a potential source of inference error, and, therefore, they should be treated with caution [111]. Despite the limitations of molecular clocks [111-112], divergence time estimates can still provide a proxy for the temporal window of evolutionary diversification in species groups of interest. Therefore and taking into account our data limitations and availability, we used BEAST v.1.6.1 [107] to estimate dates of the cladogenetic events using only ND4 and flanking tRNA-His. We used a phylogeny pruned arbitrarily to include one representative from each of the major lineages uncovered with the concatenated analysis (6 specimens in total, we excluded J. Awlime population, because of the lack of support of the branch in previous analyses). This method excludes closely related terminal taxa because the Yule tree prior (see below) does not include a model of coalescence, which can complicate rate estimation for closely related sequences [113]. Analyses were run four times for 5x107generations with a sampling frequency of 10 000. Extreme genetic diversity in Atlantolacerta 128 Models and prior specifications applied were as follows (otherwise by default): GTR+G for 12S; HKY+G for ND4 and tRNA-His; HKY+I for MC1R; HKY+I for ACM4; GTR+G+I for C-MOS; HKY+I for RAG1; GTR+I for PDC; Relaxed Uncorrelated Lognormal Clock (estimate); Yule process of speciation; random starting tree; alpha Uniform (0, 10); ucld.mean of ND4 Normal (initial value: 0.0226, mean: 0.0226, Stdev: 0.0031). Acknowledgements and funding MB is supported by the FCT grant SFRH/BD/41488/2007. This work was funded by FCT grant PTDC/BIA-BDE/74349/2006 and by grant CGL2009-11663 from the Ministerio de Educación y Ciencia, Spain to SC. Fieldwork in Morocco in 2008 and 2009 was conducted under permit decision 84° issued by Haut Commissariat aux Eaux et Forêts et à la Lutte Contre la Désertification, issued to David Donaire plus other permits issued to the latter along a 10 year period. Thanks to all colleagues from CIBIO who assisted during the hard fieldwork, especially to Anna Perera, Daniele Salvi, Fatima Jorge and Fernando Martinez-Freiria. We also want to thank to the anonymous reviewers that helped to improve this manuscript. Authors’ contributions MB carried out the molecular laboratory work, analysed the data and drafted a preliminary version of the manuscript. All authors participated in the conception and design of the study, collection of samples, writing and approval of the final manuscript. Extreme genetic diversity in Atlantolacerta 129 References 1. Arribas O, Carranza S: Morphological and genetic evidence of the full species status of Iberolacerta cyreni martinezricai (Arribas, 1996).Zootaxa 2004, 634:1-24. 2. Carretero MA: An integrated assessment of a group with complex systematics: the Iberomaghrebian lizard genus Podarcis (Squamata, Lacertidae).Integr Zool 2008, 3(4):247266. 3. Duggen S, Hoernle K, van den Bogaard P, Rupke L, Morgan JP: Deep roots of the Messinian salinity crisis.Nature 2003, 422(6932):602-606. 4. Krijgsman W, Hilgen FJ, Raffi I, Sierro FJ, Wilson DS: Chronology, causes and progression of the Messinian salinity crisis.Nature 1999, 400(6745):652-655. 5. Pinho C, Ferrand N, Harris DJ: Reexamination of the Iberian and North African Podarcis (Squamata: Lacertidae) phylogeny based on increased mitochondrial DNA sequencing. Molecular Phylogenetics and Evolution 2006, 38(1):266-273. 6. Santos X, Roca J, Pleguezuelos JM, Donaire D, Carranza S: Biogeography and evolution of the Smooth snake Coronella austriaca (Serpentes : Colubridae) in the Iberian Peninsula: evidence for Messinian refuges and Pleistocenic range expansions.Amphibia-Reptilia 2008, 29(1):35-47. 7. Schmitt T: Molecular biogeography of Europe: Pleistocene cycles and postglacial trends. Frontiers in Zoology 2007, 4(1):11. 8. Tolley KA, Chase BM, Forest F: Speciation and radiations track climate transitions since the Miocene Climatic Optimum: a case study of southern African chameleons.Journal of Biogeography 2008, 35(8):1402-1414. 9. Camargo A, Sinervo B, Sites JW: Lizards as model organisms for linking phylogeographic and speciation studies.Molecular ecology 2010, 19(16):3250-3270. 10. Carranza S, Arnold EN, Geniez P, Roca J, Mateo JA: Radiation, multiple dispersal and parallelism in the skinks, Chalcides and Sphenops (Squamata: Scincidae), with comments on Scincus and Scincopus and the age of the Sahara Desert.Molecular Phylogenetics and Evolution 2008, 46(3):1071-1094. 11. Fonseca MM, Brito JC, Paulo OS, Carretero MA, Harris DJ: Systematic and phylogeographical assessment of the Acanthodactylus erythrurus group (Reptilia: Lacertidae) based on phylogenetic analyses of mitochondrial and nuclear DNA.Molecular Phylogenetics and Evolution 2009, 51(2):131-142. 12. Fonseca MM, Brito JC, Rebelo H, Kalboussi M, Larbes S, Carretero MA, Harris DJ: Genetic variation among spiny-footed lizards in the Acanthodactylus pardalis group from North Africa.African Zoology 2008, 43(1):8-15. 13. Harris DJ, Batista V, Carretero MA: Assessment of genetic diversity within Acanthodactylus erythrurus (Reptilia: Lacertidae) in Morocco and the Iberian Peninsula using mitochondrial DNA sequence data.Amphibia-Reptilia 2004, 25(2):227-232. 14. Kaliontzopoulou A, Pinho C, Harris DJ, Carretero MA: When cryptic diversity blurs the picture: a cautionary tale from Iberian and North African Podarcis wall lizards.Biological Journal of the Linnean Society 2011, 103(4):779-800. 15. Pinho C, Harris DJ, Ferrand N: Contrasting patterns of population subdivision and historical demography in three western Mediterranean lizard species inferred from mitochondrial DNA variation.Molecular ecology 2007, 16(6):1191-1205. 16. Rato C, Harris DJ: Genetic variation within Saurodactylus and its phylogenetic relationships within the Gekkonoidea estimated from mitochondrial and nuclear DNA sequences. Amphibia-Reptilia 2008, 29(1):25-34. 17. Perera A, Harris DJ: Genetic variability within the Oudri’s fan-footed gecko Ptyodactylus oudrii in North Africa assessed using mitochondrial and nuclear DNA sequences.Molecular Phylogenetics and Evolution 2010, 54:634-639. 18. Carranza S, Romano A, Arnold EN, Sotgiu G: Biogeography and evolution of European cave salamanders, Hydromantes (Urodela: Plethodontidae), inferred from mtDNA sequences. Journal of Biogeography 2008, 35:724-738. 19. Carranza S, Arnold EN: History of West Mediterranean newts, Pleurodeles (Amphibia: Salamandridae), inferred from old and recent DNA sequences.Systematics and Biodiversity 2004, 1(3):327-337. 20. Rato C, Carranza S, Harris DJ: When selection deceives phylogeographic interpretation: The case of the Mediterranean house gecko, Hemidactylus turcicus (Linnaeus, 1758).Molecular Phylogenetics and Evolution 2011, 58:365-373. Extreme genetic diversity in Atlantolacerta 130 21. Rato C, Carranza S, Perera A, Carretero MA, Harris DJ: Conflicting patterns of nucleotide diversity between mtDNA and nDNA in the Moorish gecko, Tarentola mauritanica. Molecular Phylogenetics and Evolution 2010, 56(3):962-971. 22. Degnan JH, Rosenberg NA: Gene tree discordance, phylogenetic inference and the multispecies coalescent.Trends in ecology & evolution 2009, 24(6):332-340. 23. Edwards SV: Is a New and General Theory of Molecular Systematics Emerging? Evolution 2009, 63(1):1-19. 24. Maddison WP: Gene trees in species trees.Systematic Biology 1997, 46(3):523-536. 25. Agapow PM, Bininda-Emonds ORP, Crandall KA, Gittleman JL, Mace GM, Marshall JC, Purvis A: The impact of species concept on biodiversity studies.Q Rev Biol 2004, 79(2):161-179. 26. Mayden RL: A hierarchy of species concepts: the denouement in the saga of the species problem. In: Species: The units of diversity. Edited by Claridge MF, Dawah HA, Wilson MR. London: Chapman and Hall; 1997: 381-423. 27. Sites JW, Marshall JC: Delimiting species: a Renaissance issue in systematic biology.Trends in Ecology & Evolution 2003, 18(9):462-470. 28. Agapow M: Species: demarcation and diversity. In: Phylogeny and Conservation. Edited by Purvis A, Gittleman JL, Brooks T. Cambridge, UK: Cambridge University Press; 2005: 57–75. 29. Sattler T, Bontadina F, Hirzel AH, Arlettaz R: Ecological niche modelling of two cryptic bat species calls for a reassessment of their conservation status.Journal of Applied Ecology 2007, 44(6):1188-1199. 30. de Queiroz K, Donoghue MJ: Phylogenetic Systematics and the Species Problem.Cladistics 1988, 4(4):317-338. 31. de Queiroz K: Species Concepts and Species Delimitation.Systematic Biology 2007, 56(6):879886. 32. Dayrat B: Towards integrative taxonomy.Biological Journal of the Linnean Society 2005, 85:407–415. 33. Padial JM, Miralles A, De la Riva I, Vences M: The integrative future of taxonomy.Front Zool 2010, 7:16. 34. Schlick-Steiner BC, Steiner FM, Seifert B, Stauffer C, Christian E, Crozier RH: Integrative taxonomy: a multisource approach to exploring biodiversity.Annu Rev Entomol 2010, 55:421–438. 35. Miralles A, Vasconcelos R, Perera A, Harris DJ, Carranza S: An integrative taxonomic revision of the Cape Verdean skinks (Squamata, Scincidae).Zoologica Scripta 2010, 40:16-44. 36. Vasconcelos R, Carranza S, Harris DJ: Insight into an island radiation: the Tarentola geckos of the Cape Verde archipelago.Journal of Biogeography 2010, 37(6):1047-1060. 37. Galbreath KE, Hafner DJ, Zamudio KR: When cold is better: Climate-driven elevation shifts yield complex patterns of diversification and demography in an Alpine Specialist (American Pika, Ochotona Princeps).Evolution 2009, 63(11):2848-2863. 38. Hewitt GM: Genetic consequences of climatic oscilations in the quaternary.Philosophical Transactions of the Royal Society of London 2004, 359:183-195. 39. Knowles LL: Did the Pleistocene glaciations promote divergence? Tests of explicit refugial models in montane grasshopprers.Molecular ecology 2001, 10(3):691-701. 40. Hughes L: Climate change and Australia: Trends, projections and impacts.Austral Ecol 2003, 28(4):423-443. 41. Mouret V, Guillaumet A, Cheylan M, Pottier G, Ferchaud AL, Crochet PA: The legacy of ice ages in mountain species: post-glacial colonization of mountain tops rather than current range fragmentation determines mitochondrial genetic diversity in an endemic Pyrenean rock lizard.Journal of Biogeography 2011, 38(9):1717-1731. 42. Bons J, Geniez P: Amphibiens et reptiles du Maroc (Sahara Occidental compris) Atlas Biogéographique. Barcelone: Asociación Herpetologica Espanola; 1996. 43. Schleich HH, Kastle W, Kabisch K: Amphibians and Reptiles from North Africa. Königstein, Germany.: Koeltz Scientific Publications; 1996. 44. Arnold EN: Relationships of the Palaearctic lizards assigned to the genera Lacerta, Algyroides and Psammodromus (Reptila, Lacertidae). London: British Museum (Natural History); 1973. 45. Arnold EN: Towards a phylogeny and biogeography of the Lacertidae: relationships within an Old-World family of lizards derived from morphology. London: British Museum (Natural History); 1989. 46. Harris DJ: Molecular systematics and evolution of lacertid lizards.Natura Croatica 1999, 83(3):161-180. Extreme genetic diversity in Atlantolacerta 131 47. Mayer W, Bischoff W: Beiträge zur taxonomischen Revision der Gattung Lacerta (Reptilia: Lacertidae) Teil 1: Zootoca,Omanosaura,Timonund,Teira als eigenstandige Gattungen. Salamandra 1996, 32(3):163-170. 48. Oliverio M, Bologna MA, Mariottin P: Molecular biogeography of the Mediterranean lizards Podarcis Wagler, 1830 and Teira Gray, 1838 (Reptilia, Lacertidae) Journal of Biogeography 2000, 27:1403-1420. 49. Arnold EN, Arribas O, Carranza S: Systematics of the Palaearctic and Oriental lizard tribe Lacertini (Squamata: Lacertidae: Lacertinae), with descriptions of eight new genera. Zootaxa 2007, 1430:1-86. 50. Pavlicev M, Mayer W: Fast radiation of the subfamily Lacertinae (Reptilia: Lacertidae): History or methodical artefact? Molecular Phylogenetics and Evolution 2009, 52(3):727-734. 51. Busack SD: Notes on the biology of Lacerta andreanszkyi (Reptilia: Lacertidae).AmphibiaReptilia 1987, 8:231-236. 52. Carretero MA, Perera A, Harris DJ, Batista V, Pinho C: Spring diet and trophic partitioning in an alpine lizard community from Morocco.African Zoology 2006, 41(1):113-122. 53. Crochet PA, Chaline O, Surget-Groba Y, Debain C, Cheylan M: Speciation in mountains: phylogeography and phylogeny of the rock lizards genus Iberolacerta (Reptilia: Lacertidae). Molecular Phylogenetics and Evolution 2004, 30(3):860-866. 54. Carranza S, Arnold EN, Amat F: DNA phylogeny of Lacerta (Iberolacerta) and other lacertine lizards (Reptilia: Lacertidae): did competition cause long-term mountain restriction? Systematics and Biodiversity 2004, 2(01):57-77. 55. Funk DJ, Omland KE: Species-level paraphyly and polyphyly: Frequency, causes, and consequences, with insights from animal mitochondrial DNA.Annu Rev Ecol Evol S 2003, 34:397-423. 56. Knowles LL, Carstens BC: Delimiting species without monophyletic gene trees.Systematic Biology 2007, 56(6):887-895. 57. Hudson RR, Coyne JA: Mathematical consequences of the genealogical species concept. Evolution 2002, 56(8):1557-1565. 58. Hudson RR, Turelli M: Stochasticity overrules the "three-times rule": Genetic drift, genetic draft, and coalescence times for nuclear loci versus mitochondrial DNA.Evolution 2003, 57(1):182-190. 59. Pinho C, Harris DJ, Ferrand N: Non-equilibrium estimates of gene flow inferred from nuclear genealogies suggest that Iberian and North African wall lizards (Podarcis spp.) are an assemblage of incipient species.BMC Evolutionary Biology 2008, 8. 60. Pritchard JK, Stephens M, Donnelly P: Inference of population structure using multilocus genotype data.Genetics 2000, 155:945-959. 61. Evanno G, Regnaut S, Goudet J: Detecting the number of clusters of individuals using the software STRUCTURE: a simulation study Molecular ecology 2005, 14:2611-2620. 62. Banks MA, Eichert W: WHICHRUN (version 3.2): A computer program for population assignment of individuals based on multilocus genotype data.Journal of Heredity 2000, 91(1):87-89. 63. Corander J, Waldmann P, Sillanpaa MJ: Bayesian analysis of genetic differentiation between populations.Genetics 2003, 163(1):367-374. 64. Dawson KJ, Belkhir K: A Bayesian approach to the identification of panmictic populations and the assignment of individuals.Genet Res 2001, 78(1):59-77. 65. Manel S, Berthier P, Luikart G: Detecting wildlife poaching: Identifying the origin of individuals with Bayesian assignment tests and multilocus genotypes.Conservation Biology 2002, 16(3):650-659. 66. Pritchard JK, Donnelly P: Case-control studies of association in structured or admixed populations.Theor Popul Biol 2001, 60(3):227-237. 67. Rosenberg NA, Pritchard JK, Weber JL, Cann HM, Kidd KK, Zhivotovsky LA, Feldman MW: Genetic structure of human populations.Science 2002, 298(5602):2381-2385. 68. Turakulov R, Easteal S: Number of SNPS loci needed to detect population structure.Hum Hered 2003, 55(1):37-45. 69. Babault J, Teixell A, Arboleya ML, Charroud M: A late cenozoic age for long-wavelength surface uplift of the atlas mountains of Morocco.Terra Nova 2008, 20(2):102-107. 70. Gomez F, Beauchamp W, Barazangi M: Role of the Atlas Mountains (northwest Africa) within the African-Eurasian plate-boundary zone.Geology 2000, 28(9):775-778. 71. Brown JW, Rest JS, Garcia-Moreno J, Sorenson MD, Mindell DP: Strong mitochondrial DNA support for a Cretaceous origin of modern avian lineages.Bmc Biol 2008, 6. 138 Morphological analyses of Atlantolacerta 139 Introduction Delimiting species, despite being a controversial issue, is of major importance since species are the basic unit in areas such as ecology, biogeography and evolution, with serious implications for conservation biology (Myers et al., 2000; Sites & Marshall, 2003). However, species concepts and their delimitation are still controversial, and recently several new approaches, such as the “unified species concept” (de Queiroz, 2007) and conceptual advances in integrative taxonomy (Dayrat, 2005; Padial et al., 2010), have emerged to try to reconcile the different species concepts. The principal difficulty of determining if a population constitutes an independent evolving lineage occurs in recently separated species, which are less likely to achieve criteria such as morphological distinctiveness, reproductive isolation, ecological exclusivity and monophyly (de Queiroz, 2007). Moreover, speciation is not always accompanied with phenotypic changes, potentially leading to an underestimation of the actual levels of biodiversity. Cryptic species are an example that can be difficult to classify, particularly because morphology has been traditionally the main tool to identify and classify new species. Although in the past two decades the study of cryptic species has increased (Detwiler et al., 2010; Florio et al., 2012; Padial & de la Riva, 2009), see a review in (Bickford et al., 2007) mainly due to the advances in molecular and analytical methods, cryptic diversity remains a challenge for taxonomists. Additionally, delimitation of allopatric forms is a further challenge, as it is difficult to measure objectively some of the criteria that, usually, determine reproductive isolation. Cryptic species occur with regularity across all biogeographical regions and major metazoan taxa, and are more common than previously thought (Pfenninger & Schwenk, 2007). North Africa is a region where cryptic diversity has been described in several taxa such as plants (Abdelaziz et al., 2011), spiders (Duncan et al., 2010), mammals (Ben Faleh et al., 2012) and reptiles (Perera & Harris, 2010; Rato et al., 2012). The diverse geographical and geological features and the variety of climates exert different selective pressures that have promoted speciation processes in the region. The Atlas Mountains are especially interesting as a source of speciation. They formed at the Africa-Eurasia plate boundaries, and uplifted during the cenozoic (Gómez et al., 2000) and have been identified as refugia during the Pleistocenic climatic fluctuations (Medail & Diadema, 2009), harbouring a diversity that is still underexplored. However, there are an increasing number of examples (Brown et al., 2002; Fritz et al., 2006; Rato et al., 2010; Recuero et al., 2007) demonstrating the role of the Atlas system in species diversification. The most recent study of the Moroccan day gecko (genus Quedenfeldtia) endemic from the Atlas region, confirms once more the interest of this region as source of cryptic speciation (Barata et al., 2012b). Morphological analyses of Atlantolacerta 140 Atlantolacerta andreanskyi (Werner, 1929) is a small lacertid lizard endemic to the western and central High Atlas Mountains of Morocco. It is the lacertid found at higher altitudes in the region, being restricted to areas between 2400 m and 3800 m a.s.l. (Bons & Geniez, 1996; Schleich et al., 1996). It is often found near watercourses and in the base of cushion-like thorny plants (Bons & Geniez, 1996) that offer a buffered microclimate with humidity, food, and protection against predators and wind (Schleich et al., 1996). In a recent study, individuals from eight different geographic populations, covering the distribution range of A. andreanskyi, were compared using a multilocus approach (see Figure 1), which included two mitochondrial and five nuclear markers (Barata et al., 2012a) results revealed an extreme genetic diversity among seven of the eight populations analysed in mtDNA, showing divergence levels ranging from 1.6% to 6.6% in 12S rRNA and 5.5% to 16.5% in ND4. The nuclear markers (ACM4, MC1R, C-MOS, PDC and RAG1) were concordant with the mtDNA even if monophyly was not achieved between some populations. In view of these results, the authors suggested the possibility that A. andreanskyi might be a complex of cryptic species. Unfortunately, due to the restricted, and in many cases inaccessible, distribution range, there is a profound lack of knowledge regarding this species, and the few existing studies (Busack, 1987; Klemmer, 1969; Pasteur & Bons, 1960; Saint Girons, 1953; Stemmler, 1972; Volobouev et al., 1990; Werner, 1929, 1931, 1935) are mostly based on individuals from two geographically close populations from the High Atlas (Oukaimeden and Toubkal). The only available studies on morphology and sexual dimorphism are also based on these populations (Busack, 1987; Rykena & Bischoff, 1992; Schleich et al., 1996), although Joger & Bischoff (1989) mention the population from Jebel Ayache and the possibility that this population is distinct from the others at a subspecific or specific level, but without any detailed explanation other than its geographical isolation. In consequence, there is an urgent need to investigate the morphological variation across the distribution range of the species, in order to evaluate its real “cryptic nature”. In this study we investigate the morphological variation within A. andreanskyi in different populations, in order to compare this with the genetic variation described by Barata et al. (2012a). To achieve this, a morphological analysis of body measurements, scalation and colour characters were performed. The main objective is to assess the cryptic nature of this species, and, as necessary, to describe the new taxa observed. Morphological analyses of Atlantolacerta 141 Material and Methods Sampling and data collection The study area comprises the western and central parts of the High Atlas Mountains of Morocco, covering the whole distribution range of A. andreanskyi (Bons & Geniez, 1996). A total of 139 specimens from 6 of the 8 populations genetically characterized in Barata et al. (2012a) covering most of the species distribution range were sampled (Figure 2): Jebel Sirwa (13 males (M) and 9 females (F)), Oukaimeden (14M, 17F), Tizin Tichka (12M, 12F), Jebel Azourki (10M, 21F), Outabati (6M, 5F) and Jebel Ayache (8M, 12F). Despite considerable effort to sample in the other two populations included in the genetic study (Toubkal and Jebel Awlime), only two and three specimens could be found, respectively, and thus, although they were examined for descriptive purposes, they were not included in the multivariate analysis. Specimens were caught by hand, identified and sexed on the basis of external features (e.g. developed femoral pores, see(Schleich et al., 1996). Figure 1.Atlantolacerta andreanskyi phylogenetic trees adapted from Barata et al. (2012a). A - mtDNA tree (12S and ND4+tRNA-His) and B - nuclear DNA concatenated tree (MC1R, PDC, ACM4, C-MOS and RAG1). Populations are represented with the same colours in all the figures; Jebel Sirwa (pink), Oukaimeden (red), Tizin Tichka (dark blue), Jebel Azourki (light blue), Outabati (light green) and Jebel Ayache (dark green), Jebel Awlime (orange) and Toubkal (yellow). The last two populations were not included in this study due to the small sample size. In total, twelve linear measurements, nine pholidotic characters (Table 1) and seven colour characters were taken. Snout vent length (SVL) was measured from the tip of the snout to the cloaca opening; trunk length (TRL) was measured from the posterior edge of the forelimb insertion to the anterior edge of the hindlimb insertion. Tail width (TW) was recorded at its widest point. Head length (HL) was measured from tip of the snout to the collar, head width Morphological analyses of Atlantolacerta 142 (HW) at its widest part, usually at the level of the temporal region, and head height (HH) from occiput to jaws. Pileus length (PL) was measured from the occipital to the limit of the rostral scale. The total lengths of front (FLL) and hind (HLL) limbs were measured from the longest toe to the base of the limb. Femur length (FL) was measured from the base of the hindlimb to the knee joint, tibia length (TBL) from the knee joint to the ankle joint and fourth toe length (HTL) from the insertion of the toe to its extremity, including the claw. Also, detailed pictures of dorsal, lateral and ventral body were taken in the field, and nine pholidotic variables recorded a posteriori from them: number of ventral scales (VSN) including all the large scales counted in a midline from the collar to the anterior insertion of hindlimbs; number of gular scales (GSN) in a midline from the collar to the chin shields scales; number of collar scales (CSN), number of femoral pores (FPN) counted in males only; number of supratemporal scales (STSN); number of supra labial scales (SLSN); number of supraciliary scales (SCSN); number of supraciliary granules (SCGN), and number of enlarged side to side lamellae under the fourth toe (Lam). Regarding colour pattern, the following variables were also recorded from pictures: presence of black pigmentation (spots) in the lateral head (HPL, 1 = absent, 2 = scarce, 3 = abundant), dorsal head (HPD, 1 = absent, 2 = scarce, 3 = abundant), ventral head (HPV, 1 = absent, 2 = scarce, 3 = abundant), ventral body (VBP, 1 = absent, 2 = scarce, 3 = abundant) and cloacal region (CD, 0 = absent, 1 = one single dot in the anal plate, 2 = two dots in the anal plate; 3 = three dots, in the anal plate); presence of a central dorsal line (CBL, 1 = absent, 2 = discontinuous, 3 = continuous) and presence of light dorsolateral lines (1 = absent, 2 = present). Since morphological variation can be affected by inter-observer differences (Roitberg et al., 2011) all linear measurements were recorded in the field by the same author (MB) to the nearest 0.01 mm, using a digital calliper. Pholidotic and colour variables were also retrieved from digital pictures by the same author (MB) at least twice and the mean value was recorded. All bilateral variables were taken from the right side of the animal. In cases where data could not be collected due to member amputations or poor picture quality they were replaced by the group mean. Only adult individuals were included in this study. Individuals were released in the same place where they were caught after recording the coordinates of the location with a GPS. Tissue samples from the localities in this study were already characterized genetically in Barata et al. (2012a). Morphological analyses of Atlantolacerta 143 Figure 2. Distribution map of Atlantolacerta andreanskyi. Populations investigated in the current study: Jebel Sirwa (pink), Oukaimeden (red), Tizin Tichka (dark blue), Jebel Azourki (light blue), Outabati (light green) and Jebel Ayache (dark green). White dots represent the known distribution of the species as available in Bons & Geniez (1996). Morphological analysis Body measurements and pholidotic variables were log-transformed and checked for homoscedasticity (Levene's test) and normality (Shapiro-Wilks test) assumptions. Since linear measurements were highly correlated to body size, namely snout-vent length (SVL, Pearson correlation in all cases p < 0.01), we used an isometric correction (Somers, 1986) to estimate body-size-corrected variables that were then used to investigate the existence of possible differentiation patterns not related to body size. For this, all linear measurements (log transformed) were projected on an isometric vector, in order to obtain a multivariate representation of the isometric size of each individual (mSIZE). Each variable was then regressed on this isometric vector and the residuals were used as size-corrected variables. Thus, the multivariate representation of isometric size (mSIZE) was used as a size estimator, whereas the remaining isometric-size corrected variables were used as a representation of shape (Kaliontzopoulou et al., 2010). Correction was done using the R package (R, 2011). In order to investigate the morphological differences between sexes and populations, (multi)variate analysis of the variance (M)ANOVAs were performed on size (multivariate representation of the isometric size, mSIZE), shape (the remaining iso-corrected linear measurements) and pholidotic variables separately using POPULATION, SEX and its interaction (POPULATION*SEX), as factors. Analyses were done using the software Morphological analyses of Atlantolacerta 144 STATISTICA 7.1 (Statsoft, 2005). Since factor SEX was significant in several variables, further analyses were performed on males and females separately. Additionally, sexual dimorphism within each population was investigated using Analysis of Covariance (ANCOVAs) on linear measurements using SVL as covariate, and Student's t test on pholidotic variables. To investigate the generalized morphological relationships among the different A. andreanskyi populations, and the contribution of each dataset to the differentiation among populations, Canonical Discriminant Function Analyses (CDFA) were performed on linear measurements (mSIZE and iso-corrected variables) and pholidotic variables separately. CDFA allowed assessment at a multivariate level of which variables were the major contributors to the differentiation between populations, and to create canonical discriminant functions to calculate the probability of classifying correctly the individuals based on them. We used the leave-one-out option to cross-validate the classification results. Since this procedure (Jacknife prediction) generates individual classifications using discriminate functions based on all observations except the given case, it provides a more accurate estimate of the classification values. CDFA analyses were performed using SPSS v.20.0.0 (IBM, 2011). Finally, variation in colour pattern between populations was investigated at a multivariate level using a Multiple Correspondence Analysis (MCA) using the software R (R, 2011). In all analysis, significance level was considered at p < 0.05. Genetic diagnosis In order to identify diagnosable nucleotide positions among different populations of A. andreanskyi for descriptive purposes,sequences from Barata et al. (2012a) were investigated. Diagnosable positions were detected with the help of the software Mega 5 (Tamura et al., 2011) and confirmed by eye in BioEdit (Hall, 1999). In order to locate the exact diagnosable positions, sequences from Barata et al. (2012a) were aligned to the complete genome of Podarcis muralis (GenBank accession number NC_011607.1, (Podnar et al., 2009). The 1108 bp fragments from Barata et al. (2012a) correspond to the positions 484 to 812 (329 bp) of the P. muralis 12S rRNA and 10948 to 11539 (592 bp) of the P. muralis ND4. Morphological analyses of Atlantolacerta 145 Bayesian Species delimitation Bayesian species delimitation was used in order to support the existence of the lineages observed previously in genetic analysis (Barata et al., 2012a). Bayesian species delimitation was conducted using the program Bayesian Phylogenetics and Phylogeography, BPP v.2.0b (Rannala & Yang, 2003; Yang & Rannala, 2010) using the five nuclear loci. This method accommodates both species phylogeny and lineage sorting due to ancestral polymorphism. A gamma prior G(1, 10) was used on the population size parameters (qs). The age of the root in the species tree (t0) was assigned the gamma prior G(1, 10), while the other divergence time parameters were assigned the Dirichlet prior (Yang & Rannala, 2010): equation 2). We used algorithm 0 with the finetuning parameter = 15.0 in order to ensure adequate rjMCMC mixing. This involves specifying a reversible jump algorithm to achieve dimension matching between species delimitation models with different number of parameters. Each species delimitation model was assigned equal prior probability and each analysis was run at least twice to confirm consistency between runs. The guide tree plays an important role in the result of the species delimitation model (Leaché & Fujita, 2010); therefore we used the guide tree: ((JSi, Ouk), ((Tiz, JAz), (Out, JAy))), based on the estimate of relationships from both mtDNA and nDNA trees (Bayesian and ML analysis, see results from(Barata et al., 2012a). Results Detailed descriptive statistics for all the variables in all the lineages analysed are presented in Table 1, while Table 2 includes the detailed results regarding sexual dimorphism within each lineage. Inter-lineage variation Linear measurements The (M)ANOVA analysis showed general differences between sexes, populations and its interaction in both size (mSIZE) and shape (Table 3). Regarding size, individuals from Jebel Azourki and Jebel Ayache were larger than the other ones, although males and females had a different pattern (Table 3; Figure 3). Regarding shape, there were differences among populations in all iso-corrected linear measurements (Table 3). On the other hand, males and females also showed differences in most of them, with the exception of tail width (TW), front limbs (FLL) and hind limbs (TBL, HTL and HLL, but not FL; Table 3, Figure 3). Finally, the degree of sexual dimorphism (POPULATION*SEX interaction) was similar among populations, with the exception of trunk length (TRL; Table 3; Figure 3). Morphological analyses of Atlantolacerta 146 Table 1. Descriptive statistics of linear measurements and pholidotic variables for males and females of all the populations included in the study. For each group, mean, standard deviation (Std Dev), minimum (Min) and maximum (Max) values and sample size (n) is detailed. Males Mean Std.Dev. Min. Max. N Mean Std.Dev. Min. Max. N Mean Std.Dev. Min. Max. N Mean Std.Dev. Min. Max. N Mean Std.Dev. Min. Max. N Mean Std.Dev. Min. Max. N SVL 42,24 2,96 36,24 45,50 13 42,02 3,48 35,68 47,50 14 42,53 2,84 36,92 45,97 12 41,49 1,41 40,36 44,14 6 47,30 3,68 42,86 52,44 10 48,99 2,07 45,28 51,67 13 TRL 21,64 1,84 18,08 23,82 13 22,96 1,63 19,68 25,92 14 22,47 2,12 18,38 26,65 12 22,35 1,26 21,01 23,92 6 24,62 2,86 20,19 29,48 10 26,27 0,66 25,37 27,23 13 PL 9,85 0,75 8,73 10,74 13 9,58 0,65 8,44 10,79 14 9,50 0,60 8,39 10,30 12 9,37 0,25 9,04 9,75 6 10,95 0,53 10,32 11,80 10 11,58 0,47 10,73 12,11 13 HL 15,56 1,09 12,95 16,58 13 14,54 1,53 12,00 16,89 14 15,23 1,14 12,88 17,20 12 15,18 0,48 14,40 15,70 6 17,02 1,30 14,51 19,29 10 17,37 1,31 15,43 18,84 13 HW 6,39 0,56 5,53 7,00 13 5,77 0,29 5,32 6,40 14 5,90 0,45 5,23 6,58 12 5,47 0,14 5,26 5,67 6 6,68 0,39 6,03 7,20 10 7,18 0,37 6,65 7,78 13 HH 4,14 0,37 3,56 4,74 13 4,18 0,58 3,27 5,31 14 4,45 0,46 3,77 5,15 12 3,90 0,35 3,36 4,25 6 4,79 0,45 4,18 5,35 10 5,11 0,36 4,68 5,75 13 TW 4,92 0,17 4,64 5,30 13 3,90 0,32 3,37 4,45 14 3,88 0,46 3,20 4,31 5 3,79 0,16 3,49 3,92 6 4,59 0,42 4,11 5,21 10 4,94 0,33 4,51 5,32 13 FFL 12,67 0,81 11,40 14,01 13 12,25 0,79 11,02 13,55 14 12,18 1,00 10,12 13,23 12 12,80 0,38 12,29 13,20 6 13,85 0,98 12,05 15,34 10 14,95 0,84 13,22 15,75 13 FL 8,22 0,53 7,11 9,11 13 8,11 0,30 7,54 8,56 14 7,90 0,53 7,03 8,58 12 7,74 0,69 6,55 8,31 5 9,21 0,90 7,45 10,82 10 9,77 0,52 9,29 10,65 13 TBL 5,80 0,48 4,96 6,43 13 5,24 0,51 4,27 5,90 14 5,45 0,18 5,03 5,68 12 5,68 0,13 5,50 5,88 6 6,16 0,39 5,46 6,79 10 6,42 0,31 5,86 6,78 13 4TL 8,77 0,86 7,34 10,06 13 8,27 0,82 6,60 9,48 14 9,07 0,31 8,46 9,59 12 9,13 0,61 8,39 9,89 6 10,27 0,51 9,43 10,96 10 10,53 0,98 8,63 11,55 13 HFL 21,06 1,81 18,16 24,25 13 17,61 1,50 15,56 20,50 13 19,77 2,02 16,04 22,87 12 20,14 1,16 18,38 21,83 6 22,41 0,33 21,87 22,87 10 23,05 1,43 20,00 24,90 13 VSN 26,63 1,30 25,00 29,00 8 27,00 2,65 24,00 29,00 3 28,57 1,27 27,00 30,00 7 26,00 0,89 25,00 27,00 6 28,30 1,64 26,00 31,00 10 26,75 1,49 24,00 29,00 13 FPN 16,50 0,97 15,00 18,00 10 16,43 0,85 15,00 18,00 14 17,17 1,12 16,00 19,00 12 18,33 1,63 17,00 21,00 6 18,20 2,10 16,00 23,00 10 19,13 1,13 18,00 21,00 13 Lam 22,29 1,25 21,00 25,00 7 18,54 1,56 16,00 21,00 13 19,83 2,04 17,00 24,00 12 19,83 1,72 17,00 22,00 6 21,22 0,83 20,00 22,00 9 21,13 0,99 20,00 23,00 13 STSN 4,63 0,74 4,00 6,00 8 3,50 0,65 3,00 5,00 14 4,00 0,95 2,00 5,00 12 4,33 0,52 4,00 5,00 6 4,60 1,08 3,00 6,00 10 4,13 1,13 3,00 6,00 13 GSN 21,00 1,94 17,00 24,00 9 20,83 1,40 19,00 24,00 12 20,50 2,02 18,00 25,00 12 21,33 1,86 19,00 24,00 6 22,10 1,45 20,00 25,00 10 23,75 1,58 21,00 26,00 13 CSN 7,90 0,57 7,00 9,00 10 7,36 1,08 7,00 11,00 14 7,58 0,90 6,00 9,00 12 7,67 0,82 7,00 9,00 6 8,00 0,47 7,00 9,00 10 8,38 0,52 8,00 9,00 13 SCGN 2,14 1,22 0,00 3,00 7 4,36 1,60 2,00 7,00 14 2,40 0,52 2,00 3,00 10 3,33 1,75 2,00 6,00 6 4,20 1,69 2,00 7,00 10 3,50 1,60 1,00 6,00 13 SCSN 5,00 0,87 3,00 6,00 9 4,64 0,63 3,00 5,00 14 4,92 0,79 4,00 6,00 12 5,50 0,55 5,00 6,00 6 4,50 0,71 3,00 5,00 10 4,88 0,35 4,00 5,00 13 SLSN 6,56 0,53 6,00 7,00 9 6,00 0,39 5,00 7,00 14 6,64 0,51 6,00 7,00 11 6,67 0,82 6,00 8,00 6 6,70 0,48 6,00 7,00 10 7,38 0,92 6,00 8,00 13 Females SVL 45,51 1,12 43,88 47,70 9 41,98 2,65 37,43 46,61 17 45,03 2,48 40,80 48,30 12 44,53 1,22 43,28 46,26 5 51,79 3,47 46,50 58,91 21 44,89 5,90 33,76 52,35 12 TRL 25,74 2,92 21,27 29,35 9 24,38 2,27 20,95 28,35 17 26,75 1,86 23,30 28,70 12 26,92 2,47 24,10 29,72 5 30,87 2,28 27,52 34,80 21 26,15 3,93 18,19 33,16 12 PL 9,07 0,57 8,34 9,97 9 8,73 0,30 8,16 9,21 17 9,01 0,38 8,40 9,70 12 8,66 0,19 8,47 8,89 5 10,18 0,39 9,49 11,06 21 9,35 0,64 7,81 9,94 12 HL 14,15 1,17 12,70 15,90 9 13,81 0,87 11,98 15,45 17 13,89 0,68 12,90 15,20 12 13,73 0,68 12,89 14,76 5 15,99 0,83 14,53 17,32 21 14,11 1,32 11,50 16,31 12 HW 5,74 0,36 5,28 6,32 9 5,39 0,28 4,93 5,84 17 5,35 0,29 4,90 5,80 12 4,90 0,22 4,60 5,21 5 6,15 0,29 5,67 6,81 21 5,78 0,47 4,99 6,34 12 HH 3,96 0,51 3,31 4,64 9 3,76 0,30 3,36 4,45 17 3,96 0,39 3,20 4,50 12 3,23 0,32 2,86 3,65 5 4,39 0,33 3,69 4,95 21 4,13 0,40 3,32 4,60 12 TW 4,31 0,34 3,85 4,89 9 3,71 0,18 3,41 4,07 17 3,81 0,14 3,70 4,00 3 3,47 0,26 3,16 3,75 5 4,37 0,35 3,70 4,98 18 4,12 0,50 3,31 4,89 12 FFL 11,59 0,72 10,73 12,96 9 11,45 0,76 10,43 12,87 17 11,87 0,36 11,30 12,40 12 11,48 0,32 11,02 11,73 4 13,52 0,58 12,23 14,72 21 12,41 0,87 10,83 13,72 12 FL 7,19 0,43 6,52 7,90 9 7,10 0,71 5,71 8,14 17 7,17 0,65 5,80 8,40 12 5,66 0,59 5,06 6,40 5 8,30 0,56 7,29 9,13 21 7,29 0,59 6,16 8,17 12 TBL 5,23 0,30 4,78 5,73 9 4,95 0,29 4,34 5,30 17 5,19 0,21 4,70 5,50 12 5,30 0,29 4,88 5,55 4 5,86 0,29 5,26 6,50 21 5,42 0,28 4,95 5,92 12 4TL 8,79 0,33 8,07 9,22 9 7,76 0,70 6,30 9,14 17 8,75 0,51 7,80 9,50 12 8,34 1,10 7,03 9,36 4 10,08 0,72 8,69 11,56 21 9,05 0,64 7,67 9,93 12 HFL 19,76 1,80 16,50 21,72 9 16,69 1,01 14,94 18,56 17 18,47 1,18 16,40 19,90 12 17,55 0,91 16,80 19,01 5 21,54 1,25 19,59 25,07 21 19,23 1,44 16,47 21,17 12 VSN 30,00 2,12 27,00 32,00 5 29,50 1,00 29,00 31,00 4 30,50 1,07 29,00 32,00 8 29,20 1,92 27,00 32,00 5 31,05 1,86 27,00 35,00 21 29,67 1,83 26,00 32,00 12 Lam 20,60 1,14 19,00 22,00 5 19,00 2,34 13,00 22,00 16 20,00 0,95 18,00 21,00 12 20,20 1,79 18,00 22,00 5 20,65 1,76 18,00 24,00 20 20,58 1,62 18,00 23,00 12 STSN 5,20 1,30 4,00 7,00 5 3,88 0,99 2,00 6,00 17 4,00 0,60 3,00 5,00 12 4,80 0,45 4,00 5,00 5 3,76 0,77 2,00 5,00 21 4,83 1,12 3,00 6,00 12 GSN 21,40 2,61 17,00 23,00 5 20,93 2,02 18,00 25,00 14 19,75 2,05 17,00 24,00 12 19,80 1,64 18,00 21,00 5 21,33 1,74 19,00 24,00 21 22,08 1,78 20,00 27,00 12 CSN 7,71 0,76 7,00 9,00 7 7,13 0,99 5,00 9,00 15 6,92 0,90 6,00 8,00 12 6,80 1,10 5,00 8,00 5 8,24 0,70 7,00 9,00 21 8,50 0,67 8,00 10,00 12 SCGN 2,00 1,41 0,00 3,00 4 4,59 1,46 2,00 8,00 17 3,36 1,63 2,00 7,00 11 2,80 0,84 2,00 4,00 5 4,10 1,30 3,00 7,00 21 2,70 0,95 1,00 4,00 10 SCSN 4,00 1,00 3,00 5,00 7 4,82 0,53 4,00 6,00 17 5,17 0,84 4,00 6,00 12 5,00 1,00 4,00 6,00 5 4,86 0,48 4,00 6,00 21 4,83 0,39 4,00 5,00 12 SLSN 6,14 0,69 5,00 7,00 7 6,24 0,44 6,00 7,00 17 6,67 0,65 6,00 8,00 12 6,20 0,84 5,00 7,00 5 6,67 0,66 5,00 8,00 21 6,83 0,58 6,00 8,00 12 Jebel Sirwa Oukaimeden Tizin Tichka Outabati Jebel Azourki Jebel Ayache Morphological analyses of Atlantolacerta 147 Table 2. Sexual dimorphism in linear measurements and pholidosis within lineages. For each population, mean value of the covariate used to estimate the adjusted means, and adjusted means for each sex are shown. Significant values (p< 0.05) are in bold. Pop / Sex LogTRL LogPL LogHL LogHW LogHH LogTW LogFFL LogFL LogTBL LogHTL LogHFL VSN FPN Lam STSN GSN CSN SCGN SCSN SLSN J.Sirwa ( LogSVL 1.69) Males 1.411 1.050 1.241 0.835 0.692 0.673 1.148 0.977 0.799 1.018 1.355 26.63 16.50 22.29 4.63 21.00 7.90 2.14 5.00 6.56 Females 1.466 0.997 1.193 0.781 0.630 0.625 1.123 0.907 0.759 0.994 1.332 30.00 - 20.60 5.20 21.40 7.71 2.00 4.00 6.14 Oukaimeden (LogSVL 1.62) Males 1.357 0.978 1.159 0.758 0.616 0.589 1.085 0.908 0.715 0.915 1.245 27.00 16.43 18.54 3.50 20.83 7.36 4.36 4.64 6.00 Females 1.383 0.940 1.138 0.730 0.573 0.568 1.057 0.848 0.693 0.887 1.221 29.50 - 19.00 3.88 20.93 7.13 4.59 4.82 6.24 TizinTichka (LogSVL 1.62) Males 1.344 0.975 1.172 0.747 0.621 0.584 1.080 0.896 0.731 0.947 1.250 28.57 17.17 19.83 4.00 20.50 7.58 2.40 4.92 6.64 Females 1.381 0.950 1.138 0.717 0.542 0.582 1.075 0.837 0.698 0.925 1.227 30.50 - 20.00 4.00 19.75 6.92 3.36 5.17 6.67 J.Azourki ( LogSVL 1.64) Males 1.353 1.009 1.205 0.822 0.631 0.694 1.113 0.923 0.776 0.956 1.336 28.30 18.20 21.22 4.60 22.10 8.00 4.20 4.50 6.70 Females 1.389 0.941 1.135 0.740 0.578 0.631 1.052 0.847 0.704 0.929 1.280 31.05 - 20.65 3.76 21.33 8.24 4.10 4.86 6.67 Outabati ( LogSVL 1.62) Males 1.356 0.986 1.200 0.748 0.632 0.574 1.102 0.915 0.746 0.998 1.318 26.00 18.33 19.83 4.33 21.33 7.67 3.33 5.50 6.67 Females 1.403 0.924 1.115 0.683 0.452 0.554 1.063 0.731 0.730 0.872 1.213 29.20 - 20.20 4.80 19.80 6.80 2.80 5.00 6.20 J.Ayache ( LogSVL 1.66) Males 1.403 1.053 1.225 0.850 0.697 0.678 1.168 0.978 0.804 1.017 1.357 26.75 19.13 21.13 4.13 23.75 8.38 3.50 4.88 7.38 Females 1.429 0.980 1.161 0.766 0.624 0.626 1.099 0.872 0.737 0.959 1.288 29.67 - 20.58 4.83 22.08 8.50 2.70 4.83 6.83 Considering only males, there were differences among populations in general size and shape (Table 3), including most of the iso-corrected measurements with the exception of head length, front limb length and tibia length (HL, FLL and TBL, respectively; Table 3; Figure 3). Regarding females, general differences in size and shape were also observed, and all the linear measurements analysed were significantly different among populations with the exception of front limb length (FLL; Table 3; Figure 3). Table 3. Summary of the ANOVA/MANOVA results regarding the effect of sex, population and their interaction. Significant values are in bold. Variables mSize 50.10 <0.01 79.00 <0.01 6.80 <0.01 23.7 <0.01 37.60 <0.01 Shape 6.15 <0.01 28.75 <0.01 1.66 <0.01 3.63 <0.01 4.36 <0.01 Scales 4.14 <0.01 4.42 <0.01 0.94 0.58 2.72 <0.01 2.81 <0.01 Morphometric TrL 4.41 <0.01 237.91 <0.01 5.25 <0.01 5.69 <0.01 5.77 <0.01 PL 5.28 <0.01 19.38 <0.01 2.23 0.06 3.80 <0.01 3.68 <0.01 HL 3.09 0.01 7.23 <0.01 1.29 0.27 1.05 0.40 3.79 <0.01 HW 7.57 <0.01 17.19 <0.01 1.24 0.29 4.74 <0.01 5.02 <0.01 HH 7.49 <0.01 6.03 0.02 1.66 0.15 4.95 <0.01 4.70 <0.01 TW 19.09 <0.01 0.11 0.75 2.21 0.06 19.86 <0.01 4.39 <0.01 FFL 3.97 <0.01 1.73 0.19 0.78 0.56 2.18 0.07 2.38 0.05 FL 8.90 <0.01 46.24 <0.01 1.57 0.17 2.61 0.03 7.43 <0.01 TBL 5.63 <0.01 1.57 0.21 1.21 0.31 1.92 0.10 5.53 <0.01 HTL 5.49 <0.01 3.00 0.09 0.78 0.57 3.87 <0.01 2.98 0.02 HFL 9.52 <0.01 0.00 0.97 0.67 0.65 5.94 <0.01 3.84 <0.01 Scales VSN 4.20 <0.01 34.30 <0.01 0.40 0.82 3.23 0.02 1.41 0.24 FPN --- --- --- --- --- --- 3.30 0.02 --- --- Lam 7.18 <0.01 0.26 0.61 0.40 0.84 8.25 <0.01 2.25 0.07 STSN 2.64 0.03 0.59 0.45 2.33 0.05 1.38 0.26 4.74 <0.01 GSN 6.78 <0.01 0.13 0.72 1.15 0.34 4.66 <0.01 3.57 <0.01 CSN 6.34 <0.01 3.87 0.05 1.89 0.11 1.32 0.29 7.82 <0.01 SCGN 4.42 <0.01 0.32 0.58 0.30 0.91 1.56 0.21 3.77 <0.01 SCSN 1.44 0.22 2.69 0.11 1.66 0.16 1.96 0.12 0.61 0.69 SLSN 3.20 0.01 1.01 0.32 0.37 0.87 2.36 0.07 1.26 0.30 Total Males Females Pop Sex Pop*Sex Pop Pop