Genetic approaches to the conservation of migratory bats: a study of the eastern red bat (Lasiurus borealis)

Introduction

As concerns about anthropogenic climate change and the long-term environmental impacts of burning of fossil fuels on biological and human systems have heightened, there is increasing motivation to develop alternative sources of energy that will reduce the production of greenhouse gasses. Wind power has become an increasingly important sector of the energy industry and is one of the fastest growing sources of renewable energy (Kaldellis & Zafirakis, 2011; Leung & Yang, 2012). Despite the many positive aspects of wind power development, there have been environmental costs associated with turbine installation and operation (Morrison & Sinclair, 2004; Abbasi et al., 2014). Fatalities of bats at wind power installations have emerged as a major environmental impact of wind power development, with large mortality events being reported at a number of wind energy facilities in the United States and abroad (Erickson et al., 2001; Erickson, Johnson & Young, 2005; Kunz et al., 2007; Arnett et al., 2008). The bat species most affected by wind power in North America are migratory, tree-roosting species such as hoary bats (Lasiurus cinereus), eastern red bats (Lasiurus borealis), and silver-haired bats (Lasionycteris noctivagans), which together constitute almost three-quarters of the bat carcasses found at wind turbines (Arnett et al., 2008). Although mortalities may occur throughout April to November, most bat fatalities in North America have been reported in late summer and early autumn (reviewed by Kunz et al., 2007; Arnett et al., 2008) and appear to be concentrated during fall migration of the affected species (Cryan, 2003).

The observed high levels of mortality for these species at wind power installations raise concerns about the long-term impacts of this technology on bat populations, yet we lack the necessary information to place this mortality in context with respect to baseline population estimates and demographic trends of the affected species. For most bat species we have no knowledge of the size of populations and their demographic trends, the degree of structuring into discrete subpopulations, and whether different subpopulations use spatially segregated migratory routes. While estimates of local population sizes within particular roosts may be feasible using traditional capture-mark-recapture (CMR) methodology or survey techniques, no reliable range-wide population estimates exist for any bat species (O’Shea & Bogan, 2003; Kunz et al., 2009). Traditional demographic approaches have limitations when applied to bats, as they are nocturnal, exhibit cryptic behavior, and are difficult to follow over time during extensive seasonal movements between summer breeding areas and overwintering sites (Cryan, 2003; Rivers, Butlin & Altringham, 2006). The tree-roosting migratory bat species that are killed in high numbers at wind turbines are especially inaccessible for traditional CMR studies, given their solitary nature and restriction to forested habitats (Kunz, 1982; Shump & Shump, 1982a; Shump & Shump, 1982b). Large-scale banding studies typically experience extremely low recapture rates (e.g., Glass, 1982; reviewed in O’Shea & Bogan, 2003), and there are serious data deficiencies with respect to sex- and age-specific survival and reproductive rates that hamper our ability to widely apply demographic models to bat populations. Given these difficulties, we require other approaches to estimating population sizes and demographic trends within migratory bat populations affected by wind power development.

Genetic approaches provide an alternative to traditional demographic methods of population estimation, and allow us to estimate the degree of population structuring, demographic trends within subpopulations, and effective population size (N_e) using data on allele frequencies or the base composition of DNA sequences. Fewer individuals need to be sampled relative to CMR approaches, and individuals need only be sampled a single time for many analyses. In addition, population parameters can be estimated directly from the observed patterns of genetic variation, and age- or sex-specific demographic information may not be required. Molecular markers can also be used to examine levels of population differentiation within a species and to geographically delimit populations or groups of populations based on the observed distribution of genetic variation (Freeland, Petersen & Kirk, 2011). Importantly, such analyses can be used to define the relevant unit for population monitoring, and highlight demographic connections among populations that may not be obvious from behavioral data alone. As mating is likely to take place during migration in bats (Dodd & Adkins, 2007; Cryan, 2008; Cryan et al., 2012; Solick et al., 2012), gene flow should occur among populations that interact during migration. Therefore it is likely that any genetically distinct populations, if they exist, will be using different migratory pathways and may be subject to different mortality rates as wind turbines are concentrated heterogeneously across the landscape. The analysis of genetic population structure is therefore highly relevant to our understanding of bat—wind turbine interactions.

