Quaternary climatic fluctuations and resulting climatically suitable areas for Eurasian owlets

Abstract Aim The nested pattern in the geographical distribution of three Indian owlets, resulting in a gradient of endemicity, is hypothesized to be an impact of historical climate change. In current time, the Forest Owlet Athene blewitti is endemic to central India, and its range is encompassed within the ranges of the Jungle Owlet Glaucidium radiatum (distributed through South Asia) and Spotted Owlet Athene brama (distributed through Iran, South and Southeast Asia). Another phylogenetically close species, Little Owl Athene noctua, which is largely Palearctic in distribution, is hypothesized to have undergone severe range reduction during the Last Glacial Maximum, showing a postglacial expansion. The present study tests hypotheses on the possible role of Quaternary climatic fluctuations in shaping geographical ranges of owlets. Methods We used primary field observations, open access data, and climatic niche modeling to construct climatic niches of four owlets for four periods, the Last Interglacial (~120–140 Ka), Last Glacial Maximum (~22 Ka), Mid‐Holocene (~6 Ka), and Current (1960–1990). We performed climatic niche extent, breadth, and overlap analyses and tested if climatically suitable areas for owlets are nested in a relatively stable climate. Results Climatically suitable areas for all owlets examined underwent cycles of expansion and reduction or a gradual expansion or reduction since the Last Interglacial. The Indian owlets show significant climatic niche overlap in the current period. Climatically suitable areas for Little Owl shifted southwards during the Last Glacial Maximum and expanded northwards in the postglaciation period. For each owlet, the modeled climatic niches were nested in climatically stable areas. Main Conclusions The study highlights the impact of Quaternary climate change in shaping the present distribution of owlets. This is relevant to the current scenario of climate change and global warming and can help inform conservation strategies, especially for the extremely range‐restricted Forest Owlet.

The LGM characterized low average temperature, increased aridity, and a drop in sea levels (Clark & Huybers, 2009), leading to a change in climate, available land area, and climate-associated changes in vegetation (Anhuf et al., 2006;Bose et al., 2016). These changes possibly altered the ranges of many species. In the Holocene (~11.7 Ka to Present), a warmer climate than the LGM prevailed in the Northern Hemisphere but the tropics were colder than the present (Mayewski et al., 2004;Steig, 1999;Wanner et al., 2008).
Knowledge from the Quaternary period suggests that species responses to past climate change can provide crucial information on their current and future evolutionary and ecological trajectories.
In this paper, we examine the effect of the Quaternary climatic fluctuations on the climatic niche extents of owlets that show a gradient of endemicity and overlap in their current geographical distributions in parts of their ranges. Such comparative biogeography studies are scarce and have been recommended to comprehend community responses to global climate change (Berg et al., 2010).  (Koparde et al., 2018). Understanding how their ecological and evolutionary histories shaped their current distribution can provide vital information on how they respond to climatic changes, which will help plan their conservation strategies in the current scenario, especially in case of the Endangered Forest Owlet.
Divergence estimates from the phylogeny of Indian owlets indicate that the Plio-Pleistocene climate change may have played an important role in the speciation of Athene and Glaucidium owlets (Koparde et al., 2018) that could explain patterns of their current ranges. Pellegrino et al., (2014), suggest that the Little Owl survived in the European Southern Refugia (Iberian, Italian, and Balkan Peninsula) during the LGM, when much of its distributional range was covered in ice, and later expanded into its current range.
In the present study, we explore if Quaternary climatic fluctuations played a role in shaping the geographical distributions of owlets, using past-projected climatic niche models (CNMs) and examine if the suitable areas for the endemic and Endangered Forest Owlet were nested within climatically stable areas to a greater extent, as compared to the other relatively widespread species.

| Target species
The four target species are the Forest Owlet, Jungle Owlet, Spotted Owlet, and Little Owl and are presented in Figure 1 (Ali & Ripley, 1983;Rasmussen & Collar, 1998) and Forest Owlet distribution overlaps completely with Spotted Owlet and Jungle Owlet distributions in central India, they can be potential competitors in areas of sympatry (Ishtiaq, 2000;Jathar & Rahmani, 2004;Mehta, Kulkarni, Talmale, & Chandarana, 2018;Rasmussen & Ishtiaq, 1999). At a finer scale, however, their habitat associations are reported to be different. The Spotted Owlet and Little Owl are associated with open habitats and considered to be synanthropic while the Jungle Owlet is associated with dry to moist deciduous open forests and scrublands (Ali & Ripley, 1983).

