Demographic history, genetic structure and gene flow in a steppe-associated raptor species
BMC Evolutionary Biology
Demographic history, genetic structure and gene flow in a steppe-associated raptor species
Jesus T Garcia 2
Fernando Alda 2
Julien Terraube 1 2
Franois Mougeot 0
Audrey Sternalski 2 3
Vincent Bretagnolle 3
Beatriz Arroyo 2
0 Estacion Experimental de Zonas Aridas (EEZA-CSIC) , Ctra. de Sacramento s/ n , La Canada de San Urbano , 04120 Almeria , Spain
1 Natural Research Ltd, Banchory Business Centre , Burn O'Bennie Road, Banchory, AB31 5ZU , UK
2 Instituto de Investigacion en Recursos Cinegeticos (IREC) (CSIC-UCLM-JCCM). Ronda de Toledo s/n, E-13005 Ciudad Real , Spain
3 CEBC-CNRS , UPR1934, 79360 Beauvoir sur Niort , France
Background: Environmental preferences and past climatic changes may determine the length of time during which a species range has contracted or expanded from refugia, thereby influencing levels of genetic diversification. Connectivity among populations of steppe-associated taxa might have been maximal during the long glacial periods, and interrupted only during the shorter interglacial phases, potentially resulting in low levels of genetic differentiation among populations. We investigated this hypothesis by exploring patterns of genetic diversity, past demography and gene flow in a raptor species characteristic of steppes, the Montagu's harrier (Circus pygargus), using mitochondrial DNA data from 13 breeding populations and two wintering populations. Results: Consistent with our hypothesis, Montagu's harrier has relatively low genetic variation at the mitochondrial DNA. The highest levels of genetic diversity were found in coastal Spain, France and central Asia. These areas, which were open landscapes during the Holocene, may have acted as refugia when most of the European continent was covered by forests. We found significant genetic differentiation between two population groups, at the SW and NE parts of the species' range. Two events of past population growth were detected, and occurred ca. 7500-5500 and ca. 3500-1000 years BP in the SW and NE part of the range respectively. These events were likely associated with vegetation shifts caused by climate and human-induced changes during the Holocene. Conclusions: The relative genetic homogeneity observed across populations of this steppe raptor may be explained by a short isolation time, relatively recent population expansions and a relaxed philopatry. We highlight the importance of considering the consequence of isolation and colonization processes in order to better understand the evolutionary history of steppe species.
Contemporary patterns of genetic diversity and
population structure reflect not only current patterns of
genetic exchange but also past dispersal processes and
levels of gene flow among populations during historical
climatic events [1,2]. In many species inhabiting
temperate zones, climate-vegetation feedbacks during the
Pleistocene caused range contractions to lower latitudes
followed by range expansions during interglacial periods
[2,3], which in turn promoted much of the
diversification observed today. However, not all species responded
similarly to these past climatic events. Species-specific
responses to these changes are the result of a complex
1Instituto de Investigacin en Recursos Cinegticos (IREC) (CSIC-UCLM-JCCM).
Ronda de Toledo s/n, E-13005 Ciudad Real, Spain
Full list of author information is available at the end of the article
interplay between the behavioral, physiological and
ecological characteristics of the species, including their
biogeographic origin, habitat preferences and dispersal
capabilities . Climatic conditions prevailing at
different time periods strongly influenced the extent of each
habitat type in the past, which should in turn influence
the length of time during which a species underwent
isolation or range expansion and, consequently, the
opportunities for genetic diversification [3,5]. For
example, species inhabiting arctic or boreal areas seem to
have experienced range expansions during the long
glacial periods, but remained isolated during the short
interglacial ones [6-8], leading to a pattern of
contemporary genetic structure different from that of species
inhabiting temperate areas. The phylogeography of
temperate and arctic species is rather well studied. In
contrast, the phylogeography of steppe species, which are
biogeographically in-between the temperate and
arcticboreal fauna, remains scarcely known [4,9]. During the
glacial periods, many of these species were widely
distributed throughout the periglacial steppes of the
northern hemisphere. For these steppe species, gene flow at
large geographic scales might have been interrupted by
the postglacial retreat and reduction of steppe vegetation
during the short interglacial phases. Therefore, the
amount of time spent in isolation, and the resulting
genetic differentiation, should be smaller for steppe
species as compared with temperate ones. The distribution
of steppe species during the glacial periods is well
documented by fossil evidence [10-12], but the genetic
evidence is still relatively poor [9,13-15]. In addition,
steppe-like ecosystems (including natural and
agricultural landscapes) are probably amongst the most altered
habitats nowadays, due to human pressure and rapid
changes in land use . Consequently, the abundance
and distribution range of many steppe/farmland species
has greatly declined in recent years [17-22].
The Montagus harrier Circus pygargus is perhaps the
best example of raptor species specialized in such steppe
habitats, which has also adapted well to farmland
habitats [23,24]. This species is widely but patchily
distributed across the Palearctic region and undertakes
longdistance migrations  (Figure 1). Migration studies
based on satellite tracking and ringing recoveries show
that European and Asian breeding populations follow
different migratory pathways and overwinter in different
continents (Figure 1). Central Asian populations, where
it is accepted that the stronghold of the species occurs
(ca. 25,000-30,000 pairs ), migrate along an eastern
route and overwinter in the Indian subcontinent
(Pakistan and India ). European populations migrate
through the Mediterranean peninsulas and overwinter in
the Sahel belt from western to eastern Africa [28-30].
To date, no specific study has evaluated the genetic
structure and gene flow among populations of this
species or closely related ones.
We used molecular analyses to evaluate levels of
genetic diversity and connectivity among Montagus
harrier breeding populations. We collected genetic material
from across the whole breeding range (Spain to
Kazakhstan) and from two overwintering populations (located
in western Africa and in the Indian sub-continent).
With this information, we analyzed historical
demographic patterns, in order to better understand how past
climatic dynamics may have affected the species and
shaped its current distribution and population genetic
structuring. Specifically, our aims were: first, to evaluate
the genetic diversity of Montagus harrier populations
and analyze different scenarios that may explain the
observed genetic structure (according to physical
barriers, overwintering areas or current breeding
Figure 1 Map showing the geographic range of breeding
populations of Montagus harrier (redrawn from) and
sampling localities used in this study. Orange and red colors
correspond to breeding grounds; the large lump in the breeding
distribution in the old USSR reflects a lack of detailed knowledge on
their distribution rather than a continuous range. The figure shows
the separations between the two main distribution areas within the
range. Green colors correspond to wintering grounds. Locality
codes are as indicated in Table 1. The doted black line shows the
location of the intercontinental migratory divide proposed by
Moreau  to split the breeding populations wintering in the
Indian sub-continent from those wintering in Africa
distribution); second, to investigate past demographic
changes in the species, and patterns of gene flow
between genetically structured populations, looking for
evidence of asymmetrical exchanges between these.
Finally, we discuss the implications of these results for
the understanding of the patterns shaping the
phylogeographic structure of the fauna associated with steppe
We sequenced a total of 1063 bp of t-RNA Trp and
ND2 and 714 bp of COI genes from all samples (n =
284), which were collapsed into 36 and 13 haplotypes,
respectively (Additional file 1). The concatenated data
set (1777 bp) yielded 51 different haplotypes defined by
35 polymorphic sites. Overall, gene and nucleotide
diversities were Hd = 0.663 (SD, 0.032) and = 0.0008
(SD, 0.00007), respectively (Table 1). Most of the
observed genetic variability was in the form of single
nucleotide substitutions. The mean number of pairwise
Table 1 Sample sizes and genetic characteristics.
Summary of the number of samples (n), number of haplotypes (nH), haplotype diversity (Hd) and nucleotide diversity () for the concatenated data set in each
sampling locality. Tajimas (D), Fus (FS) measure departure of the data from neutrality. The significance (at the P = 0.05 level) of these statistics is indicated with
an asterisk (*) for each breeding population and for the SW and NE demes (see text). SD = standard deviation.
nucleotide differences for all samples (concatenated data
set) was 1.47 (SD, 0.94).
French populations showed the highest nucleotide and
haplotype diversities when considering each gene partition
separately (ND2 and COI, results not shown) or the
concatenated data set (Table 1). Haplotype diversities were
also high in the coastal regions of Spain (Galicia and
Castelln). The lowest haplotype diversities were found in
populations from Germany and Netherlands (Table 1).
Regarding the two wintering sites, levels of genetic
diversity were greater in Pakistan than in Senegal (Table 1).
The two phylogenetic methods used (Maximum
likelihood and Bayesian inference) were largely consistent in
the (lack of) relationships recovered. Despite the broad
geographic sampling (Figure 1), there was little
phylogenetic structure and low branch support (Figure 2). The
haplotype network revealed no major branching events
(Figure 2), although two groups of haplotypes could be
differentiated. The first group was distributed around
haplotype Hap2, which was observed in 57% of
individuals and had a widespread geographical distribution.
The other haplotypes were generally site-specific and
occurred at low frequencies. The second group
consisted mostly of rare haplotypes from all geographic
regions except Senegal. In general, haplotypes specific to
certain geographic regions did not form monophyletic
groups, but appeared to be randomly distributed across
the network (Figure 2).
Differences among Montagus harrier populations
accounted for a significant 5.15% of total molecular
variance (Table 2). Hierarchical analyses testing for
differences between populations in relation to their
geographic origin showed low but significant
differentiation among groups of populations located in the
southwestern (SW) vs. north-eastern (NE) parts of the
breeding range (scenario 1), with low levels of genetic
variation among populations within each group. In contrast,
the structure scenarios based on geographic barriers
(scenario 2) or based on an intercontinental migratory
divide (scenario 3) revealed no significant among-group
differences (Table 2). The FCT estimate obtained for
scenario 3 was similar to that obtained for scenario 1,
and close to significance (Table 2). This suggested some
structuring in relation to the migratory divide, although
the higher and significant value of FST indicated larger
heterogeneity among populations within groups than in
scenario 1. Including samples from wintering areas
might affect, to some extent, these results, but the
exclusion of these wintering samples also entails the loss
of important genetic information. We repeated these
analyses without including the samples from the two
overwinter sites, and results were qualitatively the same.
Figure 2 Phylogenetic relationships among mitochondrial DNA
haplotypes. A) Maximum-likelihood (ML) tree of Circus pygargus
based on mtDNA sequences. Numbers above branches indicate ML
bootstrap values (1000 replicates) and numbers below branches
indicate BI posterior probabilities. Only bipartitions with bootstrap or
posterior probability values above 50 and 0.5, respectively, are
shown. B) Median-joining haplotype network based on the mtDNA
concatenated data set. Circle size is proportional to haplotype
frequency and connecting lines are proportional to mutation steps
between haplotypes. White circles represent hypothetical
intermediate or un-sampled haplotypes. Each color refers to one
geographic area (orange: SW populations; red: NE populations; dark
green: Senegal; light green: Pakistan).
Overall, 27 out of 105 pairwise FST comparisons
between populations were significant. Comparisons
involving Pakistan were significant in all but four cases
(Table 3). Differences were smallest between Pakistan
Table 2 AMOVA summary.
and Czech Republic (FST = 0.037) and largest between
Pakistan and western Spain (PAK-SpW: FST = 0.324).
In contrast, comparisons involving Senegal were
significant only with central France (PLO) and Pakistan
(SENPLO: FST = 0.159, and SEN-PAK: FST = 0.224).
Considering only breeding populations, eastern Spain (SpE)
showed the highest number of significant pairwise
comparisons (seven out of 12), while the largest
differentiation among populations was found between western
Spain (SpW) and central France (PLO) (SpW-PLO: FST
= 0.216). Finally, the correlation between population
pairwise FST values and their geographical distances
was not significant (Mantel test, r = -0.06, p = 0.70).
The whole breeding population exhibited significantly
negative values for both Tajimas D (D = -2.28, p =
0.0023) and Fus FS (FS = -68.94, p < 0.001), suggesting
that the overall population size has fluctuated in the past.
The majority of D and FS values were negative for most
of the studied populations, but Fus FS rejected neutrality
(p < 0.02) in only two breeding populations, whereas
Tajimas D rejected neutrality (p < 0.05) in five
populations (Table 1). We further considered the two groups of
populations (from the SW and NE regions) whose
genetic structure differed according to AMOVA analysis
(see above). Each group of populations also showed
negative D and FS values that were all significant except for
the Tajimas D of the NE group (Table 1).
The effective population sizes and demographic trends
estimated by the Bayesian Skyline Plot (BSP) analysis
indicated recent population size increases in both
regions (SW and NE). However, the overall increase was
much less marked for the SW populations, which
showed a roughly linear increase during the second part
of the Holocene. Based on a range of mutation rates
from 0.02 to 0.055 s/s/Myr, population growth started
at approximately 7500 years before present (BP) in SW
and 4000 years BP in NE (Figure 3). The time to the
2. Geographic barriers
% of total variance (p-value)
Grouping Populations in group Between Between
group population within
FCT group FST
Southwestern vs. Northeastern (SpW, SpC, SpNW, SpE, PACA, PLO, PCH, AQPY, 4.23 3.21 (0.011)
CHAR) vs. (NED, GER, CZE, KZ, PAK) (0.0009)
Spain vs. Western Europe vs. (SpW, SpC, SpNW, SpE) vs. (PACA, PLO, PCH, 2.48 (0.114) 3.42 (0.007)
Eastern Europe vs. central Asian AQPY, CHAR, NED, GER) vs. (CZE) vs. (KZ, PAK)
Africa vs. India (SpW, SpC, SpNW, SpE, PACA, PLO, PCH, AQPY, 4.55 (0.087) 3.66 (0.0009)
CHAR, NED, GER, CZE, SEN) vs. (KZ, PAK)
Results of AMOVA (concatenated data, Tamura and Nei; gamma = 0.2487) testing for differences between Montagus harrier populations and geographic areas.
The analysis partitions out total molecular variance into different components, and statistical significance is obtained by randomization after 5000 permutations.
Significant differences between groups are shown in bold.
Table 3 Between-population genetic differentiation in Circus pygargus.
Below the diagonal, pairwise FST values (based on Tamura and Nei distances between haplotypes; gamma = 0.248). Above the diagonal, p-values for significant
FST values obtained after 5000 permutations (* p < 0.05; ** p < 0.01, *** p < 0.001).
most recent common ancestor (TMRCA) was estimated
at 21,000 (40,000-6000 95% highest posterior density,
HPD) and 35,000 years BP (61,000-13,000 95% HPD)
for SW and NE groups, respectively. Independent runs
of IMa gave similar results, and plots of parameter
trends indicated sufficient mixing among chains. The
estimated effective population size for the SW region
Figure 3 Bayesian skyline plot generated with BEAST v1.5.4.
and effective population size (Ne) and time estimates under
IMa. The solid red lines are the median and 95% confidence
intervals (CI) estimates for SW populations while dotted blue lines
for NE populations. The shaded box in the centre of the figure
represents the mean and 95%CI for the time of population splitting
estimated by IMa for SW and NE groups of Montagus harriers.
Boxes represent 95% CI for the Ne estimated by IMa with SW and
NE groups and the ancestral population in grey. Lines inside the
boxes represent the peak of the marginal distribution values.
(peak SW= 75.18, 90% HPD = 41.4 - 123.1; Figure 3)
was similar to that of the NE region (peak NE at 63.45,
90% HPD = 16.5 - 397.0; Figure 3), despite larger HPD
in NE than in SW group. Both descendant populations
had effective population sizes that were over one order
of magnitude larger than the ancestral population (A),
which peaked at 10.28 (90% HPD = 2.82 - 17.57). Based
on these values of , both populations appeared to have
grown substantially following divergence. The average
estimate for the scaled splitting time was t = 0.37 (90%
HPD = 0.22 - 0.50; Figure 3), suggesting that the NE
and SW groups of populations started to diverge about
5200 years BP (7100 - 3100 years BP considering the
range of mutation rates used). The gene flow estimate
from SW into NE was low (m = 0.66, 90% HPD =
0.0035 - 4.12), and close to null from NE into SW (m =
0.0035, 90% HPD = 0.0035 - 2.88). Conversion of these
values of m resulted in an estimated number of migrants
of approximately 0.004 females per generation from NE
to SW populations (one female every 250 generations),
and 22 females per generation from SW to NE
populations. However, the associated error to these estimates
was large (90% HPDSW = 0.1 - 108; 90% HPDNE = 0.1
130 female migrants per generation).
Genetic diversity of Montagus harrier
Overall, we found little mitochondrial DNA variation
among populations of Montagus harrier throughout the
breeding range. Although variability was low, we
nevertheless found significant differentiation between
southwestern and northeastern populations. In contrast to the
low among-population variability, our study revealed a
relatively high degree of genetic diversity at the mtDNA
for this species as compared with those found in other
sympatric raptor species. These other raptors often
show much less polymorphism, even when using genetic
markers that are more variable than the mitochondrial
control region (examples from vultures in ). Some of
these studies have associated these low levels of genetic
diversity with recent decreases in population sizes, like
for instance in the Spanish imperial eagle Aquila
adalberti , Red kite Milvus milvus , Bonellis eagle
Hieraaetus fasciatus , White-tailed eagle Haliaeetus
albicilla , and Bearded vulture Gypaetus barbatus
. The mtDNA diversity of the Montagus harrier is
not consistent with a recent population decline, unlike
in these other studied raptor species. However, genetic
diversity is rarely homogeneously distributed within a
species range, because it reflects the long-term effects of
historical events . Despite relatively small sample
sizes for some localities, we were able to identify two
strongholds of genetic diversity for the Montagus
harrier, one in the western part of the range (France and
coastal Spain), and another in central Asia. This
coincides with the distribution of open vegetation in the mid
Holocene, which was restricted to the Iberian Peninsula
and coastal areas of France and to central Asia [37-39]
suggesting that these areas may have acted as refuges
for the species at a time when most of Europe was
covered by dense forest . Fossil records confirm the
presence of Montagus harriers in the western region
from the Pleistocene through the Holocene .
However, our lack of sampling locations in East Asia
prevents us from pinpointing primary areas of Montagus
harrier diversity in this region. Additional sampling
would help to establish if the diversity found in
Kazakhstan reflects higher levels of diversity in the Asian
steppes, as suggested by the high diversity found in the
wintering grounds of Asian breeding populations (PAK)
compared to their west European counterparts (SEN).
The haplotype network consisted of two distinct, but
closely related, lineages (Figure 2) that could also be
indicative of two refugia or nuclei from where the
species expanded. The existence in Eurasia of a western
and an eastern or central refuge area has been proposed
for other bird species [31,35,41,42], and the highest
values of genetic diversity observed in France and in
Kazakhstan would be consistent with such scenario.
However, we found no relationship between lineages
and their geographic origin (Figure 2). This could be a
consequence of the low divergence between the two
haplotype groups, suggesting a short time of isolation
for these populations. This would be expected for a
steppe bird whose habitat was potentially reduced
during a short period within the present interglacial period,
rather than during the longer glacial periods, as
postulated for temperate species . During the Pleistocene,
the duration of the steppe-favorable period would have
allowed harriers to fully expand their range, leading to a
high connectivity between different demes. This could
have erased any phylogeographic signal, as is commonly
observed in species characterized by a high dispersal
capability and exploiting a wide ecological niche .
Population structure and phylogeography
In many animal species, the patterns of genetic
differentiation and gene flow are highly influenced by the
geographical characteristics of the places they inhabit and
by their migratory behavior. However, in the Montagus
harrier, neither the geographical barriers (mountain
ranges) nor the intercontinental migratory divide
(populations wintering in Africa or in the Indian
sub-continent) represented significant barriers against gene flow.
Interestingly, the levels of genetic differentiation
between populations were unrelated to the geographical
distances separating these populations. This could be
explained by the long-distance migration behavior of the
species. Western European breeding birds (Spain,
France) winter primarily in west Africa (Senegal,
Mauritania and Mali [30,43]), while central and
eastern-European harriers (the Netherlands, Czech Republic and
Germany in this study) winter further east (Niger,
Nigeria and Chad ). Montagus harriers may perform
winter movements following outbreaks of locusts and
grasshoppers (their main food in winter) .
Additionally, wind conditions may influence the migratory
routes, particularly in spring . Some individuals
might thus migrate through a different route upon their
return (spring) migration, and finally disperse and breed
far away from their natal place. This could explain the
lack of significant genetic structure between populations
a priori assigned to different overwintering areas. This
hypothesis of different migration routes between spring
and autumn is supported by ring recoveries  and
counts of migratory birds in the central Mediterranean
, although satellite telemetry-based studies
(conducted on a limited number of individuals) indicate that
most birds follow the same route in spring and autumn
migrations . While the correlation between genetic
and geographical distances for all populations was
nonsignificant, the AMOVA test gave statistical evidence for
a differentiation between the SW and NE populations
(although only 5% of the genetic variation was explained
by this partition). In other words, location of the
breeding areas may be important when explaining the genetic
structure of populations at a broad scale, whereas
geographical distance between populations is not. These
observations, together with the low number of
significant comparisons between pairs of populations, points
to dispersal as a major factor preventing genetic
differentiation within these two regions of the breeding range
of Montagus harriers. Such a relaxed philopatric
behavior has indeed been described for the species .
Our results strongly support a recent population
expansion as an important cause of the relative homogeneity
across populations. Such expansion is indicated by the
shallow phylogenetic tree and by the star-like haplotype
network (Figure 2). The genetic signature observed in
the two main Montagus harrier groups of populations
(SW vs. NE) is consistent with the occurrence of
postglacial demographic expansions during the second half
of the Holocene, as evidenced by the BSP analyses.
These revealed two events of population growth that
occurred first in the SW (ca. 7500-5500 years BP), and
later in the NE part of the range (ca. 3500 to 1000 years
BP; Figure 3). The effective population size increases
also appeared more pronounced in the NE than in the
SW populations (Figure 3). Furthermore, these two
groups of populations diverged around 5200 years ago,
between the two waves of population growth. This
pattern may be explained by regional differences in the
impact of both climate changes and human activity on
vegetation during the Holocene. Ample evidence
support that during the last glacial maximum
(37,00016,000 years BP), the tundra-steppe vegetation was
widespread from France to the Bering strait (e.g. ).
Therefore, according to the preference of Montagus
harrier for open steppe-like landscapes, the species
would have been widely distributed throughout Eurasia.
Then, during the following interglacial period at the
beginning of the Holocene, about 10,000 years ago, this
vast system disappeared almost completely as a
consequence of the quick expansion of temperate forests (in
Europe) and taigas (in Asia) from their ice-Age refugia
. This probably led to strong range contraction and
population declines in steppe-associated communities.
Unfortunately, based on our data set, coalescent times
go back only 10,000 years, therefore changes in
population size before that date (e.g. after the last glacial
maximum ca.18,000 to 10,000 years BP) cannot be inferred.
Although this particular issue has received little
attention so far (phylogeographic studies of steppe species
are still scarce), our data are consistent with the recent
idea of interglacial refugia, which proposes that, in
addition to the traditional high-latitude refugia of boreal
species, cryptic refugia might have existed in other areas
in the south during the interglacials [4,9]. In fact,
climate reconstructions based on pollen records have
shown that, during the Holocene, the climate was
neither stable nor uniform across Eurasia [50-52].
Therefore, the occurrence of these cryptic refugia, and
consequently the severity of climate change impacts on
species, might have been qualitatively different among
regions of Europe and Asia . During the Holocene,
steppe biomes occurred recurrently [54,55], coinciding
with major dry events leading up to glacial conditions at
different time intervals: ~ 11.000-9.500 years BP, ~
8000-7000 years BP and 4000-3000 years BP [54,56].
More recently, an overall more arid period has been
described in Eurasia during the last 4500 years (e.g.
[39,57]). This, together with the increasing
anthropogenic landscape transformations from 40003000 years
ago (e.g. clearance of forested areas, cultivation, cattle
grazing [49,58,59]), may have provided, either naturally
or artificially, new steppe-like habitats for many species
to colonize [60,61]. This temporal pattern is consistent
with our data and could explain the recent population
growths detected in both groups. The earlier and slower
population growth detected in the SW group might
have been associated with an increase in the extent of
suitable open habitats after a major Holocene climate
change dated around 8000 years BP . In contrast,
the existence of large steppe extensions in the East
together with a lack of evidence supporting the
aforementioned cooling event 8000 years ago in this region
 would explain the lack of synchronous expansions