While it is not possible to directly estimate adult census population size (N_c) using molecular data (although genetic markers can be used to identify individuals for traditional CMR analyses; Luikart et al., 2010), it is possible to estimate effective population size (N_e). N_e is defined as the number of individuals in an ideal Wright–Fisher population (a large, constant-sized, randomly-mating, hermaphroditic population with discrete generations) that would lose genetic variation through genetic drift at the same rate as the actual population (Crow & Denniston, 1988). It provides information on how quickly genetic variation is being lost, or relatedness is increasing, in a population of interest, and may be interpreted as an estimate of the number of individuals actually contributing genes to the next generation. The estimation of N_e has seen wide application in studies of threatened or isolated populations, as the magnitude of genetic drift, and hence loss of genetic variation, is inversely proportional to N_e (Leberg, 2005; Wang, 2005; Luikart et al., 2010). Current estimates of N_e can be used to assess the ‘genetic health’ of populations and their capability to respond to future environmental change or anthropogenic changes via selection (Frankham, Ballou & Briscoe, 2002). Estimation of N_e is also common in phylogeographic studies exploring past changes in population sizes in relation to changing climatic conditions or vicariant events in the evolutionary history of species (Avise, 2000; Russell et al., 2011), thus providing important insight into the demographic history of populations and species.

Here, we utilize genetic data from eastern red bats (Lasiurus borealis), one of the species most highly affected by wind power development in North America, to (1) evaluate patterns of population structure and whether different subpopulations use spatially segregated migratory routes, (2) estimate effective population size (N_e), and (3) assess signals of growth or decline in population size. This species was chosen because it is one of the three bat species of greatest concern with regard to the biodiversity impacts of wind energy, and has the highest fatality rate at a number of wind power installations in the eastern United States (Arnett et al., 2008). Although estimates of census population size would be preferable for understanding the size of bat populations and the potential impact of fatalities at wind power installations, N_e estimates may provide us with valuable information on the size of the evolutionarily relevant portion of the population (that portion contributing genes to the next generation). Further, regular monitoring of N_e might serve as a proxy for tracking changes in population size over time. Our study provides valuable data for understanding the population-level impacts of mortalities due to wind power for this migratory bat species by assessing whether there are discrete subpopulations that may represent independent management units and may undergo different migratory behavior, whether populations from different regions may be connected demographically, and the relative magnitude and historical population trends of the population or subpopulations we identify.

Methods

Sampling

Tissue samples from eastern red bats were collected by researchers capturing bats in the field or collecting carcasses at wind power developments. We asked researchers across the range to collect samples, but eastern red bats were not encountered in all areas due to regional differences in encounter rates. Therefore, we have the largest sample sizes per site and the greatest number of samples from the eastern portion of the species’ range. All researchers were required to have appropriate state and federal collecting permits. A small number of samples from Michigan were collected by one of the authors (MJV) under permit from the state of Michigan (Michigan Department of Natural Resources permit SC-1257) with appropriate Institutional Animal Care and Use Committee approval (Western Michigan University protocol 05-03-01).

We compiled a collection of tissue samples from known sample sites collected in the summer months (June to mid-August when bats are likely to be resident) primarily between 2000 and 2006, for the purpose of assessing levels of genetic population structure and estimating N_e (Table 1 and Table S1). We received tissue samples for 1–39 bats from any given site. We had sufficient sample size (N > 15) for each of 12 sites with which to carry out site-level population genetic analyses (Fig. 1 and Table 1). Unlike colonial bats roosting in buildings or trees where bats can be captured in numbers from a single point location during a single sampling session, tree-roosting bats such as eastern red bats are solitary. Sampling these bats therefore must involve the capture of foraging individuals and may encompass individuals from a wider area over a longer time scale. Therefore we define a ‘site’ as a collection of capture localities within a set of nearby counties within a single state or province. For six of our sites, bats were captured either within a single county or at a single capture location (AR, GA, MO, ON, TX, WV-Ma), while the other six sites consisted of individuals captured in several counties within a given state (IL, MD, MI, NC, TN, WV-Pe; the site label for the latter site represents one of the two counties included; Table S1). There were no consistent differences in diversity measures within sites or levels of differentiation between sites associated with sites containing samples from a single versus multiple counties (see ‘Results’).

Table 1:

Sites sampled and diversity statistics for 16-locus microsatellite genotypes and mitochondrial HV2 sequences.