| Data collection-climate data
We extracted the climate dataset available for four time periods, LIG (~120-140 Ka), LGM (~22 Ka), MDH (~6 Ka), and current  from <http://www.worldclim.org/> (Hijmans, Cameron, Parra, Jones, & Jarvis, 2005). Datasets for the LGM, MDH, and current time period were available at 2.5′ (around 5 km 2 ); this being the highest resolution for the LGM and MDH datasets. The available dataset for the LIG was 30″ (around 1 km 2 ) resolution (Otto-Bliesner et al., 2006). Therefore, we scaled the LIG dataset to 2.5′. We used the LGM and MDH dataset from the Community Climate System Model (CCSM4) (Gent et al., 2011) following Fuentes-Hurtado, Hof, and Jansson (2016). We clipped raster files of the bioclimatic variables at two extents to be used in the analysis, the Indian Subcontinent for the Indian owlets restricted to the Indian Subcontinent (5 N to 39.1 N and 55.1 E to 109.9 E) and Eurasia and parts of North Africa for Little .08 E to 134.46 E). The geographical extent under modeling is a crucial factor in determining the accuracy of species distribution models (Barve et al., 2011;VanDerWal, Shoo, Graham, & Williams, 2009). Therefore, we used two different extents to capture the predictor range better.

| Climatic niche models
We created bias files for all owlets to correct for sampling bias in modeling (Kramer-Schadt et al., 2013). We used MaxEnt v 3.4.1 (Phillips, Anderson, & Schapire, 2006) for CNMs. We performed all the Pre-and Post-MaxEnt data analyses in ArcGIS v10.1 (ESRI, 2011) and SDMToolbox (Brown, 2014) in ArcGIS. For CNMs, we followed We first performed a correlation analysis on all 19 bioclimatic predictors for the current time period to detect highly correlated (r > 0.8, r < −0.8) variables. To select the appropriate variables from pairs of highly correlated ones, all 19 variables were used in a MaxEnt run (replicate type = bootstraps, replicate runs = 50) and variables that contributed maximally in jackknifing runs were noted. For further analysis, we retained only those variables (from a correlated pair) that had high contributions in the MaxEnt output and were important considering the natural history of each owlet. This procedure has been used elsewhere (Carroll, 2010;Peterson & Robins, 2003).
The selected variables for each species are shown in Table 1. We followed this procedure to fine-tune the models and avoid overfitting (Lee-Yaw et al., 2016;Radosavljevic & Anderson, 2014). The final MaxEnt models were run with 50 bootstrap iterations. We set the regularization parameter to 1.5 to avoid overfitting of data. To determine the robustness of the model in terms of Test and Training AUC values, we randomly picked 25% of points as test points. We performed backward-time simulations by projecting CNM for the current period for each owlet at three time periods, MDH, LGM, and LIG. We used a 10th percentile logistic training presence threshold to convert continuous raster maps into binary maps to better visualize the change in the extent of climatic niche. The 10th percentile logistic threshold is a conservative estimator of predicted climatic niches and has been applied to avoid overfitting of models (Kumar & Stohlgren, 2009;Pearson, Raxworthy, Nakamura, & Townsend Peterson, 2007). The fossil data available on the study owls are scanty, mainly available for the Little Owl from Europe (Bedetti & Pavia, 2013;Mlikovsky, 2002;Pavia, Manegold, & Haarhoff, 2014); hence, validation of the past CNMs was not possible. Here, as with other CNM approaches, we assume that current species-climate relationships have been maintained in the past.

