\* Taxon: *B. sylvicola*; "X" = the bees from Alberta resembling melanopygus, but with an electrophoretic profile inconsistent with the other *melanopygus*; *B*. *melanopygus* AB = from Alberta; *B. melanopygus* OR/CA = from Oregon and California; *B. edwardsii* OR/CA = from Oregon and California

Table 6. Electrophoresis results for the five enzymes exhibiting either differences between taxa and/or variation within taxa. The other 11 loci were invariant within, and showed no differences between, all taxa. The body of the table gives the number of individual bees of each electromorph. Electromorph mobilities (mm) are standardized relative to those of B. occidentalis (= index 100, Scholl et al. 1990). Modified from Owen et al. (2010).

All bees had identical electrophoretic mobilities, and were invariant at 11 of the 16 enzyme loci examined. Five loci exhibited either differences between taxa and/or variation within taxa (Table 6). The nominate forms of *sylvicola* and *melanopygus* from Alberta clearly have different electrophoretic profiles (Table 6). The electrophoretic profiles of *melanopygus* and *edwardsii* from all locations were entirely consistent with each other. There was a very small amount of variation present, with heterozygotes being

What was really interesting was the six bees ("MEL X"), collected in Alberta, that were assigned to *melanopygus* by eye when they were collected but turned out to have an electrophoretic profile inconsistent with that of *melanopygus* but consistent with that of *sylvicola* (Table 6). Going back to the collection records it was found that these bees (plus another three that were not electrophoresed) came from high elevations in the Kananaskis Valley (Fortress Mountain and Highwood Pass) where typical *sylvicola* had been collected. These were later reassigned to *sylvicola* on the basis of the wing morphometric analysis (see

**Taxon\*** 72/82 82 93 93/100 92/96 96 95 100 102 100 100/105 105 100 105

**(***n***=18)** 1 17 18 16 2 18 1 17

**"X" (***n***=6)** 6 6 6 1 5 6

**OR/CA** *(n***=48***)* 24 2 2 46 35 13 48 48

\* Taxon: *B. sylvicola*; "X" = the bees from Alberta resembling melanopygus, but with an electrophoretic profile inconsistent with the other *melanopygus*; *B*. *melanopygus* AB = from Alberta; *B. melanopygus* OR/CA = from Oregon and California; *B. edwardsii* OR/CA = from Oregon and California

Table 6. Electrophoresis results for the five enzymes exhibiting either differences between taxa and/or variation within taxa. The other 11 loci were invariant within, and showed no differences between, all taxa. The body of the table gives the number of individual bees of each electromorph. Electromorph mobilities (mm) are standardized relative to those of B.

occidentalis (= index 100, Scholl et al. 1990). Modified from Owen et al. (2010).

*Pgm Gpi Idh* (NAD) *Hk Sdh*

16 16 13 3 16 16

23 2 1 24 15 10 25 25

 **Enzyme electromorph**

detected at a few locations (Table 6).

below).

*B. sylvicola*

*melanopygus***" AB (***n* **= 16)**

*melanopygus"*  **OR/CA (***n***=25)**

*"B. edwardsii"*

*"B.*

*"B.*

Fig. 9. Plot of the first two Canonical scores for *B. sylvicola (*S), the Alberta *B. melanopygus* (M) and the anomalous Alberta *B. melanopygus* ("Mel X"). Modified from Owen et al. (2010).

Fig. 10. Plot of the first two Canonical scores for the total data set. S = *B. sylvicola*, R = red "*melanopygus*" from Oregon and California, B = black *"edwardsii*" from Oregon and California, M = *B. melanopygus* from Alberta. Modified from Owen et al. (2010).

The discriminant functions analysis was run three times. Initially only specimens from Alberta were included. This was to verify that the technique could separate closely related species (*melanopygus* and *sylvicola*) in sympatry, and to determine the status of the aberrant *melanopygus* ("MEL X"). In addition to the six "MEL X" bees that were electrophoresed (Table 6) three other queens that were collected on the same dates and at the same locations were reassigned from *melanopygus* and included in the "MEL X" category. The plot of the first two canonical scores is shown in Figure 9. *B. melanopygus* is clearly separated from *B. sylvicola* by

Applications of Morphometrics to the Hymenoptera, Particularly Bumble Bees (*Bombus*, Apidae) 25

(*Pyrobombus*) *perplexus*. Ninety two bands using four different PCR primers were generated. Fixed differences occurred between all groups, and individual haplotypes did not occur in more than one taxonomic group, although there was overlap in haplotype components. The corrected average number of pairwise differences of between *B. moderatus* and *B. terricola* and *B. moderatus* and *B. occidentalis* was 6.98 and 5.92 respectively, and that between *B.*

**Species (n)** *B. terricola B. occidentalis B. moderatus B. perplexus*

*B. terricola* (87) 1.28 6.27 7.91 54.76 *B. occidentalis* (79) 5.07 1.11 6.77 53.87 *B. moderatus* (104*)* 6.98 5.92 0.59 53.55 *B. perplexus* (54) 53.21 52.41 52.34 1.81

Table 7. Average pairwise differences between and within bumble-bee species.Above diagonal: Average number of pairwise differences between groups (*PiXY*).Diagonal elements: Average number of pairwise differences within groups (*PiX)*. Below diagonal: Corrected average number of pairwise differences (*PiXY-( PiY)/2*). Sample sizes are given in

> *occidentalis* 1986

*terricola*1985 2 0 17 4 23 91.3% *terricola* 1986 1 3 6 15 25 84.0%

**Total** 36 57 29 29 151 85.4%

Table 8**.** Classification count (actual and predicted) of the *B. occidentalis* and *B. terricola* from

The taxa are clearly separated by both the genetic and morphological evidence. The corrected average number of pairwise differences of between *B. moderatus* and *B. terricola* and *B. moderatus* and *B. occidentalis* was 6.98 and 5.92 respectively, and that between *B. occidentalis* and *B. terricola* was 5.07. Therefore since *B. terricola* and *B. occidentalis* are

Traditional morphometric analysis was done on some specimens of *B. occidentalis* and *B. terricola* queens collected in 1985 and 1986. The left forewing was removed and measured using the methods of Owen et al. (2010) as described above, and discriminant analysis performed. The classification counts are given in Table 8, and the first and third canonical

> *terricola*  1985

16 15 4 5 40 77.5%

17 39 2 5 63 88.9%

*terricola*  1986

Total (n)

% correctly classified

*occidentalis* and *B. terricola* was 5.07 (Table 7).

parentheses.

*occidentalis* 1985

*occidentalis* 1986

scores are plotted in Fig. 11.

**Actual1** *occidentalis*

**Species Predicted**

1985

1 Reduction in classification error due to variables measured = 43.5%.

1985 and 1986 using discriminant analysis of wing venation.

the first canonical score. Similarly the "MEL X" bees are obviously distinct from *melanopygus* and are grouped with *sylvicola*. Next, the analysis was run using the complete data set (Figure 10) with the "MEL X" bees now being reclassified as *sylvicola*. Again, *B. sylvicola* is clearly separated by canonical score one, but *melanopygus* and *edwardsii* are not obviously resolved.

Enzyme electrophoresis and wing morphometrics failed to distinguish the nominate species *B. edwardsii* and *B. melanopygus*, yet clearly separated *B. sylvicola* from the latter. This, together with the colour dimorphism genetic data (Owen and Plowright 1980), and the lack of other morphological differences led Owen et al. (2010) to conclude that *melanopygus* and *edwardsii* are conspecific. If *B. melanopygus* is a "good" species, then there is a gene frequency cline for the color dimorphism (Fig. 8).

#### **4.2** *B. occidentalis/B. terricola*

Two other taxa, where the evolutionary status and taxonomic classification are also unclear, are *B. terricola* Kirby and *B. occidentalis* Greene. The basis of this confusion originates with their classification being based primarily on pile colour pattern. Greene's original description of *B. occidentalis* reads "…first four abdominal segments black… "(Franklin 1913). Given that this is the type specimen description, specimens with the first four abdominal segments being black should be considered 'typical' *B. occidentalis.* In contrast, typical *B. terricola* have TIII and TIV that are consistently and clearly defined by complete yellow bands, and lack the large amount of white to cream-coloured pile typical of *B. occidentalis* on TV and TVI. However, in some parts of its distribution including areas of overlap with *B. terricola*, *B. occidentalis* exhibits considerable pile colour variation with some specimens closely resembling *B. terricola* (Stephen, 1957; Milliron, 1971). The primary ambiguous components of these bees are the complete to incomplete yellow bands on gastral terga III and IV. Nevertheless Stephen (1957) noted that *B. terricola* was "one of the most color stable species in western America" (p. 82) showing little or no variation throughout its range, and that it could be distinguished from *B. occidentalis* in having TII always yellow and TIV black. On this basis many authors have regarded *B. occidentalis* and *B. terricola* to be separate species (Stephen 1957; Thorp et al. 1983). However, Milliron (1971) reduced *B. occidentalis* to subspecific status under *B. terricola*, citing a lack of evident reliable or constant morphological features by which to differentiate specimens in areas of overlap. Milliron (1971) also suggested that these two subspecies most probably interbreed, producing numerous perplexing subspecific hybrids. This is certainly one possible explanation for the rare occurrence of colonies headed by definite *B. occidentalis* queens which produce *B. terricola*-like offspring.

Recently Bertsch et al. (2010) sequenced part (1005 bp) of the mitochondrial cytochrome oxidase subunit I (COI) gene and found a difference of 30 nucleotides between *B. occidentalis* and *B. terricola*, which is significantly larger than that found within a species. On this basis Bertsch et al. (2010) concluded that *B. occidentalis* and *B. terricola* do represent good biological species. They also suggested that to clarify the situation these taxa should be studied in greater detail in their area of contact in British Columbia and southern Alberta.

Whidden (2002) studied sympatric populations of *B. occidentalis* and *B. terricola* in Alberta using randomly amplified polymorphic DNA (RAPD) analysis. For comparison he also analyzed one consubgeneric species, *B. moderatus*, and one non-consubgeneric species *B.*

the first canonical score. Similarly the "MEL X" bees are obviously distinct from *melanopygus* and are grouped with *sylvicola*. Next, the analysis was run using the complete data set (Figure 10) with the "MEL X" bees now being reclassified as *sylvicola*. Again, *B. sylvicola* is clearly separated by canonical score one, but *melanopygus* and *edwardsii* are not obviously resolved.

Enzyme electrophoresis and wing morphometrics failed to distinguish the nominate species *B. edwardsii* and *B. melanopygus*, yet clearly separated *B. sylvicola* from the latter. This, together with the colour dimorphism genetic data (Owen and Plowright 1980), and the lack of other morphological differences led Owen et al. (2010) to conclude that *melanopygus* and *edwardsii* are conspecific. If *B. melanopygus* is a "good" species, then there is a gene frequency

Two other taxa, where the evolutionary status and taxonomic classification are also unclear, are *B. terricola* Kirby and *B. occidentalis* Greene. The basis of this confusion originates with their classification being based primarily on pile colour pattern. Greene's original description of *B. occidentalis* reads "…first four abdominal segments black… "(Franklin 1913). Given that this is the type specimen description, specimens with the first four abdominal segments being black should be considered 'typical' *B. occidentalis.* In contrast, typical *B. terricola* have TIII and TIV that are consistently and clearly defined by complete yellow bands, and lack the large amount of white to cream-coloured pile typical of *B. occidentalis* on TV and TVI. However, in some parts of its distribution including areas of overlap with *B. terricola*, *B. occidentalis* exhibits considerable pile colour variation with some specimens closely resembling *B. terricola* (Stephen, 1957; Milliron, 1971). The primary ambiguous components of these bees are the complete to incomplete yellow bands on gastral terga III and IV. Nevertheless Stephen (1957) noted that *B. terricola* was "one of the most color stable species in western America" (p. 82) showing little or no variation throughout its range, and that it could be distinguished from *B. occidentalis* in having TII always yellow and TIV black. On this basis many authors have regarded *B. occidentalis* and *B. terricola* to be separate species (Stephen 1957; Thorp et al. 1983). However, Milliron (1971) reduced *B. occidentalis* to subspecific status under *B. terricola*, citing a lack of evident reliable or constant morphological features by which to differentiate specimens in areas of overlap. Milliron (1971) also suggested that these two subspecies most probably interbreed, producing numerous perplexing subspecific hybrids. This is certainly one possible explanation for the rare occurrence of colonies headed by definite *B. occidentalis* queens which produce *B.*

Recently Bertsch et al. (2010) sequenced part (1005 bp) of the mitochondrial cytochrome oxidase subunit I (COI) gene and found a difference of 30 nucleotides between *B. occidentalis* and *B. terricola*, which is significantly larger than that found within a species. On this basis Bertsch et al. (2010) concluded that *B. occidentalis* and *B. terricola* do represent good biological species. They also suggested that to clarify the situation these taxa should be studied in greater detail in their area of contact in British Columbia and southern Alberta. Whidden (2002) studied sympatric populations of *B. occidentalis* and *B. terricola* in Alberta using randomly amplified polymorphic DNA (RAPD) analysis. For comparison he also analyzed one consubgeneric species, *B. moderatus*, and one non-consubgeneric species *B.*

cline for the color dimorphism (Fig. 8).

**4.2** *B. occidentalis/B. terricola*

*terricola*-like offspring.

(*Pyrobombus*) *perplexus*. Ninety two bands using four different PCR primers were generated. Fixed differences occurred between all groups, and individual haplotypes did not occur in more than one taxonomic group, although there was overlap in haplotype components. The corrected average number of pairwise differences of between *B. moderatus* and *B. terricola* and *B. moderatus* and *B. occidentalis* was 6.98 and 5.92 respectively, and that between *B. occidentalis* and *B. terricola* was 5.07 (Table 7).

Table 7. Average pairwise differences between and within bumble-bee species.Above diagonal: Average number of pairwise differences between groups (*PiXY*).Diagonal elements: Average number of pairwise differences within groups (*PiX)*. Below diagonal: Corrected average number of pairwise differences (*PiXY-( PiY)/2*). Sample sizes are given in parentheses.

Traditional morphometric analysis was done on some specimens of *B. occidentalis* and *B. terricola* queens collected in 1985 and 1986. The left forewing was removed and measured using the methods of Owen et al. (2010) as described above, and discriminant analysis performed. The classification counts are given in Table 8, and the first and third canonical scores are plotted in Fig. 11.

1 Reduction in classification error due to variables measured = 43.5%.

Table 8**.** Classification count (actual and predicted) of the *B. occidentalis* and *B. terricola* from 1985 and 1986 using discriminant analysis of wing venation.