of Montagus harrier population groups. Additionally,
studies on pollen spectra have clearly indicated that
around 3000 and 1000-500 years BP steppe biomes were
relatively abundant in the eastern part of the range (e.g.
), thus agreeing with the fast population growth
observed in the northeastern populations of Montagus
harrier around that period of time (Figure 3).
Our results point to a short isolation time, relatively
recent population expansions and relaxed philopatry as
the main factors determining the relative genetic
homogeneity observed across populations of a
steppe-associated raptor species. In contrast to the traditional view
[3,63], our findings do not support an important role of
southern Mediterranean peninsulas for extensive
colonization of formerly treeless northern regions. In our case,
rather than a source of postglacial colonization, the
Iberian Peninsula would represent an area of postglacial
refuge for steppe fauna. This finding implies that the
population genetic models of glacial isolation and
postglacial colonization developed for temperate taxa might
have limited applications for steppe species. However,
there is still little evidence for the direct effect of past
climatic events on the genetic variability and
phylogeographic structure in steppe-associated fauna at a
regional or continental scale (but see [9,13]). Our study has
added new insights into the knowledge of how genetic
variation in steppe-associated taxa has been influenced
by late Pleistocene and Holocene climatic changes.
Future research should include a comparative approach,
which would allow the comparison of phylogeographic
patterns in a wider range of co-distributed species. This
would contribute to a better understanding of how
glacial cycles have sculpted the genetic variation of
steppeassociated taxa in Eurasia.
We analyzed genetic material from 284 Montagus
harrier specimens collected in 13 localities across the
species breeding range, from Spain to Kazakhstan, and
from two wintering areas (Senegal, and Pakistan) (Table
1). Samples were grouped a priori according to sampling
locality, and these groups were considered as
populations for genetic analyses (see Table 1 and Figure 1). All
samples are contemporary (collected in 1999-2009) and
consisted of blood (n = 204) or feathers (n = 80). When
nestlings were used as a source for DNA (<10% of
samples) we used only one chick per brood to avoid
pseudo-replication of mitochondrial haplotypes.
DNA isolation, polymerase chain reaction (PCR) and
Blood samples were digested (8 h) in 250 L SET buffer in
the presence of SDS (2%) and proteinase K (10 ng/L).
Feathers were processed like blood samples but increasing
proteinase K (20 ng/L) and time of digestion (16 h). Total
genomic DNA was extracted using standard NH4Ac
protocol. Purified DNA was diluted to a working concentration
of 25 ng/L. We amplified two mitochondrial regions
including partial tRNA-Trp and NADH dehydrogenase
subunit 2 (ND2) and a portion of the cytochrome oxidase
subunit I (COI) via polymerase chain reaction (PCR).
Primers L5216-H5766 and L5758-H6313 were used for the
amplification of ND2 gene , and BirdF1 and BirdR1
 for the COI fragment. PCR reactions were run using
