BACKGROUND
interspecific gene flow is an important evolutionary force.
it may enrich genetic variation of individual species, facilitate the origin of new phenotypes, lead to the origin of new species, or, conversely, to species fusion  <cit> .
it has been suggested that high speciation rates and rapid phenotypic changes observed during adaptive radiations may be facilitated by increased genetic variation and novel genetic combinations produced by introgressive hybridization among the species  <cit> .
patterns of interspecific gene flow can also provide an important insight into the genetic basis of reproductive isolation  <cit> .
in species with incomplete reproductive isolation, gene flow occurs at loci with alleles that are neutral or beneficial on the genomic background of the other species, but is limited at loci harboring alleles involved in reproductive isolation  <cit> .
most work exploring patterns and rates of interspecific gene flow focused on pairs of closely related sister species  <cit>  because, until recently  <cit> , models allowing the analysis of interspecific gene flow were limited to only two species.
several general findings emerged from these studies.
first, in the case of speciation driven by intrinsic barriers, interspecific gene flow is often reduced on the sex chromosomes   and in genomic regions with low recombination rates  <cit> , indicating that genetic incompatibilities tend to preferentially accumulate in these genomic regions.
second, in the case of speciation driven by divergent ecological selection, gene flow is correlated with ecological divergence  <cit>  and is reduced around genes underlying ecological differences  <cit> .
interspecific gene flow across larger sections of the phylogenetic tree is, however, much less studied.
such studies would be particularly interesting in instances of rapid radiations, where lineage divergence occurs on a shorter timescale than the completion of reproductive isolation .
knowledge of patterns and rates of interspecific gene flow in such systems allows to address questions, which are difficult or even impossible to study using models of only two sister species.
one can, for example, examine whether the same genes contribute to reproductive isolation in different speciation events, assess the importance of various reproductive barriers in restricting the gene flow and determine the order in which different reproductive barriers arise during lineage divergence.
this information can help us to understand why lineage divergence occurs so quickly in some taxonomic groups, while not in others.
acrocephalus warblers are small passerines that underwent an adaptive radiation during the last approximately  <dig> million years  <cit> , which gave rise to  <dig> extant species occupying mainly eurasia, africa, australia and pacific islands  <cit> .
the members of this genus are phenotypically quite uniform.
they, however, display a great diversity in ecology and behaviour, and have been used extensively as models in ecological and behavioral research  <cit> .
interspecific hybridization occasionally occurs in acrocephalus warblers.
it has been reported not only between sister species , but also between more distantly related species  .
in some of these cases, hybridization occurs between species of very different body size  or even between species belonging to different subgenera .
so far very little is known about the degree of postzygotic reproductive isolation between the species.
in some cases, e.g.
between sister species a. arundinaceus and a. stentoreus, postzygotic isolation seems to be quite strong as no backcross hybrids have been observed in a hybridizing population  <cit> .
in other cases , however, observation of fertile f <dig> hybrids  <cit>  suggests that postzygotic isolation is incomplete and interspecific gene flow might occur between the species.
here we performed the first study exploring patterns and rates of interspecific gene flow in this avian group, focusing on three species of the subgenus notiocichla: a. scirpaceus, a. palustris and a. dumetorum.
all three species are very similar in morphology, but have different habitat requirements, migration strategies and song patterns  <cit> .
a. scirpaceus breeds across europe, northern africa and western asia typically in reed beds and winters in sub-saharan africa.
its sister species, a. palustris, has similar geographical distribution, but occupies mostly damp herbaceous vegetation and winters in southeastern africa.
males of each species sing different songs, nevertheless, in a. palustris, mixed singers that incorporate song phrases from a. scirpaceus to its repertoire are known  <cit> .
both species occasionally form mixed pairs  <cit>  and are able to produce viable hybrid offspring .
a. dumetorum, the sister species to the a. scirpaceus and a. palustris clade  <cit> , breeds in asia and northeastern europe in various herbaceous or bushy vegetation and winters in india.
in eastern europe and western asia, it often co-occurs with a. palustris at the same sites and hybridization between the two species has been observed .
the exact frequency of hybridization between a. palustris and a. scirpaceus and between a. palustris and a. dumetorum is hard to estimate owing to morphological similarity of the species and difficulty to distinguish hybrids without genetic analysis.
the total number of observed mixed pairs is, however, comparable in both species pairs .
the cases of hybridization between a. dumetorum and a. scirpaceus, which are ecologically more divergent, are not known, although one putative hybrid individual has been reported  <cit> .
to get the first insight into mechanisms and selective forces driving acrocephalus radiation, we analyzed nucleotide variation at four autosomal and four loci on the z chromosome in sympatric populations of a. scirpaceus, a. palustris and a. dumetorum.
using an isolation-with-migration model for multiple populations we tested whether interspecific gene flow occurred among the studied species and if so, whether it was higher between more genetically similar or between more ecologically similar species.
in addition, using several tests of neutrality, we explored patterns of positive selection in all three species.
methods
study area and samples
we analyzed samples from  <dig> individuals of a. scirpaceus,  <dig> individuals of a. palustris and  <dig> individuals of a. dumetorum.
all birds were caught in  <dig> during the breeding season  in southern finland, where the breeding ranges of all three species overlap .
coordinates of the sites, where individual birds were captured, are shown in additional file  <dig>  from each individual, a blood sample was collected by brachial vein puncture and stored in pure ethanol for further extraction of genomic dna.
in addition, a sample from one individual of a. schoenobaenus was used as the outgroup.fig.
<dig> breeding ranges of the three studied acrocephalus warblers, a. scirpaceus, a. palustris and a. dumetorum.
sampling sites are indicated by dark grey circles
molecular sex determination and sequencing
genomic dna was purified by dneasy tissue kit  according to manufacturer’s instructions.
sex of individual birds was determined following the method of griffiths et al.
<cit>  using p <dig> and p <dig> primers.
these primers amplify a fragment of two homologous genes located on the z and w chromosomes.
the amplified fragments from the z and w chromosomes differ in length, which can be visualized on agarose gel.
we further amplified and sequenced introns of four autosomal and four z-linked loci.
the four autosomal loci are located on different chromosomes in the zebra finch  genome.
of the four z-linked genes, ppwd <dig> and adamts <dig> are close physically , and therefore may not have independent evolutionary histories.
however, the other two loci on the z chromosome are separated from each other and from ppwd <dig> and adamts <dig> by > 19 mb.
primers for pcr amplification were partially obtained from published studies and partially designed by ourselves.
all primers were designed in conserved exonic regions of the zebra finch or chicken  genome in order to amplify intronic sequences.
primer sequences and lengths of the obtained pcr products are provided in additional file  <dig>  pcr conditions are same as in storchová et al.
<cit> .
all pcr products were sequenced in both directions with the pcr primers using sanger sequencing at the dna sequencing laboratory at charles university in prague.
data analyses
sequences were manually edited using codoncode aligner software .
alignments were generated by clustalw as implemented in the program bioedit  <cit> .
all alignments were visually checked and manually adjusted.
exonic sequences as well as indel polymorphisms were excluded from the analyses.
individuals and/or positions with missing data  were eliminated from the dataset to get sequences of the same length for each locus.
diploid sequences at each locus were further separated into two haplotypes using the program phase, version  <dig> .
<dig>  <cit>  with the following parameters: number of iterations =  <dig> , thinning interval =  <dig>  burnin =  <dig>  we used the default recombination model, which is the general model for varying recombination rate.
for each dataset, we applied the algorithm five times with different random seeds, and we checked for consistency of results across independent runs.
we obtained identical haplotypes across all runs for all loci.
the final number of haploid sequences for each locus and the length of sequences are shown in table  <dig> table  <dig> polymorphism statistics for eight studied loci
17483
a.s.
a.p.
a.d.
21281
a.s.
a.p.
a.d.
24972
a.s.
a.p.
a.d.
rpl5-4
a.s.
a.p.
a.d.
adamts6
a.s.
a.p.
a.d.
ppwd1
a.s.
a.p.
a.d.
tg401
a.s.
a.p.
a.d.
tg1505
a.s.
a.p.
a.d.
aautosome , chromosome z
blength of sequence
c
a. scirpaceus , a. palustris , a. dumetorum
dnumber of haploid sequences
enumber of segregating sites
faverage number of nucleotide differences
gproportion of polymorphic sites
htajima’s d
ifu and li’s d
jdivergence to outgroup measured as average pairwise divergence, dxy
basic population genetic analyses of polymorphism, divergence, recombination, and tests of neutrality based on the allele frequency spectrum were performed with the program dnasp  <cit> .
hudson-kreitman-aguadé  tests of positive selection  <cit>  were performed using the hka program .
this program compares the ratio of polymorphism to divergence at multiple loci.
loci under positive selection are expected to show lower levels of polymorphism relative to divergence.
for each locus, we set the information whether the locus is autosomal or z-linked, based on which the hka program accounts for lower effective population size and thus the lower nucleotide variation on the z chromosome.
genealogical relationships among haplotypes for each locus were reconstructed with network software  <cit>  using median joining algorithm.
the data were fit to the isolation-with-migration model  for multiple populations  <cit>  using the program ima <dig>  the program estimates several demographic parameters based on markov chain monte carlo simulations of genealogies.
these parameters include the effective population size for the current and ancestral populations, migration rates between the populations and the population-split times.
because ima <dig> assumes no recombination within loci, we determined the longest region without observed recombination for each locus using the program imgc  <cit> .
this program removes either sites or haplotypes to produce the most information-rich contiguous dna sequence segment that passes the four-gamete test.
nonrecombinant regions represented 84 % of the length of each locus, on average, and were used as an input for the ima <dig> program.
for each locus, we provided the information whether it is autosomal or z-linked, based on which ima <dig> program accounts for different effective population sizes between autosomes and the z chromosome.
we ran the program three times with identical starting conditions, with the exception of the random number seed, to assess convergence.
to facilitate mixing of the markov chains, we used metropolis coupling with  <dig> chains and a geometric heating model.
upper bounds for the prior distributions of parameter values were set as suggested in ima <dig> documentation .
all runs began with a burn-in period of  <dig>  steps and were allowed to continue for 7– <dig> million steps.
we were able to achieve adequate mixing of the markov chains as indicated by trend line plots, effective sample size  values , and very similar parameter estimates in the first half and the second half of the run.
independent runs converged to the same result .
to perform likelihood-ratio tests of nested models, we combined results of the three independent ima <dig> runs  in a single l mode run.
all estimated parameters of the im model are scaled to the mutation rate.
to convert these parameters to biologically meaningful quantities  we calculated the neutral mutation rate for each locus  using divergence to the outgroup.
the neutral mutation rate per year was calculated for each locus by using the formula d =  <dig> μt <dig> where d is the average pairwise divergence, dxy, between the three studied species and the outgroup, and t <dig> is the divergence time.
we assumed that the divergence time between a. scirpaceus/a.
palustris/a.
dumetorum and a. schoenobaenus is  <dig>  million years.
this was estimated on the basis of cytochrome b sequence divergence, which is  <dig>  %  <cit> , and assuming approximately  <dig>  % sequence divergence per million years  <cit> .
the geometric mean of the locus-specific mutation rates was then used for the parameter conversion.
a generation time of one year  was assumed to estimate mutation rates per generation.
RESULTS
levels of intraspecific polymorphism and tests of neutrality
for each of the eight sequenced loci, we obtained  <dig> –  <dig> computationally reconstructed haploid sequences .
levels of nucleotide variation were generally high in all three studied species, suggesting a high effective population size of the species.
when averaged over all loci, π =  <dig>  % and θ =  <dig>  % in a. scirpaceus, π =  <dig>  % and θ =  <dig>  % in a. palustris, and π =  <dig>  % and θ =  <dig>  % in a. dumetorum.
we, however, observed a substantial heterogeneity in the levels of nucleotide variation among loci.
notably, within-species nucleotide variation was zero or very low in three out of four z-linked loci: adamts <dig>  tg <dig> and tg <dig>  for adamts <dig> and tg <dig>  π and θ reached 0 % in a. scirpaceus and were also quite low in the remaining two species.
for tg <dig>  π and θ reached 0 % even in all three species .
this can be also seen on haplotype networks showing zero or very low haplotype diversity within the species in the three loci .fig.
<dig> haplotype networks of four autosomal  and four z-linked  loci.
size of the circles are proportional to the number of haplotypes.
a. scirpaceus is indicated in red, a. palustris in green and a. dumetorum in grey.
in the case of  <dig> and  <dig>  haplotype networks were constructed using only sequences with no recombination witin loci obtained with program imgc
to test whether this heterogeneity in nucleotide variation could be attributed to the action of positive selection at some loci, we performed a multilocus hka test.
the test was applied separately to each species, in each case using a single sequence from a. schoenobaenus as the outgroup.
the null model was rejected in all three species , suggesting that positive selection affects levels of nucleotide variation at one or more loci in all three species.
a closer inspection of levels of polymorphism and divergence for individual loci revealed that loci tg <dig> in a. scirpaceus and a. palustris, and tg <dig> in a. scirpaceus, a. palustris and a. dumetorum showed the highest deviations from expected values .
when these two genes were excluded from the hka test, the null model was not rejected in a. palustris and a. dumetorum, but was still rejected in a. scirpaceus .
this suggests that a recent selective sweep might have occurred at these two loci independently in multiple species.table  <dig> hka test of positive selection
speciesa
p-valued
p-valuee
a.s.
<dig>
<dig>
a.p.
<dig>
<dig>
a.d.
<dig>
<dig>
a
a. scirpaceus , a. palustris , a. dumetorum
bcounted according to the formula ∑ 2/variance)
cdegree of freedom
dprobability from chi-square distribution
eprobability from simulations
we also looked for evidence of selection by comparing the distribution of allele frequencies with the expectations under a neutral equilibrium model using tajima’s d and fu-li’s d tests .
these tests were performed for all studied loci separately in each species.
there was a tendency toward negative values of tajima’s d in all three species , which may reflect mild population expansions.
nonetheless, none of the  <dig> values of tajima’s d within species were significantly different from the neutral expectation of  <dig>  no significant deviations from neutral expectations were either found in fu-li’s d tests.
these tests, however, could not be performed for loci with zero within-species nucleotide variation.
estimation of divergence times, effective population sizes and rates of interspecific gene flow
using im model for multiple populations, we estimated that the time since divergence between a. scirpaceus and a. palustris is  <dig>  mya and between the common ancestor of these species and a. dumetorum  <dig>  mya.
the estimated ne was similar for all three species,  <dig>  for a. scirpaceus,  <dig>  for a. palustris, and  <dig>  for a. dumetorum.
ne for the common ancestor of a. scirpaceus and a. palustris was substantially lower,  <dig> .
ne for the common ancestor of all three species could not be estimated accurately since the posterior probability distribution for this parameter was very flat .
the estimates of migration rates between a. scirpaceus and a. palustris were relatively high, although only in one direction .
migration rates between a. palustris and a. dumetorum were lower  and no migration was detected between a. scirpaceus and a. dumetorum .
estimates of migration rates between a. dumetorum and the common ancestor of a. scirpaceus and a. palustris were quite high, though again only in one direction ; however, these estimates might not be accurate since posterior probability distributions of these parameter were flat .
maximum-likelihood estimates  and 95 % highest posterior density  intervals of all model parameters are given in table  <dig>  and the marginal posterior probability distributions are shown in fig.
<dig> fig.
<dig> the marginal posterior probability distributions for the demographic parameters of the im model.
effective population sizes of a. scirpaceus , a. palustris , a. dumetorum , the common ancestor of a. scirpaceus and a. palustris , and the common ancestor of all three species .
divergence time between a. scirpaceus and a. palustris  and between a. dumetorum and the common ancestor of a. scirpaceus and a. palustris .
migration rates per year
ns
c
np
c
nd
c
nsp
c
ms>p
d
mp>s
d
ms>d
d
md>s
d
mp>d
d
md>p
d
md>sp
d
msp>d
d
2nsms>p
e
2nsms>d
e
2npmp>s
e
2npmp>d
e
2ndmd>s
e
2ndmd>p
e
2ndmd>sp
e
2nspmsp>d
e
t0
f
t1
f
alower 95 % hpd
bhigher 95 % hpd
ceffective population sizes of a. scirpaceus , a. palustris , a. dumetorum , and the common ancestor of a. scirpaceus and a. palustris  in numbers of individuals
dmigration rates per year
ethe population migration rate
fdivergence time between a. scirpaceus and a. palustris  and between a. dumetorum and the common ancestor of a. scirpaceus and a. palustris  in numbers of years
to test whether our estimated model with gene flow fits significantly better the data than models without gene flow between particular species, we performed log-likelihood ratio tests of nested models  <cit>  as implemented in the ima <dig> program.
we tested four nested models  where both migration parameters between  a. scirpaceus and a. palustris,  a. scirpaceus and a. dumetorum,  a. palustris and a. dumetorum, and  a. dumetorum and the common ancestor of a. scirpaceus and a. palustris were set to zero.
the full model with gene flow was a significantly better fit to the data than all nested models, except the second one .table  <dig> log-likelihood ratio tests of nested models
modela
< <dig>
<dig>
<dig>
afour nested models with zero migration rates between  a. scirpaceus and a. palustris,  a. scirpaceus and a. dumetorum,  a. palustris and a. dumetorum, and  a. dumetorum and the common ancestor of a. scirpaceus and a. palustris were compared to the estimated full model
bestimates of the posterior density function under the full model
clog-likelihood ratio statistics calculated as the difference between the highest posterior probability for the full model and the highest posterior probability for the nested model
dthe degrees of freedom.
models in which migration at least in one direction is equal to  <dig> have distributions of 2llr that are a mixture and χ <dig> mixture distribution was thus used to calculate the p-value
ethe probability of achieving the test statistics by chance under the null model
we were further interested in whether patterns of interspecific gene flow between a. scirpaceus and a. palustris, and between a. palustris and a. dumetorum are similar or not.
as ima <dig> does not allow to estimate locus-specific migration rates, we compared fst values, which reflect between-species differentiation and are correlated with levels of gene flow in hybridizing taxa .
we should, however, note that fst values might be also affected by the action of positive selection that leads to a reduced within-species nucleotide variation  <cit> .
no statistically significant correlation in fst between the two species pairs was found , suggesting that the patterns of interspecific gene flow and/or positive selection are in general different for the two species pairs.
discussion
adaptive radiations represent an outstanding system for studying mechanisms of speciation and adaptive divergence.
although divergent ecological selection seems to be an important driver of many adaptive radiations  <cit> , it has been suggested that special genomic properties, such as gene duplications, activation of transposable elements or presence of gene flow among species, can also contribute to high speciation rates and rapid phenotypic changes during adaptive radiations  <cit> .
acrocephalus warblers represent an interesting model system for studying mechanisms of radiation in birds, although their mainly continental radiation was relatively slow compared to the well-known cases of adaptive radiations on islands, such as darwin finches  <cit> .
our estimates of divergence times are approximately  <dig>  mya for the sister species a. scirpaceus and a. palustris, and  <dig>  mya for their common ancestor and a. dumetorum.
despite this relatively deep divergence, many currently existing acrocephalus species still hybridize in nature .
this hybridization might in principle lead to gene flow among the species if f <dig> hybrids are fertile.
however, since very little is known about the degree of postzygotic isolation in this avian group, it is hard to predict whether gene flow could occur among acrocephalus species and how strong it might be.
here we provided the first evidence that interspecific gene flow occurred during acrocephalus radiation and not only between sister species.
we found evidence of gene flow between a. scirpaceus and a. palustris, which are the sister species with somewhat different ecological requirements, and between a. palustris and a. dumetorum, which have similar ecological requirements and often co-occur at the same sites where their ranges overlap.
significant gene flow was detected also between a. dumetorum and the ancestor of a. scirpaceus and a. palustris.
no gene flow was, however, detected between a. scirpaceus and a. dumetorum.
although the isolation-with-migration model estimates not only contemporary gene flow, but average gene flow among populations since the time of divergence  <cit> , it is noteworthy that hybridization currently occurs between both pairs of species for which gene flow has been detected, but not between a. scirpaceus and a. palustris for which no clear cases of interspecific hybridization are known.
interestingly, the estimated levels of gene flow were in all cases asymmetric.
gene flow occurred from a. scirpaceus to a. palustris , from a. palustris to a. dumetorum , and from the a. scirpaceus/a.
palustris ancestor to a. dumetorum .
zero gene flow was detected in the opposite directions.
such unidirectional gene flow is expected when a population of one species expands to the area already occupied by a related species and interbreeding is not prevented between the two species.
introgression of neutral alleles then occurs almost exclusively from the local into the invading species  <cit> .
unfortunately, history of breeding range changes are not well known for the studied species, but historically recent breeding range expansion of a. dumetorum to eastern and northern europe where it encountered a. palustris  <cit>  is consistent with this scenario.
an important issue in speciation research is to understand evolutionary forces that are responsible for the establishment of reproductive barriers between the species.
here we tested whether gene flow is higher between genetically more similar, but ecologically somewhat different, a. scirpaceus and a. palustris, or between ecologically more similar, but genetically more distant, a. palustris and a. dumetorum.
our results show that gene flow is more than three times higher between genetically more similar a. scirpaceus and a. palustris than between a. palustris and a. dumetorum, although the latter two species have higher chance to meet in the breeding sites and very likely more often hybridize.
lower levels of gene flow between a. palustris and a. dumetorum might be caused by historically more recent contact between the species  <cit> , but might also suggest that divergent ecological selection is not the major cause of reproductive isolation and that speciation in acrocephalus warblers is more likely driven by gradual accumulation of postzygotic barriers, such as hybrid sterility.
hybrid sterility have been thought to evolve quite slowly in birds compared to other taxonomic groups; however, recent findings suggest that at least in some bird species it can evolve very rapidly  <cit> .
it is, however, also possible that very different migratory behaviors of a. dumetorum, wintering in india, and a. palustris, that winters similarly as a. scirpaceus in africa , might have contributed to their stronger reproductive isolation.
indeed, studies in swainson’s thrush  show that different seasonal migratory behavior of its two subspecies can be a source of selection against hybrids  <cit> .
analysis of more species with different genetic and ecological similarity would be needed to conclude which selective forces drive radiation of acrocephalus warblers.
although the patterns of between species differentiation were in general different for the two species pairs, we detected signatures of apparently independent instances of positive selection at the same two z-linked loci in multiple species.
selective sweeps reducing the within-species nucleotide variation occurred independently in a. scirpaceus and a. palustris at tg <dig> locus and even in all three species at tg <dig> locus.
although we do not know the physical location of the two loci in the acrocephalus genome, the two loci lie in different parts of the z chromosome in the zebra finch genome.
given the very conserved karyotypes in birds  <cit> , it is likely that the two loci are also unlinked in acrocephalus genome and that selection has led to reduced nucleotide variation independently at the two loci.
similar examples of reduced nucleotide diversity at the same loci in multiple species have been observed for example in a genomic scan of islands of divergence in two hybridizing ficedula flycatchers  <cit>  or between two closely related and hybridizing plants, silene latifolia and s. dioica  <cit> .
such independent instances of positive selection in multiple species suggest that the same  loci might have been subject to parallel selection in two or even all three acrocephalus species.
this might happen for example in genes involved in intragenomic conflicts, such as meiotic drive, which is expected to be especially common on the sex chromosomes.
whether divergence of these genes could contribute to acrocephalus speciation, however, still remains to be explored.
in addition, it should be tested whether similar patterns of adaptive evolution in the three acrocephalus species could not result from linked selection in low-recombinant regions rather than parallel adaptive evolution  <cit> .
CONCLUSIONS
our study provides the first evidence that interspecific gene flow occurred during acrocephalus radiation and not only between sister species.
interspecific gene flow could thus play an important role in shaping genetic variation in individual species and possibly accelerate adaptive evolution and speciation rates in this group.
interestingly, the gene flow was in all cases asymmetric, suggesting that it might have occurred mainly during range expansions of individual species into areas of closely related and not completely reproductively isolated species.
interspecific competition leading to ecological character displacement and habitat segregation  <cit>  and/or reinforcement of prezygotic isolation  <cit>  could then strengthen reproductive barriers between the species.
however, in contrast to darwin finches  <cit> , our results showed that gene flow was higher between genetically more closely related but ecologically more divergent species suggesting that gradual accumulation of genetic incompatibilities causing intrinsic isolation rather than divergent ecological selection are more efficient in restricting interspecific gene flow in this avian group.
thus in contrast to darwin finches  <cit> , the ecological differentiation of acrocephalus species likely arose in later stages of speciation.
data on more than three acrocephalus species would be, however, needed to generalize this conclusion to the whole acrocephalus radiation.
in addition, genome-wide data across autosomes and the whole z chromosome would be needed to understand mechanisms responsible for the parallel adaptive evolution at the same loci in multiple species and its possible role in acrocephalus speciation.
nevertheless, occurrence of interspecific hybridization among many acrocephalus species suggests that interspecific gene flow could be common across the whole radiation.
additional files
additional file 1: published cases of interspecific mating  and interspecific hybrids  among acrocephalus species.
additional file 2: list of analyzed samples including information of sex of each individual and geographic coordinates of the sites, where individual birds were captured.
additional file 3: pcr primers and length of pcr products.
additional file 4: observed and expected levels of polymorphism and divergence in hka test of positive selection.

 interspecific gene flow is an important evolutionary force. it may enrich genetic variation of individual species, facilitate the origin of new phenotypes, lead to the origin of new species, or, conversely, to species fusion. it has been suggested that high speciation rates and rapid phenotypic changes observed during adaptive radiations may be facilitated by increased genetic variation and novel genetic combinations produced by introgressive hybridization among the species. patterns of interspecific gene flow can also provide an important insight into the genetic basis of reproductive isolation. in species with incomplete reproductive isolation, gene flow occurs at loci with alleles that are neutral or beneficial on the genomic background of the other species, but is limited at loci harboring alleles involved in reproductive isolation. most work exploring patterns and rates of interspecific gene flow focused on pairs of closely related sister species because, until recently, models allowing the analysis of interspecific gene flow were limited to only two species. several general findings emerged from these studies. first, in the case of speciation driven by intrinsic barriers, interspecific gene flow is often reduced on the sex chromosomes and in genomic regions with low recombination rates, indicating that genetic incompatibilities tend to preferentially accumulate in these genomic regions. second, in the case of speciation driven by divergent ecological selection, gene flow is correlated with ecological divergence and is reduced around genes underlying ecological differences. interspecific gene flow across larger sections of the phylogenetic tree is, however, much less studied. such studies would be particularly interesting in instances of rapid radiations, where lineage divergence occurs on a shorter timescale than the completion of reproductive isolation. knowledge of patterns and rates of interspecific gene flow in such systems allows to address questions, which are difficult or even impossible to study using models of only two sister species. one can, for example, examine whether the same genes contribute to reproductive isolation in different speciation events, assess the importance of various reproductive barriers in restricting the gene flow and determine the order in which different reproductive barriers arise during lineage divergence. this information can help us to understand why lineage divergence occurs so quickly in some taxonomic groups, while not in others. acrocephalus warblers are small passerines that underwent an adaptive radiation during the last approximately million years, which gave rise to extant species occupying mainly eurasia, africa, australia and pacific islands. the members of this genus are phenotypically quite uniform. they, however, display a great diversity in ecology and behaviour, and have been used extensively as models in ecological and behavioral research. interspecific hybridization occasionally occurs in acrocephalus warblers. it has been reported not only between sister species, but also between more distantly related species. in some of these cases, hybridization occurs between species of very different body size or even between species belonging to different subgenera. so far very little is known about the degree of postzygotic reproductive isolation between the species. in some cases, e.g. between sister species a. arundinaceus and a. stentoreus, postzygotic isolation seems to be quite strong as no backcross hybrids have been observed in a hybridizing population. in other cases, however, observation of fertile f hybrids suggests that postzygotic isolation is incomplete and interspecific gene flow might occur between the species. here we performed the first study exploring patterns and rates of interspecific gene flow in this avian group, focusing on three species of the subgenus notiocichla: a. scirpaceus, a. palustris and a. dumetorum. all three species are very similar in morphology, but have different habitat requirements, migration strategies and song patterns. a. scirpaceus breeds across europe, northern africa and western asia typically in reed beds and winters in sub-saharan africa. its sister species, a. palustris, has similar geographical distribution, but occupies mostly damp herbaceous vegetation and winters in southeastern africa. males of each species sing different songs, nevertheless, in a. palustris, mixed singers that incorporate song phrases from a. scirpaceus to its repertoire are known. both species occasionally form mixed pairs and are able to produce viable hybrid offspring. a. dumetorum, the sister species to the a. scirpaceus and a. palustris clade, breeds in asia and northeastern europe in various herbaceous or bushy vegetation and winters in india. in eastern europe and western asia, it often co-occurs with a. palustris at the same sites and hybridization between the two species has been observed. the exact frequency of hybridization between a. palustris and a. scirpaceus and between a. palustris and a. dumetorum is hard to estimate owing to morphological similarity of the species and difficulty to distinguish hybrids without genetic analysis. the total number of observed mixed pairs is, however, comparable in both species pairs. the cases of hybridization between a. dumetorum and a. scirpaceus, which are ecologically more divergent, are not known, although one putative hybrid individual has been reported. to get the first insight into mechanisms and selective forces driving acrocephalus radiation, we analyzed nucleotide variation at four autosomal and four loci on the z chromosome in sympatric populations of a. scirpaceus, a. palustris and a. dumetorum. using an isolation-with-migration model for multiple populations we tested whether interspecific gene flow occurred among the studied species and if so, whether it was higher between more genetically similar or between more ecologically similar species. in addition, using several tests of neutrality, we explored patterns of positive selection in all three species. methods study area and samples we analyzed samples from individuals of a. scirpaceus, individuals of a. palustris and individuals of a. dumetorum. all birds were caught in during the breeding season in southern finland, where the breeding ranges of all three species overlap. coordinates of the sites, where individual birds were captured, are shown in additional file from each individual, a blood sample was collected by brachial vein puncture and stored in pure ethanol for further extraction of genomic dna. in addition, a sample from one individual of a. schoenobaenus was used as the outgroup.fig. breeding ranges of the three studied acrocephalus warblers, a. scirpaceus, a. palustris and a. dumetorum. sampling sites are indicated by dark grey circles molecular sex determination and sequencing genomic dna was purified by dneasy tissue kit according to manufacturer’s instructions. sex of individual birds was determined following the method of griffiths et al. using p and p primers. these primers amplify a fragment of two homologous genes located on the z and w chromosomes. the amplified fragments from the z and w chromosomes differ in length, which can be visualized on agarose gel. we further amplified and sequenced introns of four autosomal and four z-linked loci. the four autosomal loci are located on different chromosomes in the zebra finch genome. of the four z-linked genes, ppwd and adamts are close physically, and therefore may not have independent evolutionary histories. however, the other two loci on the z chromosome are separated from each other and from ppwd and adamts by > 19 mb. primers for pcr amplification were partially obtained from published studies and partially designed by ourselves. all primers were designed in conserved exonic regions of the zebra finch or chicken genome in order to amplify intronic sequences. primer sequences and lengths of the obtained pcr products are provided in additional file pcr conditions are same as in storchová et al.. all pcr products were sequenced in both directions with the pcr primers using sanger sequencing at the dna sequencing laboratory at charles university in prague. data analyses sequences were manually edited using codoncode aligner software. alignments were generated by clustalw as implemented in the program bioedit. all alignments were visually checked and manually adjusted. exonic sequences as well as indel polymorphisms were excluded from the analyses. individuals and/or positions with missing data were eliminated from the dataset to get sequences of the same length for each locus. diploid sequences at each locus were further separated into two haplotypes using the program phase, version. with the following parameters: number of iterations = , thinning interval =  burnin =  we used the default recombination model, which is the general model for varying recombination rate. for each dataset, we applied the algorithm five times with different random seeds, and we checked for consistency of results across independent runs. we obtained identical haplotypes across all runs for all loci. the final number of haploid sequences for each locus and the length of sequences are shown in table  table polymorphism statistics for eight studied loci 17483 a.s. a.p. a.d. 21281 a.s. a.p. a.d. 24972 a.s. a.p. a.d. rpl5-4 a.s. a.p. a.d. adamts6 a.s. a.p. a.d. ppwd1 a.s. a.p. a.d. tg401 a.s. a.p. a.d. tg1505 a.s. a.p. a.d. aautosome, chromosome z blength of sequence c a. scirpaceus, a. palustris, a. dumetorum dnumber of haploid sequences enumber of segregating sites faverage number of nucleotide differences gproportion of polymorphic sites htajima’s d ifu and li’s d jdivergence to outgroup measured as average pairwise divergence, dxy basic population genetic analyses of polymorphism, divergence, recombination, and tests of neutrality based on the allele frequency spectrum were performed with the program dnasp. hudson-kreitman-aguadé tests of positive selection were performed using the hka program. this program compares the ratio of polymorphism to divergence at multiple loci. loci under positive selection are expected to show lower levels of polymorphism relative to divergence. for each locus, we set the information whether the locus is autosomal or z-linked, based on which the hka program accounts for lower effective population size and thus the lower nucleotide variation on the z chromosome. genealogical relationships among haplotypes for each locus were reconstructed with network software using median joining algorithm. the data were fit to the isolation-with-migration model for multiple populations using the program ima the program estimates several demographic parameters based on markov chain monte carlo simulations of genealogies. these parameters include the effective population size for the current and ancestral populations, migration rates between the populations and the population-split times. because ima assumes no recombination within loci, we determined the longest region without observed recombination for each locus using the program imgc. this program removes either sites or haplotypes to produce the most information-rich contiguous dna sequence segment that passes the four-gamete test. nonrecombinant regions represented 84 % of the length of each locus, on average, and were used as an input for the ima program. for each locus, we provided the information whether it is autosomal or z-linked, based on which ima program accounts for different effective population sizes between autosomes and the z chromosome. we ran the program three times with identical starting conditions, with the exception of the random number seed, to assess convergence. to facilitate mixing of the markov chains, we used metropolis coupling with chains and a geometric heating model. upper bounds for the prior distributions of parameter values were set as suggested in ima documentation. all runs began with a burn-in period of steps and were allowed to continue for 7– million steps. we were able to achieve adequate mixing of the markov chains as indicated by trend line plots, effective sample size values, and very similar parameter estimates in the first half and the second half of the run. independent runs converged to the same result. to perform likelihood-ratio tests of nested models, we combined results of the three independent ima runs in a single l mode run. all estimated parameters of the im model are scaled to the mutation rate. to convert these parameters to biologically meaningful quantities we calculated the neutral mutation rate for each locus using divergence to the outgroup. the neutral mutation rate per year was calculated for each locus by using the formula d =  μt where d is the average pairwise divergence, dxy, between the three studied species and the outgroup, and t is the divergence time. we assumed that the divergence time between a. scirpaceus/a. palustris/a. dumetorum and a. schoenobaenus is million years. this was estimated on the basis of cytochrome b sequence divergence, which is  %, and assuming approximately  % sequence divergence per million years. the geometric mean of the locus-specific mutation rates was then used for the parameter conversion. a generation time of one year was assumed to estimate mutation rates per generation.  levels of intraspecific polymorphism and tests of neutrality for each of the eight sequenced loci, we obtained – computationally reconstructed haploid sequences. levels of nucleotide variation were generally high in all three studied species, suggesting a high effective population size of the species. when averaged over all loci, π =   % and θ =   % in a. scirpaceus, π =   % and θ =   % in a. palustris, and π =   % and θ =   % in a. dumetorum. we, however, observed a substantial heterogeneity in the levels of nucleotide variation among loci. notably, within-species nucleotide variation was zero or very low in three out of four z-linked loci: adamts tg and tg for adamts and tg π and θ reached 0 % in a. scirpaceus and were also quite low in the remaining two species. for tg π and θ reached 0 % even in all three species. this can be also seen on haplotype networks showing zero or very low haplotype diversity within the species in the three loci.fig. haplotype networks of four autosomal and four z-linked loci. size of the circles are proportional to the number of haplotypes. a. scirpaceus is indicated in red, a. palustris in green and a. dumetorum in grey. in the case of and haplotype networks were constructed using only sequences with no recombination witin loci obtained with program imgc to test whether this heterogeneity in nucleotide variation could be attributed to the action of positive selection at some loci, we performed a multilocus hka test. the test was applied separately to each species, in each case using a single sequence from a. schoenobaenus as the outgroup. the null model was rejected in all three species, suggesting that positive selection affects levels of nucleotide variation at one or more loci in all three species. a closer inspection of levels of polymorphism and divergence for individual loci revealed that loci tg in a. scirpaceus and a. palustris, and tg in a. scirpaceus, a. palustris and a. dumetorum showed the highest deviations from expected values. when these two genes were excluded from the hka test, the null model was not rejected in a. palustris and a. dumetorum, but was still rejected in a. scirpaceus. this suggests that a recent selective sweep might have occurred at these two loci independently in multiple species.table hka test of positive selection speciesa p-valued p-valuee a.s. a.p. a.d. a a. scirpaceus, a. palustris, a. dumetorum bcounted according to the formula ∑ 2/variance) cdegree of freedom dprobability from chi-square distribution eprobability from simulations we also looked for evidence of selection by comparing the distribution of allele frequencies with the expectations under a neutral equilibrium model using tajima’s d and fu-li’s d tests. these tests were performed for all studied loci separately in each species. there was a tendency toward negative values of tajima’s d in all three species, which may reflect mild population expansions. nonetheless, none of the values of tajima’s d within species were significantly different from the neutral expectation of no significant deviations from neutral expectations were either found in fu-li’s d tests. these tests, however, could not be performed for loci with zero within-species nucleotide variation. estimation of divergence times, effective population sizes and rates of interspecific gene flow using im model for multiple populations, we estimated that the time since divergence between a. scirpaceus and a. palustris is mya and between the common ancestor of these species and a. dumetorum mya. the estimated ne was similar for all three species, for a. scirpaceus, for a. palustris, and for a. dumetorum. ne for the common ancestor of a. scirpaceus and a. palustris was substantially lower,. ne for the common ancestor of all three species could not be estimated accurately since the posterior probability distribution for this parameter was very flat. the estimates of migration rates between a. scirpaceus and a. palustris were relatively high, although only in one direction. migration rates between a. palustris and a. dumetorum were lower and no migration was detected between a. scirpaceus and a. dumetorum. estimates of migration rates between a. dumetorum and the common ancestor of a. scirpaceus and a. palustris were quite high, though again only in one direction; however, these estimates might not be accurate since posterior probability distributions of these parameter were flat. maximum-likelihood estimates and 95 % highest posterior density intervals of all model parameters are given in table  and the marginal posterior probability distributions are shown in fig.  fig. the marginal posterior probability distributions for the demographic parameters of the im model. effective population sizes of a. scirpaceus, a. palustris, a. dumetorum, the common ancestor of a. scirpaceus and a. palustris, and the common ancestor of all three species. divergence time between a. scirpaceus and a. palustris and between a. dumetorum and the common ancestor of a. scirpaceus and a. palustris. migration rates per year ns c np c nd c nsp c ms>p d mp>s d ms>d d md>s d mp>d d md>p d md>sp d msp>d d 2nsms>p e 2nsms>d e 2npmp>s e 2npmp>d e 2ndmd>s e 2ndmd>p e 2ndmd>sp e 2nspmsp>d e t0 f t1 f alower 95 % hpd bhigher 95 % hpd ceffective population sizes of a. scirpaceus, a. palustris, a. dumetorum, and the common ancestor of a. scirpaceus and a. palustris in numbers of individuals dmigration rates per year ethe population migration rate fdivergence time between a. scirpaceus and a. palustris and between a. dumetorum and the common ancestor of a. scirpaceus and a. palustris in numbers of years to test whether our estimated model with gene flow fits significantly better the data than models without gene flow between particular species, we performed log-likelihood ratio tests of nested models as implemented in the ima program. we tested four nested models where both migration parameters between a. scirpaceus and a. palustris, a. scirpaceus and a. dumetorum, a. palustris and a. dumetorum, and a. dumetorum and the common ancestor of a. scirpaceus and a. palustris were set to zero. the full model with gene flow was a significantly better fit to the data than all nested models, except the second one.table log-likelihood ratio tests of nested models modela < afour nested models with zero migration rates between a. scirpaceus and a. palustris, a. scirpaceus and a. dumetorum, a. palustris and a. dumetorum, and a. dumetorum and the common ancestor of a. scirpaceus and a. palustris were compared to the estimated full model bestimates of the posterior density function under the full model clog-likelihood ratio statistics calculated as the difference between the highest posterior probability for the full model and the highest posterior probability for the nested model dthe degrees of freedom. models in which migration at least in one direction is equal to have distributions of 2llr that are a mixture and χ mixture distribution was thus used to calculate the p-value ethe probability of achieving the test statistics by chance under the null model we were further interested in whether patterns of interspecific gene flow between a. scirpaceus and a. palustris, and between a. palustris and a. dumetorum are similar or not. as ima does not allow to estimate locus-specific migration rates, we compared fst values, which reflect between-species differentiation and are correlated with levels of gene flow in hybridizing taxa. we should, however, note that fst values might be also affected by the action of positive selection that leads to a reduced within-species nucleotide variation. no statistically significant correlation in fst between the two species pairs was found, suggesting that the patterns of interspecific gene flow and/or positive selection are in general different for the two species pairs. discussion adaptive radiations represent an outstanding system for studying mechanisms of speciation and adaptive divergence. although divergent ecological selection seems to be an important driver of many adaptive radiations, it has been suggested that special genomic properties, such as gene duplications, activation of transposable elements or presence of gene flow among species, can also contribute to high speciation rates and rapid phenotypic changes during adaptive radiations. acrocephalus warblers represent an interesting model system for studying mechanisms of radiation in birds, although their mainly continental radiation was relatively slow compared to the well-known cases of adaptive radiations on islands, such as darwin finches. our estimates of divergence times are approximately mya for the sister species a. scirpaceus and a. palustris, and mya for their common ancestor and a. dumetorum. despite this relatively deep divergence, many currently existing acrocephalus species still hybridize in nature. this hybridization might in principle lead to gene flow among the species if f hybrids are fertile. however, since very little is known about the degree of postzygotic isolation in this avian group, it is hard to predict whether gene flow could occur among acrocephalus species and how strong it might be. here we provided the first evidence that interspecific gene flow occurred during acrocephalus radiation and not only between sister species. we found evidence of gene flow between a. scirpaceus and a. palustris, which are the sister species with somewhat different ecological requirements, and between a. palustris and a. dumetorum, which have similar ecological requirements and often co-occur at the same sites where their ranges overlap. significant gene flow was detected also between a. dumetorum and the ancestor of a. scirpaceus and a. palustris. no gene flow was, however, detected between a. scirpaceus and a. dumetorum. although the isolation-with-migration model estimates not only contemporary gene flow, but average gene flow among populations since the time of divergence, it is noteworthy that hybridization currently occurs between both pairs of species for which gene flow has been detected, but not between a. scirpaceus and a. palustris for which no clear cases of interspecific hybridization are known. interestingly, the estimated levels of gene flow were in all cases asymmetric. gene flow occurred from a. scirpaceus to a. palustris, from a. palustris to a. dumetorum, and from the a. scirpaceus/a. palustris ancestor to a. dumetorum. zero gene flow was detected in the opposite directions. such unidirectional gene flow is expected when a population of one species expands to the area already occupied by a related species and interbreeding is not prevented between the two species. introgression of neutral alleles then occurs almost exclusively from the local into the invading species. unfortunately, history of breeding range changes are not well known for the studied species, but historically recent breeding range expansion of a. dumetorum to eastern and northern europe where it encountered a. palustris is consistent with this scenario. an important issue in speciation research is to understand evolutionary forces that are responsible for the establishment of reproductive barriers between the species. here we tested whether gene flow is higher between genetically more similar, but ecologically somewhat different, a. scirpaceus and a. palustris, or between ecologically more similar, but genetically more distant, a. palustris and a. dumetorum. our results show that gene flow is more than three times higher between genetically more similar a. scirpaceus and a. palustris than between a. palustris and a. dumetorum, although the latter two species have higher chance to meet in the breeding sites and very likely more often hybridize. lower levels of gene flow between a. palustris and a. dumetorum might be caused by historically more recent contact between the species, but might also suggest that divergent ecological selection is not the major cause of reproductive isolation and that speciation in acrocephalus warblers is more likely driven by gradual accumulation of postzygotic barriers, such as hybrid sterility. hybrid sterility have been thought to evolve quite slowly in birds compared to other taxonomic groups; however, recent findings suggest that at least in some bird species it can evolve very rapidly. it is, however, also possible that very different migratory behaviors of a. dumetorum, wintering in india, and a. palustris, that winters similarly as a. scirpaceus in africa, might have contributed to their stronger reproductive isolation. indeed, studies in swainson’s thrush show that different seasonal migratory behavior of its two subspecies can be a source of selection against hybrids. analysis of more species with different genetic and ecological similarity would be needed to conclude which selective forces drive radiation of acrocephalus warblers. although the patterns of between species differentiation were in general different for the two species pairs, we detected signatures of apparently independent instances of positive selection at the same two z-linked loci in multiple species. selective sweeps reducing the within-species nucleotide variation occurred independently in a. scirpaceus and a. palustris at tg locus and even in all three species at tg locus. although we do not know the physical location of the two loci in the acrocephalus genome, the two loci lie in different parts of the z chromosome in the zebra finch genome. given the very conserved karyotypes in birds, it is likely that the two loci are also unlinked in acrocephalus genome and that selection has led to reduced nucleotide variation independently at the two loci. similar examples of reduced nucleotide diversity at the same loci in multiple species have been observed for example in a genomic scan of islands of divergence in two hybridizing ficedula flycatchers or between two closely related and hybridizing plants, silene latifolia and s. dioica. such independent instances of positive selection in multiple species suggest that the same loci might have been subject to parallel selection in two or even all three acrocephalus species. this might happen for example in genes involved in intragenomic conflicts, such as meiotic drive, which is expected to be especially common on the sex chromosomes. whether divergence of these genes could contribute to acrocephalus speciation, however, still remains to be explored. in addition, it should be tested whether similar patterns of adaptive evolution in the three acrocephalus species could not result from linked selection in low-recombinant regions rather than parallel adaptive evolution.  our study provides the first evidence that interspecific gene flow occurred during acrocephalus radiation and not only between sister species. interspecific gene flow could thus play an important role in shaping genetic variation in individual species and possibly accelerate adaptive evolution and speciation rates in this group. interestingly, the gene flow was in all cases asymmetric, suggesting that it might have occurred mainly during range expansions of individual species into areas of closely related and not completely reproductively isolated species. interspecific competition leading to ecological character displacement and habitat segregation and/or reinforcement of prezygotic isolation could then strengthen reproductive barriers between the species. however, in contrast to darwin finches, our results showed that gene flow was higher between genetically more closely related but ecologically more divergent species suggesting that gradual accumulation of genetic incompatibilities causing intrinsic isolation rather than divergent ecological selection are more efficient in restricting interspecific gene flow in this avian group. thus in contrast to darwin finches, the ecological differentiation of acrocephalus species likely arose in later stages of speciation. data on more than three acrocephalus species would be, however, needed to generalize this conclusion to the whole acrocephalus radiation. in addition, genome-wide data across autosomes and the whole z chromosome would be needed to understand mechanisms responsible for the parallel adaptive evolution at the same loci in multiple species and its possible role in acrocephalus speciation. nevertheless, occurrence of interspecific hybridization among many acrocephalus species suggests that interspecific gene flow could be common across the whole radiation. additional files additional file 1: published cases of interspecific mating and interspecific hybrids among acrocephalus species. additional file 2: list of analyzed samples including information of sex of each individual and geographic coordinates of the sites, where individual birds were captured. additional file 3: pcr primers and length of pcr products. additional file 4: observed and expected levels of polymorphism and divergence in hka test of positive selection. 