The taxa are clearly separated by both the genetic and morphological evidence. The corrected average number of pairwise differences of between *B. moderatus* and *B. terricola* and *B. moderatus* and *B. occidentalis* was 6.98 and 5.92 respectively, and that between *B. occidentalis* and *B. terricola* was 5.07. Therefore since *B. terricola* and *B. occidentalis* are

Applications of Morphometrics to the Hymenoptera, Particularly Bumble Bees (*Bombus*, Apidae) 27

Aytekin A., M., Terzo, M., Rasmont P. & Çağatay, N., (2007) Landmark Based Geometric

*Bombus* Latreille). Annales de la Société Entomologique de France, 43: 95-102. Aytekin,A. M., Rasmont, P. & Çagatay, N. (2003) Molecular and morphometric variation in

Baylac, M., Villemant, C. & Simbolott, G. (2003) Combining geometric morphometrics with

Bertsch, A., Hrabe de Angelis, M. & Przemeck, G.K.H. (2010) A phylogenetic framework for

Cameron, S.A., Lozier, J.D., Strange, J.P., Koch , J.B., Cordes, N., Solter, L.F. & Griswold, T.L.

Clarke, G.M. (1997) The genetic basis of developmental stability. III. Haplo-diploidy: Are

Colla, S.R. & Packer, L. (2008) Evidence for decline in eastern North American bumblebees

Crespi, B.J. & Vanderkist, B.A. (1997) Fluctuating asymmetry in vestigial and functional

Danaraddi, C.S. & Viraktamath, S. (2009) Morphometric studies on the stingless bee, *Trigona iridipennis* Smith. Karnataka Journal of Agricultural Science, 22: 796-797. Eickwort, K.R. (1969) Differential variation of males and females in *Polistes exclamans*.

Evans, E., Thorpe, R., Jepson, S. & Black, S.H. (2008) Status review of three formerly

Fisher, R.A. (1936) The use of multiple measurements in taxonomic problems. Annals of

Franklin, H.J. (1913) The Bombidae of the New World. Transactions of the American

Francoy, T.M., Wittmann, D., Steinhage, V., Drauschke, M., Müller, S., Cunha, D.R.,

Gelin, L.F.F., Da Cruz, J.D., Noll, F.B., Giiannotti, E., Dem Santos, G.M., & Bichara-Filho,

common species of bumble bee in the subgenus *Bombus*. 63 pp. The Xerces Society,

Nascimento, A.M., Figueiredo, V.L.C. Simões ,Z.L.P., De Jong, D., Arias,M.C. & Gonçalves, L.S. (2009) Morphometric and genetic changes in a population of *Apis mellifera* after 34 years of Africanization. Genetics and Molecular Research, 8: 709-717. Fruciano, C., Tigano, C. & Ferrito, V. (2011) Traditional and geometric morphometrics detect

morphological variation of lower pharyngeal jaw in *Coris julis* (Teleostei, Labridae).

C.C. (2008) Morphological Caste Studies In The Neotropical Swarm-Founding

males more unstable than females? Evolution 51: 2012-2028.

Cresson (Hymenoptera: Apidae). Canadian Journal of Zoology, 66: 1172-1178.

traits of a haplodiploid insect. Heredity, 79: 624-630.

Torre (Hymenoptera: Apidae) Mellifera 3:34-40

markers. Beitrage zur Entomologie, 60, 229-242.

the Linnean Society 80: 89-98.

*January 3, 2011*

1379-1391

Evolution, 23: 391-405.

Portland, Oregon.

Eugenics, 7: 179-188.

Entomological Society 38:177-486.

Italian Journal of Zoology, 78: 320-327.