Site labels represent two-letter state codes as in Fig. 1.

Site	State or province	N Gen	NA	HO	HE	AR	ARPriv	F IS	N Seq	NH	h	π
AR	Arkansas	39	18.25	0.84	0.88	13.14	0.50	0.044	25	21	0.987	0.016
GA	Georgia	30	16.75	0.81	0.87	13.12	1.16	0.064	17	13	0.963	0.009
IL	Illinois	26	15.31	0.80	0.87	12.88	0.56	0.084	26	22	0.985	0.013
MD	Maryland	21	13.31	0.81	0.86	12.19	0.80	0.057	15	15	1.000	0.012
MI	Michigan	17	12.69	0.82	0.88	12.69	0.84	0.073	16	16	1.000	0.013
MO	Missouri	27	16.25	0.84	0.89	13.20	0.80	0.056	34	21	0.961	0.009
NC	North Carolina	18	13.19	0.81	0.88	12.87	0.76	0.079
ON	Ontario	19	14.13	0.87	0.88	13.43	1.05	0.021	19	17	0.983	0.012
TN	Tennessee	22	14.50	0.82	0.87	12.98	0.79	0.065	26	23	0.991	0.010
TX	Texas	20	14.19	0.79	0.88	13.14	0.81	0.105	21	20	0.995	0.011
WV-Pe	West Virginia	20	13.25	0.83	0.87	12.35	0.79	0.050	19	18	0.994	0.010
WV-Ma	West Virginia	25	15.44	0.85	0.88	13.02	0.45	0.036
Overall		284	14.77	0.82	0.88	12.92	0.78	0.061	218	18.6	0.986	0.011

DOI: 10.7717/peerj.983/table-1

Notes:

N _Gen

number of individuals genotyped

N_A

number of alleles

H_O

observed heterozygosity

H_E

expected heterozygosity

A_R

allelic richness

AR_Priv

private allelic richness

F _IS

inbreeding coefficient

N _Seq

number of individuals sequenced at mitochondrial HV2 locus

N_H

number of haplotypes

h

haplotype diversity

π

nucleotide diversity

Overall values represent means for all measures except NGen and NSeq, which represent sums.

Results

N_e estimation

We used three coalescent methods to estimate N_e for eastern red bats: IMa2, Lamarc, and msvar. These methods utilize different suites of data (microsatellites only for msvar, nuclear and mitochondrial sequence data only for Lamarc, all three data types for IMa2), and therefore were expected to provide complementary estimates based on differences in the mutation rates of the markers used and differences in the underlying models assumed.

IMa2

This analysis converged on an unambiguous, unimodal posterior distribution for the single population parameter θ (=4N_eμ) for the panmictic eastern red bat population. The most probable value of θ was estimated to be 37.95 (95% CI [32.15–45.55]). We used Pesole et al.’s (1999) estimate of mammalian mitochondrial mutation rates (=2.740 × 10⁻⁸ substitutions per site per year) to calculate locus-specific mutation rates for our data. The geometric mean of these rates (=8.03 × 10⁻⁶ substitutions per locus per year =1.61 × 10⁻⁵ substitutions per locus per generation; Table S5) was used to convert coalescent-scaled estimates of θ into estimates of N_e. Our analysis thus supports an effective population size of approximately 5.91 × 10⁵ individuals (95% CI [5.00–7.09] ×10⁵; Fig. 3).

Lamarc

We used coalescent-based analyses in Lamarc to provide estimates of θ and population growth independently for the nuclear CHY and mitochondrial HV2 loci. Analyses of both markers provided unambiguous, unimodal posterior probability distributions for both parameters. Utilizing the relative mutation rates estimated from IMa2, estimates of N_e using HV2 across three runs in Lamarc were 5.18 × 10⁵ (95% CI [4.25–7.22] ×10⁵; Table 4). The estimate of N_e using CHY (males and females) was significantly larger, with a mean of 1.52 × 10⁶ (95% CI [1.05–2.18] ×10⁶; Table 4). There was a clear signal of historical population growth recovered from both loci (Table 4); however, the time scale over which this growth occurred is not estimated in the Lamarc model.

Table 4:

Estimates of θ, N_e, and population growth (g) based on Lamarc analyses.

	θ (95% CI)	Ne (95% CI)	g (95% CI)