independent instances of positive selection at the same loci in multiple species indicate an interesting possibility that the same loci might have contributed to reproductive isolation in several speciation events. gene flow was higher between genetically more closely related a. scirpaceus and a. palustris than between ecologically more similar a. palustris and a. dumetorum, suggesting that gradual accumulation of intrinsic barriers rather than divergent ecological selection are more efficient in restricting interspecific gene flow in acrocephalus warblers. here we studied patterns of interspecific gene flow and selection across three european acrocephalus warblers to get a first insight into mechanisms of radiation of this avian group. using an isolation-with-migration model for multiple populations, we found evidence for unidirectional gene flow from a. scirpaceus to a. palustris and from a. palustris to a. dumetorum. our study brings the first evidence that gene flow occurred during acrocephalus radiation and not only between sister species. interspecific gene flow could thus be an important source of genetic variation in individual acrocephalus species and could have accelerated adaptive evolution and speciation rate in this avian group by creating novel genetic combinations and new phenotypes. adaptive radiationspeciationgene flowparallel adaptive evolutionz chromosomeacrocephalus warblersczech science foundation 15-10884yissue-copyright-statement© the author 2016 acrocephalus warblers are small passerines that underwent an adaptive radiation in the last approximately million years that gave rise to extant species, many of which still hybridize in nature. acrocephalus warblers have served as model organisms for a wide variety of ecological and behavioral studies, yet our knowledge of mechanisms and selective forces driving their radiation is limited. although levels of genetic differentiation between different species pairs were in general not correlated, we found signatures of apparently independent instances of positive selection at the same two z-linked loci in multiple species. understanding the mechanisms and selective forces leading to adaptive radiations and origin of biodiversity is a major goal of evolutionary biology. we analyzed nucleotide variation at eight nuclear loci in three hybridizing acrocephalus species with overlapping breeding ranges in europe.
 understanding the mechanisms and selective forces leading to adaptive radiations and origin of biodiversity is a major goal of evolutionary biology. acrocephalus warblers are small passerines that underwent an adaptive radiation in the last approximately million years that gave rise to extant species, many of which still hybridize in nature. acrocephalus warblers have served as model organisms for a wide variety of ecological and behavioral studies, yet our knowledge of mechanisms and selective forces driving their radiation is limited. here we studied patterns of interspecific gene flow and selection across three european acrocephalus warblers to get a first insight into mechanisms of radiation of this avian group.  we analyzed nucleotide variation at eight nuclear loci in three hybridizing acrocephalus species with overlapping breeding ranges in europe. using an isolation-with-migration model for multiple populations, we found evidence for unidirectional gene flow from a. scirpaceus to a. palustris and from a. palustris to a. dumetorum. gene flow was higher between genetically more closely related a. scirpaceus and a. palustris than between ecologically more similar a. palustris and a. dumetorum, suggesting that gradual accumulation of intrinsic barriers rather than divergent ecological selection are more efficient in restricting interspecific gene flow in acrocephalus warblers. although levels of genetic differentiation between different species pairs were in general not correlated, we found signatures of apparently independent instances of positive selection at the same two z-linked loci in multiple species.  our study brings the first evidence that gene flow occurred during acrocephalus radiation and not only between sister species. interspecific gene flow could thus be an important source of genetic variation in individual acrocephalus species and could have accelerated adaptive evolution and speciation rate in this avian group by creating novel genetic combinations and new phenotypes. independent instances of positive selection at the same loci in multiple species indicate an interesting possibility that the same loci might have contributed to reproductive isolation in several speciation events. electronic supplementary material the online version of this article contains supplementary material, which is available to authorized users. keywords adaptive radiationspeciationgene flowparallel adaptive evolutionz chromosomeacrocephalus warblersczech science foundation 15-10884yissue-copyright-statement© the author 2016 