| Post-CNM analysis
We performed intersection and stability analyses on climatically suitable areas by superimposing suitable area polygons for a specific time period with another time period. We treated areas common to both the polygons (intersection analysis) as conserved areas (niche stable) and nonoverlapping areas as shift (contraction/expansion/displacement) in the climatic niche. We performed niche overlap analysis to compute I statistic (Warren, Glor, & Turelli, 2008) and niche breadth analysis to compute the B2 statistic (uncertainty index) (Nakazato, Warren, & Moyle, 2010) using ENMTools v1.4.4 (Warren, Glor, & Turelli, 2010) to explore niche overlap across time periods and species. The niche overlap index (I statistic) varies between 0 (no overlap) and 1 (complete overlap). The higher values in case of B2 index represent broader niche. Finally, we created a climatic heterogeneity layer for each time period using SDMToolbox in ArcGIS. The climatic heterogeneity information is in percentage (0-100), 0 signifying highly homogenous (stable) climate and a value of 100 indicating highly heterogeneous climate. Climatic heterogeneity information was extracted based on 1,000 random points generated for each suitable area polygon to test if the predicted climatically suitable areas of owlets fall in areas with higher climatic stability, expecting that suitable areas for the endemic Forest Owlet will be nested in climatically stable zones as compared to the widespread owlets. We performed one-way ANOVA test on the climatic heterogeneity data to test for variation within a species across time periods and between species for each time period.

| Climatically suitable areas and niche breadth of owlets
The models had a low false positive rate (model summaries in TA B L E 1 Summary of the best-fit climatic niche models (CNMs) for the current time period. Variables used are the same as for the past-projections underwent a cyclic reduction and expansion throughout the four time periods (Figure 2a-d, Table 2, Supporting Information Appendix S1: Figure S1 Table 2, Supporting Information Appendix S1: Figures S1.4. and S1.5.).
We detected a southward shift during the LGM and northward expansion post-LGM in climatically suitable areas for the Little Owl ( Figure 2m, Table 2, Supporting Information Appendix S1: Figure S1.5.).  (Figures 3 and 4, Supporting Information Appendix S1: Figures S1.6. and S1.7., Appendix S2: Table S2.2.). In the current time period, all owlets occupied areas with higher climatic stability than in the past.

| Quaternary climatic fluctuations and climatically suitable areas for owlets
Assuming that the owlets have tracked the climatically suitable areas predicted by our models, we detected variable responses of the four owlets to the Quaternary climatic fluctuations. The climatically suitable areas for the currently severely range-restricted Endangered Forest Owlet showed distinct cycles of reduction and expansion; whereas suitable areas for other currently widespread owlets showed either an overall progressive expansion or reduction ( Figure 2, Table 2). The change in climatically suitable areas for the Indian owlets might be a function of climate and climate-mediated change in habitat, prey, and interactions among these species. Currently, the Forest Owlet is sympatric with Jungle Owlet and Spotted Owlet whereas the Jungle Owlet and Spotted Owlet overlap in occurrence in parts of their overall distributional ranges.
Interspecific interactions and resource use when factored into niche models could improve predictions (Araújo & Luoto, 2007;Wisz et al., 2013). Incorporating data on recent as well as paleoclimatic fluctuations generated for time periods not covered in this study are recommended, to obtain a more comprehensive picture of species responses to climate change.

| Climatically suitable areas for owlets and climatic heterogeneity
Following the hypothesis of endemic species occupying climatically stable areas (Dynesius & Jansson, 2000;Jansson, 2003), over a gradient of range extents, we expected the following pattern: the

| Caveats and conclusion
There is no definitive way to empirically validate the past distribution models constructed for the study species and hence interpretations are presumptive. Nonavailability of fossil data for focal species makes it difficult to examine the accuracy of the past distribution models. The two central assumptions of the study are (a) climatically suitable area is a proxy for the geographical area occupied by a species and (b) the current species-climate relationships have been maintained in the past. We recommend validating the past distribution models with the help of fossil occurrence data whenever available. Apart from these major issues, the quality and accuracy of the predictor dataset and its projections are of concern (reviewed in Nogués-Bravo, 2009).
Our results suggest that Quaternary climatic fluctuations might have played a significant role in shaping the present distribution of owlets. Such information can help in deciphering the biogeography of species, with varying habitat associations albeit with overlapping geographical distributions. Future research in this area should focus on more substantial datasets incorporating information on interspecific interactions and regional climate to understand the effect of Holocene climatic fluctuations on species. Pradesh, Gujarat, and Chhattisgarh) for providing permits to collect data on owls. Pankaj Koparde (P.K.) would like to thank the Ecology Lab members at IISER Tirupati for discussions.

CO N FLI C T O F I NTE R E S T
None declared.

AUTH O R CO NTR I B UTI O N S
P.K., V.V.R., and S.M conceived the ideas; P.K. and P.M. collected the data; P.K. analyzed the data; and P.K., V.V.R., and S.M. led the writing.

DATA ACCE SS I B I LIT Y
The climatic niche model output raster files can be accessed from