HV2
Run 1	0.022	5.0 × 105	964.25
	(0.018, 0.031)	(4.16–7.26 × 105)	(361.03, 1007.18)
Run 2	0.024	5.52 × 105	965.75
	(0.019, 0.029)	(4.33–6.78 × 105)	(358.34, 1007.50)
Run 3	0.022	5.0 × 105	965.95
	(0.018, 0.033)	(4.25–7.61 × 105)	(382.04, 1006.35)
Overall	0.022	5.18 × 105	965.32
	(0.018, 0.031)	(4.25–7.22 × 105)	(367.14, 1007.01)
CHY
Run 1	0.048	1.54 × 106	958.85
	(0.033, 0.067)	(1.07–2.15 × 106)	(496.01, 1002.27)
Run 2	0.046	1.50 × 106	957.10
	(0.032, 0.067)	(1.03–2.17 × 106)	(486.19, 1002.01)
Run 3	0.047	1.52 × 106	952.73
	(0.033, 0.069)	(1.06–2.21 × 106)	(479.76, 1001.04)
Overall	0.047	1.52 × 106	956.23
	(0.033, 0.068)	(1.05–2.18 × 106)	(487.32, 1001.77)

DOI: 10.7717/peerj.983/table-4

msvar

Although we found considerable variation from run to run, there were some clear patterns that emerged from these analyses. Importantly, we found no consistent difference between parameter estimates from the 2002 vs. 2010 time points (Fig. 2 and Figs. S1–S2). We also found no consistent difference between independent subsamples (A vs. B, each run 2–3 times) of the full dataset (runs A1–A3 vs. B1 and B3 for 2002; runs A1–A3 vs. B1–B3 for 2010). For the current effective population size N₁, we recovered generally consistent estimates on the order of 10⁴–10⁵ (average N₁ ≈ 74,500). Estimates of ancestral effective population size N_A were less consistent among runs, but did result in estimates ranging in the same order of magnitude as N₁ (average N_A ≈ 194,300; Fig. S1). These analyses yielded differing signals of population growth vs. decline between runs (Table 5), although a majority of runs (8 of 11) support a model of population decline rather than growth. The time of this population size change (t) was also variable among runs, but generally was on the order of 10³–10⁴ years (average t ≈ 21,600 years; Table 5 and Fig. S2). While the time of population size change is difficult to pinpoint with great accuracy, these analyses clearly are not informative regarding very recent population size change.

Table 5:

Estimates of current and ancestral N_e, time of growth and population trend based on msvar analyses.

Year	Current Ne (mode ± variance)	Ancestral Ne (mode ± variance)	Time of growth (mode ± variance)	Trend
2002_A1	125,786 ± 4.9	24,191 ± 4.5	5,353 ± 1.9	Growth
2002_A2	21,120 ± 6.2	57,497 ± 3.1	6,924 ± 1.9	Decline
2002_A3	106,925 ± 3.3	22,460 ± 6.0	27,256 ± 5.8	Growth
2002_B1	137,848 ± 6.1	14,626 ± 4.6	11,710 ± 3.3	Growth
2002_B3	195,164 ± 4.4	651,754 ± 3.1	21,915 ± 1.1	Decline
2010_A1	46,279 ± 3.5	59,872 ± 2.6	88,776 ± 1.9	Decline
2010_A2	36,766 ± 5.6	427,688 ± 3.8	16,088 ± 4.3	Decline
2010_A3	24,733 ± 2.4	44,036 ± 5.5	32,866 ± 1.6	Decline
2010_B1	22,845 ± 3.8	81,332 ± 7.6	5,161 ± 4.2	Decline
2010_B2	12,670 ± 5.4	29,191 ± 3.8	9,978 ± 1.3	Decline
2010_B3	89,050 ± 10.0	724,656 ± 2.8	11,552 ± 7.2	Decline

DOI: 10.7717/peerj.983/table-5

Discussion