the following parameters: denaturation at 95 C for 3 min,
followed by 35 cycles of 94 C for 60 s, 54 C for 60 s, and
72 C for 60 s, and a final extension at 72C for 5 min.
PCRs contained approximately 25 ng of template DNA, 1
PCR buffer (Biotools), 0.25 mM of each dNTP, 0.3 M of
each primer, 2 mM MgCl2, and 0.5 U of Taq DNA
polymerase (Biotools) in a total volume of 10 L.
PCR-products were purified with Exonuclease I and
Shrimp Alkaline Phosphatase enzymatic reactions (United
States Biochemical). Purified reactions were sequenced in
an ABI 3130 automated sequencer (Applied Biosystems)
using dye-terminator chemistry (BigDye kit 3.1, Applied
Biosystems) with the same primers used for PCR. All
sequences are accessible at GenBank (Additional File 1).
Genetic diversity and population structure analyses
We edited and aligned DNA sequences using Bioedit
 and Clustal W . Arlequin 184.108.40.206  was
used to determine the number of haplotypes and
variable sites, and to calculate genetic diversity in each
population (at the haplotype and nucleotide levels).
Genetic differentiation between populations was tested
using pairwise FST comparisons for each mtDNA
region and for the concatenated data set, using a
Tamura-Nei evolutionary model , as this is the
closest model to the one inferred for our data set in
Modeltest 3.7 . The significance of pairwise FST
comparisons was given by a P value calculated using
10,000 random permutation tests; p-values were
further adjusted according to sequential Bonferroni
corrections for multiple tests . Evidence for
population genetic structure was assessed using an analysis
of molecular variance (AMOVA) as implemented in
Arlequin 220.127.116.11. Tamura-Nei distances plus gamma
correction (a = 0.2487) were selected for the
concatenated data set. We first examined overall differences
among populations (i.e. sampling localities without
grouping) and then we determined the contributions
of different grouping scenarios to the partitioning of
genetic variation in the dataset. For this purpose, we
tested three hypothetical scenarios (Table 2): (1)
differentiation explained by the species current geographic
distribution. We compared south-western (SW) vs.
north-eastern populations (NE), which represent
currently samples from a large and continuous breeding
area (SW) and smaller populations with a more overall
patchy distribution, both separated by a gap in the
overall breeding distribution range (; see Figure 1);
(2) differentiation between populations separated by
geographic barriers (i.e. mountain ranges: Pyrenees/
Alps/Urals) as potential barriers to gene flow. For this,
we compared Spanish, Western-central European,
Eastern European, and Asian populations; (3)
differentiation in relation to the intercontinental migratory divide
(i.e. differences between populations wintering in
Africa vs. those wintering in the Indian subcontinent).
Since all Montagus harrier overwintering in the
Indian subcontinent come from Asian breeding
populations , winter samples from the Indian sub-continent
(Pakistan) were included into the NE population group
(in scenario 1) or into the Asian population group (in
scenarios 2 and 3). Likewise, in scenario 3, samples
from Senegal were pooled into the group of Montagus
harrier populations wintering in Africa. In all cases, we
repeated these analyses without including samples from
the two overwinter sites.
A pattern of isolation-by-distance was explicitly tested
using Mantel tests to compare pairwise geographic and
genetic distances between populations. These were
statistically tested using linearized pairwise differentiation
indexes (FST /(1-FST) in Arlequin 18.104.22.168. The statistical
significance of correlations between distance matrices
Phylogenetic relationships of Montagus harrier were
reconstructed by examining mtDNA sequence variation
in all samples. The best-fit evolutionary model for the
concatenated data set was determined using the Akaike
information criterion implemented in Modeltest 3.7. A
Maximum-Likelihood (ML) tree was built using a
heuristic search starting from a neighbour-joining tree and a
tree bisection reconnection (TBR) algorithm for branch
swapping, with random addition of sequences in PAUP*
4.0 . The statistical support for internal branches of
the tree was estimated by 1,000 bootstrap replicates.
This model was also used to carry out Bayesian
inference (BI) of phylogeny as implemented in MrBayes
v3.1.2 , simulating four simultaneous Monte Carlo
Markov Chains (MCMC) for 5 106 generations each.
The first 250,000 generations were discarded as burn-in.
Bayesian posterior probabilities were obtained to assess
the robustness of the BI trees. Trees were rooted with
one sequence of C. macrourus, which was used as
outgroup (Additional File 1).
We also represented the genealogical relationships of
all the analyzed samples of C. pygargus with a haplotype
network calculated using the median-joining algorithm
 in NETWORK 4.5 http://www.fluxus-engineering.
Signatures of demographic changes or selection in the
recent history of C. pygargus, considering a model of
mutation-drift equilibrium, were addressed using
analyses based on different coalescent approaches. Firstly,
we tested the data against a neutral Wright-Fisher
model using Fus Fs  and Tajimas D , which
aim to identify an excess of recent single nucleotide
substitution caused by population growth, bottleneck, or
background selection. We performed this test in
Arlequin 22.214.171.124. Significance of the statistics was determined
by 1000 coalescent simulations of the neutral model,
where P must be less than 0.02 to be significant due to
the non-normal distribution of the Fs statistic .
Secondly, because departures from neutrality are often
caused by changes in effective population size, we
generated a Bayesian Skyline Plot (BSP) to explore changes in
genetic diversity occurring at a certain time period
within a given genealogy using MCMC based sampling,
and to generate the posterior distribution of the effective
population size at that time (Ne) . We used a strict
molecular clock and a range of substitution rates
estimated for other bird species (0.02-0.055 substitutions/
site/Myr; [13,80,81]). Four independent analyses were
performed in BEAST v1.5.4  and ran for 4 107
generations with a sampling frequency of 1000 steps.
Convergence was assessed using Tracer v1.5  and
uncertainty in parameter estimates reflected in values of
the 95% highest posterior density (HPD). We also used
BEAST v1.5.4 to estimate the time to the most recent
common ancestor (TMRCA) for each group of
sequences analyzed as well as for the complete data set.
Thirdly, because AMOVA tests revealed genetic
breaks between geographical regions (SW vs. NE,
scenario 1), we used the program IMa  to test the
hypothesis of a shared-history scenario of isolation
with migration for the SW and NE population groups.
This model assumes that an ancestral population of
constant size and population parameter A separated
into two populations (SW and NE) at time T, to
simultaneously determine (1) time since divergence (t), (2)
effective population sizes of each population (SW and
NE) and the ancestral population (A) at time of split,
and (3) immigration rates (mSW and mNE). We ran
three replicate runs with a random seed to initiate each
run. In all analyses, we used at least 20
Markovcoupled chains with a geometric heating scheme, a
burn-in of 200,000 steps, and run until the effective
sample sizes (ESS; see ) for each parameter were at
least 500. To ensure proper chain mixing and
parameter convergence, all parameter trend lines were
visually inspected and three independent runs, which
differed only in starting random seed, were compared.
To convert IMa parameter estimates to biologically
meaningful values, the parameters were scaled to a
substitution rate of = 4 10-8 substitutions per site
per year (s/s/y). We also used a lower ( = 2 10-8 s/
s/y) and an upper ( = 5.5 10-8 s/s/y) limit for this
conversion, as indicated for the BSP analysis. To
estimate generation time (g) we used the equation g = a +
(s/(1- s)), where a is the age of first reproduction in
females and s is the expected adult survival rate .
Age at which Montagus harrier females reach maturity
was set at 2 years and estimated adult survival rate at
0.67, according to published data [24,85]. Therefore,
we considered a generation time (g) of 4 years to
express the output parameters (SW, NW, A, mSW,
mNW, t) in demographic units as follows: effective
population size: N = /4g; number of migrants per
generation: M = m/2g; divergence time in years: T
Two anonymous reviewers contributed to improving the MS. Samples were
kindly provided by R. Limiana, X. Vzquez, C. Alonso, A. Guerrero, and E. de
Prada (Spain), C. Prckhauer (Senegal), C. Trierweiller (The Netherlands), I.
Kunstmuller (Czech Republic), A. Baqri (Pakistan), and B. Van Hecke, T.
Printemps (PLO, France), F. Arrias, P. Maigre (PACA, France), S. Paris, JL.
Bourrioux, F. Burda (CHAR, France), MF Canevet (AQPY, France), S. Augiron,
JFB Blanc (PCH, France). We also thank the staff of the collection of DNA
and Tissues from the Museo Nacional de Ciencias Naturales (CSIC). Sally
Bach reviewed the English text. For research permits we thank the Ministerio
de Medio Ambiente, Rural y Marino. This work was supported by JCCM
project (Ref: PAC06-0137) and CSIC-MICINN project (Ref: PIE 201030I019), and
by Natural Research (fieldwork in Kazakhstan). We performed all the
laboratory work, including sequencing, at the Genetics Laboratory of IREC
JTG conceived and designed the study, collected samples, helped in
molecular genetic work, performed analysis and drafted the manuscript. FA
participated in the study design, carried out molecular genetic work,
performed analysis and helped to draft the manuscript. FM, JT, AS, VB
participated in the study design, collected samples and helped to draft the
manuscript. BA conceived and designed the study, collected samples and
helped to draft the manuscript. All authors read and approved the final
1. Avise JC : Phylogeography: the history and formation of species Harvard University Press; 2000 .
2. Hewitt G : The genetic legacy of the Quaternary ice ages . Nature 2000 , 405 ( 6789 ): 907 - 913 .
3. Hewitt GM : Some genetic consequences of ice ages, and their role, in divergence and speciation . Biol J Linn Soc 1996 , 58 ( 3 ): 247 - 276 .
4. Stewart JR , Lister AM , Barnes I , Daln L : Refugia revisited: individualistic responses of species in space and time . Proc R Soc Lond B 2010 , 277 ( 1682 ): 661 .
5. Hewitt GM : Post glacial re-colonization of European biota . Biol J Linn Soc 1999 , 68 ( 1 2): 87 - 112 .
6. Fedorov V , Goropashnaya A , Jaarola MM , Cook A : Phylogeography of lemmings (Lemmus): no evidence for postglacial colonization of Arctic from the Beringian refugium . Mol Ecol 2003 , 12 : 725 - 731 .
7. Hamill R , Doyle D , Duke E : Spatial patterns of genetic diversity across European subspecies of the mountain hare , Lepus timidus L. Heredity 2006 , 97 ( 5 ): 355 - 365 .
8. Hewitt G : Genetic consequences of climatic oscillations in the Quaternary . Philos Trans Roy Soc Lond B 2004 , 359 ( 1442 ): 183 .
9. Garca JT , Maosa S , Morales MB , Ponjoan A , Garcia de la Morena EL , Bota G , Bretagnolle V , Dvila JA : Genetic consequences of interglacial isolation in a steppe bird . Mol Phylogenet Evol 2011 , 61 : 671 - 676 .
10. Li Y-X , Yun-Xianga Z , Xiang-Xu X : The composition of three mammal faunas and environmental evolution in the last glacial maximum , Guanzhong area, Shaanxi Province, China. Quatern Int 2011 .
11. Snchez-Marco A : Implications of the avian fauna for paleoecology in the Early Pleistocene of the Iberian Peninsula . J Hum Evol 1999 , 37 ( 3- 4 ): 375 - 388 .
12. Snchez-Marco A : Avian zoogeographical patterns during the Quaternary in the Mediterranean region and paleoclimatic interpretation . Ardeola 2004 , 51 ( 1 ): 91 - 132 .
13. Garca JT , Surez F , Garza V , Calero-Riestra M , Hernndez J , Prez-Tris J : Genetic and phenotypic variation among geographically isolated populations of the globally threatened Dupont's lark Chersophilus duponti . Mol Phylogenet Evol 2008 , 46 ( 1 ): 237 - 251 .
14. Warmuth V , Eriksson A , Bower M , Caon J , Cothran G , Distl O , GlowatzkiMullis M-L , Hunt H , Lus C , Oom M , et al: European Domestic Horses Originated in Two Holocene Refugia . Plos One 2011 , 6 ( 3 ): e18194 .
15. Melo-Ferreira J , Boursot P , Randi E , Kryukov A , Suchentrunk F , Ferrand N , Alves C : The rise and fall of the mountain hare (Lepus timidus) during Pleistocene glaciations: expansion and retreat with hybridization in the Iberian Peninsula . Mol Ecol 2007 , 16 : 605 - 618 .
16. Moreno V , Morales M , Traba J : Avoiding over-implementation of agrienvironmental schemes for steppe bird conservation: A species-focused proposal based on expert criteria . J Environ Manage 2010 , 91 ( 8 ): 1802 - 1809 .
17. Fuller R , Gregory R , Gibbons D , Marchant J , Wilson J , Baillie S , Carter N : Population declines and range contractions among lowland farmland birds in Britain . Conserv Biol 1995 , 9 ( 6 ): 1425 - 1441 .
18. Chamberlain D , Fuller R , Bunce R , Duckworth J , Shrubb M : Changes in the abundance of farmland birds in relation to the timing of agricultural intensification in England and Wales . J Appl Ecol 2000 , 37 ( 5 ): 771 - 788 .
19. Siriwardena GM , Crick HQP , Baillic SR , Wilson JD : Agricultural land use and the spatial distribution of granivorous lowland farmland birds . Ecography 2000 , 23 ( 6 ): 702 - 719 .
20. Donald PF , Green RE , Heath MF : Agricultural intensification and the collapse of Europe's farmland bird populations . Proc R Soc Lond B 2001 , 268 ( 1462 ): 25 - 29 .
21. Ormerod S , Watkinson A: Editors' introduction: birds and agriculture . J Appl Ecol 2000 , 37 ( 5 ): 699 - 705 .
22. Sderstrm B , Prt T : Influence of Landscape Scale on Farmland Birds Breeding in Semi Natural Pastures . Conserv Biol 2000 , 14 ( 2 ): 522 - 533 .
23. Clarke R : Montagu's harrier Arlequin Press; 1996 .
24. Arroyo B , Garca JT , Bretagnolle V : Conservation of the Montagu's harrier (Circus pygargus) in agricultural areas . Anim Cons 2002 , 5 ( 04 ): 283 - 290 .
25. Arroyo BE , Garca JT , Bretagnolle V : Circus pygargus Montagu's harrier . BWP Update 2004 , 6 ( 1 ).
26. Ferguson-Lees J , Christie DA , Franklin K , Mead D , Burton P : Raptors of the world . CSIRO ; 2001 .
27. Moreau RE : The Palaearctic-African bird migration systems . Academic Press; 1972 .
28. Garcia JT , Arroyo B : Migratory movements of western European Montagu's Harrier Circus pygargus: a review . Bird Study 1998 , 45 : 188 - 194 .
29. Trierweiler C , Koks B , Bairlein F , K-M E , Komdeur J , Dijkstra C : Migratory routes and wintering behavior of NW European Montagu's Harriers revealed by satellite telemetry . J Ornithol 2006 , 147 ( suppl 1 ): 265 .
30. Limiana R , Soutullo V , Urios V , Reig-Ferrer A : Migration and wintering areas of adult Montagus harriers (Circus pygargus) breeding in Spain . J Ornithol 2011 .
31. Poulakakis N , Antoniou A , Mantziou G , Parmakelis A , Skartsi T , Vasilakis D , Elorriaga J , de la Puente J , Gavashelishvili A , Ghasabyan M , et al: Population structure, diversity, and phylogeography in the nearthreatened Eurasian black vultures Aegypius monachus (Falconiformes; Accipitridae) in Europe: insights from microsatellite and mitochondrial DNA. Biol J Linn Soc 2008 , 95 : 859 - 872 .
32. Martnez-Cruz J , Godoy J , Negro J : Population genetics after fragmentation: the case of the endangered Spanish imperial eagle Aquila adalberti . Mol Ecol 2004 , 13 : 2243 - 2255 .
33. Roques S , Negro J : MtDNA genetic diversity and population history of a dwindling raptorial bird, the red kite (Milvus milvus) . Biol Cons 2005 , 126 : 41 - 50 .
34. Cadaha L , Urios V , Negro J : Low mitochondrial DNA diversity in the endangered Bonelli's eagle (Hieraaetus fasciatus) from SW Europe (Iberia) and NW Africa. J Ornithol 2007 , 148 ( 1 ): 99 - 104 .
35. Hailer F , Helander B , Folkestad A , Ganusevich S , Garstad S , Hauff P , Koren C , Masterov V , Nygrd T , Rudnick J , et al: Phylogeography of the white-tailed eagle, a generalist with large dispersal capacity . J Biogeogr 2007 , 34 : 1193 - 1206 .
36. Godoy J , Negro J , Hiraldo F , Donzar J : Phylogeography, genetic structure and diversity in the endangered bearded vulture (Gypaetus barbatus L.) as revealed by mitochondrial DNA . Mol Ecol 2004 , 13 : 371 - 390 .
37. Pantalon-Cano J , Errikarta-Imanol Y , Prez-Obiol R , Roure J : Palynological evidence for vegetational history in semi-arid areas of the western Mediterranean (Almeria, Spain) . The Holocene 2003 , 13 ( 1 ): 109 - 119 .
38. Gallimore R , Jacob R , Kutzbach J : Coupled atmosphere-ocean-vegetation simulations for modern and mid-Holocene climates: role of extratropical vegetation cover feedbacks . Clim Dynam 2005 , 25 : 755 - 776 .
39. Jalut G , Esteban Amat A , Bonnet L , Gauquelin T , Fontugne M : Holocene climatic changes in the western Mediterranean, from south-east France to south-east Spain . Palaeogeogr Palaeocl 2000 , 160 : 255 - 290 .
40. Mitchell F : How open were European primeval forests? Hypothesis testing using palaeoecological data . J Ecol 2005 , 93 : 168 - 177 .
41. Tiedemann R , Paulus K , Scheer M , al e: Mitochondrial DNA and microsatellite variation in the eider duck (Somateria mollissima) indicate stepwise postglacial colonization of Europe and limited current longdistance dispersal . Mol Ecol 2004 , 13 : 1481 - 1494 .
42. Liebers D , De Knijff P , Helbig A : The herring gull complex is not a ring species . Proc R Soc Lond B 2004 , 271 : 893 - 901 .
43. Limiana R , Soutullo A , Urios V : Autumn migration of Montagu's Harriers tracked by satellite telemetry . J Ornithol 2007 , 148 : 517 - 523 .
44. Klaassen R , Hake M , Strandberg R , Alerstam T : Geographical and temporal flexibility in the response to crosswinds by migrating raptors . Proc R Soc Lond B 2011 , 278 : 1339 - 1346 .
45. Panuccio M , Agostini N : Spring migration of Pallid (Circus macrourus) and Montagu's harriers (Circus pygargus) in relation to sex and age classes at two watchsites of the central mediterranean . Buteo 2006 , 15 : 3 - 10 .
46. Trierweiler C , Koks B : Montagu's Harrier Circus pygargus . In Living on the edge: Wetlands and birds in a changing Sahel . Edited by: Zwarts L, Bijlsma R , van der Kamp J , Wymenga E. Zeist: KNNV Publishing ; 2009 : 312 - 327 .
47. Limiana R , Garca JT , Guerrero A , Gonzlez JM , Lavedn J , Moreno JD , Muoz A-R , Palomares LE , Pinilla A , Ros G , et al: Phylopatry and natal dispersal of Montagu's harriers Circus pygargus breeding in Spain: a review of existing data . Eur J Wildl Res 2011 .
48. Peyron O , Guiot J , Cheddadi R , Tarasov P , Reille R , de Beaulieu J-L , Bottema S , Andrieu V : Climatic reconstruction in Europe for 18,000 yr BP from pollen data . Quaternary Res 1998 , 49 : 183 - 196 .
49. Novenko E , Volkova E , Nosova N , Zuganova I : Late Glacial and Holocene landscape dynamics in the southern taiga zone of East European Plain according to pollen and macrofossil records from the Central Forest State Reserve (Valdai Hills, Russia) . Quatern Int 2009 , 207 : 93 - 103 .
50. Klimanov V : Paleoclimatic reconstructions based on the information statistical method . In Late Quaternary Environments of the Soviet Union. Edited by: Velichko A. Minneapolis: University of Minnesota Press ; 1984 : 297 - 304 .
51. Savina S , Khotinsky N : Holocence paleoclimatic reconstructions based on the zonal method . In Late Quaternary Environments of Soviet Union. Edited by: Velichko A. Minneapolis: University of Minnesota Press ; 1984 : 287 - 296 .
52. Zubakov V , Borzenkova I : Global Paleoclimate of Late Cenozoic Amsterdam: Elsevier ; 1990 .
53. Mayewski P , Rohling E , Stager J , Karln W , Maasch K , L D , Meyerson E , Gasse F , van Kreveld S , Holmgren K , et al: Holocene climate variability . Quaternary Res 2004 , 62 : 243 - 255 .
54. Tarasov P , Bezrukova E , Karabanov E , Nakagawa T , Wagner M , Kulagina N , Letunova P , Abzaeva A , Granoszewski W , Riedel F : Vegetation and climate dynamics during the Holocene and Eemian interglacials derived from Lake Baikal pollen records . Palaeogeogr Palaeocl 2007 , 252 : 440 - 457 .
55. Sorrel P , Popescu S-M , Klotz S , Suc J-P , Oberhnsli H : Climate variability in the Aral Sea basin (Central Asia) during the late Holocene based on vegetation changes . Quaternary Res 2007 , 67 : 357 - 370 .
56. Gasse F , van Campo E : Abrupt post-glacial climatic events in West Asia and North Africa mlonsoon domlains . Earth Planet Sc Lett 1994 , 126 : 435 - 456 .
57. Prez-Obiol R , Jalut G , Juli R , Plachs A , Iriarte M , Otto T , HernndezBeloqui B : Mid-Holocene vegetation and climatic history of the Iberian Peninsula . The Holocene 2011 , 21 : 75 - 93 .
58. Zohary D , Hopf M : Domestication of plants in the old world Oxford : Oxford University Press ; 1993 .
59. Ruddiman W : The anthropogenic greenhouse era began thousands of years ago . Climate Change 2003 , 61 : 261 - 293 .
60. O'Connor R , Shrubb M : Farming and Birds Cambridge : Cambridge University Press ; 1986 .
61. Bouma J , Varralyay G , Batjes N : Principal land use changes anticipated in Europe . Agriculture, Ecosystems and Environment 1998 , 67 : 103 - 119 .
62. Alley R , gstsdttir A : The 8 k event: cause and consequences of a major Holocene abrupt climate change . Quaternary Sci Rev 2005 , 24 : 1123 - 1149 .
63. Taberlet P , Fumagalli L , Wust-Saucy A , Cosson J : Comparative phylogeography and postglacial colonization routes in Europe . Mol Ecol 1998 , 7 : 453 - 464 .
64. Sorenson M , Ast J , Dimcheff D , Yuri T , Mindell D : Primers for a PCR-based approach to mitochondrial genome sequencing in birds and other vertebrates . Mol Phylogenet Evol 1999 , 12 ( 2 ): 105 - 114 .
65. Hebert PDN , Stoeckle MY , Zemlak TS , Francis CM : Identification of birds through DNA barcodes . PLoS Biology 2004 , 2 ( 10 ): e312 .
66. Hall T : BioEdit: a user-friendly biological sequence alignment editor and analysis program for Windows 95 /98/NT. 1999 , 1999 : 95 - 98 .
67. Thompson JD , Higgins DG , Gibson TJ : CLUSTAL W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice . Nucl Acids Res 1994 , 22 ( 22 ): 4673 - 4680 .
68. Excoffier L , Laval G , Schneider S : Arlequin ( version 3. 0 ) : an integrated software package for population genetics data analysis . Evol Bioinf Online 2005 , 1 : 47 .
69. Tamura K , Nei M : Estimation of the number of nucleotide substitutions in the control region of mitochondrial DNA in humans and chimpanzees . Mol Biol Evol 1993 , 10 ( 3 ): 512 .
70. Posada D , Crandall K : MODELTEST: testing the model of DNA substitution . Bioinformatics 1998 , 14 ( 9 ): 817 .
71. Holm S : A simple sequentially rejective multiple test procedure . Scand J Stat 1979 , 6 ( 2 ): 65 - 70 .
72. Del Hoyo J , Elliott A , Sargatal J , Cabot J : Handbook of the birds of the world: New world vultures to guineafowl . Barcelona: Lynx Editions; 1994 .
73. Swofford D : PAUP*. Phylogenetic analysis using parsimony (* and other methods) . Version 4. Sinauer Associates , Sunderland, Massachusetts; 2002 .
74. Ronquist F , Huelsenbeck JP : MrBayes 3: Bayesian phylogenetic inference under mixed models . Bioinformatics 2003 , 19 ( 12 ): 1572 - 1574 .
75. Bandelt HJ , Forster P , Rhl A : Median-joining networks for inferring intraspecific phylogenies . Mol Biol Evol 1999 , 16 ( 1 ): 37 .
76. Fu YX : Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection . Genetics 1997 , 147 ( 2 ): 915 .
77. Tajima F : Statistical method for testing the neutral mutation hypothesis by DNA polymorphism . Genetics 1989 , 123 ( 3 ): 585 .
78. Drummond A , Rambaut A : BEAST v1 .4. 2006 [http://beast.bio.ed.ac.uk/].
79. Arbogast B , Drovetski S , Curry R , Boag P , Seutin G , Grant P , Grant R , Anderson D : The origing and diversification of Galapagos Mockingbirds . Evolution 2006 , 60 ( 2 ): 370 - 382 .
80. Peck D , Congdon B : Reconciling historical processes, population structure in the sooty tern Sterna fuscata . J Avian Biol 2004 , 35 : 327 - 335 .
81. Rambaut A , Drummond A : Tracer v1.4. 2007 [http://beast.bio.ed.ac.uk/ Tracer].
82. Hey J , Nielsen R : Multilocus Methods for Estimating Population Sizes, Migration Rates and Divergence Time , With Applications to the Divergence of Drosophila pseudoobscura and D. persimilis. Genetics 2004 , 167 : 747 - 760 .
83. Saether B , Engen S , Mller A , Visser M , Matthysen E , Fiedler W , Lambrechts M , Becker P , Brommer J , Dickinson J : Time to extinction of bird populations . Ecology 2005 , 86 ( 3 ): 693 - 700 .
84. Millon A , Bretagnolle V : Predator population dynamics under a cyclic prey regime: numerical responses, demographic parameters and growth rates . Oikos 2008 , 117 : 1500 - 1510 .
85. del Hoyo J , Elliott A , Sargatal J : Handbook of the birds of the world: New world vultures to guineafowl . Barcelona: Lynx Editions; 19942 .