morphometric analysis of wing shape in *Sibiricobombus* Vogt (Hymenoptera: Apidae:

Bombus terrestris Lucoformis Krüger and Bombus terrestris Dalmatinus Dalla

pattern recognition for the investigation of species complexes. Biological Journal of

the North American species of the subgenus *Bombus* sensu stricto (*Bombus affinis, B. franklini, B. moderatus, B. occidentalis* & *B. terricola*) based on mitochondrial DNA

(2011) Patterns of widespread decline in North American bumble bees. Published online before print January 3, 2011, doi: 10.1073/pnas.1014743108, *PNAS, U.S.A.*

(Hymenoptera:Apidae), with special focus on *Bombus affinis*. Biodiv and Cons 17:

differentiated from each other to the same extent as they are from *B. moderatus*, they should regarded as distinct taxa. Discriminant function analysis of wing morphometric data correctly classified over 85% of the specimens of *B. occidentalis* and *B. terricola*, indicating significant morphological divergence.

Fig. 11. Plot of the first and third Canonical scores for the 1985 and 1986 specimens of *B. occidentalis* (O) and *B. terricola* (T).

#### **5. Conclusions**

Morphometric analysis has been applied in a number of different ways to problems in the Hymenoptera and has proved to have an important and useful set of techniques for answering interesting questions. It is particularly useful for species identification and classification. The more traditional approaches appear to be as sensitive as geometric morphometrics for many problems. A powerful approach is to combine morphometric genetic methods, particularly to help answer questions of systematic and taxonomy.

#### **6. Acknowledgements**

Funding for some of this research was provided by a Mount Royal University Internal Research grant. I thank Daniel Lee-Owe for measuring the bee wings.

#### **7. References**

Adams, D.C., Rohlf, F.J. & Slice, D.E. (2004) Geometric morphometrics: Ten years of progress following the 'revolution'. Italian Journal of Zoology, 71: 5-16.

Alexander, B. (1991) A cladistic analysis of the genus *Apis*. In: *Diversity in the genus* Apis*,*  Smith, D.R. (ed.) pp. 1-28, Westview Press, Boulder, Colorado.

differentiated from each other to the same extent as they are from *B. moderatus*, they should regarded as distinct taxa. Discriminant function analysis of wing morphometric data correctly classified over 85% of the specimens of *B. occidentalis* and *B. terricola*, indicating

Fig. 11. Plot of the first and third Canonical scores for the 1985 and 1986 specimens of *B.*

Morphometric analysis has been applied in a number of different ways to problems in the Hymenoptera and has proved to have an important and useful set of techniques for answering interesting questions. It is particularly useful for species identification and classification. The more traditional approaches appear to be as sensitive as geometric morphometrics for many problems. A powerful approach is to combine morphometric

Funding for some of this research was provided by a Mount Royal University Internal

Adams, D.C., Rohlf, F.J. & Slice, D.E. (2004) Geometric morphometrics: Ten years of progress following the 'revolution'. Italian Journal of Zoology, 71: 5-16. Alexander, B. (1991) A cladistic analysis of the genus *Apis*. In: *Diversity in the genus* Apis*,*

Smith, D.R. (ed.) pp. 1-28, Westview Press, Boulder, Colorado.

genetic methods, particularly to help answer questions of systematic and taxonomy.

Research grant. I thank Daniel Lee-Owe for measuring the bee wings.

significant morphological divergence.

*occidentalis* (O) and *B. terricola* (T).

**5. Conclusions**

**6. Acknowledgements**

**7. References**

Applications of Morphometrics to the Hymenoptera, Particularly Bumble Bees (*Bombus*, Apidae) 29

Owen, R.E. & Plowright, R.C. (1980) Abdominal pile color dimorphism in the bumble bee

Owen, R.E. & Plowright, R.C. (1988) Inheritance of metasomal pile colour variation in the

Owen, R.E. (1986) Gene frequency clines at X-linked or haplodiploid loci. *Heredity*, 57:209-219. Owen, R.E. (1989) Differential size variation of male and female bumble bees (Hymenoptera,

Owen, R.E. 1988 Body size variation and optimal body size of bumble bee queens

Owen, R.E., T. L. Whidden & Plowright, R.C. (2010) Genetic and morphometric evidence for

Palmer, A. R. (1994) Fluctuating asymmetry analyses: A primer. In: T. A. Markow (ed.) pp.

Palmer, A. R., & Strobeck, C. (1986) Fluctuating asymmetry: measurement, analysis,

Parsons, P.A. (1990) Fluctuating asymmetry: an epigenetic measure of stress. Biological

Pekkarinen, A. (1979) Morphometric, colour and enzyme variation in bumblebees

Plowright, R.C. & Owen, R.E. (1980) The evolutionary significance of bumble bee color

Plowright, R.C. & Pallett, M.J. (1979) A morphometric study of the taxonomic status of *Bombus sandersoni* (Hymenoptera: Apidae) Canadian Entomologist, 110: 647-654. Plowright, R.C., & Stephen, W.P. (1973) A numerical taxonomic analysis of the evolutionary

Plowright, R.C., & Stephen, W.P. (1980) The taxonomic status of *Bombus franklini*

Pretorius, E. (2005) Using geometric morphometrics to investigate wing dimorphism in males

Rohlf F. J. & Marcus L. F. (1993) A revolution in morphometrics. Trends Ecol. Evol., 8: 129-132 Rohlf, F.J. & Slice, D. (1990) Extensions of the procrustes method for the optimal

Rosenberg, M.S. (1997) Evolution of shape: differences between the major and minor

chelipeds of *Uca pugnax* (Decapoda: Ocypodidae). Journal of Crustacean Biology,

patterns. Annual Review of Ecology and Systematics*,* 17: 391-421.

Reviews of the Cambridge Philosophical Society, 65: 131-145.

patterns: A mimetic interpretation. Evolution, 34, 622-637.

(Hymenoptera: Apidae). Canadian Entomologist, 112: 475-479.

superimposition of landmarks, Systematic Zoology, 39: 40-59.

the conspecific status of the bumble bees *Bombus melanopygus* and *B. edwardsii* (Hymenoptera, Apidae). Journal of Insect Science, 10:109, available online:

335-364, *Developmental Instability: Its Origins and Evolutionary Implications.* Kluwer,

(Hymenoptera, Apidae, *Bombus*) in Fennoscandia and Denmark. Acta Zoologica

relationships of *Bombus* and *Psithyrus* (Apidae: Hymenoptera). Canadian

and females of Hymenoptera – a case study based on the genus *Tachysphex* Kohl (Hymenoptera: Specidae: Larinae). Australian Journal of Entomology, 44: 113-121. Prŷs-Jones, O.E. & Corbet, S.A. (1987) *Bumblebees*. Cambridge University Press, Cambridge. Richards, K.W. (1984) Alfalfa leafcutter bee management in Western Canada. Agriculture

*Bombus melanopygus*. Journal of Heredity, 71: 241-247.

(Hymenoptera: Apidae). Canadian Entomologist, 120:19-27.

Apidae, *Bombus*). Journal of Heredity, 80:39-43.

bumble bee *Bombus rufocinctus*

insectscience.org/10.109.

Dordrecht, Netherlands.

Fennica 158: 1-60.

Entomologist, 105: 733-743.

Canada Publication 1495E: 1-53.

17: 52-59.

Polistinae Wasp *Angiopolybia pallens* (Lepeletier) (Hymenoptera: Vespidae) *Neotropical Entomology 37: 691-701*

Goulson, D., Lye, & B. Darvill, B. (2008) Decline and conservation of bumble bees. Ann. Rev.

Grixti, J.C., Wonga, L.T., Cameron, S.A. & Favreta, C. (2009) Decline of bumble bees

Harder, L.D. (1982) Measurement and estimation of functional proboscis length in

Harder, L.D. (1985) Morphology as a predictor of flower choice by bumble bees. Ecology 66:

Hastings, J. M., C.W. Holliday and J. R. Coelho (2008) Body size relationship between

Hintze, J.L. (1996) *NCSS 6.0.21-2 Statistical System for Windows. User's Manual II*. Hintze,

Huber, J.T. (2009) Biodiversity of Hymenoptera. In: Foottit, R.G. and Adler, P.H. (eds.), pp. 303-323, *Insect Biodiversity; Science and Society*, Wiley-Blackwell, Chichester, UK. Huxley, J. (1972) *Problems of Relative Growth 2nd edition,* Dover Publications Inc., New York. Hastings, J. M., C.W. Holliday and J. R. Coelho. (2010) Size-specific provisioning by cicada

Klingenberg, C.P., Badyaev, A.V., Sowry, S.M. & Beckwith, N.J. (2001) Inferring

Kozmus, P., Virant-Doberlet, M., Meglič & Dovč, P. (2011) Identification of *Bombus* species

Marcus, J.M. (2001) The development and evolution of crossveins in insect wings. Journal of

Medler, J.T. (1962) Morphometric studies on bumble bees. Annals of the Entomological

Michener, C.D. (2000) *The Bees of the World*. Johns Hopkins University Press, Baltimore and

Milankov, V., Francuski, L., Ludoški, J., Stăhls, G. & Vujié, A. (2010) Estimating genetic and

Milliron, H.E. (1971) A monograph of the western hemisphere bumblebees (Hymenoptera:

O'Neill, K.M., ,Pearce, A.M., O'Neill, R.P, & Miller, R.S. (2010) Offspring Size and Sex Ratio

phenotypic diversity in a northern hoverfly reveals lack of heterozygosity correlated with significant fluctuating asymmetry of wing traits. Journal of Insect

Apidae; Bombinae) - The Genera Bombus, Megabombus Subgenus Bombias.

Variation in a Feral Population of Alfalfa Leafcutting Bee*s* (Hymenoptera: Megachilidae) Annals of the Entomological Society of America, 103: 775- 784 Owen, R.E. & McCorquodale, D.B. (1994) Quantitive variation and heritability of

postdiapause development time and body size in the Alfalfa leafcutting bee (Hymenoptera: Megachilidae). Annals of the Entomological Society of America,

*Sphecius speciosus* (Hymenoptera: Crabronidae), and their prey: Prey size

killers, *Sphecius speciosus*, (Hymenoptera : Crabronidae) in north Florida. Florida

developmental modularity from morphological integration: Analysis of individual variation and asymmetry in bumblebee wings. The American Naturalist, 157: 11-23.

(*Bombus*) in the North American Midwest. Biol. Cons. 142: 75-84

bumblebees. Canadian Journal of Zoology, 60: 1073-1079.

determines wasp size. Florida Entomologist 91: 657-663.

based on wing venation structure. Apidologie 42: 472-480.

*Memoirs of the Entomological Society of Canada,* 82, 1-80.

*Neotropical Entomology 37: 691-701*

Ent*.* 53:191–208

Kaysville, Utah.

Entomologist 93: 412-421.

Anatomy, 199: 211-216.

Conservation, 14: 77-88.

London.

87:922-927.

Society of America, 55: 212-218.

198-210.

Polistinae Wasp *Angiopolybia pallens* (Lepeletier) (Hymenoptera: Vespidae)

**0**

**2**

<sup>1</sup>*France* 2,3*USA*

**The Mosquito Fauna: From Metric**

<sup>1</sup>*UMR 5090 MIVEGEC, Avenue Agropolis, IRD, Montpellier*

Jean-Pierre Dujardin1, P. Thongsripong2 and Amy B. Henry<sup>3</sup>

<sup>2</sup>*Department of Tropical Medicine, Medical Microbiology, and Pharmacology,*

<sup>3</sup>*Department of Microbiology, University of Hawaii at Manoa, Honolulu, Hawaii*

Biodiversity, encompassing diversity of genes, species, and ecosystems, is fundamental to biology (Gaston & Spicer, 2004), yet tools to monitor it are insufficient. Biodiversity can be estimated by using the number of species (species richness) in a community, and/or by this number together with the proportion of species (species evenness), and/or by other more

Studies on morphological and biological diversity have highlighted the complexity of the diversity structure, that is, the relationship between morphological and taxonomic diversity (Foote, 1992). Most studies acknowledge a certain level of dissociation between morphological diversity and species richness, suggesting that taxonomic and morphological diversity patterns are distinct ones (Foote, 1993; Moyne & Neige, 2007; Roy et al., 2001; Roy & Foote, 1997; Vasil'ev et al., 2010). Occasionally, the use of metric diversity as a proxy for species

However, appropriate morphological disparity metrics and sets of morphological characters are not clearly defined (Navarro, 2003; Roy & Foote, 1997; Wills et al., 1994). Although most studies have used one or two measures of disparity to quantify and characterize the occupation of morphospace, multiple measures might be necessary to fully detect changes in

Moreover, organisms present an indefinitely large number of potentially quantifiable traits, and in practice only a small number of features can be studied. Therefore, we cannot strictly measure morphological diversity, but diversity with respect to some set of traits. The common practice is to seek broad coverage of morphology, using as many characters as possible. The present study is disputing this common practice, showing that the opposite strategy could be

Our hypothesis is that when grouping many traits morphological variation becomes too complex to validly reflect a single factor like species richness. Organismal morphology is the result of many biological causes that are not only evolutionary factors but also environmental

indirect estimators among which is morphological variation.

patterns of morphospace occupation (Ciampaglio et al., 2001).

more informative: to seek for elementary coverage of morphology.

richness has been suggested (Dolan et al., 2006).

**1. Introduction**

**Disparity to Species Diversity**

*University of Hawaii at Manoa, Honolulu, Hawaii*

## **The Mosquito Fauna: From Metric Disparity to Species Diversity**

Jean-Pierre Dujardin1, P. Thongsripong2 and Amy B. Henry<sup>3</sup> <sup>1</sup>*UMR 5090 MIVEGEC, Avenue Agropolis, IRD, Montpellier* <sup>2</sup>*Department of Tropical Medicine, Medical Microbiology, and Pharmacology, University of Hawaii at Manoa, Honolulu, Hawaii* <sup>3</sup>*Department of Microbiology, University of Hawaii at Manoa, Honolulu, Hawaii* <sup>1</sup>*France* 2,3*USA*

#### **1. Introduction**

30 Morphometrics

Rothschild, M. (1979) Factors influencing size and sex ratio in *Megachile rotundata*

Ruttner, F. (1986) Geographical variability and classification. In: *Bee Genetics and Breeding*,

Scholl, A., E. Obrecht & Owen, R.E. (1990) The genetic relationship between *Bombus*

Scholl, A., Thorp, R.W., Owen, R.E. & Obrecht, E. (1992) Specific distinctiveness of *Bombus*

Sharkey, M.J. (2007) Phylogeny and classification of Hymenoptera. Zootaxa, 1668, 521–548 Shih, C., Feng, H., Liu, C., Zhao, Y., & Ren, D. (2010) Morphology, phylogeny, evolution,

Slice, D.E., Bookstein, F.L., Marcus, L. F. & Rohlf, F.J. (Revised Feb. 12, 2009) A Glossary for

Smith, D.R., Crespi, B.J. & Bookstein, F.L. (1997) Fluctuating asymmetry in the honey bee,

Sneath, P.H.A. & Sokal, R.R. (1973) *Numerical Taxonomy*. W.H. Freeman & Co., San

Stephen, W.P. (1957) Bumble Bees of Western America - (Hymenoptera: Apoidea). *Technical*

Thorp, R.W., Horning, D.S., and Dunning, L.L. (1983) Bumble bees and cuckoo bumblebees of California (Hymenoptera: Apidae). *Bulletin of the California Insect Survey* 23:1-79. Tofilski, A. (2004) DrawWing, aprogram for numerical description of insect wings. 5pp. *Journal of Insect Science*, 4:17, Available online: http://www.insectscience.org/4.17 Tüzün, A. (2009) Significance of wing morphometry in distinguishing some of the hymenoptera species. African Journal of Biotechnology 8: 3353-3363. Whidden, T.L. (2002) Applications of randomly amplified and microsatellite DNA to

Williams, P.H. (1998) An annotated checklist of bumble bees with an analysis of patterns of

Wilson, E.O. (1971) *The Insect Societies*. The Belknap Press of Harvard University Press,

years. Annals of the Entomological Society of America, 103: 875-885. Silva, M.C., Lomônaco, C., Augusto, S.C. & Kerr, W.E. (2009) Climatic and anthropic

30.09.2011, Available from: http://life.bio.sunysb.edu/morph/

Rinderer, T.E. (ed.) pp. 23-56 Academic Press, inc., New York.

Journal of Zoology, 68, 2264-2268.

Research, 8: 730-737.

10: 551-574.

Franscisco.

Calgary.

*of London* 67, 79-152.

Cambridge, Mass.

*Bulletin,* 40, 2-163.

relationships. Zoologica Scripta, 31, 57–66

401.

(Hymenoptera, Megachilide). Journal of the Kansas Entomological Society, 52: 392-

*moderatus* Cresson and the *Bombus Iucorum* Auct. species complex. Canadian

*nevadensis* Cresson and *B.auricomus* (Robertson)(Hymenoptera: Apidae) - enzyme electrophoretic data. Journal of the Kansas Entomological Society, 65, 134-140. Sharkey, M.J. & Roy, A. (2002) Phylogeny of the Hymenoptera: a reanalysis of the Ronquist

*et al.* (1999) reanalysis, with an emphasis on wing venation and apocritan

and dispersal of pelecinid wasps (Hymenoptera: Pelecinidae) over 165 million

influence on size and fluctuating asymmetry of Euglossine bees (Hymenoptera, Apidae) in a semideciduous seasonal forest reserve. Genetics and Molecular

Geometric Morphometrics In: *Morphometrics at SUNY Stony Broo,k* Date of access,

*Apis mellifera*: effects of ploidy and hybridization. Journal of Evolutionary Biology,

problems in bumble bee biology. 155 pp. Doctoral dissertation, University of

description (Hymenoptea: Apidae, Bombini). *Bulletin of the National History Museum*

Biodiversity, encompassing diversity of genes, species, and ecosystems, is fundamental to biology (Gaston & Spicer, 2004), yet tools to monitor it are insufficient. Biodiversity can be estimated by using the number of species (species richness) in a community, and/or by this number together with the proportion of species (species evenness), and/or by other more indirect estimators among which is morphological variation.

Studies on morphological and biological diversity have highlighted the complexity of the diversity structure, that is, the relationship between morphological and taxonomic diversity (Foote, 1992). Most studies acknowledge a certain level of dissociation between morphological diversity and species richness, suggesting that taxonomic and morphological diversity patterns are distinct ones (Foote, 1993; Moyne & Neige, 2007; Roy et al., 2001; Roy & Foote, 1997; Vasil'ev et al., 2010). Occasionally, the use of metric diversity as a proxy for species richness has been suggested (Dolan et al., 2006).

However, appropriate morphological disparity metrics and sets of morphological characters are not clearly defined (Navarro, 2003; Roy & Foote, 1997; Wills et al., 1994). Although most studies have used one or two measures of disparity to quantify and characterize the occupation of morphospace, multiple measures might be necessary to fully detect changes in patterns of morphospace occupation (Ciampaglio et al., 2001).

Moreover, organisms present an indefinitely large number of potentially quantifiable traits, and in practice only a small number of features can be studied. Therefore, we cannot strictly measure morphological diversity, but diversity with respect to some set of traits. The common practice is to seek broad coverage of morphology, using as many characters as possible. The present study is disputing this common practice, showing that the opposite strategy could be more informative: to seek for elementary coverage of morphology.

Our hypothesis is that when grouping many traits morphological variation becomes too complex to validly reflect a single factor like species richness. Organismal morphology is the result of many biological causes that are not only evolutionary factors but also environmental

80,000 mosquitoes were collected. For the morphometric study, only a subset of the female

The Mosquito Fauna: From Metric Disparity to Species Diversity 33

Genus Species *n* F FF R RF SU U *Aedes aegypti 12* . . 13 5 3

*Aedomyia catasticta 1* . 1. . . . *Anopheles baezai 1* .... 1.

*Armigeres dentatus 1* 1. . . . .

*Coquillettidia crassipes 8* .. .8. .

*Culex bitaeniorhynchus 22* 54 14 7 1

*Ficalbia minima 3* . 3 .. *Heizmania* unknown *7* 7 . .. *Mansonia annulifera 5* . 3. . 2

*albopictus 7* . 4 21 . . *lynetopennis 14* . 4 10 . . . *mediolineatus 12* . 6 6. . . *vexans 25* . 7 16 . 2 . unknown *20* 15 72 5 .

*barbirostris 10* 23. 1 3 1 *kochi 6* . 6. . . . *minimus 2* . 2. . . . *peditaeniatus 5* . . 12 2 . *phillippines 5* . 3 2. . . *tessellatus 9* . 4. . 5 . *vagus 37* . 21 12 . 1 3 unknown *7* . 7.  1. . . . *subalbatus 10* 24 1. 1 2 unknown *2* . 2. . . .

unknown *5* . 5. . . .

*brevipalpis 7* . 2 22 . 1 *fuscocephala 6* . 5. . . 1 *gelidus 18* . 4 15 4 4 *mocthogenes 1* . 1. . . . *nigropunctatus 12* 10 . . . 1 1 *quinquefasciatus 11* . 1. . 4 6 *sinensis 28* 3 11 1 7 6 . *tritaeniorhynchus 1* . 1. . . . *vishnui 53* 1 29 8 1 11 3 unknown *24* 2 14 1 3 2 2

*indiana 6* . .15. *uniformis 11* . 7. 1 3 . *chamberlainai 5* . 31 1 . *hybrida 7* 22 21 . . *luzonensis 7* 21 1. . 3 *metallica 2* .... 2. unknown *1* . 1. . .

mosquitoes (Table 1) were examined.

and historical. Instead of considering a conglomerate of traits to " cover " the organismal morphology, we suggest to decompose the global morphology into more elementary units and test their variation against species richness.

A set of anatomical landmarks in the wing represents a suitable tool to explore that hypothesis since it can be decomposed into subsets of different landmarks. In our approach, the complete set of anatomical landmarks of the wing would represent the broadest morphological coverage, and its decomposition into smaller configurations of landmarks would provide more elementary units of morphology. Do different subsets (combinations) of landmarks reflect biodiversity in the same way as the total set? Landmark-based geometric morphometrics, which is applied here, provides a convenient way for measuring morphological variety, requiring only the recognition of homologous landmarks in all individuals under comparison. This condition applies very well to mosquito wings, because their venation pattern is almost identical among different species and higher taxa, including different tribes.

#### **2. Materials & methods**

#### **2.1 The insects**

We used a total sample of 480 individuals (one wing per individual). They were tentatively identified using available morphological keys (Rattanarithikul et al., 2005; Rattanarithikul, Harrison, Panthusiri, Peyton & Coleman, 2006; Rattanarithikul, Harrison, Harbach, Panthusiri & Coleman, 2006.; Rattanarithikul et al., 2010). A total of 10 genera and 43 species of Culicidae (mosquitoes) was found, with unidentified species pooled into one putative "species" (Table 1). This collection was a representative set of higher taxa within Culicidae: it contained indeed the two subfamilies (Anophelinae and Culicinae) and, within the Culicinae, 6 tribes out of 11.

#### **2.2 Mosquito collection**

The mosquito collection was done during the rainy season of 2008 (June-August) along a forest-agro-urban landscape gradient within Nakhon Nayok province, central Thailand. Six habitat types: forest, fragmented forest, rice field, rural, suburban, and urban were identified and characterized across the landscape gradient (Table 1). For each habitat type, four sites were picked as representative of similar habitat range. Forest sites were situated along the border of the pristine Khao Yai National Park. Fragmented forest sites were on the edge of disturbed vegetation patch not far from the National Park, where human settlements were sparse and traditional small-scale agricultures were practiced. Rice field sites were further away from the National park and situated in the lowland closer to the main river hence the big-scale and irrigated rice agricultures were possible. The rural, suburban, and urban sites distributed based on the distance from the centre of town. The rural sites were more than 7 km from town; the suburban were within 5 km from town; and the urban sites were in the center of the town.

Four types of adult mosquito trap: BG sentinel, Mosquito Magnet, CDC UV light traps, and CDC backpack aspirator, were used in order to maximize the variety of mosquito samples collected. In each site, mosquitoes were collected for 24 hours using 8 BG traps, 2 Mosquito Magnet traps, and 3 of 3 to 10 minute-long aspirations for the day trapping, and 8 UV light traps, 8 BG traps, and 2 Mosquito Magnet traps for the night trapping. A total of over 2 Will-be-set-by-IN-TECH

and historical. Instead of considering a conglomerate of traits to " cover " the organismal morphology, we suggest to decompose the global morphology into more elementary units

A set of anatomical landmarks in the wing represents a suitable tool to explore that hypothesis since it can be decomposed into subsets of different landmarks. In our approach, the complete set of anatomical landmarks of the wing would represent the broadest morphological coverage, and its decomposition into smaller configurations of landmarks would provide more elementary units of morphology. Do different subsets (combinations) of landmarks reflect biodiversity in the same way as the total set? Landmark-based geometric morphometrics, which is applied here, provides a convenient way for measuring morphological variety, requiring only the recognition of homologous landmarks in all individuals under comparison. This condition applies very well to mosquito wings, because their venation pattern is almost identical among different species and higher taxa, including

We used a total sample of 480 individuals (one wing per individual). They were tentatively identified using available morphological keys (Rattanarithikul et al., 2005; Rattanarithikul, Harrison, Panthusiri, Peyton & Coleman, 2006; Rattanarithikul, Harrison, Harbach, Panthusiri & Coleman, 2006.; Rattanarithikul et al., 2010). A total of 10 genera and 43 species of Culicidae (mosquitoes) was found, with unidentified species pooled into one putative "species" (Table 1). This collection was a representative set of higher taxa within Culicidae: it contained indeed the two subfamilies (Anophelinae and Culicinae) and, within the Culicinae, 6 tribes out of 11.

The mosquito collection was done during the rainy season of 2008 (June-August) along a forest-agro-urban landscape gradient within Nakhon Nayok province, central Thailand. Six habitat types: forest, fragmented forest, rice field, rural, suburban, and urban were identified and characterized across the landscape gradient (Table 1). For each habitat type, four sites were picked as representative of similar habitat range. Forest sites were situated along the border of the pristine Khao Yai National Park. Fragmented forest sites were on the edge of disturbed vegetation patch not far from the National Park, where human settlements were sparse and traditional small-scale agricultures were practiced. Rice field sites were further away from the National park and situated in the lowland closer to the main river hence the big-scale and irrigated rice agricultures were possible. The rural, suburban, and urban sites distributed based on the distance from the centre of town. The rural sites were more than 7 km from town; the suburban were within 5 km from town; and the urban sites were in the

Four types of adult mosquito trap: BG sentinel, Mosquito Magnet, CDC UV light traps, and CDC backpack aspirator, were used in order to maximize the variety of mosquito samples collected. In each site, mosquitoes were collected for 24 hours using 8 BG traps, 2 Mosquito Magnet traps, and 3 of 3 to 10 minute-long aspirations for the day trapping, and 8 UV light traps, 8 BG traps, and 2 Mosquito Magnet traps for the night trapping. A total of over

and test their variation against species richness.

different tribes.

**2.1 The insects**

**2. Materials & methods**

**2.2 Mosquito collection**

center of the town.

80,000 mosquitoes were collected. For the morphometric study, only a subset of the female mosquitoes (Table 1) were examined.

(4) (5) 0 (10) (11)

00 . 0 0 0

taxa where groups represent communities of various taxonomic richness.

of high evenness is generally not the one found in natural conditions.

landmarks, and the total number of landmarks (13).

**2.5 Species assemblage**

**2.6 Sampling individuals**

**2.7 Metric disparity**

1 8-9-11-13 1-4-6-8-11 0 1-2-5-7-8-10-11-12-13 1-3-4-5-7-8-10-11-12-13 2 4-5-7-13 3-5-9-10-11 0 1-3-4-5-6-7-9-10-11-12 1-2-3-4-5-6-9-10-11-12-13 3 2-6-7-13 5-9-10-11-12 0 1-2-3-4-5-7-8-10-11-12 1-2-3-4-7-8-9-10-11-12-13 4 2-3-8-13 4-6-7-10-13 0 2-3-5-6-7-8-9-10-12-13 2-3-4-5-7-8-9-10-11-12-13 5 1-2-3-10 1-4-6-9-13 0 3-4-5-6-7-8-9-10-11-13 1-2-3-6-7-8-9-10-11-12-13 6 4-5-12-13 1-2-3-4-5 0 3-4-5-6-8-9-10-11-12-13 1-2-4-6-7-8-9-10-11-12-13 7 5-8-10-13 2-3-7-9-11 0 1-2-4-6-7-8-9-10-11-12 2-3-4-5-6-7-8-10-11-12-13 8 2-9-12-13 2-3-5-7-12 0 1-2-6-7-8-9-10-11-12-13 1-2-3-4-5-6-7-9-11-12-13 9 1-3-8-12 1-3-7-8-10 0 1-2-3-5-6-7-9-11-12-13 2-3-4-6-7-8-9-10-11-12-13 10 1-4-5-9 1-5-7-12-13 0 1-2-4-5-6-9-10-11-12-13 1-2-3-4-5-6-7-9-10-11-12 11 6-8-10-13 2-4-7-8-10 0 2-3-4-6-7-8-9-10-11-13 1-2-3-5-6-7-8-9-10-11-12 12 4-7-8-9 3-4-8-11-12 0 1-2-3-4-5-6-7-8-12-13 1-3-4-5-6-7-8-9-10-11-12 13 1-3-6-9 1-5-8-9-10 0 1-2-3-5-7-8-10-11-12-13 1-2-3-4-5-7-8-9-10-11-13 14 1-5-8-10 3-5-9-11-12 0 2-3-4-5-6-7-8-9-10-13 1-2-4-5-6-7-8-9-10-12-13 15 1-4-8-9 3-8-9-10-11 0 1-2-4-5-6-7-8-10-12-13 1-2-3-4-5-6-7-8-9-10-12

The Mosquito Fauna: From Metric Disparity to Species Diversity 35

30 3-5-8-12 2-5-8-9-13 0 1-2-3-4-5-7-9-11-12-13 1-2-3-4-5-6-7-8-9-12-13 Table 2. A partial representation of the sets of landmarks used to examine a sequence of increasing taxonomic richness (see Table 3). Landmarks were numbered from 1 to 13 (Fig. 1).

In the building of sample sets of increasing taxonomic richness (a sequence), species were not " added " to previous ones, they were resampled at each step from the total available. This can be seen in the sequence shown Table 3. The program was not simulating a temporal variation where species progressively accumulate in a given environment, but a spatial sampling of

To measure metric disparity (MD), a total of 50 specimens by assemblage was typically used. However, for assemblages of low SR (2 to 8 species), the sample size could be less than 50 (28.6 +- 9.6). On average, the sample size was 45.8 +- 8.7. The abundance of each species within an assemblage was approximately the same. For instance, an assemblage of 10 species contained approximately 5 specimens per species, while an assemblage of 24 species could contain for instance 20 species with 2 individuals and 4 species containing 3 individuals. This situation

To estimate morphological diversity, we considered only the geometric, landmark-based approach. Morphological diversity was estimated by the metric disparity (MD) index. For

The first row refers to the number of landmarks (between brackets). The first column enumerates the 30 combinations of landmarks randomly generated for each number of landmarks. The second column partially shows the 30 combinations to represent 4 landmarks. To save space, only some combinations for 4, 5, 10 and 11 landmarks are represented here. To these 240 (30\*8) combinations, we added the 13 combinations of 12

Table 1. List of genera and species of Culicidae identified on morphological ground, with their repartition according to the habitat. F, forest; FF, fragmented forest; RF, rice field; R, rural; SU, semi-urban and U, urban. For statistical tests, unknown species have been pooled into one single taxon (92 specimens). *n*, number of specimens submitted to morphometric analyses.

#### **2.3 Shape of the wing**

The shape of the mosquito wings was described by 13 landmarks (see Fig. 1).

Fig. 1. Mosquito wing. Landmarks are labelled according to the order of digitization.

The decomposition of shape used a total of 254 different LM configurations out of the 7814 possible ones, i.e. the totality (13) of landmarks (LM), the 13 combinations of 12 LM, and 240 combinations of decreasing numbers of LM. For each subset of LM, going from 11 LM to 4 LM, 30 different combinations were tested (see Table 2).

The above testing of 254 LM configurations was replicated on ten different species sequences (see Table 3).

#### **2.4 Species sequence**

A sequence refers here to a succession of species assemblages with increasing richness. The mosquitoes were randomly sampled with replacement into 22 assemblages of increasing species richness. These 22 assemblages constituted a sequence where each unit represented a community with different species richness, starting from 2 species, 4 species, 6 species, and so on until 44 species. A partial representation of a sequence is shown Table 3; ten such sequences were constituted. The range of categorical units (44 species) was not modified.

4 Will-be-set-by-IN-TECH

Genus Species *n* F FF R RF SU U *Uranotaenia campestris 2* . 2. . . .

Table 1. List of genera and species of Culicidae identified on morphological ground, with their repartition according to the habitat. F, forest; FF, fragmented forest; RF, rice field; R, rural; SU, semi-urban and U, urban. For statistical tests, unknown species have been pooled into one single taxon (92 specimens). *n*, number of specimens submitted to morphometric

The shape of the mosquito wings was described by 13 landmarks (see Fig. 1).

Fig. 1. Mosquito wing. Landmarks are labelled according to the order of digitization.

LM, 30 different combinations were tested (see Table 2).

The decomposition of shape used a total of 254 different LM configurations out of the 7814 possible ones, i.e. the totality (13) of landmarks (LM), the 13 combinations of 12 LM, and 240 combinations of decreasing numbers of LM. For each subset of LM, going from 11 LM to 4

The above testing of 254 LM configurations was replicated on ten different species sequences

A sequence refers here to a succession of species assemblages with increasing richness. The mosquitoes were randomly sampled with replacement into 22 assemblages of increasing species richness. These 22 assemblages constituted a sequence where each unit represented a community with different species richness, starting from 2 species, 4 species, 6 species, and so on until 44 species. A partial representation of a sequence is shown Table 3; ten such sequences

were constituted. The range of categorical units (44 species) was not modified.

analyses.

(see Table 3).

**2.4 Species sequence**

**2.3 Shape of the wing**

*lateralis 1* .... 1. *lutescens 4* .. .13. *micans 7* . 1. 3 3 . *nivipleuraura 3* 2. . . 1 . *subnormalis 2* .... 2. unknown *22* 17 2 1 1 1

Table 2. A partial representation of the sets of landmarks used to examine a sequence of increasing taxonomic richness (see Table 3). Landmarks were numbered from 1 to 13 (Fig. 1). The first row refers to the number of landmarks (between brackets). The first column enumerates the 30 combinations of landmarks randomly generated for each number of landmarks. The second column partially shows the 30 combinations to represent 4 landmarks. To save space, only some combinations for 4, 5, 10 and 11 landmarks are represented here. To these 240 (30\*8) combinations, we added the 13 combinations of 12 landmarks, and the total number of landmarks (13).

#### **2.5 Species assemblage**

In the building of sample sets of increasing taxonomic richness (a sequence), species were not " added " to previous ones, they were resampled at each step from the total available. This can be seen in the sequence shown Table 3. The program was not simulating a temporal variation where species progressively accumulate in a given environment, but a spatial sampling of taxa where groups represent communities of various taxonomic richness.

#### **2.6 Sampling individuals**

To measure metric disparity (MD), a total of 50 specimens by assemblage was typically used. However, for assemblages of low SR (2 to 8 species), the sample size could be less than 50 (28.6 +- 9.6). On average, the sample size was 45.8 +- 8.7. The abundance of each species within an assemblage was approximately the same. For instance, an assemblage of 10 species contained approximately 5 specimens per species, while an assemblage of 24 species could contain for instance 20 species with 2 individuals and 4 species containing 3 individuals. This situation of high evenness is generally not the one found in natural conditions.

#### **2.7 Metric disparity**

To estimate morphological diversity, we considered only the geometric, landmark-based approach. Morphological diversity was estimated by the metric disparity (MD) index. For

to compute MD (see above). Five ANOVA were performed, one using the totality of LM to compute the PW, and four using selected combinations of LM. The latter were chosen according to their relationship to species diversity: they were two configurations of LM that produced MD highly correlated to the species richness (SR), and two others not related to SR

The Mosquito Fauna: From Metric Disparity to Species Diversity 37

To estimate the biodiversity, we used the species richness index (SR, or the number of categorical units) which is the total number of species found in the community (an assemblage). Since evenness was maintained as high as possible, biodiversity indexes like the Shannon index (Shannon & Weaver, 1949) or the Simpson one (Simpson, 1949) were not

Diversity structure (DS) is described here as the relationship between MD and SR (Foote, 1992). It was estimated by the determination coefficients (squared linear correlation coefficients) (see Table 5), and illustrated graphically (Fig. 2). An average estimate of MD was used during all correlations. For each specific combination of LM, at each level of SR, MD was an average value derived from the 30 different assemblages (Table 3). For one specific number of LM, for example 4, this average was performed taking into account also 30 combinations of

A special TclTk (http://www.tcl.tk/software/tcltk/) script was written where Procrustes superposition (GPA), partial warps (PW) computations, as well as metric diversity

(http://www.mpl.ird.fr/morphometrics/clic/index.html). Table 3 was

Globally, the effects of SR on various LM configurations of the same wing confirmed our initial hypothesis: different aspects of shape did not vary the same way in response to the same factor. The relationship between species richness (SR) and metric disparity (MD) was estimated by the determination coefficient (see DS, for "diversity structure" in Tables 4, 5 and 6). This coefficient does not inform about the sign of the correlation between SR and MD, which was not necessary since these correlation coefficients were all positive. The species richness contributed to the metric disparity according to the number of landmarks used to compute MD (Tables 4 and 5), and also according to the identity of landmarks involved (Table 6). We subdivide our results according to these two aspects of shape composition: the number (see paragraph 3.1) and the identity (see paragraph 3.2) of LM. We then show the ANOVA output for some remarquable configurations of LM (high DS, low DS): it estimates the possible

(MD) estimations, made use of procedures extracted from the CLIC package

computed using STATA (**?**). Figure 2 used the GNUMERIC spreadsheet

role of habitat heterogeneity on the metric disparity (see paragraph 3.3).

(http://projects.gnome.org/gnumeric/).

(see Table 7).

examined.

4 LM (Table 2).

**2.11 Software**

**3. Results**

**2.10 Diversity structure**

**2.9 Biodiversity**

Table 3. A partial representation of a sequence of increasing taxonomic richness (numbers between brackets) as used by our simulation programme. Each number represents a mosquito species. There were 44 different taxa. The first set of rows shows two assemblages of two randomly selected species (species "13", species "21" and species "3", species "14") among the 30 species pairs randomly generated. The second set of rows shows two such assemblages of 4 species among the 30 random combinations of 4 species, and so on till reaching an assemblage of 44 species. This makes a total of 22 assemblages of increasing species richness, each one sampled 30 times among species. The last assemblage, containing the totality of the species, was sampled 30 times among individuals (*1-30*). Ten such sequences were generated, and each one was explored by using 254 different landmarks configurations (Table 2).

each of N wings, after Procrustes superposition using the Generalized Procrustes Algorithm (GPA) (Rohlf, 1990) and partial warps (PW) computation as in Rohlf (1996), the sum of squared PW was obtained. This sum was divided by the degrees of freedom (N-1) to compute MD (Zelditch et al., 2004). To estimate MD for an assemblage with, for instance, 2 species, 30 random pairs of species were used (see Table 3) and the average MD value was considered.

#### **2.8 Habitat heterogeneity**

The sample composition did not allow valid statistics using the habitat as a categorical unit (instead of the species). To evaluate the importance of the environment on the metric disparity (MD), we performed a simple two-way ANOVA with taxa and habitat as effects. The variable used was the individual sum of squared PW because this sum is the term directly used to compute MD (see above). Five ANOVA were performed, one using the totality of LM to compute the PW, and four using selected combinations of LM. The latter were chosen according to their relationship to species diversity: they were two configurations of LM that produced MD highly correlated to the species richness (SR), and two others not related to SR (see Table 7).

#### **2.9 Biodiversity**

6 Will-be-set-by-IN-TECH

(*SR*) (*2*) (*4*) (*6*) (*8*) (*etc*.) (*40*) (*44*)

*1-30* 1 2 3 4 5 6 7 8 etc. 42 43 44 Table 3. A partial representation of a sequence of increasing taxonomic richness (numbers between brackets) as used by our simulation programme. Each number represents a

mosquito species. There were 44 different taxa. The first set of rows shows two assemblages of two randomly selected species (species "13", species "21" and species "3", species "14") among the 30 species pairs randomly generated. The second set of rows shows two such assemblages of 4 species among the 30 random combinations of 4 species, and so on till reaching an assemblage of 44 species. This makes a total of 22 assemblages of increasing species richness, each one sampled 30 times among species. The last assemblage, containing

each of N wings, after Procrustes superposition using the Generalized Procrustes Algorithm (GPA) (Rohlf, 1990) and partial warps (PW) computation as in Rohlf (1996), the sum of squared PW was obtained. This sum was divided by the degrees of freedom (N-1) to compute MD (Zelditch et al., 2004). To estimate MD for an assemblage with, for instance, 2 species, 30 random pairs of species were used (see Table 3) and the average MD value was considered.

The sample composition did not allow valid statistics using the habitat as a categorical unit (instead of the species). To evaluate the importance of the environment on the metric disparity (MD), we performed a simple two-way ANOVA with taxa and habitat as effects. The variable used was the individual sum of squared PW because this sum is the term directly used

the totality of the species, was sampled 30 times among individuals (*1-30*). Ten such sequences were generated, and each one was explored by using 254 different landmarks

*1* 13 21 *etc*.. . *30* 3 14

configurations (Table 2).

**2.8 Habitat heterogeneity**

*1* 15 18 23 34 *etc*.. . . . *30* 3 12 14 38

*1* 7 17 22 35 37 41 *etc*.. . . . . . *30* 3 5 6 10 26 38

*1* 11 12 16 22 24 35 38 42 *etc*.. . . . . . . . *30* 3 6 7 13 18 20 26 37 *1* 6 9 18 19 22 28 34 38 etc. *etc*. . . . . . . . . etc. *30* 3 6 9 13 15 24 27 38 etc. *etc*. etc. etc. etc. etc. etc. etc. etc. etc. etc. *1* 1 2 5 6 8 11 12 13 etc. 41 *etc*. . . . . . . . . etc. . *30* 1 3 4 6 7 8 10 11 etc. 44 To estimate the biodiversity, we used the species richness index (SR, or the number of categorical units) which is the total number of species found in the community (an assemblage). Since evenness was maintained as high as possible, biodiversity indexes like the Shannon index (Shannon & Weaver, 1949) or the Simpson one (Simpson, 1949) were not examined.

#### **2.10 Diversity structure**

Diversity structure (DS) is described here as the relationship between MD and SR (Foote, 1992). It was estimated by the determination coefficients (squared linear correlation coefficients) (see Table 5), and illustrated graphically (Fig. 2). An average estimate of MD was used during all correlations. For each specific combination of LM, at each level of SR, MD was an average value derived from the 30 different assemblages (Table 3). For one specific number of LM, for example 4, this average was performed taking into account also 30 combinations of 4 LM (Table 2).

#### **2.11 Software**

A special TclTk (http://www.tcl.tk/software/tcltk/) script was written where Procrustes superposition (GPA), partial warps (PW) computations, as well as metric diversity (MD) estimations, made use of procedures extracted from the CLIC package (http://www.mpl.ird.fr/morphometrics/clic/index.html). Table 3 was computed using STATA (**?**). Figure 2 used the GNUMERIC spreadsheet (http://projects.gnome.org/gnumeric/).

#### **3. Results**

Globally, the effects of SR on various LM configurations of the same wing confirmed our initial hypothesis: different aspects of shape did not vary the same way in response to the same factor. The relationship between species richness (SR) and metric disparity (MD) was estimated by the determination coefficient (see DS, for "diversity structure" in Tables 4, 5 and 6). This coefficient does not inform about the sign of the correlation between SR and MD, which was not necessary since these correlation coefficients were all positive. The species richness contributed to the metric disparity according to the number of landmarks used to compute MD (Tables 4 and 5), and also according to the identity of landmarks involved (Table 6). We subdivide our results according to these two aspects of shape composition: the number (see paragraph 3.1) and the identity (see paragraph 3.2) of LM. We then show the ANOVA output for some remarquable configurations of LM (high DS, low DS): it estimates the possible role of habitat heterogeneity on the metric disparity (see paragraph 3.3).

NLM DS SD MIN MAX STABILITY

The Mosquito Fauna: From Metric Disparity to Species Diversity 39

Table 5. Determination coefficients (column DS) as percentages representing the contribution of species richness (SR) to metric disparity (MD) according to the number of landmarks (NLM) used to compute MD. These coefficients allow comparisons of the diversity structures (DS), i.e. the relationship between MD and SR according to NLM. Columns DS lists the determination coefficients and SD their standard deviation; MIN is the minimum value of DS and MAX its maximum. The table, except last column, represents an average DS derived from an average sequence of groups having increased SR; it was computed from the ten replicated sequences used in the study (see Table 3). The last column (STABILITY) is the standard deviation of the DS mean scores obtained from one sequence to another; it indicates how stable were the DS for the same LM configurations across 10 different random sequences of species. For each number of landmarks (each row), DS and SD values were averaged over 30 combinations of different landmarks (except for 12 LM which had only 13 possible configurations and of course for 13LM, see Table 2) and over ten replicated sequences, each one providing an average estimate from 30 random assemblages of species (Table 3)

For the same specimens and the same species arrangements, some shape components (i.e. landmarks configurations) of the wing varied in accordance with the number of species and others did not. Furthermore, a LM configuration highly predictive of SR, like the set of landmarks 2, 3, 8, 13 or the landmarks 2, 6, 7, 13, remained predictive regardless of the species sequences. The same observation applied for non-predictive sets of landmarks, like for instance the set of landmarks 1, 4, 5, 6 or 1, 3, 6, 9. This stability was verified across the ten

The highly predictive LNLM configurations, like 2-3-8-13 and 2-6-7-13 used in the ANOVA, were both affected by species richness only (P < 0.0001), not by the habitat heterogeneity (P > 0.0500), while the TNLM and two poorly predictive LNLM (1-4-5-6 and 1-3-6-9) were affected

To explore the diversity structure (DS), i.e. the relationship between metric and biological diversity, our model tested the effect of species richness (SR) on the metric disparity (MD) computed from 254 possible combinations of landmarks. Our data showed that the DS

replicates (see last column of Table 5).

by both species and habitat (Table 7).

**3.3 The habitat heterogeneity**

**4. Discussion**

4 64.7% 13.8 27 94 5.8 5 63.4% 7.7 35 85 6.0 6 65.4% 8.3 40 89 5.7 7 65.3% 5.8 41 86 6.4 8 64.1% 6.5 40 83 6.6 9 66.4% 3.6 49 83 6.1 10 64.5% 3.9 47 79 6.8 11 65.8% 3.4 48 80 6.3 12 67.5% 2.5 59 71 2.2 13 64.0% 0.6 57 62 3.5

#### **3.1 Species richness and the number of landmarks**

We observed that the configurations involving a low number of landmarks varied more frequently in very close accordance with species richness. Actually, for the same species arrangements, some configurations of the wing involving a very low number of landmarks (LNLM configurations) had either very high (column "DS > 75", see Table 4) or very low (column "DS < 50", see Table 4) prediction power on species diversity, while more complex anatomical configurations of the same wing showed a more stable but less predictive relationship with SR.

Table 4. For each number of landmarks (LM\_nb), the column "total" gives the total number of possible configurations among the possible landmark positions on the wing (a total of 13), the column "sampled" indicates the percentage of such configurations that have been studied; except for 12 and 13 LM, we always examined 30 random configurations for each number of landmark. Thus, 4% for instance is 30 out of 715. DS, or the "diversity structure" was estimated here by the determination coefficient (expressed in percentages) between metric disparity and species richness. The three last columns refers to the frequency at which a given DS was observed. It can be seen that a determination coefficient lower than 55% or higher than 75% between species richness and metric disparity was observed more frequently when using a low number of LM

Both the best and worst DS scores between SR and MD were obtained with configurations made from a low number of landmarks (LNLM). The range of scores progressively decreased with the addition of more LM (see column SD of Table 5). With more numerous LM however, the best predictive values did not reach such high levels as with fewer LM. Thus, more LM meant a more stable assessment of diversity structure (DS, or the relationship between MD and species richness), with no occurrence of very high values (Table 5).

#### **3.2 Species richness and the identity of landmarks**

Not only were we able to disclose different diversity structures according to the number of landmarks, but also according to specific configurations of landmarks. A partial output is presented Table 6. The highest determination coefficient (94%) observed was obtained with a specific combination of 4 LM (see 2-3-8-13, Table 4); other combinations involving the same number of LM gave much lower predictability (see 1-3-8-12, Table 4). It could be as low as 27% (see the columns MAX and MIN of Table 5).

8 Will-be-set-by-IN-TECH

We observed that the configurations involving a low number of landmarks varied more frequently in very close accordance with species richness. Actually, for the same species arrangements, some configurations of the wing involving a very low number of landmarks (LNLM configurations) had either very high (column "DS > 75", see Table 4) or very low (column "DS < 50", see Table 4) prediction power on species diversity, while more complex anatomical configurations of the same wing showed a more stable but less predictive

> LM\_nb total sampled DS ≤ 55 55 < DS < 75 DS ≥ 75 4 715 (4%) 26% 45% 29% 5 1287 (2%) 21% 71% 8% 6 1716 (2%) 17% 71% 12% 7 1716 (2%) 12% 82% 9% 8 1287 (2%) 16% 77% 7% 9 715 (4%) 6% 87% 7% 10 286 (10%) 10% 83% 6% 11 78 (38%) 5% 86% 9% 12 13 (100%) 0% 98% 2% 13 1 (100%) 0% 100% 0%

Table 4. For each number of landmarks (LM\_nb), the column "total" gives the total number of possible configurations among the possible landmark positions on the wing (a total of 13), the column "sampled" indicates the percentage of such configurations that have been studied; except for 12 and 13 LM, we always examined 30 random configurations for each number of landmark. Thus, 4% for instance is 30 out of 715. DS, or the "diversity structure" was estimated here by the determination coefficient (expressed in percentages) between metric disparity and species richness. The three last columns refers to the frequency at which a given DS was observed. It can be seen that a determination coefficient lower than 55% or higher than 75% between species richness and metric disparity was observed more

Both the best and worst DS scores between SR and MD were obtained with configurations made from a low number of landmarks (LNLM). The range of scores progressively decreased with the addition of more LM (see column SD of Table 5). With more numerous LM however, the best predictive values did not reach such high levels as with fewer LM. Thus, more LM meant a more stable assessment of diversity structure (DS, or the relationship between MD

Not only were we able to disclose different diversity structures according to the number of landmarks, but also according to specific configurations of landmarks. A partial output is presented Table 6. The highest determination coefficient (94%) observed was obtained with a specific combination of 4 LM (see 2-3-8-13, Table 4); other combinations involving the same number of LM gave much lower predictability (see 1-3-8-12, Table 4). It could be as low as

and species richness), with no occurrence of very high values (Table 5).

**3.1 Species richness and the number of landmarks**

frequently when using a low number of LM

**3.2 Species richness and the identity of landmarks**

27% (see the columns MAX and MIN of Table 5).

relationship with SR.

Table 5. Determination coefficients (column DS) as percentages representing the contribution of species richness (SR) to metric disparity (MD) according to the number of landmarks (NLM) used to compute MD. These coefficients allow comparisons of the diversity structures (DS), i.e. the relationship between MD and SR according to NLM. Columns DS lists the determination coefficients and SD their standard deviation; MIN is the minimum value of DS and MAX its maximum. The table, except last column, represents an average DS derived from an average sequence of groups having increased SR; it was computed from the ten replicated sequences used in the study (see Table 3). The last column (STABILITY) is the standard deviation of the DS mean scores obtained from one sequence to another; it indicates how stable were the DS for the same LM configurations across 10 different random sequences of species. For each number of landmarks (each row), DS and SD values were averaged over 30 combinations of different landmarks (except for 12 LM which had only 13 possible configurations and of course for 13LM, see Table 2) and over ten replicated sequences, each one providing an average estimate from 30 random assemblages of species (Table 3)

For the same specimens and the same species arrangements, some shape components (i.e. landmarks configurations) of the wing varied in accordance with the number of species and others did not. Furthermore, a LM configuration highly predictive of SR, like the set of landmarks 2, 3, 8, 13 or the landmarks 2, 6, 7, 13, remained predictive regardless of the species sequences. The same observation applied for non-predictive sets of landmarks, like for instance the set of landmarks 1, 4, 5, 6 or 1, 3, 6, 9. This stability was verified across the ten replicates (see last column of Table 5).

#### **3.3 The habitat heterogeneity**

The highly predictive LNLM configurations, like 2-3-8-13 and 2-6-7-13 used in the ANOVA, were both affected by species richness only (P < 0.0001), not by the habitat heterogeneity (P > 0.0500), while the TNLM and two poorly predictive LNLM (1-4-5-6 and 1-3-6-9) were affected by both species and habitat (Table 7).

#### **4. Discussion**

To explore the diversity structure (DS), i.e. the relationship between metric and biological diversity, our model tested the effect of species richness (SR) on the metric disparity (MD) computed from 254 possible combinations of landmarks.  Two-way ANOVA showing the contribution of species richness and habitat heterogeneity to the individual sum of squared PW: this sum per individual is the value on which metric disparity is directly estimated (as an average). The ANOVA was performed on values computed from the total number of landmarks (TNLM configuration, top) or from the following configuration of four landmarks: 2-3-8-13 (bottom). For the latter (and also for LNLM highly predictive configurations like 2-6-7-13, not shown here), only species richness contributed to metric variation. Both species richness and habitat heterogeneity contributed

to the variation of the total set of landmarks (as well as of LNLM poorly predictive

likely to affect the diversity structure of some landmark configurations.

The model did not randomize the total number of species, which was always fixed to 44 taxa, nor did it randomize the sequence of species richness (2, 4, 6,..., 44), two parameters that are

There is however no *a priori* reason to think that these shortcomings would reduce the interest to use smaller configurations of shape to simplify the interpretation of the diversity structure. More critical may be that our model did not take into account variation in species evenness. In each assemblage of species for a given sequence, the model was designed to get an approximately equal number of specimens within each species. In the natural conditions, there is generally no such evenness, and, because the MD index is an average value, this might affect the correlation between MD and SR. Rare taxa are often assumed to exhibit unusual morphologies because of specialized life habits and could thus contribute disproportionately to the disparity of an assemblage (Deline, 2009). But wings of mosquitoes, because they are generally very similar, are unlikely to have such an effect. Their general similarity may cause the opposite effect to take place: rare taxa may fill the mophospace among the common species

Finally, in our simulation the successive assemblages contained an increasing number of species randomly selected from a total pool of mosquitoes, they were not obtained by " adding " species to previously selected ones. Because of this specific design, our simulation is not relevant to temporal follow-up of metric diversity, however it is specific to the spatial

configurations like 1-4-6-9 or 1-3-5-6, not shown here)

and would likely lower MD (Deline, 2009).

comparison of biodiversity.

Table 6. Minimum (MIN) and maximum (MAX) determination coefficients as percentages representing the contribution of species richness to metric disparity according to the configurations of landmarks used. The 8-9-11-13 formula means the configuration of four landmarks using landmarks 8, 9, 11 and 13 as represented in Fig. 1. As in Table 2 the LM configurations are classified according to the number of LM (number between brackets). To save sapece, we present only a partial output of the data for 4LM and 11LM. For each configuration of landmarks (each row), values were obtained from ten replicated sequences, each one providing an average estimate from 30 random assemblages of species (Table 3). It can be seen that the contribution of species richness to metric disparity (MD) can be very high when MD is computed from a low number of LM (column 4LM), which does not seem to be the case for configurations involving more landmarks (column 11LM). Table 5 provides the average DS over all the tested configurations for each number of LM (not only 4 and 11).

depended on the aspects of shape that were considered (the number of LM, the configuration of LM).

#### **4.1 Our model**

The ability to detect morphological trends and occupation patterns within morphospace might depend on using the appropriate measure(s) of disparity. Since there is no clear indication about which index is best, some authors used a variety of measurements (Ciampaglio et al., 2001; Navarro, 2003). Our model presented results regarding only the MD index: it is the most commonly applied measurement of morphological disparity in landmark-based geometric studies (Zelditch et al., 2004).

In trying to include 50 individuals per group within our model, we were unable to consider the effect of sampling variation. However, metric disparity, as measured here, is relatively insensitive to sample size variation and has been shown to be a stable estimate when using a sample of 50 individuals (Navarro, 2003).

Table 6. Minimum (MIN) and maximum (MAX) determination coefficients as percentages representing the contribution of species richness to metric disparity according to the configurations of landmarks used. The 8-9-11-13 formula means the configuration of four landmarks using landmarks 8, 9, 11 and 13 as represented in Fig. 1. As in Table 2 the LM configurations are classified according to the number of LM (number between brackets). To save sapece, we present only a partial output of the data for 4LM and 11LM. For each configuration of landmarks (each row), values were obtained from ten replicated sequences, each one providing an average estimate from 30 random assemblages of species (Table 3). It can be seen that the contribution of species richness to metric disparity (MD) can be very high when MD is computed from a low number of LM (column 4LM), which does not seem to be the case for configurations involving more landmarks (column 11LM). Table 5 provides the average DS over all the tested configurations for each number of LM (not only 4 and 11).

depended on the aspects of shape that were considered (the number of LM, the configuration

The ability to detect morphological trends and occupation patterns within morphospace might depend on using the appropriate measure(s) of disparity. Since there is no clear indication about which index is best, some authors used a variety of measurements (Ciampaglio et al., 2001; Navarro, 2003). Our model presented results regarding only the MD index: it is the most commonly applied measurement of morphological disparity in

In trying to include 50 individuals per group within our model, we were unable to consider the effect of sampling variation. However, metric disparity, as measured here, is relatively insensitive to sample size variation and has been shown to be a stable estimate when using a

landmark-based geometric studies (Zelditch et al., 2004).

sample of 50 individuals (Navarro, 2003).

of LM).

**4.1 Our model**

(4) MIN MAX (11) MIN MAX 8-9-11-13 71 93 1-3-4-5-7-8-10-11-12-13 56 75 4-5-7-13 44 62 1-2-3-4-5-6-9-10-11-12-13 55 67 2-6-7-13 71 90 1-2-3-4-7-8-9-10-11-12-13 58 70 2-3-8-13 77 94 2-3-4-5-7-8-9-10-11-12-13 58 76 1-2-3-10 44 72 1-2-3-6-7-8-9-10-11-12-13 57 76 4-5-12-13 54 74 1-2-4-6-7-8-9-10-11-12-13 57 76 5-8-10-13 52 72 2-3-4-5-6-7-8-10-11-12-13 58 75 2-9-12-13 62 80 1-2-3-4-5-6-7-9-11-12-13 56 70 1-3-8-12 37 48 2-3-4-6-7-8-9-10-11-12-13 61 77 1-4-5-9 39 62 1-2-3-4-5-6-7-9-10-11-12 49 63 6-8-10-13 51 71 1-2-3-5-6-7-8-9-10-11-12 49 68 4-7-8-9 74 89 1-3-4-5-6-7-8-9-10-11-12 50 68 1-3-6-9 64 78 1-2-3-4-5-7-8-9-10-11-13 56 77 1-5-8-10 47 65 1-2-4-5-6-7-8-9-10-12-13 52 73 1-4-8-9 65 79 1-2-3-4-5-6-7-8-9-10-12 48 67 .. . . . . . 3-5-8-12 42 63 1-2-3-4-5-6-7-8-9-12-13 54 76

Table 7. Two-way ANOVA showing the contribution of species richness and habitat heterogeneity to the individual sum of squared PW: this sum per individual is the value on which metric disparity is directly estimated (as an average). The ANOVA was performed on values computed from the total number of landmarks (TNLM configuration, top) or from the following configuration of four landmarks: 2-3-8-13 (bottom). For the latter (and also for LNLM highly predictive configurations like 2-6-7-13, not shown here), only species richness contributed to metric variation. Both species richness and habitat heterogeneity contributed to the variation of the total set of landmarks (as well as of LNLM poorly predictive configurations like 1-4-6-9 or 1-3-5-6, not shown here)

The model did not randomize the total number of species, which was always fixed to 44 taxa, nor did it randomize the sequence of species richness (2, 4, 6,..., 44), two parameters that are likely to affect the diversity structure of some landmark configurations.

There is however no *a priori* reason to think that these shortcomings would reduce the interest to use smaller configurations of shape to simplify the interpretation of the diversity structure.

More critical may be that our model did not take into account variation in species evenness. In each assemblage of species for a given sequence, the model was designed to get an approximately equal number of specimens within each species. In the natural conditions, there is generally no such evenness, and, because the MD index is an average value, this might affect the correlation between MD and SR. Rare taxa are often assumed to exhibit unusual morphologies because of specialized life habits and could thus contribute disproportionately to the disparity of an assemblage (Deline, 2009). But wings of mosquitoes, because they are generally very similar, are unlikely to have such an effect. Their general similarity may cause the opposite effect to take place: rare taxa may fill the mophospace among the common species and would likely lower MD (Deline, 2009).

Finally, in our simulation the successive assemblages contained an increasing number of species randomly selected from a total pool of mosquitoes, they were not obtained by " adding " species to previously selected ones. Because of this specific design, our simulation is not relevant to temporal follow-up of metric diversity, however it is specific to the spatial comparison of biodiversity.

**4.3 Species diversity: The identity of landmarks**

observation raises two questions:

richness estimations?

still to be investigated.

For a given number of landmarks, the same configurations produced high (or low) correlation with SR regardless of the taxonomic composition of the species sequence. For instance, the 2-3-8-13 configuration always produced high correlation with SR regardless of the sequence of species considered, and the 1-3-6-9 configuration always produced one of the worst correlations. Thus, there seemed to be specific configurations of LM that could reflect the number of species, and others that could not, regardless of the species sequence used. This

The Mosquito Fauna: From Metric Disparity to Species Diversity 43

• Are there any LM configurations that could be used in mosquitoes as a proxy for species

• What are the sources (other than randomness) that generated variation for the LM

The first question is about the effect of taxonomic diversity on wing shape variation: can the

Before starting to answer that question, one could ask what the interests are in making an indirect estimation of species richness. The answer is about time. Metric disparity extracted from wing venation is much faster to obtain than metrics requiring a taxonomic identification of each collected individual. Mosquitoes are often distinguished on the basis of labile characters which may be lost on damaged specimens, making species diagnostic difficult in a group where species are numerous (3500 species) and taxonomists are unfortunately not. An

Our results showed that some landmarks configurations were highly predictive for SR. If these landmark configurations were known in advance, their variation could be used indeed as a proxy for SR. In our study, their performances seemed stable regardless of species arrangement in different groups (Table 5, last column). However, whether the high scoring of these LM configurations would be maintained in other studies on other mosquito species is

The second question refers to other possible meanings of MD in our data set. Other than SR, what could be the cause of shape variation? It probably has many causes, certainly among which is species evenness (Deline, 2009), but also functional and ecological attributes (Roy & Foote, 1997), environmental conditions (Vasil'ev et al., 2010), founder effects (David, 1999; Whitlock & Fowler, 1999), endemicity (Magniez-Jannin et al., 2000) or reproduction mode (Baltanas et al., 2002). The only factor we could discuss with our data was the habitat, an environmental parameter to which insect metric properties are known to be sensible (Benítez et al., 2011; Tantowijoyo & Hoffmann, 2011). Insects were collected in the forest, in rural and urban areas. It was difficult to address the question of habitat influence for a given species since the within species sample sizes did not allow valid statistics inference (see Table 1). Thus, the habitat as explored here was not free from the possible effects due to different species compositions. Our ANOVA analyzes could however provide some indications. When both species richness and habitat heterogeneity were significantly contributing to the MD, a good relationship between MD and SR was not observed. A broad morphological covering (TNLM), or some unpredictive elementary units of morphology (LNLM configurations like 1-4-6-9 or 1-3-5-6), were under the influence of both species and habitat heterogeneities (Table

configurations that were not influenced by the species richness?

latter be used as a proxy for indirect species richness estimation?

indirect but fast comparative estimation of SR would be welcome.

#### **4.2 Species diversity: The number of landmarks**

Although correlations between MD and SR could be weak or not significant, they were positive. A positive correlation has an intuitive explanation: more diversity means more forms. Weak and not significant correlation between SR and MD obtained from some combinations of LM could be the consequence of a few species immediately occupying the extremes of the morphospace. In that situation, by filling in the morphospace between disparate taxa, increasing taxonomic diversity was not able to have a significant effect on MD estimates (Roy & Foote, 1997).

The intensity of the correlation between SR and MD depended on the number of LM. Unexpectedly, the low number of LM (LNLM) configurations gave the best predictive power on SR (Fig. 2, left side). This is counter-intuitive since it is the current belief that more shape would have more taxonomic contents and hence would be more useful to distinguish species. We verified that by using one of the most predictive configuration of 4 LM (2-3-8-13) only 30% of individuals could be correctly assigned to their species, while more than 75% could be correctly attributed with 13 LM (the TNLM configuration). How could it be that a better discriminating configuration of LM could give a lower taxonomic prediction than a poorly discriminating one? The answer could lie in the observation that a LNLM configuration also could have no relationship at all with SR, and that the TNLM configuration is gathering both the predictive and non-predictive subsets of LM. Thus, the mixing of better predictive with poorly predictive components of shape into a larger set of LM was likely to mix opposite trends and blur the relationship between MD and SR (Fig. 2, right side). In other words, when a larger number of LM is used as a proxy for SR, a larger number of influences is also allowed which is not limited to SR. This explanation raises the question whether in our material there were other influences affecting the variation of the LM configurations that were not related to SR.

Fig. 2. Two graphs representing the diversity structure (DS), i.e. the relationship between species richness (SR, on the horizontal axis) and metric diversity (MD, on vertical axis). The left graph shows the DS appearing when using a set of four landmarks (namely the landmarks 2, 3, 8 and 13; see Fig. 1). The right graph shows the DS using the total number (13) of LM. When using the totality of LM (right graph) there is a first list of increasing MD values up to a SR of 30 species. This first list is common to the two graphs. In the right one however, after 30 species, there is a " plateau " breaking the global correlation between SR and DM

#### **4.3 Species diversity: The identity of landmarks**

12 Will-be-set-by-IN-TECH

Although correlations between MD and SR could be weak or not significant, they were positive. A positive correlation has an intuitive explanation: more diversity means more forms. Weak and not significant correlation between SR and MD obtained from some combinations of LM could be the consequence of a few species immediately occupying the extremes of the morphospace. In that situation, by filling in the morphospace between disparate taxa, increasing taxonomic diversity was not able to have a significant effect on

The intensity of the correlation between SR and MD depended on the number of LM. Unexpectedly, the low number of LM (LNLM) configurations gave the best predictive power on SR (Fig. 2, left side). This is counter-intuitive since it is the current belief that more shape would have more taxonomic contents and hence would be more useful to distinguish species. We verified that by using one of the most predictive configuration of 4 LM (2-3-8-13) only 30% of individuals could be correctly assigned to their species, while more than 75% could be correctly attributed with 13 LM (the TNLM configuration). How could it be that a better discriminating configuration of LM could give a lower taxonomic prediction than a poorly discriminating one? The answer could lie in the observation that a LNLM configuration also could have no relationship at all with SR, and that the TNLM configuration is gathering both the predictive and non-predictive subsets of LM. Thus, the mixing of better predictive with poorly predictive components of shape into a larger set of LM was likely to mix opposite trends and blur the relationship between MD and SR (Fig. 2, right side). In other words, when a larger number of LM is used as a proxy for SR, a larger number of influences is also allowed which is not limited to SR. This explanation raises the question whether in our material there were other influences affecting the variation of the LM configurations that were not related to

> 0.006 0.007 0.008 0.009 0.01 0.011 0.012

Fig. 2. Two graphs representing the diversity structure (DS), i.e. the relationship between species richness (SR, on the horizontal axis) and metric diversity (MD, on vertical axis). The

left graph shows the DS appearing when using a set of four landmarks (namely the landmarks 2, 3, 8 and 13; see Fig. 1). The right graph shows the DS using the total number (13) of LM. When using the totality of LM (right graph) there is a first list of increasing MD values up to a SR of 30 species. This first list is common to the two graphs. In the right one however, after 30 species, there is a " plateau " breaking the global correlation between SR

0 10 20 30 40 50

**4.2 Species diversity: The number of landmarks**

0 10 20 30 40 50

MD estimates (Roy & Foote, 1997).

SR.

0.002 0.0022 0.0024 0.0026 0.0028 0.003 0.0032

and DM

For a given number of landmarks, the same configurations produced high (or low) correlation with SR regardless of the taxonomic composition of the species sequence. For instance, the 2-3-8-13 configuration always produced high correlation with SR regardless of the sequence of species considered, and the 1-3-6-9 configuration always produced one of the worst correlations. Thus, there seemed to be specific configurations of LM that could reflect the number of species, and others that could not, regardless of the species sequence used. This observation raises two questions:

The first question is about the effect of taxonomic diversity on wing shape variation: can the latter be used as a proxy for indirect species richness estimation?

Before starting to answer that question, one could ask what the interests are in making an indirect estimation of species richness. The answer is about time. Metric disparity extracted from wing venation is much faster to obtain than metrics requiring a taxonomic identification of each collected individual. Mosquitoes are often distinguished on the basis of labile characters which may be lost on damaged specimens, making species diagnostic difficult in a group where species are numerous (3500 species) and taxonomists are unfortunately not. An indirect but fast comparative estimation of SR would be welcome.

Our results showed that some landmarks configurations were highly predictive for SR. If these landmark configurations were known in advance, their variation could be used indeed as a proxy for SR. In our study, their performances seemed stable regardless of species arrangement in different groups (Table 5, last column). However, whether the high scoring of these LM configurations would be maintained in other studies on other mosquito species is still to be investigated.

The second question refers to other possible meanings of MD in our data set. Other than SR, what could be the cause of shape variation? It probably has many causes, certainly among which is species evenness (Deline, 2009), but also functional and ecological attributes (Roy & Foote, 1997), environmental conditions (Vasil'ev et al., 2010), founder effects (David, 1999; Whitlock & Fowler, 1999), endemicity (Magniez-Jannin et al., 2000) or reproduction mode (Baltanas et al., 2002). The only factor we could discuss with our data was the habitat, an environmental parameter to which insect metric properties are known to be sensible (Benítez et al., 2011; Tantowijoyo & Hoffmann, 2011). Insects were collected in the forest, in rural and urban areas. It was difficult to address the question of habitat influence for a given species since the within species sample sizes did not allow valid statistics inference (see Table 1). Thus, the habitat as explored here was not free from the possible effects due to different species compositions. Our ANOVA analyzes could however provide some indications. When both species richness and habitat heterogeneity were significantly contributing to the MD, a good relationship between MD and SR was not observed. A broad morphological covering (TNLM), or some unpredictive elementary units of morphology (LNLM configurations like 1-4-6-9 or 1-3-5-6), were under the influence of both species and habitat heterogeneities (Table

**8. References**

Baltanas, A., Alcorlo, P. & Danielopol, D. L. (2002). Morphological disparity in populations

The Mosquito Fauna: From Metric Disparity to Species Diversity 45

Benítez, H. A., Briones, R. & Jerez, V. (2011). Intra and inter-population morphological

the Baker River basin, Chilean Patagonia, *Journal of Insect Science, 11(94): 1–9*. Ciampaglio, C. N., Kemp, M. & McShea, D. W. (2001). Detecting changes in morphospace

David, P. (1999). A quantitative model of the relationship between phenotypic variance and heterozygosity at marker loci under partial selfing, *Genetics, 153: 1463–1474*. Deline, B. (2009). The effects of rarity and abundance distributions on measurements of local

Dolan, J. R., Jacquet, S. & Torréton, J.-P. (2006). Comparing taxonomic and morphological

Foote, M. (1992). Rarefaction Analysis Of Morphological And Taxonomic Diversity,

Foote, M. (1993). Discordance and Concordance Between Morphological and Taxonomic

Gaston, K. J. & Spicer, J. I. (2004). Biodiversity: An introduction, *Oxford: Blackwell Publishing*. Magniez-Jannin, F., David, B. & Dommergues, J. L. (2000). Analysing disparity by applying

Moyne, S. & Neige, P. (2007). The space-time relationship of taxonomic diversity and

Navarro, N. (2003). MDA: a MATLAB-based program for morphospace-disparity analysis,

Rattanarithikul, R., Harrison, B. A., Panthusiri, P. & Coleman, R. E. (2005). Illustrated keys to

Rattanarithikul, R., Harrison, B. A., Panthusiri, P., Peyton, E. L. & Coleman, R. E. (2006).

Rattanarithikul, R., Harrison, B., Harbach, R., Panthusiri, P. & Coleman, R. (2006.). Illustrated

Rattanarithikul, R., Harbach, R. E., Harrison, B. A., Panthurisi, P., Coleman, R. E. &

Rohlf, F. J. (1990). Rotational fit (Procrustes) methods, *in* F. Rohlf & F. Bookstein (eds),

*Southeast Asian J Trop Med Public health, Vol. 41 (Suppl.1): pp. 225*.

beetles, *Biological Journal of the Linnean Society 71(2): 343–358*.

Ostracoda), *Biological Journal of the Linnean Society 75(1): 9–19*.

disparity, *Paleobiology* 27(4): 695–715.

*Oceanography, 51 (2)* pp. 950–958.

Diversity, *Paleobiology* 19(2): 185–204.

*Palaeoclimatol. Palaeoecol., 248: 82–95*.

*Computers & Geosciences, 29: 655–664*.

*Med Public Health, 37(Suppl 1): 1–85*.

*Health, 36(Suppl 1): 1–80*.

*Health 37(Suppl.2): 1–128*.

*Paleobiology* 18(1): 1–16.

morphological disparity, *Paleobiology* 35(2): 175–189.

with and without sexual reproduction: a case study in *Eucypris virens* (Crustacea :

variation of shape and size of the chilean magnificent beetle, *Ceroglossus chilensis* in

occupation patterns in the fossil record: characterization and analysis of measures of

biodiversity of Tintinnids (planktonic ciliates) of New Caledonia, *Limnology and*

combined morphological and molecular approaches to French and Japanese carabid

morphological disparity in the Middle Jurassic Ammonite radiation, *Palaeogeogr.*

the mosquitoes of Thailand I. Background; geographic distribution; lists of genera, subgenera, and species; and a key to the genera., *Southeast Asian J Trop Med Public*

Illustrated keys to the mosquitoes of Thailand III. Genera *Aedeomyia*, *Ficalbia*, *Mimomyia*, *Hodgesia*, *Coquillettidia*, *Mansonia*, and *Uranotaenia*., *Southeast Asian J Trop*

keys to the mosquitoes of Thailand IV. *Anopheles*., *Southeast Asian J Trop Med Public*

Richardson, J. (2010). Illustrated keys to the mosquitoes of Thailand. VI. Tribe Aedini,

*Proceedings of the Michigan Morphometrics Workshop. Special Publiation Number 2.*

7). On the contrary, the highly predictive LNLM configurations (like 2-3-8-13 or 2-6-7-13) were apparently under the influence of SR only (Table 7).

#### **5. Perpectives**

Almost two decades ago, Foote (1993) claimed that " *discordances between morphological and taxonomic diversity demand to be interpreted biologically, not explained away as artifact of taxonomic practice* ". Morphological disparity is obviously under the influence of factors other than the mere species number, and we showed here for instance the likely influence of the mosquito habitat. However, this does not preclude the possibility to use shape variation as an indicator of species richness. We considered this possibility through decomposing wing shape into elementary components and comparing their respective relationships with taxonomic richness.

We suggest that the use of elementary units of shape (LNLM) could allow one to focus on a single factor, like species richness, and that in this regard the use of many characters (TNLM) has the inconvenience of mixing various effects, making a clear interpretation difficult. In a recent study showing good parallelism between SR and metric diversity, a single character was used (Dolan et al., 2006).

If the objective of the morphometric analysis was to accurately reflect one factor, then the use of LNLM is recommended, although not any LNLM configuration. The remaining question is: which LNLM configuration to use?

The answer to such a question certainly implies to explore the relationship of the LNLM configurations with known factors other than possible species richness, like the habitats, the localities of collection, etc. Any combination of landmarks which would vary under the influence of such parameters would be less likely to reliably reflect species richness alone. However, a more definitive answer cannot be provided through the use of a model which did not reflect natural conditions closely enough. As explained above, an investigation is still needed to evaluate the interference of species evenness on the relationship between MD and SR (Deline, 2009), an issue not contemplated in our model.

#### **6. Conclusion**

At this stage, our data confirmed that metric properties of a given community contain hidden but accurate information about species richness. Our model suggests that this information is likely to be found through the examination of some elementary shape configurations rather than of a global multivariate projection of many morphological traits.

#### **7. Acknowledgment**

We thank R. Rattanarithikul for her kind help in taxonomic identification, as well as P. Kittayapong and C. Apiwhatnasorn (Mahidol University, Bangkok, Thailand) who encouraged this work. This study has been supported by the TRF/BIOTEC Special Program for Biodiversity Research and Training Grant BRT R\_352052 and the Commission on Higher Education, RMU5080060 (Thailand), as well as by the IRD grants number HC3165-3R165-GABI-ENT2 and HC3165-3R165-NV00-THA1.

#### **8. References**

14 Will-be-set-by-IN-TECH

7). On the contrary, the highly predictive LNLM configurations (like 2-3-8-13 or 2-6-7-13) were

Almost two decades ago, Foote (1993) claimed that " *discordances between morphological and taxonomic diversity demand to be interpreted biologically, not explained away as artifact of taxonomic practice* ". Morphological disparity is obviously under the influence of factors other than the mere species number, and we showed here for instance the likely influence of the mosquito habitat. However, this does not preclude the possibility to use shape variation as an indicator of species richness. We considered this possibility through decomposing wing shape into elementary components and comparing their respective relationships with taxonomic

We suggest that the use of elementary units of shape (LNLM) could allow one to focus on a single factor, like species richness, and that in this regard the use of many characters (TNLM) has the inconvenience of mixing various effects, making a clear interpretation difficult. In a recent study showing good parallelism between SR and metric diversity, a single character

If the objective of the morphometric analysis was to accurately reflect one factor, then the use of LNLM is recommended, although not any LNLM configuration. The remaining question

The answer to such a question certainly implies to explore the relationship of the LNLM configurations with known factors other than possible species richness, like the habitats, the localities of collection, etc. Any combination of landmarks which would vary under the influence of such parameters would be less likely to reliably reflect species richness alone. However, a more definitive answer cannot be provided through the use of a model which did not reflect natural conditions closely enough. As explained above, an investigation is still needed to evaluate the interference of species evenness on the relationship between MD and

At this stage, our data confirmed that metric properties of a given community contain hidden but accurate information about species richness. Our model suggests that this information is likely to be found through the examination of some elementary shape configurations rather

We thank R. Rattanarithikul for her kind help in taxonomic identification, as well as P. Kittayapong and C. Apiwhatnasorn (Mahidol University, Bangkok, Thailand) who encouraged this work. This study has been supported by the TRF/BIOTEC Special Program for Biodiversity Research and Training Grant BRT R\_352052 and the Commission on Higher Education, RMU5080060 (Thailand), as well as by the IRD grants number

apparently under the influence of SR only (Table 7).

**5. Perpectives**

richness.

**6. Conclusion**

**7. Acknowledgment**

was used (Dolan et al., 2006).

is: which LNLM configuration to use?

SR (Deline, 2009), an issue not contemplated in our model.

than of a global multivariate projection of many morphological traits.

HC3165-3R165-GABI-ENT2 and HC3165-3R165-NV00-THA1.

**3**

Nikhil Whitaker,

*Tamil Nadu,*

*India*

**Reproduction and Morphohlogy**

**(***Indotestudo travancorica***) Boulenger, 1907**

The most threatened chelonians occur in Asia, where virtually all species are heavily harvested for food and traditional medicinal trades (van Dijk *et al.* 2000). This includes the two endemic chelonians, *Indotestudo travancorica*, and the sympatric *Vijayachelys silvatica* in the Region of the Western Ghats, South India. *Indotestudo travincorica* is listed as Vulnerable under the IUCN Red

Phylogenetic relationships between the three species of *Indotestudo*, these being *forstenii, elongata*, and *travancorica*, have been in flux. Iverson *et al* (2001) examined these clades, and found that *I.travancorica*, the species in question here, was found to be more closely related to *I. elongata*. They also place *Indotestudo forstenii* from Sulawesi and Halmahera, as a distinct species unrelated to *I.travancorica* and *I. elongata*, refuting information that *I. forstenii* were

In 1982, 14 sub-adult *Indotestudo travancorica* were collected from Kerala (10° 5, 76°4), by the then MCBT (Madras Crocodile Bank Trust) researcher J. Vijaya, who also worked on *Vijayachelys silvatica* (Whitaker & Jaganathan, 2009). Between 1988 – 1995, the captive group of Travancore tortoises (*Indotestudo travancorica*) at the MCBT have produced 21 clutches of eggs (Whitaker and Andrews, 1997). However, not much information resulted from these clutches, as most of them were infertile (but see below). Apart from this publication, there are only a small number of other publications on the biology of this species, by Appukuttan (1991), Ramesh (2002, 2003, 2007, 2008), Bhupathy and Choudhury (1995), Vasudevan et al (2010), and Vijaya (1983). Here I examine factors related to the reproductive biology, morphology, and temperature selection in the species, in between the years 1999 – 2001, and

The adult *Indotestudo travancorica* breeding pen at the MCBT houses two males, fifteen females , and two sub-adults. It consists of two interconnected exclosures, measuring 105 sq meters. One of the two exclosures was meshed off for a *Geochelone elegans* enclosure, on 18th

List, and comes under Schedule IV of the Indian Wildlife (Protection) Act.

**1. Introduction**

2008 – 2011.

**2. Methods**

introduced from India to Indonesia.

**of the Travancore Tortoise**

*Madras Crocodile Bank Trust, Mamallapuram,*

*The University of Michigan Museum of Zoology. Ann Arbor, MI, pp380*, University of Michigan Museums, Ann Arbor, pp. 227–236.

## **Reproduction and Morphohlogy of the Travancore Tortoise (***Indotestudo travancorica***) Boulenger, 1907**

Nikhil Whitaker, *Madras Crocodile Bank Trust, Mamallapuram, Tamil Nadu, India*

#### **1. Introduction**

16 Will-be-set-by-IN-TECH

46 Morphometrics

Rohlf, F. J. (1996). Morphometric spaces, shape components and the effects of linear

Roy, K., Balch, D. P. & Hellberg, M. E. (2001). Spatial patterns of morphological diversity

Roy, K. & Foote, M. (1997). Morphological approaches to measuring biodiversity, *Trends In*

Shannon, C. E. & Weaver, W. (1949). The mathematical theory of communication, *University*

Tantowijoyo, W. & Hoffmann, A. A. (2011). Variation in morphological characters of two

Vasil'ev, A. G., Vasil'eva, I., Gorodilova, Y. V. & Chibiryak, M. V. (2010). Morphological

Wills, M. A., Briggs, D. E. G. & Fortey, R. A. (1994). Disparity As An Evolutionary Index - A Comparison Of Cambrian And Recent Arthropods, *Paleobiology* 20(2): 93–130. Zelditch, M. L., Swiderski, D. L., Sheets, H. D. & Fink, W. L. (2004). *Geometric morphometrics*

*for biologists: a primer*, Elsevier, Academic Press. New-York.

Michigan Museums, Ann Arbor, pp. 227–236.

Simpson, E. H. (1949). Measurement of diversity, *Nature, 163: 688*.

gradient, *Journal of Insect Science, 11(69): 1–16*.

*Ecology & Evolution* 12(7): 277–281.

*of Illinois Press, Urbana*.

*Evolution, 53(4): 1143–1156*.

*2503–2508*.

Plenum Publ. NATO ASI, ser. A, Life Sciences, pp. 117–129.

*The University of Michigan Museum of Zoology. Ann Arbor, MI, pp380*, University of

transformations, *in* L. F. Marcus, M. Corti, A. Loy, G. Naylor & D. Slice (eds), *Advances in Morphometrics. Proceedings of the 1993 NATO-ASI on Morphometrics*, New York:

across the Indo-Pacific: analyses using strombid gastropods, *Proc. R. Soc. Lond., 268:*

invasive leafminers, *Liriomyza huidobrensis* and *L. sativae*, across a tropical elevation

disparity in populations with and without sexual reproduction: a case study in *Eucypris virens* (Crustacea : Ostracoda), *Russian Journal of Ecology, 41(2): 153–158*. Whitlock, M. C. & Fowler, K. (1999). The distribution of phenotypic variance with inbreeding, The most threatened chelonians occur in Asia, where virtually all species are heavily harvested for food and traditional medicinal trades (van Dijk *et al.* 2000). This includes the two endemic chelonians, *Indotestudo travancorica*, and the sympatric *Vijayachelys silvatica* in the Region of the Western Ghats, South India. *Indotestudo travincorica* is listed as Vulnerable under the IUCN Red List, and comes under Schedule IV of the Indian Wildlife (Protection) Act.

Phylogenetic relationships between the three species of *Indotestudo*, these being *forstenii, elongata*, and *travancorica*, have been in flux. Iverson *et al* (2001) examined these clades, and found that *I.travancorica*, the species in question here, was found to be more closely related to *I. elongata*. They also place *Indotestudo forstenii* from Sulawesi and Halmahera, as a distinct species unrelated to *I.travancorica* and *I. elongata*, refuting information that *I. forstenii* were introduced from India to Indonesia.

In 1982, 14 sub-adult *Indotestudo travancorica* were collected from Kerala (10° 5, 76°4), by the then MCBT (Madras Crocodile Bank Trust) researcher J. Vijaya, who also worked on *Vijayachelys silvatica* (Whitaker & Jaganathan, 2009). Between 1988 – 1995, the captive group of Travancore tortoises (*Indotestudo travancorica*) at the MCBT have produced 21 clutches of eggs (Whitaker and Andrews, 1997). However, not much information resulted from these clutches, as most of them were infertile (but see below). Apart from this publication, there are only a small number of other publications on the biology of this species, by Appukuttan (1991), Ramesh (2002, 2003, 2007, 2008), Bhupathy and Choudhury (1995), Vasudevan et al (2010), and Vijaya (1983). Here I examine factors related to the reproductive biology, morphology, and temperature selection in the species, in between the years 1999 – 2001, and 2008 – 2011.

#### **2. Methods**

The adult *Indotestudo travancorica* breeding pen at the MCBT houses two males, fifteen females , and two sub-adults. It consists of two interconnected exclosures, measuring 105 sq meters. One of the two exclosures was meshed off for a *Geochelone elegans* enclosure, on 18th

Reproduction and Morphohlogy of the Travancore

from nests following compaction of the nest by the female.

abnormally large yolk sacs, with these hatchings not surviving.

**3.1 Clutch and egg sizes of** *Indotestudo travancorica*

all analyses was *P <0.05*.

**3. Results & discussion**

Tortoise (*Indotestudo travancorica*) Boulenger, 1907 49

Incubation temperature was monitored utilizing a HOBO XT ™ automatic temperature logger (Onset Instruments, P.O. Box 3450, Pocasset, MA 02559, USA) set to record temperature 4 times a day throughout the incubation period, at 800, 1200, 1600, and 2000 hrs for a group of 11 clutches laid in 1999. To get mean temperatures of all of the incubation boxes, the thermocouple from the logger was positioned in the middle of all boxes. Eggs were placed in 10 centimeter by 8 centimeter Tupperware boxes, within a plywood cupboard, to avoid predation of eggs by *Rattus bengalensis*. Incubation substrate was fine grain vermiculite. In addition, a bottle of beer with ca. ¼ of the contents left, and a piece of small beef within, was used in 2009 – 2011 to provide additional protection of eggs from fruit flies, as described by Wolff (2007) When oviposition was observed, eggs were removed

Eggs that have rotted or cracked were discarded to avoid attracting flies and ants to other eggs in the same box. This was clearly evident from loss in egg weight and/or a sulphur-like smell. Morphometric measurements were recorded once hatching occurred, these being carapace length (millimeters), carapace width (millimeters), plastron length (millimeters), plastron width (millimeters), and weight (grams). A "successful hatching" was defined as an event where a hatchling managed to "pip", emerge from the egg, and survive >7 – 10 days post hatch This was the number of days typically required to absorb external yolk from pre-mature hatchlings. Hatchlings were measured and weighed within 12 hours of hatching, residual yolk outside the abdominal cavity varied from small, normal quantities, to

To record tortoise temperature, an I-button TM was lodged on top of the carapace, and affixed with Gorilla glue ™. At between 30 – 90 days, animals were captured, and temperatures for the preceding days were recorded via the USB interface provided with the I-button software/I-buttons onto a Lenovo ™ laptop. Resolution of the I-buttons was at ± 0.5 °C. One male and four female tortoises was were measured for a period of 400+ days, and a sample is given here for 157 days, resulting in 943 observations (temperature loggers were programmed to record temperature every six hours), between 30th May – 3rd December 2011. The number of observations per individual varies, as loggers fell off of animals. Searches in the leave litter had

Statistical analysis follows that of McGuinness (1999) and Fowler et al*.* (1998). In addition, SPSS 10.0 was used to confirm manual calculation of statistics. The level of significance for

Clutch size of *Indotestudo travancorica* clutches collected ranged from one to six eggs. Clutch sizes and morphology of eggs is presented in Table 1. One particularly large egg from a clutch laid on the of 9th September 2011 measured 152.39 mm (egg length), 123.41 mm (egg width), and weighed 85 grams. Sane & Sane (1988) hatched a single *Indotestudo travancorica*, incubation period for this one egg was 139 days. Das (1995) puts the average incubation period of this species at 146-149 days. Incubation period for the seven hatched *Indotestudo travancorica* in 1999 ranged from 128 – 159 days (*X* = 141.5 days), with temperatures ranging between 22.4 º C - 28.7

to be done to retrieve lost loggers. Temperatures for these periods were eliminated.

º C. All clutches with more than one egg had hatchlings emerging on different days.

December 2010, which left a single exclosure for *Indotestudo travancorica*, reducing area to 87 sq meters*.* The exclosures contain *Pongamia* and *Bambusa*, shade trees. Leaf litter is removed and replaced every six months. A circular pond 100 centimeters in diameter varies from 5 – 10 cm in depth, and this provides drinking water and heat sinks. Hatchlings are housed in a separate outdoor terrarium measuring 2.54 x 0.55 x 0.28 m. The substratum is beach sand covered with dry bamboo leaf litter, and both substrata are changed on a regular basis. A bowl of water is provided filled with pebbles leaving 2 – 4 cm exposed, to prevent accidental drowning. V-shaped roof tiles are used as shelters.

Juvenile and adult tortoises are fed on tomato, carrots, beans, and various types of spinach, pumpkin, grasses, and beef*.* Food is placed at three feeding stations, on granite slabs measuring 1.2 x 1.2 meters, and 5 centimeters off of the exclosure floor. This occurs between 1000 – 1030 hrs, and 1530 – 1600 hrs, and any remnant feed is removed the following morning/evening Precise quantities of feed offered are not recorded. Feed is provided every day.

Chelonian embryos are at the gastrula stage (at a presomite stage of development) at oviposition ( Ewert, 1979), as compared to crocodilians, with 16 – 18 somites at oviposition (Ferguson, 1985). The most advanced eggs of the Reptilia at egg laying are the lepidosaurs with 20-30 pairs of somites at oviposition (Muthukkaruppan *et. al.* 1970).Eggs are candled utilizing the technique used to candle crocodilian eggs as described by Hutton & Webb (1990), under a focused light, to determine the presence or absence of sub-embryonic fluid. This is found to be inaccurate in *I. travancorica*, and eggs deemed non viable hatch and vice versa (see section 3.3) Part of the chalking process in turtles, the development of an opaque band, involves adhesion of the vitelline membrane to the shell membrane, and this determines to some degree subsequent embryonic orientation (Ewert 1985). In addition, Andrews and Mathies (2000), note that adhesion of the embryo to the egg shell following oviposition may have a respiratory purpose "since chalking (drying) of the shell that is associated with adhesion increases the conductance to gases". However, in the travancore tortoise, viability of eggs is best confirmed by the presence of vascularisation under a candling lamp around a month into incubation, and by embryonic movement later in incubation.

Eggs were collected from 22 nests laid between December 1999 and September 2011 (Table 1). Eggs and hatchling tortoises were measured with vernier calipers (±0.1 mm) and weighed with an ACCULAB ™electronic weighing scale (±1 gram resolution). Eggs were then candled for presence or absence of sub-embryonic fluid. Clutches were segregated into different boxes, half immersed in vermiculite media, and allowed to develop at ambient room temperature. When egg laying was observed, female morphometrics, namely carapace length (mm), carapace width (mm), plastron length (mm), and weight (grams) were measured after females had compacted the nest. Where sexual size dimorphism in adults was done, four male and eighteen females were used from both the current MCBT population and animals measured in-situ by the late MCBT researcher J. Vijaya. Two females were included from a study by Bhupathy & Choudhury (1995), and nineteen females were from a dataset collected by Arun Kalagaven, from Eastern Kerala. The online database of the Global Diversity Information Facility revealed 26 specimens of unknown sex/size residing at the Florida Museum of Natural History (GBIF, 2011). Locations and GPS coordinates are not given on purpose, to protect in-situ populations.

December 2010, which left a single exclosure for *Indotestudo travancorica*, reducing area to 87 sq meters*.* The exclosures contain *Pongamia* and *Bambusa*, shade trees. Leaf litter is removed and replaced every six months. A circular pond 100 centimeters in diameter varies from 5 – 10 cm in depth, and this provides drinking water and heat sinks. Hatchlings are housed in a separate outdoor terrarium measuring 2.54 x 0.55 x 0.28 m. The substratum is beach sand covered with dry bamboo leaf litter, and both substrata are changed on a regular basis. A bowl of water is provided filled with pebbles leaving 2 – 4 cm exposed, to prevent accidental

Juvenile and adult tortoises are fed on tomato, carrots, beans, and various types of spinach, pumpkin, grasses, and beef*.* Food is placed at three feeding stations, on granite slabs measuring 1.2 x 1.2 meters, and 5 centimeters off of the exclosure floor. This occurs between 1000 – 1030 hrs, and 1530 – 1600 hrs, and any remnant feed is removed the following morning/evening Precise quantities of feed offered are not recorded. Feed is provided every

Chelonian embryos are at the gastrula stage (at a presomite stage of development) at oviposition ( Ewert, 1979), as compared to crocodilians, with 16 – 18 somites at oviposition (Ferguson, 1985). The most advanced eggs of the Reptilia at egg laying are the lepidosaurs with 20-30 pairs of somites at oviposition (Muthukkaruppan *et. al.* 1970).Eggs are candled utilizing the technique used to candle crocodilian eggs as described by Hutton & Webb (1990), under a focused light, to determine the presence or absence of sub-embryonic fluid. This is found to be inaccurate in *I. travancorica*, and eggs deemed non viable hatch and vice versa (see section 3.3) Part of the chalking process in turtles, the development of an opaque band, involves adhesion of the vitelline membrane to the shell membrane, and this determines to some degree subsequent embryonic orientation (Ewert 1985). In addition, Andrews and Mathies (2000), note that adhesion of the embryo to the egg shell following oviposition may have a respiratory purpose "since chalking (drying) of the shell that is associated with adhesion increases the conductance to gases". However, in the travancore tortoise, viability of eggs is best confirmed by the presence of vascularisation under a candling lamp around a month into incubation, and by embryonic movement later in

Eggs were collected from 22 nests laid between December 1999 and September 2011 (Table 1). Eggs and hatchling tortoises were measured with vernier calipers (±0.1 mm) and weighed with an ACCULAB ™electronic weighing scale (±1 gram resolution). Eggs were then candled for presence or absence of sub-embryonic fluid. Clutches were segregated into different boxes, half immersed in vermiculite media, and allowed to develop at ambient room temperature. When egg laying was observed, female morphometrics, namely carapace length (mm), carapace width (mm), plastron length (mm), and weight (grams) were measured after females had compacted the nest. Where sexual size dimorphism in adults was done, four male and eighteen females were used from both the current MCBT population and animals measured in-situ by the late MCBT researcher J. Vijaya. Two females were included from a study by Bhupathy & Choudhury (1995), and nineteen females were from a dataset collected by Arun Kalagaven, from Eastern Kerala. The online database of the Global Diversity Information Facility revealed 26 specimens of unknown sex/size residing at the Florida Museum of Natural History (GBIF, 2011). Locations and GPS

coordinates are not given on purpose, to protect in-situ populations.

drowning. V-shaped roof tiles are used as shelters.

day.

incubation.

Incubation temperature was monitored utilizing a HOBO XT ™ automatic temperature logger (Onset Instruments, P.O. Box 3450, Pocasset, MA 02559, USA) set to record temperature 4 times a day throughout the incubation period, at 800, 1200, 1600, and 2000 hrs for a group of 11 clutches laid in 1999. To get mean temperatures of all of the incubation boxes, the thermocouple from the logger was positioned in the middle of all boxes. Eggs were placed in 10 centimeter by 8 centimeter Tupperware boxes, within a plywood cupboard, to avoid predation of eggs by *Rattus bengalensis*. Incubation substrate was fine grain vermiculite. In addition, a bottle of beer with ca. ¼ of the contents left, and a piece of small beef within, was used in 2009 – 2011 to provide additional protection of eggs from fruit flies, as described by Wolff (2007) When oviposition was observed, eggs were removed from nests following compaction of the nest by the female.

Eggs that have rotted or cracked were discarded to avoid attracting flies and ants to other eggs in the same box. This was clearly evident from loss in egg weight and/or a sulphur-like smell. Morphometric measurements were recorded once hatching occurred, these being carapace length (millimeters), carapace width (millimeters), plastron length (millimeters), plastron width (millimeters), and weight (grams). A "successful hatching" was defined as an event where a hatchling managed to "pip", emerge from the egg, and survive >7 – 10 days post hatch This was the number of days typically required to absorb external yolk from pre-mature hatchlings. Hatchlings were measured and weighed within 12 hours of hatching, residual yolk outside the abdominal cavity varied from small, normal quantities, to abnormally large yolk sacs, with these hatchings not surviving.