We observed extremely low levels of population structure and effective panmixia across the sampled sites for eastern red bats using both nuclear and mitochondrial DNA markers. Furthermore, there is no evidence for the historical use of different migratory pathways and no evidence for any barriers to gene flow among any of the sampled localities. Few geographic barriers to the movement of vagile organisms such as bats exist east of the Rocky Mountains, and therefore there are likely few impediments to the movement of individuals across the landscape. Phylogeographic studies of widespread bats and birds have shown low levels of genetic differentiation among eastern North American populations (Gibbs, Dawson & Hobson, 2000; Kimura et al., 2002; Jones et al., 2005; Turmelle, Kunz & Sorenson, 2011; Irwin, Irwin & Smith, 2011; but see Miller-Butterworth et al., 2014). When present, genetic structure in these species is often restricted to broad-scale differentiation between eastern and western populations on either side of the Rocky Mountains. In the case of eastern red bats, evidence from museum records indicates that they most likely migrate from northern parts of their range to the southeastern United States (Cryan, 2003) where they roost in trees during warmer periods and may hibernate beneath leaf litter for short durations during colder temperatures (Saugey et al., 1998; Moorman et al., 1999; Mormann & Robbins, 2007). However, there are summer resident populations in the southeastern United States that likely do not migrate, and it is possible that there is variation in migratory tendency across the range of eastern red bats, much like tricolored bats (Perimyotis subflavus; Fraser et al., 2012). Mating likely takes place before or during migration in eastern red bats (Dodd & Adkins, 2007; Cryan, 2008; Cryan et al., 2012; Solick et al., 2012), and can take place before bats hibernate or during warm periods on the wintering grounds. Thus, the potential for mating, and hence gene flow, among individuals that spent their summers in geographically disparate areas during migration or on the wintering grounds is likely very high.

In most colonial temperate bat species, females are philopatric to natal nursery colonies or undergo short dispersal distances to nearby colonies while mating takes place during swarming and/or hibernation at distant sites that act as hotspots of gene flow between bats occupying distant roosts during the summer (Kerth et al., 2003; Veith et al., 2004; Furmankiewicz & Altringham, 2007). As a consequence, levels of mitochondrial differentiation (indicative of female movements) are often quite high among summer maternity colonies while levels of nuclear differentiation (indicative of gene flow through mating) are typically low (Castella, Ruedi & Excoffier, 2001; Bilgin et al., 2008; Kerth et al., 2008; Vonhof, Strobeck & Fenton, 2008; Bryja et al., 2009; Lack, Wilkinson & van den Bussche, 2010; Turmelle, Kunz & Sorenson, 2011). Eastern red bats and other members of the genus Lasiurus roost solitarily in foliage during the summer (Shump & Shump, 1982a; Shump & Shump, 1982b), and if they exhibited philopatry it would likely occur within broader landscape units such as forest patches or stands rather than a single roost. The absence of significant mitochondrial differentiation among samples of eastern red bats suggests that females may be exhibiting high levels of dispersal, and that gene flow likely takes place via both male and female movements and mating (e.g., Russell, Medellín & McCracken, 2005; Vonhof, Strobeck & Fenton, 2008).

Before undertaking our study, we had no prior knowledge of whether the eastern red bat was divided into a series of discrete subpopulations, possibly undertaking migration along different pathways and possibly varying in size, or whether it functioned as a single, panmictic population of unknown size. Our estimates of N_e varied considerably (over an order of magnitude) among the different approaches we used, ranging from 7.45 × 10⁴ based on microsatellite genotypes only (msvar), to 1.52 × 10⁶ for sequence data only (CHY in Lamarc), with intermediate estimates of 5.18 × 10⁵ for HV2 (Lamarc) and 5.91 × 10⁵ using all markers combined (IMa2). This variation is the result of methodological differences among the approaches we used, which all utilize different aspects of the data and make varying assumptions about the underlying historical population processes that may have occurred. Further, the analyses each used different marker data, which vary in their mutation rates, and so are providing estimates across varying time scales. Nevertheless, in combination with the results of population structure analyses, our data indicate that eastern red bats form a single, large, panmictic population across their range and that minimum effective population sizes are likely in the hundreds of thousands.

The parameter most relevant to management of this species, the actual number of individuals in the population (N_c), is not obtainable from our estimates of N_e. A variety of factors may reduce N_e relative to N_c, including fluctuations in population size over time, overlapping generations, and variation among individuals in reproductive success. Attempts have been made to compare estimates of N_e to N_c, and across a wide range of organisms the average N_e/N_c ratio is 0.11–0.14 (Frankham, 1995; Palstra & Ruzzante, 2008); for mammals alone, the average ratio is 0.34 (Frankham, 1995). If we applied this latter mean ratio (0.34) to our point estimates of N_e, we would obtain N_c estimates of 2.19 × 10⁵ to 4.5 × 10⁶ individuals. However, there are a number of serious problems with the use of our coalescent estimates in this way. N_e is a theoretical concept that relates the genetic characteristics of a population to those expected of an ideal population under a Wright–Fisher model. We can evaluate N_e as a measure of the evolutionary potential of populations, but there is no clear relationship between current demography and changes in genetic variation that influence coalescent estimates of N_e. Further, there are a number of methodological concerns. First, N_e has most often been estimated for very small populations of less than 1,000 individuals, and we do not know how the N_e/N_c ratio may vary with the magnitude of N_c. Second, the majority of the ratios provided by Frankham (1995) utilize demographic, rather then genetic, estimates of N_e, and demographic estimates may differ substantially from genetic estimates even when population sizes are small (Luikart et al., 2010). Third, the majority of estimates in Frankham (1995) come from organisms with very different life histories than bats, and we do not know to what extent the N_e/N_c ratio might vary from the overall mean for bats (or most other organisms). Fourth, the calculation of N_e using coalescent-based methods requires division of estimates of θ by the mutation rate (μ) to obtain values of N_e, but mutation rates are extremely difficult to estimate and few good estimates exist for any gene (Ho et al., 2006; Montooth & Rand, 2008; Nabholz, Glémin & Galtier, 2009), much less for any bat species. As a result, any inaccuracy in the mutation rate estimate is amplified arithmetically in the subsequent calculation of N_e (Ovenden et al., 2007; Luikart et al., 2010). Therefore, applying a standard conversion to convert N_e to N_c is highly problematic, and it is best to use our estimates to indicate relative orders of magnitude of bat population sizes rather than to provide any specific population size estimates.

The potential value of our estimates of N_e is that they may be used as a baseline for future monitoring. Assuming fatality rates at wind turbines remain high and continue to grow as wind energy development continues, it is possible that regular estimates of N_e could be utilized to document population trends of affected species (Antao, Perez-Figueroa & Luikart, 2011). Regional projections of bat fatalities predict annual fatality rates numbering in the tens of thousands (Kunz et al., 2007), and the total number of fatalities is likely to continue to rise as wind power development expands. However, the loss of genetic variation from populations and declines in N_e estimates based on linkage disequilibrium are only apparent when population sizes are very small (e.g., Waples & Do, 2010), suggesting that cumulative population declines may have to be very severe before they affect genetic estimates. Had our estimates of N_e been considerably smaller, or had we detected numerous subpopulations among which gene flow was restricted, then there may have been greater potential to document population size changes using genetic approaches. Given our results supporting a large, panmictic population, simulation studies are required to assess the sensitivity of coalescent-based estimates of N_e to population decline and to assess the utility of this approach for eastern red bats.

Our genetic data indicating panmixia and a lack of evidence for the use of different migratory pathways in different parts of the range highlights the need to consider the global implications of current and future fatalities associated with wind power. Despite growing conservation concern, current monitoring of bat fatalities at wind power developments is performed on an ad-hoc, site-by-site basis and may vary tremendously in scope according to local regulations. While such monitoring can provide valuable insights leading to site-level mitigation strategies or changes in turbine placement in some cases, biologists lack the necessary broader context within which to assess the long term, population-level impacts of observed fatality rates and management strategies at specific sites. For instance, site-specific, per-turbine thresholds to limit fatalities through curtailment (reducing turbine blade speed and operating time on low-wind nights in summer and fall to decrease fatalities; Baerwald et al., 2009; Arnett et al., 2011) ignore the fact that the demographic consequences of mortality extend well beyond any particular jurisdiction. Evidence from stable isotopes indicates that bats killed at wind power developments may originate from wide geographic areas (Voigt et al., 2012; Baerwald, Patterson & Barclay, 2014), and thus mortality at any given site can impact bat populations using geographically widespread catchment areas. Given that observed bat fatality rates at wind power facilities vary considerably among sites and regions (Arnett et al., 2008), our findings underscore the need for better data integration across jurisdictions and monitoring programs to adequately assess the cumulative demographic and genetic impacts of continued fatalities.

Supplemental Information