Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

The influence of spatial and temporal climate variation on species' distributions, phenologies and interactions Kharouba, Heather Marie 2013

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata

Download

Media
24-ubc_2013_fall_kharouba_heather.pdf [ 5.72MB ]
Metadata
JSON: 24-1.0058479.json
JSON-LD: 24-1.0058479-ld.json
RDF/XML (Pretty): 24-1.0058479-rdf.xml
RDF/JSON: 24-1.0058479-rdf.json
Turtle: 24-1.0058479-turtle.txt
N-Triples: 24-1.0058479-rdf-ntriples.txt
Original Record: 24-1.0058479-source.json
Full Text
24-1.0058479-fulltext.txt
Citation
24-1.0058479.ris

Full Text

The influence of spatial and temporalclimate variation on species?distributions, phenologies andinteractionsbyHeather Marie KharoubaB.Sc., University of Ottawa 2005M.Sc., University of Ottawa 2007A THESIS SUBMITTED IN PARTIAL FULFILLMENT OFTHE REQUIREMENTS FOR THE DEGREE OFDOCTOR OF PHILOSOPHYinThe Faculty of Graduate Studies(Zoology)THE UNIVERSITY OF BRITISH COLUMBIA(Vancouver)August 2013c? Heather Marie Kharouba 2013AbstractWith climate change, species are shifting their distributions polewards andupwards, and advancing their phenologies. However, there is substantialinterspecific variation in these responses and ecologists are having difficultyexplaining why. Understanding this variation is critical as it is likely tolead to widespread consequences for trophic interactions and ecological com-munities. In this thesis, I test several hypotheses concerning the causesand ecological consequences of interspecific and intertaxonomic variation inclimate-distribution and phenology-temperature relationships.First, I tested and found support for the hypothesis that species fromdifferent taxonomic groups vary in the strength of climate-distribution rela-tionships, likely because of differences in their life history strategies. How-ever, my results suggest that dispersal ability is unlikely to be the key traitaffecting species? geographic distributions at broad scales.Across broad-scales, I found that butterfly and plant phenologies werestrongly affected by temperature, suggesting that phenological shifts in re-sponse to future climate change are likely to be widespread. Flight seasontiming of early-season butterfly species and those with lower dispersal abil-ity was more sensitive to temperature than later-season species and thosewith greater dispersal ability, suggesting that ecological traits can accountiifor some of the interspecific variation in phenological sensitivity to temper-ature. Differences in phenological sensitivities of butterflies and plants totemperature imply that shifts in phenological synchrony are likely and couldbe substantial for interacting species, potentially resulting in important fit-ness consequences.Finally, experimentally warming the egg masses and larvae of the westerntent caterpillar (Malacosoma californicum pluviale) placed on the branchesof its host plant, red alder (Alnus rubra), in the field led to opposing directand indirect effects on larval development. Warming significantly advancedlarval but not leaf emergence, which initially prolonged larval development.However, once leaves were present, warming accelerated larval development,resulting in no overall effects on larval development.Taken together, this thesis demonstrates that to understand the fullimplications of climate change for species and communities, accounting forspecies? life history strategies and interactions will be essential. However,without more quantitative estimates of the fitness consequences of shifts inphenological synchrony, this understanding will be limited.iiiPrefaceChapter 2 is published in Ecography. I designed the study and compiledthe data in collaboration with Jenny McCune. I analyzed the data andwrote the first version of the manuscript. Jenny McCune also contributedto interpreting the results and writing the manuscript. W. Thuiller and B.Huntley contributed data and edited the manuscript.Kharouba, H.M., McCune, J.L., Thuiller, W., and Huntley, B. 2013. Doecological differences between taxonomic groups influence the relationshipbetween species? distributions and climate? A global meta-analysis usingspecies distribution models. Ecography 36:657-664.Chapter 3 - Mark Vellend and I conceived the idea for this chapter. Icompiled and analyzed the data, and wrote the manuscript. Mark Vellendcontributed to the interpretation of the results and editing the manuscript.Sebastien Rioux Paquette constructed the phylogenetic tree and with JeremyKerr, edited the manuscript.Chapter 4 - Mark Vellend and I conceived the idea for this chapter. Icompiled and analyzed the data, and wrote the manuscript. Mark Vellendcontributed to the interpretation of the results and editing the manuscript.Chapter 5 - Judy Myers suggested this experimental system and MarkVellend and I designed the study. I collected and analyzed the data, andwrote the manuscript. Rana Sarfraz conducted the bioassay experiment.ivMark Vellend and Judy Myers also helped with interpreting the results andediting the manuscript.vTable of ContentsAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiPreface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ivTable of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . viList of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiList of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xivAcknowledgements . . . . . . . . . . . . . . . . . . . . . . . . . . . xvi1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.1 Structure of this thesis . . . . . . . . . . . . . . . . . . . . . 72 Taxonomic comparison in the strength of species? climate-distribution relationships . . . . . . . . . . . . . . . . . . . . . 112.1 Synopsis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112.2 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . 122.3 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 152.3.1 Data compilation . . . . . . . . . . . . . . . . . . . . 152.3.2 Statistical analysis . . . . . . . . . . . . . . . . . . . . 18vi2.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 222.4.1 Relationship between discrimination ability and taxo-nomic group . . . . . . . . . . . . . . . . . . . . . . . 222.4.2 Relationship between discrimination ability, covari-ates and taxonomic group . . . . . . . . . . . . . . . 242.5 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . 252.5.1 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . 303 Explaining variation in butterfly phenological sensitivity totemperature . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 323.1 Synopsis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 323.2 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . 333.3 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 373.3.1 Butterfly data . . . . . . . . . . . . . . . . . . . . . . 373.3.2 Climate data . . . . . . . . . . . . . . . . . . . . . . . 383.3.3 Analyses . . . . . . . . . . . . . . . . . . . . . . . . . 383.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 453.4.1 Testing phenological sensitivity to spatial and tempo-ral temperature . . . . . . . . . . . . . . . . . . . . . 453.4.2 Using ecological traits to predict phenological sensi-tivity to temperature . . . . . . . . . . . . . . . . . . 463.4.3 Evaluating trends in phenology over time . . . . . . . 463.5 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . 513.5.1 Conclusions . . . . . . . . . . . . . . . . . . . . . . . 59vii4 Butterfly-plant comparison in the sensitivity of phenologyto temperature . . . . . . . . . . . . . . . . . . . . . . . . . . . 604.1 Synopsis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 604.2 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . 614.3 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 664.3.1 Butterfly data . . . . . . . . . . . . . . . . . . . . . . 664.3.2 Plant data . . . . . . . . . . . . . . . . . . . . . . . . 684.3.3 Climate data . . . . . . . . . . . . . . . . . . . . . . . 694.3.4 Statistical analysis . . . . . . . . . . . . . . . . . . . . 704.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 764.4.1 Testing phenological sensitivity to spatial and tempo-ral temperature . . . . . . . . . . . . . . . . . . . . . 764.4.2 Evaluating phenological shifts over time . . . . . . . . 774.4.3 Butterfly-plant comparison . . . . . . . . . . . . . . . 774.5 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . 874.5.1 Conclusions . . . . . . . . . . . . . . . . . . . . . . . 915 The effects of experimental warming on the timing of aplant-insect interaction . . . . . . . . . . . . . . . . . . . . . . 925.1 Synopsis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 925.2 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . 935.3 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 975.3.1 Study system . . . . . . . . . . . . . . . . . . . . . . 975.3.2 Study design overview . . . . . . . . . . . . . . . . . . 985.3.3 Warming experiments . . . . . . . . . . . . . . . . . . 100viii5.3.4 Leaf quality experiment . . . . . . . . . . . . . . . . . 1025.3.5 Statistical analyses . . . . . . . . . . . . . . . . . . . 1035.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1055.4.1 Treatment temperature effects . . . . . . . . . . . . . 1055.4.2 Testing the effect of branch-level warming on budburst 1055.4.3 Treatment differences in phenology and fitness proxies 1065.4.4 Within-treatment relationships . . . . . . . . . . . . . 1075.4.5 Leaf quality experiment . . . . . . . . . . . . . . . . . 1095.5 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1145.5.1 Conclusions . . . . . . . . . . . . . . . . . . . . . . . 1216 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1226.0.2 General conclusions . . . . . . . . . . . . . . . . . . . 1256.0.3 Limitations and future directions . . . . . . . . . . . 128Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 134Appendices . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 172Appendix A Methodological attributes of studies . . . . . . . 172Appendix B Taxonomic attributes of studies . . . . . . . . . . 175Appendix C Attributes of studies with dispersal distances . 177Appendix D Collinearity correlations . . . . . . . . . . . . . . . 179Appendix E Analysis of deviance results for ?large studies? 180ixAppendix F Description of data preparation . . . . . . . . . . 182Appendix G Description of main statistical analyses . . . . . 190Appendix H Phylogenetic data and analyses . . . . . . . . . . 193Appendix I Additional analyses . . . . . . . . . . . . . . . . . . 211Appendix J Additional methods, statistical analyses and re-sults . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 222Appendix K Description of main statistical analyses . . . . . 224xList of Tables2.1 Analysis of deviance table for the relationship between dis-crimination ability, covariates and taxonomic group . . . . . . 213.1 The relationship between species? ecological traits and phe-nological sensitivity to temperature . . . . . . . . . . . . . . . 484.1 Differences in the phenological sensitivity of butterflies andplants, and butterfly specialists and generalists, to spatialand temporal temperature. . . . . . . . . . . . . . . . . . . . 804.2 Comparison between butterflies and plants in their phenolog-ical sensitivity to temporal and spatial temperature . . . . . . 825.1 Differences in phenology and insect fitness proxies betweenwarming and control treatments at Totem Field and UBCFarm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1105.2 Within-treatment relationships of insect fitness proxies withmean daily maximum temperature and degree of mismatchat both sites . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111xiA.1 Methodological attributes used to build SDMs for each studyused in the analysis. . . . . . . . . . . . . . . . . . . . . . . . 173B.1 Taxonomic attributes of studies . . . . . . . . . . . . . . . . . 176C.1 Attributes of studies that contained dispersal distances . . . . 178D.1 Collinearity correlations between all continuous covariates . . 179E.1 Analysis of deviance table for the relationship between modelaccuracy, covariates and taxonomic group when studies thatcontributed more than half of the total number of species inone taxonomic group were removed. . . . . . . . . . . . . . . 181F.1 Relationship between seasonal temperature and timing offlight season across species . . . . . . . . . . . . . . . . . . . . 186F.2 Predicted relationships between relevant ecological traits andthe sensitivity of flight season timing to temperature. . . . . . 187H.1 Model comparison between phylogenetic and non-phylogenetic gls models for the relationship betweenspecies? traits and phenological sensitivity to overall, spatialand temporal temperature for only those species withnegative sensitivities . . . . . . . . . . . . . . . . . . . . . . . 198H.2 Model comparison between phylogenetic and non-phylogenetic gls models for phenological sensitivity tooverall temperature with each species? traits. . . . . . . . . . 201xiiH.3 Model comparison between phylogenetic and non-phylogenetic gls models for phenological sensitivity totemporal temperature with each species? traits. . . . . . . . . 204H.4 Model comparison between phylogenetic and non-phylogenetic gls models for phenological sensitivity tospatial temperature with each species? traits. . . . . . . . . . 207I.1 Relationship between the timing of flowering season and tem-perature . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 213I.2 Differences in plant phenological sensitivity to temperaturebased on two different metrics of temperature. . . . . . . . . . 215I.3 Relationship between the timing of flight season and temper-ature. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 216I.4 A comparison of mixed-effects models that evaluate the phe-nological sensitivity of butterflies and plants to temperaturebased on the inclusion of a smoothing term to assess the po-tential for non-linearity. . . . . . . . . . . . . . . . . . . . . . 221K.1 Model comparison of different temperature variables and fit-ness proxies at each site . . . . . . . . . . . . . . . . . . . . . 228xiiiList of Figures2.1 Taxonomic differences in discrimination ability (AUC) acrossall studies. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 232.2 Taxonomic differences in maximum dispersal distances. . . . 263.1 Distribution of slope values across 204 butterfly species forphenological sensitivity to temporal and spatial temperature,temporal phenological shifts and the degree of temperaturechange. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 523.2 Relationship between phenological sensitivity to spatial andtemporal temperature across species . . . . . . . . . . . . . . 533.3 The sensitivity of flight season timing to spatial variation intemperature a species with low and high mobility and acrossspecies as a function of the timing of flight season and mobilityindex . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 544.1 Plant and butterfly phenological sensitivity to temporal andspatial temperature. . . . . . . . . . . . . . . . . . . . . . . . 854.2 Butterfly-plant differences in phenological sensitivity to tem-poral temperature taking into account pair-wise interactions . 86xiv5.1 Conceptual diagram of direct and indirect effects of warming,through shifts in phenological mismatch and changes in leafquality, on insect fitness proxies . . . . . . . . . . . . . . . . . 1135.2 Comparison of the timing of budburst on branches that werecut to ones that remained attached for both control andwarmed treatments. . . . . . . . . . . . . . . . . . . . . . . . 1155.3 Development of western tent caterpillar larvae at Totem Fieldand UBC Farm represented by the cumulative appearance ofeach instar across trees for families from warmed and controltreatments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1165.4 Bivariate relationships of development time with mismatchand mean daily maximum temperature for Totem Field andUBC Farm. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117F.1 Pattern of butterfly collecting efforts in database . . . . . . . 189H.1 Phylogenetic tree of all Lepidoptera species included in thisstudy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 210xvAcknowledgementsDuring my time at UBC, I have been lucky to have three mentors, MarkVellend, Judy Myers and Chris Harley, who have provided valuable guidance.First, I thank my primary advisor, Mark Vellend, for his unfalteringacademic and emotional support. He allowed me the freedom to pursue myown research interests, while at the same time providing valuable advice andguidance. I found his energy and positivity motivating in challenging times.By being a tough editor, he taught me a lot about effectively communicatingmy ideas. Finally, he was always around for a chat- in his office or on skype.I thank Judy Myers for her generosity, patience with me in the field,research and academic advice. I am also extremely grateful for all the timeshe spent with me in the field and in the lab- helping to set up experimentsand collect data, regardless of weather.Finally, thanks to Chris Harley for graciously accepting a terrestrialecologist into his lab. Thanks for providing helpful research and academicadvice. He was always generous with his time and would always make sureI was doing okay.Thanks to all my thesis committee members, Mark Vellend, Judy Myers,Chris Harley, Dolph Schluter and Jessica Hellmann for contributing theirunique perspectives to my work. My thesis is much stronger as a result.xviAdditionally, I would like to thank Dolph for statistical advice and helpwith R programming and Jessica for her expert advice on butterfly biologyand field-related issues.I thank my extended labs (Vellend, Myers and Harley?s Herd of HQP) foremotional and academic support over the years. In particular, I am grate-ful to Anne Bjorkman, Emily Drummond, Tanis Gieselman, Jamie Leathem,Patrick Lilley, Janet Maclean, Jenny McCune and Andrea Stephens. Thanksalso to Nathan Kraft and Lizzie Wolkovich for valuable research and aca-demic advice. My work has also benefited from fruitful discussions in FLO-RUM and the Biodiversity Discussion Group over the years.Finally, thanks to my family and non-academic friends for all theirsupport and encouragement over the years. Lastly, I would like to thankmy husband Paul for all his love, support, and patience during my longhours in the field and at the computer, and for moving to Vancouver withme. I could not have finished this PhD without him!For Chapter 2, I thank Catherine Graham and three anonymous review-ers who gave valuable comments on the manuscript. Thanks to O. Allouche,M. Bakkenes, J. Elith, R. Kadmon, J. McPherson, R. Pearson, J. Pyry, andA. Tsoar for contributing supplemental data. M. Vellend, D. Srivastava, andJ. Goheen provided helpful comments on previous versions of the manuscriptand D. Schluter, B. Bolker, and D. Fournier gave statistical advice.For Chapter 3, I am grateful to Guy Baillargeon, Cris Guppy, NorbertKondla, Ross Layberry and Karen Needham for providing data. Ross Lay-berry helped with georeferencing and Laura Tremblay-Boyer helped with Rxviicode.For Chapter 4, I am grateful to Krysta Banack, Steve Healy, JessicaHenry, Ariel Kettle, Christine Lester, Russell Prentice, Jordan Pryce, CharlStafleu, Connie Tam, and Isabelle Teasdale for help in compiling and geo-referencing the herbarium specimens, and Amber Saundry for access to theUBC Herbarium database. Thanks to Ross Layberry for help with georef-erencing the butterfly collection records and Laura Tremblay-Boyer helpedwith R code.For Chapter 5, I thank Krysta Banack, Aliya Kanani, Christine Lesterand Charl Stafleu for help conducting the research, Abe Miller-Rushing forproviding supplemental data, Jenny S. Corey for technical advice, and AbeMiller-Rushing and Chris Harley for helpful feedback on this manuscript.Thanks also to UBC Farm for site access.This work was funded by a NSERC Canada Graduate Scholarship,NSERC BRITE Graduate Fellowship through the Beaty Biodiversity Re-search Centre, Department of Zoology Graduate Fellowship and a IODECanada War Memorial Scholarship.xviiiChapter 1IntroductionEcological responses to climate changeClimatic changes are expected to have a major impact on global biodi-versity (Bellard et al., 2012; Sala et al., 2000; Thomas et al., 2004). Thepast century has seen a nearly 1?C rise in global average temperature, withup to 7?C of average warming predicted by 2100 (IPCC, 2007). The biolog-ical impacts of climate change are likely to be multifaceted, involving be-havioural change, evolutionary change, and local extinction. Consequently,these changes are expected to accelerate the ongoing mass extinction (Car-penter et al., 2008; Maclean & Wilson, 2011; Thomas et al., 2004; Thuilleret al., 2005). As a result, there is a need to translate our understandingof ecology into predictions of the ecological consequences of these climaticchanges.Climate change has already had substantial effects on biodiversity.Since species? distributions and phenologies are strongly linked to climate,the two most commonly observed and reported responses to climate changethus far are shifts in species? ranges and phenologies (reviewed in Parmesan,2006; Parmesan & Yohe, 2003; Walther et al., 2002). In general, species areshifting their distributions to higher elevations and latitudes in response to1warming temperatures (e.g. Chen et al., 2011; Hickling et al., 2006; Parme-san, 2006). However, some species have been unable to track the changingclimate and have experienced range contractions or local extinctions (e.g.McLaughlin et al., 2002; Moritz et al., 2008; Sinervo et al., 2010) and othershave shifted their ranges counter to expectations (i.e. southward and lowerin elevation Hickling et al., 2006). The majority of species are also advanc-ing their phenologies, the timing of life history events (e.g. Bradley et al.,1999; Menzel et al., 2006; Parmesan, 2007; Roy & Sparks, 2000). However,some species have not shifted or even delayed their phenology, in somecases with negative fitness consequences (Lane et al., 2012; Willis et al.,2008). Therefore, while species are generally shifting polewards and springis arriving earlier, there appears to be significant taxonomic, spatial andtemporal variability in the magnitude of this change (Diez et al., 2012; Fitter& Fitter, 2002; Hickling et al., 2006; Primack et al., 2009) and it remainsdifficult for ecologists to explain this variation (Buckley & Kingsolver, 2012).Understanding variation in ecological responses to climate changeInterspecific variation in recent range shifts is likely partly a function ofthe variation in the degree to which species? distributions are determined byclimate. Indeed, the extent to which current climate predicts species? distri-butions differs considerably (e.g. Elith et al., 2006, Chapter 2). Moreover,which environmental factors are important at limiting species? range limitsand by which mechanisms are still uncertain (Chuine, 2010). While climatehas long been considered the most important factor in determining species?distribution limits at broad spatial scales (e.g. Gaston, 2003; Good, 1931;2Merriam, 1894), other factors may affect whether climate plays a primaryor secondary role. For example, other abiotic (e.g. soil properties) or bioticfactors (e.g. species interactions) may prevent a species from persisting evenwhere the climate is suitable (e.g. Hampe, 2004; Marsico & Hellmann, 2009;Pearson & Dawson, 2003). Disentangling the direct and indirect effects ofclimatic factors on species? distributions still poses considerable challenges(Sexton et al., 2009).Similarly, variation in species? recent phenological responses to climatechange is most likely related to the sensitivity of phenology to temperature.Temperature is thought to be the most consistent and dominant controllerof spring insect and plant phenology (Bale et al., 2002; Hodgson et al., 2011;Pau et al., 2011; Wolkovich et al., 2012). However, depending on the typeof study (e.g. genetic vs. physiological), species and region, other cues (e.g.photoperiod, frost, precipitation) can have greater relative importance (Diezet al., 2012; Hodgson et al., 2011; Pau et al., 2011). Recent studies havedocumented substantial variation among species in phenological sensitivityto temperature, even for single phenological phases (Diez et al., 2012; Hodg-son et al., 2011; Primack et al., 2009; Wolkovich et al., 2012), yet little isknown about the causes of this variation.A trait-based approach, the characterization of organisms in terms oftheir multiple biological attributes such as physiological, morphological orlife history traits (Webb et al., 2010), has recently been applied to attemptto explain the variation in distributional and phenological shifts due toclimate change, with some notable successes (Altermatt, 2010; Angert et al.,2011; Buckley & Kingsolver, 2012; Diamond et al., 2011). Since traits can3often be related directly to the environmental conditions experienced byan individual, they provide a bridge between processes at the communitylevel and global change predictions (Pau et al., 2011; Suding et al., 2008;Webb et al., 2010). Despite theoretical support for the effects of species?traits on variation in the rate of spatial spread of populations (Clark, 1998)and species? phenology (e.g. Bolmgren & Cowan, 2008; Bolmgren et al.,2003), the full potential of traits to accurately predict variation in the rateand extent of species? responses to climate change remains unclear (Angertet al., 2011; Buckley & Kingsolver, 2012).Consequences of variation in ecological responses to climate changeInterspecific variation in climate change responses has widespread con-sequences for trophic interactions, ecological communities, and ecosystemservices (Diez et al., 2012; Pau et al., 2011). Species interact in complexfood webs and species within or between the different trophic levels can re-spond to climate change in different ways (e.g. Gilman et al., 2010; Schweigeret al., 2008; Voigt et al., 2003). Differences in the magnitude and direction ofthese responses can decouple species? interactions, both spatially and tem-porally. Spatial mismatches, through differential shifting in the geographicranges of species, and temporal mismatches, through differential changesin the phenology of interacting species, are likely. For example, terrestrialplants might not track climate as fast as their insect herbivores: earlieremergence and range shifts in some plants has occurred at a slower ratethan those of their insect herbivores (Gilman et al., 2010; Harrington et al.,1999; Hellmann, 2002; Schweiger et al., 2008). Such mismatches could in4turn alter the dynamics, functioning and composition of many ecosystems(e.g. Edwards & Richardson, 2004; Liu et al., 2011; Winder & Schindler,2004). Therefore, some of the most profound effects of climate change arepredicted to be driven by changes in the timing of biotic interactions (Yang& Rudolf, 2010).Differences among species in distributional and phenological responsesto changing climates also have the potential to affect patterns of speciescoexistence and community assembly (Lavergne et al., 2010; Wolkovich& Cleland, 2010). For example, differences in flowering responses amongnative and exotic species in a community can affect interspecific compe-tition and play an important role in the success of invasive plant species(Wolkovich & Cleland, 2010). The heterogeneity in species responses toclimate change is also likely to change the composition of local communitiesby altering dominances of species and leading to the formation of non-analogcommunities, where extant species co-occur in historically unprecedentedcombinations (Hobbs et al., 2009; Walther, 2010; Webb III, 1986; Williams& Jackson, 2007). For example, the present species? composition of com-munities at higher altitudes on sub-Antarctic Marion Island is not simplyan analogue of past community composition of the native vascular flora atlower altitudes but rather constitutes a new combination of extant species(Roux & McGeoch, 2008). Consequently, in addition to the direct effectson individuals and populations, climate change can also lead to changes inother members of a community via species interactions (i.e. indirect effects).Indirect vs. direct effects of climate change5Predicting the impacts of climate change routinely focuses on the directeffects on individual species; yet increasing evidence suggests that speciesinteractions can strongly influence species? responses to changing climates(Barton et al., 2009; Davis et al., 1998; Harley, 2011; Harmon et al., 2009;Pateman et al., 2012; Suttle et al., 2007). For example, spatial mismatchesbetween a butterfly and its host plant can limit the poleward movement ofbutterfly populations (Hellmann et al., 2012; Pelini et al., 2009). The mostcommon proximate cause of local extinctions related to climate change, whenidentified, is changes in species? interactions rather than physiological lim-itations (Cahill et al., 2013). A review of >600 studies found that climatechange influenced virtually every type of species interaction (Tylianakiset al., 2008), yet most studies focused on predicting the responses of speciesto climate change ignore biotic interactions (Gilman et al., 2010; Harring-ton et al., 1999; Walther, 2010). Indirect effects can be stronger than thedirect effects of temperature and the two may interact, either amplifyingor countering one another (e.g. Barton et al., 2009; Boggs & Inouye, 2012;Harley, 2011; Suttle et al., 2007). Understanding the full implications ofclimate change for a focal species will thus depend on incorporating the di-rect effects of climate change on the species, as well as the sensitivity of theinteracting species and the strength of the interaction.Insect herbivores and pollinators are especially likely to be affected byclimatic changes both directly, given their physiological and behavioral re-quirements, and also indirectly, through interactions with their host plants(Bale et al., 2002; Forrest & Thomson, 2011). Insects are poikilothermic, sotheir development, activity, reproductive physiology, population dynamics,6and dispersal are heavily influenced by environmental conditions, such astemperature and moisture (Bale et al., 2002). Host plants impact insectherbivores directly through phenological and nutritional conditions (Feeny,1970; Hunter & Price, 1992), and any variation in resources may influence thedynamics and abundances of insect herbivores (Sanders & Gordon, 2003).Climate change-driven shifts in fitness, abundance and geographic distribu-tions of insect herbivores and pollinators will have implications for plant-insect interactions and host plant species. In particular, rates of herbivory(DeLucia et al., 2012), insect pest and pathogen dynamics (Bentz et al.,2010; Dukes et al., 2009) and pollination are likely to be affected (Heglandet al., 2009; Memmott et al., 2007). Therefore, understanding how climatechange influences insect herbivores and pollinators is critical to comprehend-ing some of the broader ecological consequences of climate change.1.1 Structure of this thesisThe goal of my thesis is to test several hypotheses related to the causesand ecological consequences of interspecific and intertaxonomic variation inclimate-distribution and phenology-temperature relationships. This thesisfocuses on the direct effects of temperature on insect and plant phenology,rather than other abiotic effects (e.g. photoperiod), since temperature isthought to be the most consistent and dominant controller of the phenol-ogy of these taxa (Bale et al., 2002; Hodgson et al., 2011; Pau et al., 2011;Wolkovich et al., 2012). While the direct influence of climate on species? dis-tributions and phenologies has been well established, it remains uncertain as7to why there is so much variation in species? distributional and phenologicalshifts in response to recent climate change. Moreover, the consequences ofinterspecific variation in phenological shifts for species interactions and thefitness of focal species have rarely been considered.To address these questions, I combine a broad-scale approach that usesa meta-analysis and collection records, with a fine-scale approach that usesmanipulative experiments. These approaches have complementary strengthsand weaknesses, such that their combined use is valuable to advancingclimate-change research. Whereas broad-scale approaches allow for the de-tection of general patterns and trends, they are limited in elucidating themechanisms that drive these patterns, the main strength of experiments.However, experiments are often imperfect simulations of the real world.0In Chapter 2, I test the hypothesis that species from different taxonomicgroups vary in the strength of their climate-distribution relationships be-cause of differences in life history strategies, in particular dispersal ability.Previous studies suggest that species that produce many propagules thattravel long distances are more likely to be able to cross unsuitable habi-tat, and thus should be more likely to be found everywhere the climate issuitable. Therefore, I predicted that taxonomic groups with lower disper-sal ability would have weaker species? climate-distribution relationships. Iused a meta-analysis approach and combined the discrimination ability met-rics from species distribution models using only climatic variables. I firstdetermined whether species varied predictably in their climate-distributionrelationships based on taxonomic affinities. I then compiled dispersal dis-tances for a subset of these species to determine whether dispersal ability8directly influenced the strength of species? climate-distribution relationships.In Chapter 3, I test the usefulness of museum butterfly collectionrecords, an under-utilized data source in the phenology literature, to de-tect broad-scale phenology-temperature relationships for butterfly speciesacross Canada. Since temperature has been shown to have different effectson phenology in space (warm vs. cold sites) and time (warm vs. cold years),I compare this relationship across both space and time. Finally, given thatspecies vary widely in their ecological and life history strategies, I test theability of ecological traits, such as dispersal ability, to predict interspecificvariation in the phenology-temperature relationship.Following the finding that collection records can be used to detectphenology-temperature relationships across broad-scales (Chapter 3), I com-pile herbarium records and compare the phenological sensitivity of butterfliesand their host plants to temperature across the province of British Columbiain Chapter 4. Given discrepancies among studies in the magnitude of species?recent phenological responses to climate change across taxonomic groups, itis difficult to make general predictions about the magnitude and directionof shifts in phenological synchrony (relative timing of life cycle events ofclosely interacting species) under climate change. However, differences inphenological sensitivity to temperature between plants and butterflies canbe used to quantify potential shifts in butterfly-plant phenological synchronyover time. Since specialists are thought to be more sensitive to climatechange (Gilman et al., 2010), I then evaluate the influence of butterfly eco-logical specialization on butterfly phenological sensitivities to temperatureand butterfly-plant differences in sensitivity to determine whether shifts in9synchrony are likely to differ between butterfly specialists and generalists.In Chapter 5, I test the influence of warming on the fitness consequencesof altered timing of species interactions. Few studies have experimentallytested the influence of warming on the timing of species interactions and inparticular, the potential for different direct and indirect (via phenologicalchange) effects of temperature on each species. As an experimental system,I use an insect herbivore (the western tent caterpillar (Malacosoma califor-nicum pluviale)) and its host plant (red alder (Alnus rubra)) subjected toa field warming treatment. I compare the direct and indirect effects (i.e.phenological mismatch) of warming on insect fitness proxies (larval devel-opment time, larval weight and family survival). A bioassay in a controlledenvironment was also conducted to determine the effects of warming-drivenchanges in leaf quality, the other primary indirect effect in this system, oninsect fitness.Finally, in Chapter 6, I summarize conclusions from the previous chap-ters and make general conclusions that can be drawn from this study aboutunderstanding species distributions and phenology in the context of climatechange. I also examine the limitations of this study and provide recommen-dations for future research priorities.10Chapter 2Do ecological differencesbetween taxonomic groupsinfluence the relationshipbetween species?distributions and climate? Aglobal meta-analysis usingspecies distribution models2.1 SynopsisUnderstanding whether and how ecological traits affect species? geographicdistributions is a fundamental issue that bridges ecology and biogeography.While climate is thought to be the major determinant of species? distri-11butions, there is considerable variation in the strength of species? climate-distribution relationships. One potential explanation is that species withrelatively low dispersal ability cannot reach all geographic areas where cli-matic conditions are suitable. We tested the hypothesis that species fromdifferent taxonomic groups varied in their climate-distribution relationshipsbecause of differences in life history strategies, in particular dispersal ability.We conducted a meta-analysis by combining the discrimination ability (AUCvalues) from 4317 species distribution models (SDMs) using fit as an indi-cation of the strength of the species? climate-distribution relationship. Wefound significant differences in the strength of species? climate-distributionrelationships across taxonomic groups, however we did not find support forthe dispersal hypothesis. Our results suggest that relevant ecological traitvariation among broad taxonomic groups may be related to differences inspecies? climate-distribution relationships but which ecological traits are im-portant remains unclear.2.2 IntroductionUnderstanding whether and how ecological traits affect species? geographicdistributions is a fundamental issue that bridges ecology and biogeography(Brown, 1995; Wiens, 2011). This issue has become even more relevant asecologists and biogeographers struggle to understand the variation in speciesresponses to climatic change. For example, recent studies have examined therelationship between species? ecological traits, such as dispersal ability andecological generalization, and changes in their distributions and phenology12with recent climatic changes (Angert et al., 2011; Diamond et al., 2011).Identifying characteristics of organisms that determine their sensitivity toenvironmental change is crucial to ecological forecasting and conservationplanning.Central to this work is the theory of the niche: the set of abiotic andbiotic conditions within which a species can persist (Hutchinson, 1957). Aspecies? distribution is limited to geographic areas where all these conditionsmeet the species niche requirements. At broad spatial scales, climate haslong been considered the most important factor in determining species? dis-tribution limits (e.g. Gaston, 2003; Good, 1931; Merriam, 1894). However,there seems to be considerable variation in the degree to which species? dis-tributions are predicted by climate. There are three potential reasons forthis variation. First, other abiotic or biotic factors may prevent a speciesfrom persisting even where the climate is suitable (Luoto et al., 2007). Al-ternatively, regions of suitable climate may be separated by areas that arenot suitable which the species does not have sufficient dispersal ability tocross (Blach-Overgaard et al., 2010; Graham et al., 2010). Finally, if thespecies is relatively new and/or the climate has only recently become suit-able, the species may not have had enough time to reach all suitable areas(Blach-Overgaard et al., 2010; Paul et al., 2009).Dispersal ability is thought by some to determine how closely a species?current distribution matches the geographic distribution where all abioticand biotic conditions meet its niche requirements. Species that producemany propagules that travel long distances are more likely to be able tocross any unsuitable habitat, and thus should be more likely to be found ev-13erywhere the climate is suitable. Therefore, dispersal ability may determinethe strength of the species? climate-distribution relationship. Indeed, somestudies have found evidence that dispersal ability can strongly affect species?distributions (e.g. Po?yry et al., 2008; Thuiller et al., 2004b). However, oth-ers suggest that the dispersal of individuals happens over such small timescales relative to the formation of species? geographic distributions that ithas little importance (Lester et al., 2007).Many have hypothesized that species in different taxonomic groupsshould vary in their climate-distribution relationship because of their dif-ferent life history strategies, in particular dispersal ability (e.g. Araujo &Pearson, 2005; Wisz et al., 2008). The fit of species distribution models(SDMs) has often been used to test this hypothesis (e.g. Araujo & Pear-son, 2005; Tsoar et al., 2007). SDMs use various statistical techniques todescribe the relationship between observed environmental variables, suchas mean annual temperature, and the recorded spatial occurrence (pres-ence/absence) of a species (see e.g. Guisan & Zimmermann, 2000). Theability of an SDM based only on climatic factors to predict the presence orabsence of a species can be considered an indication of the strength of thespecies? climate-distribution relationship: the greater the success of a SDMat predicting the species? presence/absence in a given location, the strongerthe correlation between climatic variables and the presence/absence of thespecies. Some studies have found species? climate-distribution relationshipdifferences between taxonomic groups (Araujo & Pearson, 2005; Tsoar et al.,2007), whereas others have not (Pearce & Ferrier, 2000; Wisz et al., 2008).It is unclear whether these varying results are due to the different geographic14regions, groupings of species, or modeling techniques of each study. Despitethe availability of SDMs for thousands of species, a comprehensive compar-ison of the fit of SDMs between different taxonomic groups has not beenmade.Here, we tested the hypothesis that taxonomic groups varied in thestrength of their species? climate-distribution relationships. We predictedthat taxonomic groups with lower dispersal ability would have weakerspecies? climate-distribution relationships. We used a meta-analysis ap-proach and combined the discrimination ability metrics that were reportedfrom 4317 SDMs in twenty studies using only climatic variables to deter-mine whether species varied predictably in their climate-distribution rela-tionships based on taxonomic affinities. We also compiled dispersal distancesfor a subset of these species to determine whether dispersal ability directlyinfluenced the strength of species? climate-distribution relationships. To fa-cilitate a quantitative comparison we used a standardized discriminationability measure and accounted statistically for methodological differencesamong studies.2.3 Methods2.3.1 Data compilationWe conducted a literature search using Web of Science for studies (publishedbefore March 2009) that reported statistical measures of goodness-of-fit forSDMs constructed for individual species based on climatic variables only. Wesearched for studies using the terms ?ecological niche model? and ?climat*?,15?species distribution model? and ?climat*?, and ?climate envelope model?and ?climat*?. Studies were excluded if: (1) one or more non-climaticvariables, such as soil fertility, land use or land cover, were included inthe SDM; (2) model fit was measured only qualitatively or not reported; or(3) model fit was reported only as averages across species. In cases wheremodel fit was not reported for all individual species modeled, we requestedthese data from the authors. Due to the small number of studies modelingaquatic species, we limited our analysis to terrestrial species.We needed a metric of model fit that was comparable across studies.We found AUC (area under a receiver operating characteristic curve) tobe the most common metric (other metrics: Cohen?s kappa, sensitivity,specificity, range filling rates), therefore our analysis was limited to stud-ies that reported AUC. AUC measures the ability of a SDM to discriminatesites where a species is present from sites where it is absent, rather thangoodness-of-fit per se. It considers the relationship between false-positivesand true-positives and ranges from zero to one, where perfect discriminationgives a value of one (Fielding & Bell, 1997). Hereafter, we use the term SDM?fit? to indicate ?discrimination ability? as measured by AUC. When studiesreported AUC for both training and test data, test AUC values were used.Although this metric has been criticised (e.g. Lobo et al., 2008), it was theonly measure in common across most of the studies.Some species? distributions were modeled several times, either by thesame study (using multiple modeling techniques (n=9) or resolutions (n=1))or by several studies (most such species were modeled by only two studies).In all cases, we randomly selected one SDM per species and used the as-16sociated AUC value and methodology. This produced a dataset of 4317species and their SDMs from twenty studies (See Appendix A,B). Thesestudies modeled species in Europe (10 studies, 2301 spp.), North America(2 studies, 67 spp.), South America (2 studies, 32 species) and Africa (6studies, 1917 spp.). We classified each species into one of five broad taxo-nomic groups: mammals (483 spp.), butterflies (116 spp.), herptiles (reptilesand amphibians; 114 spp.), birds (2099 spp.), and plants (1505 spp.).SDM fit can be affected by the type of model used (e.g. Elith et al.,2006), the number of climatic variables used (e.g. Pearce & Ferrier, 2000),the resolution or grain size used (e.g Guisan et al., 2007), the total extentover which the species? range was modeled (e.g. Luoto et al., 2005), andlatitude (Brown et al., 1996; Luoto et al., 2005). Therefore, for each SDM wenoted the modeling technique, number of distinct climatic variables used inthe model, resolution (km2), total spatial extent (km2) and average absolutelatitude and then included these as covariates in our statistical analysis.Another factor which may lead to differences in SDM fit between speciesis prevalence (McPherson et al., 2004; Santika, 2011), the number of gridcells from which a species is recorded as present expressed as a proportion ofthe total number of grid cells from which data are available. We were ableto obtain prevalence values for almost all of the SDMs (n=4089), allowingus to explore any effects of prevalence on SDM fit.Finally, we scanned the literature to find dispersal distances for as manyof our species as possible to assess whether there were significant differencesin measured dispersal ability among our taxonomic groups. True disper-sal distances are very difficult to measure due to phenomena such as very17rare long-distance dispersal events. Therefore, we used the directly mea-sured ability of an organism or its propagules to move (i.e. its mobility) asan estimate of a species? dispersal distance. We considered both maximumand mean measured dispersal distances but excluded migratory distances tostandardize measures of dispersal distances across taxonomic groups. Wheremore than one distance was reported per species or study we used the meanof mean distances, and the maximum of maximum distances. We foundmean dispersal distances for 241 species for which we also had AUC values(birds=103, butterflies=22, mammals=22, plants=94). For maximum dis-persal distance, we found 105 species that also had AUC values (birds=27,butterflies=18, mammals=30, plants=18). For further details, see AppendixC.2.3.2 Statistical analysisThere were two parts to the analysis. The first was to determine whetherthere were any significant differences in SDM fit between taxonomic groupsand whether those differences were robust to potential confounding factors(covariates). The second was to explore the relationship between SDMfit, taxonomic group and the other covariates. We used generalized linearmixed-effects models (GLMM, glmmadmb function in the ?glmmADMB?package (Skaug et al., 2012) in R (Team, 2012) with a Beta error distributionwith AUC as our response variable and ?study? as a random factor. AUCvalues of exactly one, which are not allowed with the beta distribution,were converted to 0.99 instead (eight significant digits were used to ensurea unique value and to match the maximum precision of the data, n=117).18To allow for model estimation, we collapsed the six rarest modeling typesinto one category to reduce the number of types (from 18 to 12; these sixtechniques were used for only 0.35% of all SDMs). We took the logarithmof spatial extent to improve normality (except in the collinearity test), butall other covariates were used without transformation. Taxonomic groupand model type were categorical, and all other covariates were continuous.Relationship between discrimination ability and taxonomic groupTo test whether taxonomic group explained significantly more deviancein AUC than expected at random, we compared a model with only an in-tercept to a model with only taxonomic group. We then tested whetherdifferences in discrimination ability across taxonomic groups explained sig-nificant additional deviance after accounting for the combined effect of thedifferences in the methodological approach of studies (i.e. the covariates:model type, resolution, number of climatic variables, spatial extent and lat-itude). For all model comparisons, we used a likelihood ratio test. We alsocalculated AIC for all models to evaluate the relative effects of individualcovariates.We first inspected bivariate plots of all continuous covariates beforeconstructing pairwise correlations to identify potential problems withmulti-collinearity among covariates (Appendix D). Latitude was highlycorrelated with spatial extent and resolution (Spearman?s r = - 0.903,-0.589 respectively, n = 4317, Appendix D) and explained less deviancein AUC than spatial extent or resolution (Table 2.1), therefore the ?fullmodel? included taxonomic group, model type, spatial extent and number19of climatic variables. We considered the effect of ?study? by including itas a random factor and by testing the influence of individual studies thatcontributed more than half of the total number of species in one taxonomicgroup (?large studies?) by comparing results obtained with and withouteach of these studies (Araujo & Pearson, 2005; Huntley et al., 2006; Luotoet al., 2005).20Table 2.1: Analysis of deviance table for the relationship between discrimination ability, covariates and taxonomicgroup. Presented are the differences in degrees of freedom, AIC and deviance between full and reduced models aswell as the associated p value. Models are compared for all species (n=4317) and for the subset of species withprevalence values (n=4089). Depending on the model comparison and term of interest, the full model includes allother covariates (number of variables, log(spatial extent), model type, resolution and taxonomic group). ??? meansthat no solution was found and ??? means a model solution could only be found if number of climatic variables wasnot included.Model for comparison Data Model terms ?d.f. ?AIC ?Deviance pJust intercept All species Intercept+ taxonomic group 4 38.98 46.98 <0.0001+ model type 10 100.58 120.58 <0.0001+ log(spatial extent) 1 2.58 4.58 0.03235+ resolution 1 1.38 3.38 0.0660+ number of climatic variables 1 1.20 0.8 0.3711+ latitude 1 0.58 1.42 0.2334Subset Intercept+ prevalence 1 335.36 337.36 <0.0001Full model All species Full model+ taxonomic group 4 38.64 46.64 <0.0001+ model type 10 101.52 120.14 <0.0001+ log(spatial extent) 1 1.12 3.12 0.0773+ resolution 1 -1.38 0.62 0.431+ number of climatic variables? NA NA NA NASubset Full model+ prevalence? 1 445.62 447.62 <0.000121Relationship between SDM fit, covariates and taxonomic groupWe tested whether individual covariates (including prevalence) explainedsignificantly more deviance in AUC than under random expectation andafter accounting for all other covariates (including taxonomic group) bycomparing each model to a reduced one. Finally, to test whether therewere significant differences in dispersal distance (both mean and maximum)across taxonomic groups, we used a Kruskal-Wallis rank sum test. We thentested whether dispersal distance explained significantly more deviance inAUC by comparing a model with and without dispersal distance. Dispersaldistance was log-transformed to improve normality. Lastly, to test for thepossibility that an interaction between dispersal distance and taxonomicgroup explained deviance in AUC, we compared a model with and withoutthis two-way interaction. All statistical analyses were performed using R2.14.1 (Team, 2012).2.4 Results2.4.1 Relationship between discrimination ability andtaxonomic groupMean AUC across all species was 0.941 (?0.00104 SE, n=4317). Birds hadthe highest mean AUC (0.954 ?0.00145 SE, n=2099) and butterflies had thelowest mean AUC (0.856 ?0.0114 SE, n=116; Figure 2.1a). However, theranking and pair-wise comparison of taxonomic groups changed dependingon which ?large study? was removed (Figure 2.1).Taxonomic group explained significant deviance in AUC (LRT7,3=46.98,22Figure 2.1: Taxonomic differences in discrimination ability (AUC) acrossall studies (based on 4317 species from twenty published studies (numberof species: birds n=2099; herptiles n=114; butterflies n=116; mammalsn=483; plants n=1505)) (a), without Huntley et al. 2006 (based on 2860species from nineteen published studies (number of species: birds n=642;herptiles n=114; butterflies n=116; mammals n=483; plants n=1505)) (b),without Araujo et al. 2005 (based on 2539 species from nineteen pub-lished studies (number of species: birds n=1942; herptiles n=11; butter-flies n=116; mammals n=331; plants n=139)) (c), and without Luoto etal. 2005 (based on 4238 species from nineteen published studies (number ofspecies: birds n=2099; herptiles n=114; butterflies n=37; mammals n=483;plants n=1505)) (d). Taxonomic groups represented are: ?BIRD?= birds,?HER?= herptiles, ?INV?= butterflies, ?MAM?= mammals, ?P?= plants.Taxonomic groups with different letters above them are significantly differ-ent according to pair-wise comparisons. Outliers were removed to improvevisual contrasts between taxonomic groups.23p<0.0001; Table 2.1), even after accounting for all covariates(LRT20,16=46.64, p<0.0001; Table 2.1). The effect of taxonomic group wasalso robust to the exclusion of each of the ?large studies? (AppendixE).2.4.2 Relationship between discrimination ability,covariates and taxonomic groupSDM model type explained significant deviance in AUC (LRT3,13=120.58,p<0.0001; Table 2.1), even after accounting for all the other covariates(LRT20,10=120.14, p<0.0001; Table 2.1). For the subset of species forwhich we had prevalence data, prevalence also explained significant deviancein AUC after accounting for all covariates (including taxonomic group;LRT12,11=447.62, p<0.0001; Table 2.1). SDMs with greater prevalence hadlower AUC (Spearman?s r= -0.4937).In our subset of species with dispersal distances, mean dispersal distancewas greatest for mammals (175 km) while birds had the greatest maximumdispersal distance (1305 km; Figure 2.2). Butterflies had the shortest meanand maximum dispersal distance (0.441 km and 2.25 km, respectively; Fig-ure 2.2). The ranking of groups closely matched the ranking of groups ofthe entire dataset in terms of AUC for both dispersal measures (Figure2.1a, Figure 2.2). There was also a significant difference between taxo-nomic groups in dispersal distance (mean: df=3, ?2=181.006, p<0.0001;max: df=4, ?2=291.557, p<0.0001). Taxonomic group explained significantdeviance in AUC (mean: LRT6,3=10.386, p=0.01555; max: LRT7,3=13.022,p=0.01117). However, dispersal distance did not explain significant deviancein AUC (mean: LRT4,3=2.068, p=0.1504; max: LRT4,3=0.144, p=0.7043).24There was no significant interaction between taxonomic group and dis-persal distance (mean: LRT10,7=4.508, p=0.2116; max: LRT12,8=4.506,p=0.3418).2.5 DiscussionWe found support for taxonomic differences in SDM fit suggesting a rolefor ecological traits in affecting species? geographic distributions at broadscales. However, prevalence and methodological issues, such as model type,also influenced SDM fit. Indeed, both factors have been shown previouslyto influence SDM fit (e.g. Elith et al., 2006; Santika, 2011). We also foundthat ?large studies? influenced the relationship among taxonomic groups andAUC, for example the taxonomic group with the highest mean AUC variedwith the subset of species considered (Figure 2.1). Therefore, species? tax-onomic affinities, prevalence and methodological issues, such as the modeltype, are all important in influencing species? climate-distribution relation-ships as measured by SDMs.There are a number of potential explanations for the difference in thestrength of species? climate-distribution relationships between taxonomicgroups. First, taxonomic differences may reflect differences in dispersal abil-ity among groups. Certainly, we found differences in measured dispersaldistances between broad taxonomic groups that were consistent with thedispersal hypothesis (Figure 2.1a, 2.2). However, there were inconsistenciesin the ranking and pair-wise comparisons of taxonomic groups in SDM fitdepending on the subset of species considered (Figure 2.1). Moreover, there25Figure 2.2: Taxonomic differences in log (base 10) maximum dispersaldistances (km) for 105 species (birds=27, butterflies=18, mammals=30,plants=18).26was no significant relationship between AUC and dispersal distance. There-fore, our results indicate that greater dispersal ability, at least in terms ofmeasurable differences in mobility, may not result in stronger overall species?climate-distribution relationships at broad scales. However, dispersal dis-tance is inherently difficult to measure and our estimate of dispersal abilitymay not have been the most appropriate for all species. For example, wedid not take into account migratory or rare long-distance dispersal events.Consequently, we may have underestimated the role of dispersal ability forcertain species.Alternatively, dispersal may not be an important trait in determiningspecies? climate-distribution relationships. The majority of species had lowprevalence (77% species had <0.1 prevalence) and species with lower preva-lence were more likely to have higher AUC values. If these low preva-lence species are mainly specialists (i.e. restricted range endemics) thatare adapted to uncommon climatic conditions found in small, contiguousareas, they could have strong climate-distribution relationships regardlessof dispersal ability.Third, other life history traits, for example, body size, generation time ordiet breadth, may influence the strength of species? climate-distribution re-lationships between taxonomic groups. However, determining their relativeimportance may be difficult across the broad taxonomic groups considered.Lower-order taxonomic groups, or functional groups of species within oracross taxonomic groups, might be more effective in dividing species ac-cording to relevant traits. Nevertheless, while some recent studies divid-ing species into finer taxonomic or functional divisions have found signifi-27cant differences in species? climate-distribution relationships (e.g. Syphard& Franklin, 2010), others have not (e.g. Huntley et al., 2004).On the other hand, taxonomic differences in SDM fit may be a functionof the sample unbalance (across studies and taxonomic groups; Appendix B)and the high average discrimination ability. Both of these factors could re-flect issues related to fitting, testing and publishing SDMs. SDMs have beencriticized for not using independent data to test their models (e.g. Hampe,2004; Segurado et al., 2006). Without independent test occurrence points,well-fitting models could reflect spatial autocorrelation between training andtesting points rather than relationships between species? presence/absenceand climatic variables. Moreover, SDMs may be overfitted by fitting com-plex response curves and re-fitting models until a high AUC is achieved(Araujo et al., 2005a; Guisan & Thuiller, 2005). We also suggest that therecould be a ?file-drawer? problem, whereby species that do not achieve a highenough AUC value based on the literature standard (Swets, 1988) are notpublished. In particular, when the objective of fitting the SDM is to predictspecies? potential distribution shifts under various climate change scenarios,authors (rightly) do not use SDMs with very low discrimination ability. Forexample, of the 453 species that Huntley et al. (2008) modeled, 13 nativespecies that did not yield ?useful? models (sensu Swets, 1988) were excludedfrom the synthesis. Taken together, these issues could inflate AUC valuesand reduce overall variation, making it difficult to detect the true relation-ship between taxonomic groups. While we acknowledge these limitations ofSDMs, to our knowledge, there are no other comparable published metrics toevaluate individual species? climate-distribution relationships at such large28scales. Moreover, SDMs are still being used to better understand the re-lationship between species? distributions and climate (e.g. Blach-Overgaardet al., 2010; Graham et al., 2010).Lastly, because SDMs are fitted to species? current distributions theyreflect both direct and indirect influences of climate on those distributions.Non-climatic factors that limit a species to certain broad areas (such as bioticinteractions or other abiotic factors) are generally modulated by climaticconditions. For example, since its introduction to Hawaii, avian malarianow restricts native bird species to higher elevations, where temperaturehalts development of the malaria pathogen inside its mosquito vector (vanRiper III et al., 1986). Differences among taxonomic groups in the ability ofclimate to directly limit species? distributions thus cannot be revealed by ourdata, given that the SDMs we used cannot differentiate direct from indirectclimatic effects. However, we have no a priori reason to expect cases whereclimate acts principally indirectly to occur more frequently in one taxonomicgroup than another. In addition, even if a species? distribution is indirectlylimited by climate due to the climatic tolerances of a competitor, predator,or disease, at broad scales, climate is still the ultimate determinant of thespecies? distribution.There are a number of steps to be taken in the future to clarify how eco-logical traits influence species? climate-distribution relationships. Firstly,more SDMs are needed for some taxonomic groups, particularly inverte-brates and herptiles. Secondly, we should strive to eliminate issues relatedto species distribution modeling by using spatially/temporally independenttraining and test datasets where possible (e.g.Beerling et al. 1995, Randin29et al. 2006). Third, analyzing SDM prediction errors might help to shedlight on the mechanism driving the variation in species? climate-distributionrelationships, especially in cases of poor fit (e.g. Hanspach et al., 2011). Forexample, SDMs with more false negatives overall than false positives couldsuggest that source-sink dynamics are important: even where conditions arenot favourable, individuals may still persist owing to a rescue effect, or tem-poral variation in conditions (Gaston, 2003; Pulliam, 2000). Alternatively,models with greater rates of false positives might suggest that dispersallimitation or interspecific interactions, such as competition, are limiting aspecies? distribution (Graham et al., 2010; Pulliam, 2000). Finally, exploringspatial variation in model behaviour, for example testing model performancein climatically heterogeneous regions or through patterns of spatial predic-tion errors (Hanspach et al., 2011), could also improve our understanding ofmodel performance and thus species? climate-distribution relationships.2.5.1 ConclusionWe found a statistically significant effect of membership in broad taxonomicgroups on SDM fit even after accounting for methodological issues, sug-gesting a role for ecological traits in determining the strength of species?climate-distribution relationships. However, the study itself, the model typeused to build the SDM and species? prevalence all had significant effects ondiscrimination ability. Our results did not the support the hypothesis thatdispersal ability affects the strength of species? climate-distribution relation-ships. However, more work is needed to determine which ecological traitsare important in determining the strength of this relationship, and at what30spatial scale and taxonomic level they are manifested.31Chapter 3Predicting the sensitivity ofbutterfly phenology totemperature over the pastcentury using collectionrecords and ecological traits3.1 SynopsisStudies to date have documented substantial variation among species in thedegree to which phenology responds to temperature and changes over time,but we have a limited understanding of the causes of such variation. Here,we use a spatially and temporally extensive dataset to evaluate the useful-ness of collection records in detecting broad-scale phenology-temperaturerelationships and to test for systematic differences in the sensitivity of phe-nology to temperature (days/?C) of Canadian butterfly species according torelevant ecological traits. We showed that the timing of flight season pre-32dictably responded to temperature both across space (variation in averagetemperature from site to site in Canada) and across time (variation fromyear to year within each individual site). This reveals that collection records,a relatively unexploited resource, can be used to quantify broad-scale rela-tionships between species? phenology and temperature. The timing of theflight season of earlier fliers and less mobile species was more sensitive totemperature than later fliers and more mobile species, suggesting that eco-logical traits can account for some of the interspecific variation in species?phenological sensitivity to temperature. Finally, we found that phenologicalsensitivity to temperature differed across time and space implying that bothdimensions of temperature will be needed to translate species? phenologicalsensitivity to temperature into accurate predictions of species? future pheno-logical shifts. Given the widespread temperature sensitivity of flight seasontiming, we can expect long-term temporal shifts with increased warming(?2.4 days/?C (0.18SE)) for many if not most butterfly species.3.2 IntroductionThe annual timing of vegetative and reproductive stages has been frequentlyobserved to shift in response to recent climate change (e.g. Menzel et al.,2006; Parmesan, 2007; Thackeray et al., 2010). While many of these phe-nological events now occur earlier due to warmer temperatures, the direc-tion and magnitude of these responses varies considerably (e.g. Both et al.,2009; Primack et al., 2009; Thackeray et al., 2010). Several hypotheses mayexplain interspecific variation in recent phenological shifts (Hodgson et al.332011;Diez et al. 2012). Most broadly, this variation could be a function of thedegree of temperature change experienced by the species and/or the strengthof the relationship between species phenology and temperature (i.e. theirsensitivity). Recent studies have documented substantial variation amongspecies in phenological sensitivity to temperature, even for single pheno-logical phases (e.g. first flowering; Diez et al., 2012; Hodgson et al., 2011;Wolkovich et al., 2012), yet little is known about the causes of this variation.The sensitivity of a species? phenology to temperature may be influencedby the relative importance of other environmental cues (e.g. precipitation,photoperiod, resource availability) and ecological traits linked to phenology.Species vary widely in their ecological and life history strategies, althoughlittle is known about the degree to which such traits might permit gen-eral predictions about the responsiveness of species? phenology to changesin temperature. For example, dispersal ability (ability of an organism orits propagules to move among areas of suitable habitat) determines the po-tential to escape adverse consequences of temperature changes (Berg et al.,2010; Watkinson & Gill, 2002), such that the phenology of species withgreater dispersal ability might be less sensitive to temperature. While eco-logical traits have recently been used to explore interspecific variation inrecent climate change-related shifts in phenology (e.g. Altermatt, 2010; Dia-mond et al., 2011), few have used ecological traits to predict the temperaturesensitivity of phenology (Huey et al., 2002; Hurlbert & Liang, 2012; Mous-sus et al., 2011). Here we test for systematic differences in the sensitivityof butterfly species? phenology (specifically the timing of their flight season)to temperature across Canada in relation to relevant ecological strategies.34Once quantified, phenological sensitivity to temperature can be used topredict how species? phenology will shift over time given continuing climatechange (Diez et al., 2012; Hodgson et al., 2011). This sensitivity to tem-perature can be quantified by relating phenology and temperature acrossspace (cold vs. warm sites) or across time at particular sites (cold vs. warmyears). In the context of climate change, species? phenological shifts throughtime are of particular interest, especially given potential impacts on speciesinteractions (e.g. Liu et al., 2011; Post et al., 2008; Tylianakis et al., 2008).However, recent studies show that temperature can have different effects onphenology across space vs. over time (e.g. Doi & Takahashi, 2008; Forkneret al., 2008; Hodgson et al., 2011), likely because temperature and day length(the two main cues of spring phenology) both vary geographically, but onlytemperature varies from year to year. Moreover, the relative importance oftemperature is thought to vary with latitude, which may lead to variationamong populations in their phenological sensitivity to climate fluctuations(Doi & Takahashi, 2008; Pau et al., 2011; Primack et al., 2009). Relativelyfew studies have considered both components of temperature sensitivity inphenology (Doi & Takahashi, 2008; Ellwood et al., 2012; Hodgson et al.,2011; Rock et al., 1993). Isolating the temporal component of temperaturesensitivity from the spatial component should allow better predictions offuture phenological shifts.Understanding interspecific variation in phenological sensitivity to tem-perature across both space and time requires long-term and broad-spatialscale data, yet such data can be difficult to acquire. Detecting phenology-temperature relationships in short-term or local datasets can be challenging35given inter-annual variations in weather or spatial variability in climate, re-spectively (Brown et al., 2011; Diez et al., 2012; Parmesan et al., 2011; Rob-birt et al., 2011). An underutilized source of data with immense potential isfound in museum collections, where millions of dated and spatially referencedspecimens are available for species across the globe and over time periodsof decades to centuries. For plants, herbarium records have recently beeneffectively used to document phenological changes over upwards of a centuryin plant communities (e.g. Gallagher et al., 2009; Lavoie & Lachance, 2006;Primack et al., 2004). However, relatively few studies have used collectionrecords to document similar phenological changes in animal taxa (but seeBartomeus et al., 2011; Polgar et al., 2013) and none at a near-continentalscale. While collection records have their limitations (e.g. non-systematiccollecting, relatively less information in individual locations), they have thepotential to quantify species? phenological responses to recent environmentalchanges given their extensive coverage (Lavoie, 2013; Vellend et al., 2013).In this study, we used the Canadian National Collection of Butter-flies (Layberry et al., 1998), a spatially and temporally extensive dataset(?91,500 georeferenced records for 297 species for the past 139 years), topursue three main objectives: (i) to test whether phenological sensitivityto temperature can be detected using collection records; and if so, (ii) tocompare this sensitivity across space (variation in average temperature fromsite to site in Canada) and time (variation from year to year within eachindividual site), and (iii) to test the ability of ecological traits, such as dis-persal ability, to predict interspecific variation in the sensitivity of phenologyto temperature. We also evaluated whether the sensitivity of phenology to36temperature has led to detectable, directional trends in phenology over time.We focused on the timing of the flight season as the phenophase of interest,for which median collection date (for each species at each collection site) isa rough proxy. First collection date was not used as its estimation is biasedby the intensity of collection efforts, which varied considerably across sitesand years.3.3 Methods3.3.1 Butterfly dataOur main data source was the Canadian National Collection of Butterfliesdatabase (updated as of January 2011; Layberry et al., 1998), which in-cludes ? 91,500 georeferenced (>80% of records to within 1 km) collectionrecords for 297 species dating from 1873 to present. Each collection recordincludes a specimen preserved in one of 40 Canadian natural history collec-tions: specimens were collected and identified initially by lepidopterists andre-verified by Lepidopteran systematists (see Layberry et al., 1998). We sup-plemented the database with additional data for British Columbia, Alberta,Yukon and Ontario from the Spencer Entomological Collection (Universityof British Columbia) and the personal and professional collections of Cana-dian butterfly experts Syd Cannings, Cris Guppy, Ross Layberry, NorbertKondla and Jon Sheppard (all pers. comm.). Supplemental records withoutassociated geographic coordinates were georeferenced by R. Layberry usingtheir locality descriptions, GPS software (QuoVadis, http://www.quovadis-gps.de), Google Earth and Google Maps. Only locations accurate to within372km were used. Nomenclature was standardized based on Pelham (2011).Details on the combined database are presented in Appendix F.To isolate the effect of spring temperature on the same years adult flightseason, we excluded all non-resident species in Canada (migratory, rarestrays etc.). We also restricted analyses to species for which there wereat least 10 records spread across a range of greater than 30 years. There-fore, our analysis included ?48,000 georeferenced records for 204 species forthe past 139 years. Species varied in the total number of records (11-1475),total number of years with data (4-113) and range of years (30-135) withdata.3.3.2 Climate dataDaily temperature data was extracted from the National Cli-mate Data and Information Archive (Environment Canada;http://climate.weatheroffice.gc.ca) for weather stations across Canadausing only quality-controlled data (i.e. the most reliable data available).For each butterfly collection record, weather data from the closest weatherstation within 10km in any direction was taken. We used mean maximumdaily temperature for April 1 to June 30 as our estimate of temperature asthis was the time period that best predicted the timing of flight season (seeAppendix F for details on selection of seasonal block).3.3.3 AnalysesThe day of year of collection records was used to estimate the timing offlight season for each species-site-year combination. We restricted analyses38to the 204 species for which there were at least 10 records spread across arange of ?30 years. The analysis was divided into three sections. First,we tested the sensitivity of species? phenology (i.e. flight season timing)to temporal and spatial dimensions of temperature. Next, we examinedwhether ecological traits (Table F.2) could predict interspecific variation inphenological sensitivity to spatial, temporal, and overall (based on actualsite-year temperature assigned to each specimen) temperature. Finally, wetested for trends in flight season timing shifts over time. All statisticalanalyses (see Appendix G for details) were performed using R 2.14.1 (Team,2012).As an estimate of flight season timing and to avoid pseudo-replicationwithin year, we calculated the median collection date for each species ineach location for each year. We assumed that the specimens collectedrepresent unbiased samples from the adults flying in a given location andyear. While there are substantially more collection records in recent years(since 1970; Figure F.1), this would only influence our estimate of shiftsin timing of flight seasons if people collected systematically earlier inthe year over the past forty years. However, there is no reason to thinkthat collectors would have systematically changed their sampling strategy,especially across such broad scales. If anything, we might expect anunderrepresentation of specimens collected during peak flight in especiallyearly or late years, thus making our tests conservative.i) Testing phenological sensitivity to spatial and temporal temperatureOur objective was to quantify responses in the timing of flight season to39spatial (variation in average temperature from site to site in Canada) andtemporal (variation from year to year within each individual site) dimensionsof temperature. To do so, we constructed a mixed-effects model for eachspecies with day of year as a function of two temperature variables (spatialand temporal dimensions of temperature), and the year and weather stationassociated with the specimen were included as random effects. The spatialdimension of temperature was characterized by calculating the mean tem-perature across all years of data available for each weather station (i.e. site)associated with a collection record. The temporal dimension of temperaturewas characterized by calculating the difference between the temperature forthe site and year of a given specimen and the mean temperature for that site,effectively estimating an inter-annual temperature differential at that site.The two regression coefficients from this model were used to define pheno-logical sensitivity to temperature with units of days/?C (hereafter referredto as ?phenological sensitivity?); hereafter we use the terms ?phenologicalsensitivity to spatial temperature? and ?phenological sensitivity to temporaltemperature? to refer to each coefficient. Therefore, we define sensitivity asa measure of the potential of a species? phenology to respond to tempera-ture changes (a slope in days/?C) rather than as a measure of goodness offit (e.g., R2) or statistical significance of the relationship between phenologyand temperature. We tested for an interaction between the two tempera-ture variables, potential nonlinear relationships in phenological sensitivityand the presence of temporal autocorrelation in the mixed-effects model butfound no evidence of any of these potential issues for the majority of species(>50%; see Appendix G for details).40We tested the prediction that average phenological sensitivity to temper-ature across species is negative (i.e. warmer temperatures leads to earlierflight seasons) for most species using one-tailed one-sample t-tests, sepa-rately for spatial and temporal temperature. Where distributions of species?phenological sensitivities did not meet assumptions of normality, we reportmedian values.We then compared phenological sensitivity to spatial and temporaltemperature using two approaches. First, we tested whether there was asignificant difference in the mean slope of these two sensitivities using apaired t-test. We then tested whether these sensitivities were correlatedacross species using the Pearson correlation coefficient.ii) Using ecological traits to predict phenological sensitivity to temperatureWe tested the influence of individual ecological traits on phenologicalsensitivity to spatial and temporal temperature. We first asked whethertraits could predict phenological sensitivity - regardless of direction (i.e.warmer temperatures lead to earlier or later flight seasons) - by calculatingthe absolute value of phenological sensitivity. Second, we restricted thesame analysis using only those species where warmer temperatures lead toearlier flight seasons (i.e. negative sensitivity; 82-89% of species dependingon whether spatial or temporal temperature was considered). Since resultsdid not qualitatively change between these two approaches, results basedon all species are reported in the main text, and results for species withnegative sensitivities are in Table H.1. For each trait, we ran three separatemodels predicting phenological sensitivity to spatial, temporal, and actual41temperature (the actual site-year temperature assigned to each specimen;hereafter referred to as ?overall?). For each model, the absolute value ofphenological sensitivity was calculated and then square root transformedto meet the assumption of a normal error distribution. We recognize thepotential for committing a type I error with the large number of modelswe tested. As such, our interpretations are based on the full set of resultsrather than results from individual models.We evaluated the predictive ability of eight ecological traits that de-scribed characteristics of species? flight seasons, overwintering strategy, hostplant specialization, dispersal ability, and range size (Table F.2). Thesetraits were chosen because they may influence the sensitivity of the tim-ing of flight season to temperature (Table F.2), and given data availability.Specifically, we considered the following flight season attributes as ?traits?:average number of generations across the species? range, average length offlight season across the species? range and timing of flight season. Maximumnumber of generations and maximum flight season length were not analyzedas they were strongly correlated with average number of generations andaverage flight season length (rho=0.9294, df=186; rho=0.9663, df=186, re-spectively). The number of generations and length of flight season werecompiled from Layberry et al. 1998. The timing of flight season initiationwas determined by taking the 5th quantile median collection date across allrecords for each species. We also considered overwintering strategy as atrait, with four categories: egg, larvae, chrysalis, adult (Klinkenberg, 2012;Layberry et al., 1998; Opler et al., 2012).Larval host plant breadth, mobility, wing length and range size were42taken from Burke et al. (2011). Larval host plant breadth was classified asmonophagous (one host species), oligophagous (congeneric host species), orpolyphagous (host species in more than one genus). We used two estimatesof dispersal ability, mobility and wing length. Mobility estimates were basedon a survey of naturalists who ranked species from 0-10 (Burke et al., 2011).Wing length is often used as a proxy for species-level dispersal rates inbiogeographic research. Average wingspan was measured as the distancetaken at the widest spread of forewings in mm. Range size (km2) was basedon species? North American distributions. Where differences in traits weregiven for subspecies, we took the average trait value (this only occurred for3 species complexes: Anthocharis sara, Callophrys gryneus, Coenonymphatullia).Since traits of related taxa may be similar due to common ancestryand therefore not statistically independent (Felsenstein, 1985; Harvey &Pagel, 1991), we assessed whether a phylogenetic analysis was necessary bycomparing residuals from linear models to phylogenetically-adjusted linearmodels (Revell, 2010, TableH.2,H.3,H.4). In order to account for potentialphylogenetic non-independence in our analyses, a molecular phylogenetictree of all species included in our study was constructed (Figure F.1; seeAppendix H for details). We report the results from the model (standardlinear or phylogenetically-corrected model) that best fit the data based onAIC (Appendix H).iii) Evaluating trends in phenology over timeTo test for trends in flight season timing shifts over time, we evaluated43two additional relationships for each species: (1) temperature as a function ofyear for the set of locations and times at which there were collection recordsfor that species and (2) timing of flight season as a function of year (here-after referred to as ?temporal phenological shift? with units of days/year).For both models, we included the identity of the closest weather station asa random effect. We tested for nonlinear relationships and temporal auto-correlation, and then accounted for temporal autocorrelation appropriately(see Appendix G for details). We tested the predictions that temperaturehas increased and that the timing of flight season has advanced over timeusing one tailed one sample t-tests. Where these distributions did not meetassumptions of normality, we report median values.In addition, using the regression coefficients from each analysis, we eval-uated whether a species? temporal phenological shift could be predicted by i)the degree of temperature change for the set of locations and times at whichthere were collection records for that species, ii) its phenological sensitivityto spatial temperature and iii) its phenological sensitivity to temporal tem-perature. We used generalized least squares to test the significance of theserelationships (for details on analysis see Appendix G). Finally, to determinewhether the precision of our estimate of flight season timing influenced ourability to detect phenological trends through time, we limited the analysisof temporal phenological shifts to those species with relatively shorter flightseasons (<2 months) and single generations each year (n=122 species).443.4 Results3.4.1 Testing phenological sensitivity to spatial andtemporal temperatureThe timing of Canadian butterfly species? flight seasons responded pre-dictably to temperature. On average, the timing of flight seasons was sig-nificantly earlier across species in warmer years (-2.38 days/?C (0.18SE),t203=-13.60, p<0.0001; Figure 3.1a) and warmer locations (-1.50 days/?C(0.19SE), s203=36, p<0.0001; Figure 3.1b). Negative slopes were found for89% and 82% of species (n=204), and 58 % and 42% of all slopes weresignificantly negative for phenological sensitivity to temporal and spatialtemperature, respectively (Figure 3.1ab).Flight season timing responded differently to spatial and temporaldimensions of temperature. There was a significant difference betweenphenological sensitivity to spatial and temporal temperature (t203=3.37,p=0.00092), where mean sensitivity to temporal temperature was greaterthan mean sensitivity to spatial temperature by 0.84 days/?C (0.25SE).Moreover, species? phenological sensitivities to spatial and temporal temper-ature were not correlated (r=0.053, t202=0.75, p=0.45; Figure 3.2). There-fore, a species? flight season timing might have been strongly related totemperature through space but not time, and vice versa (Figure 3.2).453.4.2 Using ecological traits to predict phenologicalsensitivity to temperatureAs predicted, earlier-season and less mobile species had stronger phenolog-ical sensitivity to temperature than later-season and more mobile species(Table 3.1; Figure 3.3cd). However, these traits explained only a modestamount of variation (Figure 3.3cd) and were only significant when pheno-logical sensitivity was evaluated with spatial temperature (Table 3.1). Noneof the other traits had any impact on the sensitivity of phenology to tem-perature (Table 3.1).3.4.3 Evaluating trends in phenology over timeTemporal phenological shifts were much weaker than phenological sensitiv-ity to temperature. For species with shorter flight seasons and a singlegeneration per year, we detected phenological trends through time: themean change across species in the timing of flight seasons through timewas marginally significantly less than zero (-0.19 days/decade (0.12SE),t120=-1.62, p=0.054). However, we detected a much weaker advancementin the timing of flight season across all species over the past century (-0.048days/decade (0.12SE), t203=-0.40, p=0.35; Figure 3.1c).While spring temperatures have increased on average for the sets of lo-cations and times for which there were collection records (0.0090 ?C/year(0.0012SE), t203=7.69, p<0.0001), species have experienced widely differ-ent magnitudes and even directions of temperature change (Figure 3.1d).Across species, greater temporal phenological shifts were associated with46sites where temperature changes have been greatest (-3.55 days/?C (0.58SE),LRT3,2=34.66, p<0.0001). Temporal phenological shifts were not signifi-cantly related to phenological sensitivity to spatial (LRT3,2=0.011, p=0.92)or temporal (LRT3,2=1.88, p=0.17) temperature.47Table 3.1: The relationship between species? ecological traits and phenological sensitivity to temperature.Phenological sensitivity to overall, spatial and temporal temperature were analyzed separately. Shown is thebest model for each trait based on phylogenetic and non-phylogenetic models (see Table H.2,H.3,H.4 for modelcomparison). ?Brownian? and ?Pagel? represent common branch length transformations in phylogenetic modelsand ?GLS? represents a non-phylogenetic model. The response variable is the square root of the absolute valueof phenological sensitivity to temperature as defined in the main text. Range size was square-root transformed,and average flight season length and wingspan were log-transformed.?*? means F value.Trait Dimension of phe-nological sensitiv-ityBestmodelTransformationparameterCoefficient (SE) df LRT pvalueAverage numberof generations(residual df=202)Overall GLS NA -0.028 (0.044) 1 0.41* 0.52Temporal BrownianNA 0.061 (0.0097) 1 1.80* 0.18Spatial Pagel ?=0.15 -0.087 (0.053) 1 2.72* 0.1Average lengthof flight season(residual df=202)Overall GLS NA -0.054 (0.066) 3,2 0.68 0.41Temporal GLS NA -0.021 (0.070) 3,2 0.093 0.76Spatial Pagel ?=0.14 -0.034 (0.082) 4,3 0.18 0.67Timing of flightseason (residualdf=202)Overall Pagel ?=0.102 -0.0037(0.0014) 5,4 6.72 0.0095Temporal Pagel ?=-0.00024 -0.002 (0.0013) 5,4 2.18 0.14Continued on next page48Trait Dimension of phe-nological sensitiv-ityBestmodelTransformationparameterCoefficient (SE) df LRT pvalueSpatial Pagel ?=0.17 -0.0039 (0.0018) 5,4 4.7 0.03Mobility (residualdf=198)Overall GLS NA -0.064 (0.025) 4,3 6.47 0.011Temporal GLS NA -0.039 (0.028) 4,3 1.95 0.16Spatial GLS NA -0.073 (0.031) 4,3 5.48 0.019Larval hostbreadth (residualdf=188)Overall GLS NA 0.040 (0.063) 4,3 0.405 0.53Temporal GLS NA 0.023 (0.062) 4,3 0.13 0.71Spatial GLS NA -0.079 (0.076) 4,3 1.09 0.3Range size (resid-ual df=166)Overall GLS NA -3.04e-5 (7.4e-5) 3,2 0.17 0.68Temporal GLS NA -6.7e10-5(7.89e10-5)3,2 0.73 0.39Spatial GLS NA -1.33e10-4(9.0e10-5)3,2 2.18 0.14Wingspan (resid-ual df=195)Overall GLS NA -0.16 (0.101) 3,2 2.42 0.12Temporal GLS NA -0.089 (0.11) 3,2 0.69 0.41Spatial GLS NA -0.14 (0.12) 3,2 1.35 0.25Overwintering(residual df=187)Overall GLS NA -0.0019 (0.049) 3,2 1.60E-050.99Continued on next page49Trait Dimension of phe-nological sensitiv-ityBestmodelTransformationparameterCoefficient (SE) df LRT pvalueTemporal GLS NA 0.039 (0.052) 3,2 0.56 0.45Spatial GLS NA -0.030 (0.060) 3,2 0.26 0.61503.5 DiscussionMuseum collections for animal taxa have rarely been used to test forphenology-climate relationships (but see Bartomeus et al., 2011; Polgaret al., 2013). Our analysis revealed that collection records for butterflies canbe used to detect broad-scale relationships between phenology and climate,an important goal for global change biology given potential consequencesof phenological shifts for trophic interactions and ecosystem services (Pauet al., 2011). We showed that the timing of flight season predictably re-sponded to temporal and spatial dimensions of temperature across Canada(Figure 3.1). The cross-species average degree of phenological sensitivity totemporal temperature (-2.38 days/?C (0.18SE)) is within the range reportedfor other butterfly species based on the dates of first appearance and peakflight (-11.8 to 8.5 days/?C; Dell et al., 2005; Ellwood et al., 2012; Polgaret al., 2013; Roy & Sparks, 2000; Sparks & Yates, 1997; Stefanescu et al.,2003). Therefore, collection records can provide critical historical long-termseries and broad spatial coverage for detecting species? phenological sensi-tivity to temperature.Studies to date have documented considerable variation among speciesin the degree to which phenology responds to temperature and changesover time (e.g. Diez et al., 2012; Hodgson et al., 2011; Wolkovich et al.,2012), but we have a limited understanding of the causes of such variation.Here we show that ecological traits can explain significant variation amongspecies in their phenological sensitivity to temperature (also Hurlbert &Liang, 2012; Moussus et al., 2011). Specifically, species with earlier flight51Figure 3.1: Distribution of slope values across 204 butterfly species for phe-nological sensitivity to temporal (a) (-2.38 days/?C (0.18SE), t203=-13.60,p<0.0001) and spatial (b) (-1.50 days/?C (0.19SE), s203=36, p<0.0001)temperature, temporal phenological shifts (c) (-0.0048 days/year (0.012SE),t203=-0.40, p=0.35), and the degree of temperature change for the set of lo-cations and times at which there were collection records (d) (0.0090 ?C/year(0.0012SE), t203=7.69, p<0.0001). The arrow represents the mean slopeacross species and a slope of zero is represented by a solid line (superim-posed in c).52-20-10010-20 -10 0 10Sensitivity to temporal temperature (days/?C)Sensitivity to spatial temperature (days/?C)Figure 3.2: Relationship between phenological sensitivity to spatial and tem-poral temperature (days/?C) across species. The solid line represents thePearson correlation (r=0.053, t202=0.75, p=0.45).53a)# b)#c)# d)#13014015016017015 16 17 18Temperature (?C)Flight season timing (doy)10015020025030012 14 16 18 20 22Temperature (?C)01234100 120 140 160 180 200Timing of flight season (doy)Phenological sensitivity (days/?C)012342 4 6 8Mobility indexFigure 3.3: The sensitivity of flight season timing (doy: day of year) tospatial variation in temperature (?C) for (a)Erora laeta, a species with lowmobility (mean mobility index=2.5; lme: -6.82 days/?C (1.99SE), n=12), (b)Vanessa annabella, a species with high mobility (mean mobility index=7.67;lme: -0.15 days/?C (3.14.SE), n=73) and across species (days/?C) as afunction of (c) the timing of flight season (-0.0039 (0.0018SE), n=203; PGLS:LRT5,4=4.70, p=0.030) and (d) mobility index (-0.073 (0.031SE), n=199;GLS: LRT5,4=5.48, p=0.019). Shown is the predicted slope with a 95%confidence interval.54seasons and lower dispersal ability appear more sensitive to temperaturethan species with later flight seasons and greater dispersal ability (Figure3.3). Earlier fliers may be more sensitive to temperature (e.g. Anthocharissara) because spring temperature is generally more variable across yearsthan summer temperature and there is a greater likelihood of spring frostin temperate regions than other areas (Cook et al., 2012; Inouye, 2000),leading to potentially heavier fitness costs of mis-timing flight seasons thanit would be for later fliers (e.g. Hesperia leonardus) (Pau et al., 2011).The timing of flight seasons for species with greater dispersal ability (e.g.Vanessa annabella) are likely to be less sensitive to temperature than weakdispersers (e.g. Erora laeta) because they are better able to track suitableclimatic conditions (e.g. finding microsites where food is available) and/orhave reduced local temperature adaptation, and are thus less responsiveto broad-scale climate (Figure 3.3; Diamond et al., 2011; Doligez & Pa?rt,2008). These results are consistent with other studies that have found thatinsect species that emerge earlier in the year have advanced their phenologyto a greater degree (e.g. Altermatt, 2010; Bartomeus et al., 2011; Diamondet al., 2011; Hassall et al., 2007). However, other studies have not found thatmobility is predictive of phenological shifts through time (Diamond et al.,2011; Sparks et al., 2006). These discrepancies could be a result of thedifficulties in measuring mobility or that previous studies related mobilityto temporal phenological shifts (Diamond et al., 2011) rather than directlyto temperature sensitivity as we have done here.Most generally, our results reinforce the conclusions from recent studiesshowing that traits can predict species? responses to climate change (Al-55termatt 2010;Angert et al. 2011; Diamond et al. 2011) with two potentialprovisos. First, traits in our analysis explained only a modest amount ofvariation in phenological sensitivity to temperature (Figure 3.3) indicatinga limited ability of traits to predict species responses to climate change(Angert et al., 2011; Buckley & Kingsolver, 2012). In some cases, evalu-ating the usefulness of traits will require better quantitative estimates ofinherently difficult traits to measure, such as dispersal ability and ecologicalspecialization. Second, these traits only influenced phenological sensitivityto spatial temperature implying that these traits are likely mediating therelationship between phenology, multiple cues, and local adaptation. It re-mains unclear whether the low predictive power of traits we observed here isthe result of poor knowledge of species life histories, drawbacks of the traitsthemselves or a lack of understanding of the processes underlying these re-sponses (Angert et al., 2011; Buckley & Kingsolver, 2012).Interspecific variation in phenological sensitivity to temperature mayhave important consequences for understanding and predicting variation infuture population trends. This sensitivity has recently been used as an in-dicator to predict species? vulnerability to future climate change (Clelandet al., 2012; Willis et al., 2008). For example, plant species in New Eng-land whose flowering-time was not responsive to temperature have greatlydecreased in abundance over the past 150 years (Willis et al., 2008). Onemight also expect negative fitness consequences for those species that arehighly sensitive to temperature, regardless of direction (delayed or advancedflight seasons with warming temperatures over time). A substantial advanceor delay in the timing of a flight season could shorten the length of an in-56dividual?s flight season, affecting their ability to acquire resources, lead topotential phenological mismatches with their host plants, and increase expo-sure to stressful abiotic conditions (Lane et al., 2012; Miller-Rushing et al.,2010).Predicting these ecological consequences will require translating our un-derstanding of phenological sensitivity to temperature into accurate fore-casts of species? phenological shifts through time. In this study, the sensi-tivity of butterfly phenology to temperature differed across time and space(Figure 3.2). Therefore, a species? flight season timing might have beenstrongly related to temperature through time but not space, and vice versa.This finding has two implications for predicting temporal phenological shifts.First, accurately predicting species? phenological shifts through time will re-quire isolating temporal from spatial dimensions of temperature since theyhave different effects on phenology. Therefore, space-for-time substitutionsare unlikely to work in predicting species? phenological responses to climatechange (Illan et al., 2012) and thus historical data are required to make fu-ture predictions. Second, butterflies are likely not simply responding phys-iologically to temperature. Instead, the timing of flight seasons for manybutterfly species is also likely being influenced by local adaptation and/orspatially fixed cues such as day length (Hodgson et al., 2011). These resultsare consistent with recent studies that show a combination of cues may de-termine the flight periods of different species of Lepidoptera (Hodgson et al.,2011; Valtonen et al., 2011). Therefore, quantifying the relative importanceof different cues and local adaptation on phenology, and taking into accountthe sensitivity of phenology to spatio-temporal dimensions of temperature57will be critical for making better predictions of phenological responses toclimate change.The sensitivity of species? flight season timing to temperature only trans-lated into shifts in flight season phenology through time for a subset ofspecies: those with shorter flight seasons and a single generation (59% ofspecies). While most evidence suggests that species have been recently ad-vancing their reproductive or vegetative phenologies (e.g. Parmesan, 2007),such a signal was unlikely to be found here given the substantial noise com-mon in collection records (Robbirt et al., 2011) and in this database (e.g.single records as estimates of phenology, variation in the degree of temper-ature change and phenological sensitivity to temperature). Our estimate ofthe timing of flight season was imprecise as we only had a single observationin each year for the majority of species-site combinations, thus contributingto inter-annual variation in phenology and making it difficult to detect a con-sistent advancement in phenology over time. Temporal phenological shiftswere greater not only in species for which the precision of the estimate ofthe flight season timing was greatest (short flight season, single generation),but also when the degree of temperature change across a species? range wasgreater. Therefore, on average, we suspect that past increases in tempera-ture were not large or consistent enough, relative to inter-annual variation,to allow detection of directional shifts in the timing of flight seasons for manyspecies, and that some species simply have not shifted the timing of theirflight seasons because this phenological phase is not responsive to temper-ature. However, the clear sensitivity of flight season timing to temperaturesuggests that with increased warming, future temporal shifts are likely for58many if not most species. Such shifts could lead to fitness consequences(Lane et al., 2012; Willis et al., 2008) and to shifts in the relative timing oflife cycle events of closely interacting species (Both et al., 2006; Post et al.,2008).3.5.1 ConclusionsMuseum collection records, an under-exploited source of phenological data,can provide a critical resource to explore broad-scale relationships betweenphenology and temperature. However, collection records are likely to beless applicable to i) making precise estimates of species-specific phenology-climate relationships, unless there is exceptionally good spatial-temporalcoverage for a given species, and ii) to detecting temporal trends in species?phenology for species with extended flight seasons and multiple generations.Our results suggest that ecological traits can help account for interspecificvariation in species? phenological sensitivity to temperature. Finally, isolat-ing the temporal-spatial dimensions of temperature sensitivity of phenologywill be critical in accurately predicting species? phenological responses tofuture climate change. It remains uncertain, however, to what degree lo-cal adaptation and day length, among other factors, interact to modify thetiming of butterfly flight seasons.59Chapter 4Evaluating differences inphenological sensitivity totemperature betweenbutterflies and their nectarfood plants using museumand herbarium specimens4.1 SynopsisVariation among species in their phenological responses to climate changesuggests that shifts in the relative timing of interacting species are likely tooccur. However, the expected magnitude and direction of these shifts formost interactions is unknown. Here, we use collection records of butter-flies and their nectar food plants to compare their phenological sensitivity60to temperature (days/?C) and estimate the magnitude of potential pheno-logical shifts that could influence the synchrony between these groups oforganisms. We also evaluate whether butterfly ecological specialization af-fects the difference in phenological sensitivities between these plants andbutterflies to temperature. On average, the phenology of both butterfliesand plants advanced in response to warmer temperatures. For pair-wisebutterfly-nectar food plant interactions, we found that plant phenology wasmore sensitive to temperature than butterfly phenology, suggesting thatshifts in phenological synchrony for plant-insect interactions with climatewarming are likely. Substantial interspecific variation in phenological sen-sitivity to temperature occurred for plants and for most plant species, thissensitivity differed across space and time (i.e. species? timing of floweringmight have been strongly related to temperature through space (cold vs.warm sites) but not through time at particular sites (cold vs. warm years)).This suggests the additional importance of non-temperature cues and/orlocal adaptation for many species. Finally, we found no significant differ-ences in phenological sensitivity to temperature between butterfly specialistsand generalists. Nevertheless, specialists are likely to be more vulnerable tothose shifts in synchrony than generalists since they are dependent on thephenology of fewer plant species.4.2 IntroductionPhenological shifts in response to recent climate change have been reportedfrom many different taxonomic groups (e.g. Menzel et al. 2006; Parmesan612007; Roy & Sparks 2000). While the most commonly reported phenologi-cal response has been advancement in seasonal timing, substantial variationwithin and across taxonomic groups has occurred (Both et al., 2009; Parme-san, 2007; Thackeray et al., 2010). One of the potential consequences ofthis variation is a shift in the relative timing of life cycle events of closelyinteracting species (phenological synchrony). Altered timing of ecologicalinteractions could be disruptive given the potential for fitness consequences(Both et al., 2006; Klapwijk et al., 2010; Liu et al., 2011; Post et al., 2008).For example, changes to plant-insect interactions could influence pest out-breaks and mediate pollination services (Bentz et al., 2010; Burkle et al.,2013; Dukes et al., 2009; Hegland et al., 2009). However, little is knownabout the magnitude and direction of shifts in phenological synchrony un-der climate change and how sensitive ecological interactions are to theseshifts.Here, we used ? 100 years of butterfly collection records and plantherbarium specimens to compare the sensitivity of phenology to temper-ature of two interacting trophic groups across British Columbia, Canada,a climatically heterogeneous region. Despite their limitations (e.g. non-systematic collecting, relatively less information in individual locations), col-lection records have recently been used to document phenological responsesto climate change (e.g. Bartomeus et al., 2011; Gallagher et al., 2009; Polgaret al., 2013; Primack et al., 2004, Chapter 3) and they facilitate a broad-scale comparative approach to phenological research (Vellend et al., 2013).Understanding whether interacting taxa have different sensitivities to thesame abiotic cue can help estimate the magnitude of potential shifts in phe-62nological synchrony in response to climate change for interacting species.Therefore, we used differences in phenological sensitivity to temperature,as indicated by the date of collection, between butterflies and plants usedas nectar food plants to quantify potential shifts in synchrony over time.Since temperature can have different effects on phenology across space vs.over time (e.g. Doi & Takahashi, 2008; Forkner et al., 2008; Hodgson et al.,2011) and we are particularly interested in predicting species? phenologicalshifts through time to estimate shifts in phenological synchrony, we quan-tified temperature sensitivity by relating phenology and temperature acrossboth space (variation in average temperature from site to site in Canada)and across time at particular sites (variation from year to year within eachindividual site).Temperature is thought to be the most consistent and dominant con-troller of spring plant and insect phenology (Bale et al., 2002; Hodgsonet al., 2011; Pau et al., 2011; Wolkovich et al., 2012). However, within taxo-nomic groups, there can be substantial interspecific variation in phenologicalsensitivities to temperature (Chapter 3; Diez et al., 2012; Hodgson et al.,2011; Primack et al., 2009; Wolkovich et al., 2012), suggesting there is likelyto be variation across taxonomic or trophic groups, thus increasing the po-tential for shifts in phenological synchrony. While a relatively small numberof studies have compared the temperature sensitivity of phenology for dif-ferent trophic levels (Cook et al., 2008; Doi et al., 2008; Forrest & Thomson,2011; Gordo & Sanz, 2005; Huey et al., 2002; Primack et al., 2009), to ourknowledge, no study has compared the sensitivity across many interactingspecies.63At present, it is difficult to make general predictions about the differencesin butterfly-plant phenological sensitivities to temperature, and thus, themagnitude of shifts in phenological synchrony. Most studies have found thatplants have greater phenological sensitivity to temperature than animal taxa(Cook et al. 2008; Doi et al. 2008; Forrest & Thomson 2011; Huey et al. 2002;Primack et al. 2009, but see Gordo & Sanz 2005), suggesting that mobilityconstraints might be an important mechanism. Plant phenology could bemore sensitive to temperature than animal phenology (for mobile life historystages) because individual plants are sessile, such that they cannot modulatethe warming they experience by moving among microhabitats, thus increas-ing the adaptive advantage of phenological plasticity (Huey et al., 2002;Thackeray et al., 2010). However, some studies have documented greaterphenological shifts over time for insects relative to plants (Parmesan, 2007;Visser & Holleman, 2001) and others have found no difference (Bartomeuset al., 2011; Polgar et al., 2013). While these studies did not measure temper-ature sensitivity of phenology directly, their results suggest that metabolicconstraints, rather than mobility, could be an important mechanism. As aresult of differences in metabolic complexes (respiration-limited metabolismvs. photosynthesis-limited metabolism), the metabolism and abundance ofheterotrophs appears to be more sensitive to nonlethal temperature shiftsthan autotrophs (Allen et al., 2005; O?Connor et al., 2011). Therefore, in-sect phenology could be more sensitive to temperature than plants due tometabolic differences.Interaction strength, the degree to which one species? survival or repro-duction is limited by the other species or factors involved in the interaction,64should be a predictor of the extent to which a species? fitness is affected bya shift in phenological synchrony (EncinasViso et al., 2012; Miller-Rushinget al., 2010). Since data on interaction strengths are quite rare, ecologicalspecialization can be used as a rough proxy. For example, Lepidopteranspecies with a broader diet breadth are less dependent on individual hostplant species and are thus better able to spatially track suitable abiotic con-ditions and resource availability (Altermatt, 2010; Diamond et al., 2011).As such, the phenology of generalists is likely to be less dependent on thephenology of particular plant species, provided there is adequate variationin the phenology of their nectar food plants and that abundances and geo-graphic distributions of nectar food plants are similar among specialists andgeneralists, and less sensitive to temperature (Miller-Rushing et al., 2010;Visser & Both, 2005). Therefore, there should be a difference in phenologi-cal sensitivity to temperature between generalists and specialists, leading todifferences in potential shifts in phenological synchrony. Indeed, butterflygeneralists have recently experienced relatively smaller temporal advances inphenology than species with a narrower larval diet breadth (Altermatt, 2010;Diamond et al., 2011), suggesting this trait might be useful in identifyingspecies that are especially sensitive to shifts in phenological synchrony.In this study of butterflies and plants in British Columbia, we focusedon the timing of the flowering season and the timing of the flight season asthe two phenological phases of interest given our use of collection records.We had three main objectives: (i) to compare the phenological sensitivity totemperature of adult butterflies and plants used as nectar food plants; ii) todetermine whether this sensitivity differed across space and time for these65taxa; and (iii) to evaluate the influence of butterfly ecological specializationon butterfly phenological sensitivities to temperature and butterfly-plant dif-ferences in sensitivity. We also evaluated whether the sensitivity of species?phenology to temperature has been translated into detectable, directionaltrends in phenology over time.4.3 Methods4.3.1 Butterfly dataOur main data source for butterflies was the Canadian National Collectionof Butterflies database (Layberry et al., 1998, updated as of January 2011),which includes ?16,000 georeferenced (>80% of records to within 1 km)collection records for 187 species in British Columbia dating from the late19th century to the present. Each collection record includes a specimenpreserved in one of 40 Canadian natural history collections: specimens werecollected and identified initially by lepidopterists and re-verified by Lepi-dopteran systematists (see Layberry et al., 1998). We supplemented thedatabase with additional data (?59,500 records) from the Spencer Ento-mological Collection (University of British Columbia) and the personal andprofessional collections of Canadian butterfly experts Syd Cannings, CrisGuppy, Ross Layberry, Norbert Kondla and Jon Sheppard (all pers. comm.).Supplemental records without associated geographic coordinates were geo-referenced by R. Layberry using their locality descriptions, GPS software(QuoVadis, http://www.quovadis-gps.de/), Google Earth and Google Maps.Only locations accurate to within 2km were used. Nomenclature was stan-66dardized based on Pelham (2011).To avoid pseudo-replication, duplicate records of a species on a particu-lar day at a specific location (within ?2km) were considered a single record.We worked at the species level so records identified only to genus level orhybrids were not used and sub-species were grouped to the species level.Based on previous tests of sampling bias (Bedford et al., 2012), we removedspecies suspected of having large geographical collection biases (large ar-eas of severe undersampling): Plebejus glandon, Boloria chariclae, Boloriafrigga, Lycaena phlaeas, and Oeneis melissa. To avoid including erroneousrecords, we eliminated records representing unique location-year combina-tions with recorded collection dates in December, January, or February ondays that had a daily maximum temperature less than 10?C (there were 12such records). We assumed that these dates were incorrect as none of thespecies in our database are known to be in flight during those months (Lay-berry et al., 1998). To isolate the effect of spring temperature on the sameyear?s adult flight season, we excluded all non-resident species in BritishColumbia (migratory, rare strays etc.). Only collection records with match-ing weather station data (see description below) were included. Finally,we restricted analyses to species for which there were at least 10 collectionrecords spread across a range of ?30 years. In total, there were 14,629 but-terfly records specimens of 122 species that met our criteria. Sample size perspecies ranged from 12 to 475 records (median 89.5) that covered a medianof 102.5 years.674.3.2 Plant dataFor plants, we focused on species that are known to be or are potentialnectar sources for adult butterflies found in British Columbia. Plant specieswere chosen by examining lists of adult nectar food plants in butterfly atlases(Guppy & Shepard, 2001; Opler et al., 2012; Pyle, 2002, EH Strickland Ento-mological Collection: http://entomology.museums.ualberta.ca/index.html,Insects of Alberta: http://www.insectsofalberta.com). Where only agenus or species assemblage (e.g. thistles) was specified, E-Flora BC(http://www.geog.ubc.ca/biodiversity/eflora/) was consulted to determinewhich species of the genus were present in British Columbia. Both exotic andnative species were included. We compiled an initial list of 578 plant speciesin British Columbia that could be used by butterflies for nectar sources.Herbarium records were compiled from the University of BritishColumbia (UBC) herbarium, the largest collection of plant specimens inWestern Canada. We checked each specimen for the presence of flowers andused only these specimens in our analyses. As for butterflies, specimensidentified only to genus level or hybrids were not used and sub-species weregrouped to the species level. For 330 of 578 species, we found at least onespecimen that was in flower at the time of collection. Specimens with nodate of collection or locality description were not included. To avoid pseudo-replication, duplicate records of a species on a particular day at a specificlocation were considered a single record. About half of the specimens hadno geographic coordinates associated with them. Therefore, these recordswere georeferenced based on locality descriptions associated with the speci-68men using Biogeomancer (http://bg.berkeley.edu/latest/) and Google Maps.We excluded records with >10km uncertainty (i.e. those coordinates onlyaccurate to a single decimal degree or where uncertainty was explicitly cal-culated based on the detail included in the locality description). This dataset included 4689 georeferenced records from 289 plant species. However,we restricted records to only those with matching weather station data (seesection below) and analyses to species for which there were at least 10 col-lection records spread across a range of ?30 years. In the final dataset, therewere 2100 records of 60 plant species with 11-58 records each (median=16)covering an average of 75.3 years.Of the plant species, 87% (52/60) were herbaceous and 13% were shrubs.Twenty-four families were represented with the most common being Aster-aceae, Fabaceae and Apiaceae (25%, 13%, 8%, respectively). Most of thespecies were perennial (73%), 12% were annual and the duration of life cyclefor the rest reflects a combination of strategies (e.g. annual and biennial).4.3.3 Climate dataDaily temperature data was extracted from the National Cli-mate Data and Information Archive (Environment Canada;http://climate.weatheroffice.gc.ca) for weather stations across BritishColumbia using only data that was quality-controlled (i.e., the mostreliable data available). For each butterfly collection record and herbariumspecimen, weather data from the closest weather station within 10km inany direction was taken. As our estimate of temperature for plants, weused mean daily temperature for March 1 to May 31. For butterflies we69used mean daily temperature for May 1 to July 31. These variables werethe best predictors of flowering season and flight season timing, respectively(see Appendix I for details on selection of seasonal block).4.3.4 Statistical analysisThe day of year of collection records was used to estimate the flowering andflight season timing for each species-site-year combination. The analysis wasdivided into three sections. First, we tested the sensitivity of phenology totemporal and spatial dimensions of temperature. Next, we tested for trendsin phenology over time. Finally, we compared the phenological sensitivityto temperature of butterflies and plants by taking two approaches: a broadcomparison between the two sets of species and one focusing on pair-wisebutterfly-plant interactions.As an estimate of flowering and flight season timing and to avoid pseudo-replication within years, we calculated the median collection date for eachspecies in each location for each year. The majority of locations had a sin-gle specimen per year per species (71% of butterfly specimens and >97%of herbarium specimens). We assumed that the specimens collected repre-sent unbiased samples from the individuals flowering and flying in a givenlocation and year. While there are substantially more collection records inrecent years (since 1970), this would only influence our estimate of shiftsin the timing of flowering and flight seasons if people collected systemat-ically earlier in the year over the past forty years. However, there is noreason to think that collectors systematically changed their sampling strat-egy, especially across such broad scales. If anything, we might expect an70underrepresentation of specimens collected during peak flowering or flightin especially early or late years, thus making our tests conservative. Meanoverall day of flowering was 175.9 (25 June; 0.93SE) and mean overall dayof flight was 178.3 (27 June; 0.32SE).We used linear mixed-effects models (?nlme? package in R (Pinheiro et al.,2012) and lmer; ?lme4? package in R (Bates et al., 2011)) to test phenologicalsensitivity to temperature, and the degree of temperature change and phe-nological shifts over time. All models were fitted using restricted maximumlikelihood (REML) which handles correlated error structures, unbalanceddesigns and gives less biased variance estimates than the method of maxi-mum likelihood (Pinheiro & Bates, 2000). Likelihood ratio tests were usedfor model comparison and AIC was used for model selection following Burn-ham & Anderson (2002). If sample size was less than 40, AICc (AIC cor-rected for small sample size) was used. We judged significant improvementin model fit when AIC or AICc improved by greater than two (Burnham &Anderson, 2002).We evaluated the assumptions of all parametric tests and used non-parametric tests where assumptions were not met. All statistical analyseswere performed using R 2.14.1 (Team, 2012).i) Testing phenological sensitivity to spatial and temporal temperatureOur objective was to quantify responses in the timing of flowering andflight seasons to spatial (cold vs. warm sites) and temporal (across timeat particular sites (cold vs. warm years)) dimensions of temperature. Todo so, we constructed a mixed-effects model for each species with day of71year as a function of two temperature variables (spatial and temporal di-mensions of temperature), and the year and weather station associated withthe specimen were included as random effects. The spatial dimension oftemperature was characterized by calculating the mean temperature acrossall years of data available for each weather station (i.e., site). The temporaldimension of temperature was characterized by calculating the differencebetween the temperature for the site and year of a given specimen and themean temperature for that site, effectively estimating an inter-annual tem-perature differential at that site. The two regression coefficients from thismodel were used to define phenological sensitivity to temperature with unitsof days/?C (hereafter referred to as ?phenological sensitivity?); hereafter weuse the terms ?phenological sensitivity to spatial temperature? and ?phe-nological sensitivity to temporal temperature? to refer to each coefficient.Therefore, we define sensitivity as a measure of the potential of a species?phenology to respond to temperature changes (a slope in days/?C) ratherthan as a measure of goodness of fit (e.g., R2) or statistical significance ofthe relationship between phenology and temperature. We tested for poten-tial nonlinear relationships and the presence of temporal autocorrelation butfound no evidence of either in the majority of species (>50%; see AppendixI for details).We tested the prediction that average phenological sensitivity to temper-ature across species is negative (i.e. warmer temperatures leads to earlierflight seasons) for most species using one-tailed one-sample t-tests, sepa-rately for spatial and temporal temperature. Where distributions of species?phenological sensitivities did not meet assumptions of normality, we report72median values.We then compared phenological sensitivity to spatial and temporal tem-perature using two approaches. First, we tested whether there was a signif-icant difference in the mean slope of these two sensitivities using a pairedt-test. We then tested whether these sensitivities were correlated acrossspecies using the Pearson correlation coefficient.Finally, to evaluate the effect of butterfly ecological specialization on thetemperature sensitivity of the timing of flight season, we compared butterflyspecialists and generalists in their phenological sensitivity to both dimen-sions of temperature using a two-tailed two-sample t-test. Butterflies werecategorized as specialists if they use a single plant family for nectar as adults(n=25) and generalists if they use more than one plant family for nectar(n=56). Only families and genera found in British Columbia were counted.Larval host plant preference was not considered because phenological infor-mation was not available at this scale for the larval stages of the butterfly orfor the plant (e.g. budburst). Nectar plant use was not known for all adultbutterflies in the database and some adult butterflies do not use nectar as afood resource. Therefore, this analysis was restricted to 81 butterfly species.ii) Evaluating phenological shifts over timeTo test for trends in flowering and flight season timing shifts over time,we evaluated two additional relationships for each species: (1) temperatureas a function of year for the set of locations and times at which therewere collection records for that species and (2) timing of flowering or flightseason as a function of year (hereafter referred to as ?temporal phenological73shift? with units of days/year). For both models, we included the identityof the closest weather station as a random effect. We tested for potentialnonlinear relationships and the presence of temporal autocorrelation butfound no evidence of either in the majority of species (>50%; see AppendixI for details). We tested the predictions that temperature has increasedand that the timing of flowering and flight season has advanced over timeusing one tailed one sample t-tests. Where these distributions did not meetassumptions of normality, we report median values.iii) Butterfly-plant comparisonsTo compare the phenological sensitivity to temperature of butterfliesvs. plants, we took two approaches. First, we compared the central ten-dency of both groups using a two-tailed two-sample t-test. All butterflieswere included in this comparison, regardless of whether their adult nectarplant preference was identified. We also considered the effect of any outliersgiven preliminary analyses, which suggested the presence and influence ofoutlier species on butterfly-plant differences in phenological sensitivity totemperature. Second, we evaluated the butterfly-plant difference in phe-nological sensitivity to temperature based on pair-wise species interactions.Pairs (hereafter referred to as ?butterfly-plant interactions?) were identifiedbased on adult nectar food plants or food groups (e.g. thistles) listed foreach butterfly species in the butterfly atlases. Where only a genus or group(e.g. thistles) was specified, we assumed any species of the genus present inBritish Columbia could be a potential nectar source. Therefore, both knownand potential nectar sources for adult butterflies were included. We calcu-74lated the difference in butterfly and plant phenological sensitivity for eachspecies pair. In cases where a butterfly species was paired with multipleplant species, we took the median sensitivity across the plant species. Wethen tested whether the mean difference in sensitivity across butterfly-plantinteractions differed from zero using a paired t-test. In total, there were 58butterfly-plant interactions.To determine whether the butterfly-plant difference in phenological sen-sitivity to temperature differed for butterfly specialists and generalists, wegrouped the butterfly-plant interactions based on butterfly specialization(specialists: 9 interactions; generalists: 42 interactions). For both groups,we tested whether the mean difference in sensitivity across butterfly-plant in-teractions differed from zero using a paired t-test. To test whether butterfly-plant differences in phenological sensitivity differed between specialists andgeneralists, we used a two-sample two-tailed test. There were fewer special-ists and generalists in these comparisons than the phenological sensitivityanalysis (previous section) because butterfly-plant interactions were con-structed at the species level, rather than the family and genus level, thusreducing the number of pair-wise interactions given the relatively small plantdataset.754.4 Results4.4.1 Testing phenological sensitivity to spatial andtemporal temperatureThe timing of flowering and flight seasons predictably responded to spatialand temporal dimensions of temperature. In warmer years and locations,the timing of flowering and flight season occurred significantly earlier (Table4.1; Figure 4.1).Butterfly and plant phenology responded differently to spatial and tem-poral dimensions of temperature. For both taxa, mean phenological sensi-tivity to temporal temperature was significantly different from phenologicalsensitivity to spatial temperature (plants: Sign test; s=22, p=0.052; butter-flies: t121=-3.21, p=0.0017, Table 4.1; Figure 4.1). On average, phenologicalsensitivity to temporal temperature was greater than phenological sensitivityto spatial temperature for plant and butterfly phenology by 1.70 days/?Cand 1.44 days/?C, respectively (Table 4.1; Figure 4.1). Butterfly species?phenological sensitivity to spatial and temporal temperature was correlated(r120=0.27, p=0.0031) i.e. species? flight season timing varied similarly withtemperature through space and time. There was no significant correlationbetween phenological sensitivity to spatial and temporal temperature forplants (r=0.17, p=0.20).No significant difference in the phenological sensitivity to temperatureoccurred between butterfly specialists and generalists (temporal: t79=-0.49,p=0.63; spatial: t79=0.16, p=0.87).764.4.2 Evaluating phenological shifts over timeOn average, neither butterflies nor plants showed significant shifts in phenol-ogy through time (Table 4.1). The mean change in phenology did not differsignificantly from zero for either group (Table 4.1). While mean daily maxi-mum temperatures for May-July have increased on average for the sets of lo-cations and times for which collection records were available for each species,mean temperatures did not significantly increase for March-May (Table 4.1).Moreover, species experienced widely different magnitudes and even direc-tions of temperature change over time (range in temperature change: plants:-0.080 to +0.16 ?C/year, butterflies: -0.14 to +0.082 ?C/year).4.4.3 Butterfly-plant comparisonThe difference in phenological sensitivity to temporal temperature betweenbutterflies and plants used as nectar sources was influenced by the subsetof species and approach used. When all species were considered, the meanbutterfly-plant difference in phenological sensitivity to temporal tempera-ture was not significantly different from zero (Welch t-test: t64.75=-0.47,p=0.64; Table 4.2, Figure 4.1ac). However, interspecific variation in thesensitivity of the timing of flowering season to temporal temperature wasmore than twofold greater than the interspecific variation in the sensitivityof the timing of flight season to temporal temperature, even with outliersexcluded (plant: n=59, sd=10.10, butterfly: n=121, sd=4.17; Figure 4.1ac).When only those species with identified butterfly-plant interactions were in-cluded (i.e. this resulted in a non-random subset of the full dataset), the77mean butterfly-plant difference in phenological sensitivity to temperatureincreased from 0.77 to 2.57 days/?C (Welch t-test: t39.83=1.088, p=0.28),with plants demonstrating greater phenological sensitivity than butterflies(Table 4.2). However, a single plant species outlier had a large effect onthis difference: removing this species increased the mean butterfly-plantdifference to 4.015 days/?C (t=20.93, p=0.043, Table4.2).Taking into account pair-wise interactions between butterflies and theirnectar food plants further increased the mean butterfly-plant difference inphenological sensitivity to temporal temperature. The difference betweenthe sensitivity of the timing of flowering season and flight season to temporaltemperature was significant across all butterfly-plant interactions (t57=5.64,p<0.0001; Table 4.2). The mean temperature sensitivity of the timing ofplant flowering season was greater than the mean sensitivity of the timingof butterfly flight season by 6.87 days/?C(Figure 4.2) and in 77.6% of theinteractions, plant phenological sensitivity to temperature was greater thanbutterfly sensitivity. Plant and butterfly species in 27.6% of the pair-wiseinteractions had sensitivities that differed in sign (e.g. plant species hadnegative sensitivity and butterfly species had positive sensitivity). Acrossbutterfly-plant interactions, there was no correlation in phenological sensi-tivity to temporal temperature (r=-0.13, t56=-1.0, p=0.32; Figure 4.2b).Grouping the butterfly-plant interactions based on the specialization ofbutterflies did not substantially change mean butterfly-plant difference insensitivity (Table4.2). Nectar plant specialists and generalists did not sig-nificantly differ in mean butterfly-plant differences in phenological sensitivityto temperature (Table4.2).78There were no butterfly-plant differences in phenological sensitivity tospatial temperature regardless of species grouping or approach (Table4.2).However, interspecific variation in the sensitivity of the timing of floweringseason to spatial temperature was greater than the interspecific variation inthe spatial temperature sensitivity of the timing of flight season, even withoutliers excluded (plant: n=57, sd=5.57, butterfly: n=122, sd=3.18; Figure4.1bd).79Table 4.1: Differences in the phenological sensitivity of butterflies and plants, and butterfly specialists andgeneralists, to spatial and temporal temperature. Also shown is the mean change across species in phenologyand temperature over time for the set of locations and times at which there were observations. Where a signtest instead of a one-sample t-test was used, the s test statistic and median slope value is shown. Significant pvalues are in bold.Taxa Temperaturevariable or groupof species (num-ber of species)Model (units of coefficient) Mean coeffi-cient (SE)T teststatisticP valuePlants Mean for March-May (n=60)Sensitivity to temporal temperature(days/?C)-3.69 (1.57) -2.35 0.011Sensitivity to spatial temperature(days/?C)-1.82 (1.31) S=22 0.026Phenological change (days/year) 0.13 (0.045) 2.82 0.99Temperature change (?C /year) 0.0050(0.0050)1.01 0.16Butterflies Mean for May-July (n=122)Sensitivity to temporal temperature(days/?C)-2.92 (0.43) -6.84 <0.0001Sensitivity to spatial temperature(days/?C)-1.48 (0.29) -5.14 <0.0001Phenological change (days/year) 0.0036 (0.016) S=65 0.79Temperature change (?C /year) 0.0075(0.0020)3.8 0.00012Continued on next page80Taxa Temperaturevariable or groupof species (num-ber of species)Model (units of coefficient) Mean coeffi-cient (SE)T teststatisticP valueSpecialists(n=25)Sensitivity to temporal temperature(days/?C)-3.13 (0.60) S=1 <0.0001Sensitivity to spatial temperature(days/?CC)-1.24 (0.72) -1.72 0.049Phenological change (days/year) -0.0015 -0.049 0.48Temperature change (?C /year) 0.0081(0.0035)2.3 0.015Generalists(n=56)Sensitivity to temporal temperature(days/?C)-2.90 (0.47) -6.2 <0.0001Sensitivity to spatial temperature(days/?C)-1.36 (0.36) -3.81 0.00018Phenological change (days/year) 0.01 S=34 0.96Temperature change (?C/year) 0.0064(0.0021)S=43 <0.000181Table 4.2: Comparison between butterflies and plants in their phenological sensitivity to temporal and spa-tial temperature. Also shown is whether any outliers were removed in the analysis, the statistical approachused, mean phenological sensitivity to temperature (days/?C) of each group and mean butterfly-plant differ-ences in phenological sensitivity to temperature (days/?C). ?Interactions? here refer to pair-wise butterfly-plantinteractions. Significant p values are in bold.Group ofspeciesTemperaturevariationSample size(butterflies,plants)Outliersre-moved?StatisticalapproachMean sen-sitivity(days/?C)(SE)(butterfly,plant)Differencein sen-sitivity(days/?C)Teststatistic,p valueAll species Temporal 122,60 No Welch twosample t-test-2.92 (0.43), -3.69 (1.57)0.77 T=-0.47,p=0.64121,59 Yes Welch twosample t-test-3.12 (0.38), -4.57 (1.31)1.45 T=-1.062,p=0.29Spatial 122, 60 No Welch twosample t-test-1.48 (0.29), -0.83 (1.31)0.65 T=0.48,p=0.63122,58 Yes Two sam-ple t-test-1.48 (0.29), -1.31 (0.73)0.17 T=0.26,p=0.80Continued on next page82Group ofspeciesTemperaturevariationSample size(butterflies,plants)Outliersre-moved?StatisticalapproachMean sen-sitivity(days/?C)(SE)(butterfly,plant)Differencein sen-sitivity(days/?C)Teststatistic,p valueSpecieswith iden-tifiedinteractionsTemporal 58,38 No Welch twosample t-test-2.38 (0.45), -4.95 (2.32)2.57 T=1.088,p=0.2858,37 Yes Welch twosample t-test-2.38 (0.45), -6.40 (1.86)4.015 T=20.93,p=0.043Spatial 58,38 No Welch twosample t-test-1.36 (0.34),0.22 (1.84)1.14 T=-0.84,p=0.4158,37 Yes Welch twosample t-test-1.36 (0.34), -1.32 (1.04)0.04 T=-0.031,p=0.98Temporal 58 interac-tionsNo Paired t-test-2.38 (0.45), -9.25(1.07)6.87 T=5.64,<0.0001Spatial 58 interac-tionsNo SIGN test -1.36 (0.34), -0.83 (0.66)-0.48 S=25,p=0.36Specialists Temporal 9 interac-tionsNo Paired t-test-2.23 (1.24), -8.97 (3.90)6.7 T=1.47,p=0.18Continued on next page83Group ofspeciesTemperaturevariationSample size(butterflies,plants)Outliersre-moved?StatisticalapproachMean sen-sitivity(days/?C)(SE)(butterfly,plant)Differencein sen-sitivity(days/?C)Teststatistic,p valueSpatial 9 interac-tionsNo Paired t-test-1.65 (0.43),0.10 (1.71)-1.75 T=-1.011,p=0.30Generalists Temporal 42 interac-tionsNo Paired t-test-2.43 (0.54), -9.52 (1.18)7.091 T=5.17,<0.0001Spatial 42 interac-tionsNo SIGN test -1.35 (0.46),-0.87 (0.84)0.028 S=21,0.99840246810-40 -20 0 20 40Number of speciesTemporal05101520-20 0 20 40 60Spatial010203040-40 -20 0 20 40Sensitivity (days/?C)Number of species0204060-20 0 20 40 60Sensitivity (days/?C)a)# b)#c)# d)#Figure 4.1: Plant (a,b) and butterfly (c,d) phenological sensitivity to tem-poral (a,c) and spatial (b,d) temperature. Mean phenological sensitivityis a:-3.69, b:-1.82, c:-2.92, d:-1.48 days/?C. Zero temperature sensitivity isrepresented by a solid line85Figure 4.2: Butterfly-plant differences in phenological sensitivity to temporal temperature taking into accountpair-wise interactions. (a) An example of the difference in the phenology-temperature relationship for a pair-wiseinteraction between Glaucopsyche piasus (Arrowhead blue butterfly; circles) and Lupinus sericeus (silky lupine(plant); triangles). Shown are best-fit lines from linear models (solid=G. piasus, -0.49 days/?C (2.56SE); t72=-0.49, p=0.85; dashed= L. sericeus, -3.62 days/?C (6.41SE); t 9=-0.56, p=0.59). (b) Butterfly-plant differences inphenological sensitivity (days/?C) for 58 species-pair interactions with a dotted line representing the 1:1 relation-ship, showing that mean plant phenological sensitivity is greater and more variable than butterfly phenologicalsensitivity864.5 DiscussionWhile there are likely to be important ecological and evolutionary conse-quences of shifts in phenological synchrony due to climate change (Klapwijket al., 2010; Liu et al., 2011; Post et al., 2008), it is uncertain how extensivethese shifts are likely to be and which species are likely to be more sensitive tothese shifts. This study contributes three main findings about the likelihoodand magnitude of shifts in phenological synchrony for butterflies and plants.First, while there was no general tendency for a butterfly-plant difference inphenological sensitivity to temperature, incorporating pair-wise interactionsrevealed significant butterfly-plant differences in sensitivity. Thus, focusingon the central tendency of taxonomic groups, a common approach in thephenology literature (e.g. Both et al., 2009; Parmesan, 2007), likely under-estimates the differences between pairs of interacting species in phenologicalsensitivity to temperature. This is to be expected given the substantial in-terspecific variation in sensitivity of both butterflies and their nectar plants(Figure 4.1). Our results suggest that shifts in phenological synchrony arelikely to be widespread for interacting plant-insect species. However, syn-chrony could increase or decrease depending on which species? phenologyis earliest in the year, and on the direction and magnitude of shift of eachspecies. For example, butterfly-nectar food plant synchrony could decrease(timing between events widens) if the nectar plant?s peak flowering initiallyoccurs before the butterfly?s peak flight season, and subsequently advancesmore with warming than the butterfly?s peak flight season.Second, plants had greater overall phenological sensitivity to temporal87temperature than butterflies based on pair-wise interactions (Figure 4.2).These results are consistent with studies that have compared plant pheno-logical sensitivity to temperature to animal taxa (Cook et al. 2008; Doi et al.2008; Forrest & Thomson 2011; Huey et al. 2002; Primack et al. 2009,but seeGordo & Sanz 2005), and with some trends in phenological advances due torecent climate change (Huey et al., 2002; Thackeray et al., 2010). However,other studies have found that butterflies have experienced greater pheno-logical shifts over time than plants (Both et al., 2009; Gordo & Sanz, 2005;Parmesan, 2007). Discrepancies here are likely as a result of not directly re-lating phenology to temperature as we have done here. Phenological trendsover time do not necessarily indicate responsiveness of phenology to tem-perature as they are dependent on the degree and pattern of temperaturechange (Chapter 3; Cook et al., 2008), timeframe of available data (Badecket al., 2004; Diez et al., 2012; Iler et al., 2013) and can be influenced bychanging population sizes (Miller-Rushing et al., 2008). The greater sen-sitivity of plant phenology to temperature than adult butterfly phenologysuggests that mobility is an important trait that can influence the abilityof individuals to modulate the local temperatures they experience (Chapter3; Huey et al., 2002; Thackeray et al., 2010). Given these butterfly-plantdifferences in phenological sensitivity to temperature, shifts in phenologicalsynchrony are likely to be driven by more pronounced responses of plantsthan butterflies to climate warming.Third, substantial interspecific variation in the sensitivity of the timingof flowering season to temperature, and the lack of spatio-temporal correla-tion across plant species in phenological sensitivity to temperature, suggests88that making general predictions for shifts in phenological synchrony (e.g.direction of shifts) based on plants could be more challenging than for but-terflies. Our results suggest that the phenology of some plant species maybe more influenced by cues other than temperature (e.g. precipitation, frost;Crimmins et al., 2010; Inouye, 2008) and that spatially fixed cues such asday length and local adaptation are also influencing the timing of flowering(Doi & Takahashi, 2008; Hodgson et al., 2011). Therefore, quantifying therelative importance of different cues and local adaptation on plant phenol-ogy, and taking into account the sensitivity of plant phenology to spatial andtemporal dimensions of temperature will be critical for making general pre-dictions about shifts in phenological synchrony between insects and plantsin response to climate change.Butterfly ecological specialization did not influence the direct effect oftemperature on the timing of flight season nor the butterfly-plant differencesin phenological sensitivity to temperature. These results suggest that dietbreadth does not influence phenological sensitivity to temperature, which issupported by results based on larval diet breadth from Chapter 3. Althoughrecent studies have found that Lepidopteran specialists have advanced theirphenologies to a greater degree in response to recent climate change thangeneralists (Altermatt, 2010; Diamond et al., 2011), it is unknown how wellthis trait directly predicted phenological sensitivity to temperature for thosespecies. Our results suggest that the likelihood and magnitude of shiftsin phenological synchrony is unlikely to differ among butterfly specialistsand generalists but given our rough classification of specialization (usingplant families rather than genera or species) and relatively small sample89size, further studies are required. Nevertheless, provided the nectar plantsused by both groups have similar abundances and geographic distributions,butterfly specialists are likely to be more vulnerable to shifts in phenologicalsynchrony than generalists since they are dependent on the phenology offewer plant species.The sensitivity of species? phenology to temperature did not trans-late into shifts in phenology through time. While most evidence suggeststhat species have been recently advancing their phenologies (e.g. Parmesan,2007), a strong signal was unlikely here given the substantial noise in thedata (e.g. single records as estimates of phenology, variation in the degreeof phenological sensitivity to temperature; Chapter 3). Indeed, results fromChapter 3 suggest that the strength of temporal phenological shifts increasesi) with the degree of temperature change across a species? range and ii) theprecision of the estimate of timing. Therefore, on average, we suspect thatpast increases in temperature were not large or consistent enough, relative tointer-annual variation, to allow detection of directional shifts in phenology,and that some species simply have not shifted the timing of their seasons be-cause these phenological phases are not responsive to temperature (Chapter3). Also, estimates of butterfly-plant phenology were imprecise as we onlyhad a single observation each year for the majority of species-site combina-tions, thus contributing to inter-annual variation in phenology and makingit difficult to detect a consistent advancement in phenology over time (Chap-ter 3). Nevertheless, the clear temperature sensitivity of phenology of bothtrophic groups suggests that with increased warming future temporal shiftsare likely for many if not most species.904.5.1 ConclusionsFor pair-wise butterfly-plant interactions, we found that plant phenologywas more sensitive to temperature than butterfly phenology, suggesting thatsubstantial shifts in phenological synchrony for plant-insect interactions withclimate warming are likely already occurring and likely to continue. Our re-sults also suggest that the consideration and identification of ecological in-teractions will be critical for making accurate predictions about the effects offuture climate change for ecological communities. While phenological shiftsin response to future climate change are likely to be widespread for plantsand insects, accurately predicting the direction and magnitude of shifts inphenological synchrony with plants will require a better understanding ofthe relative importance of different abiotic cues influencing plant phenology.Finally, whereas butterfly specialists may not necessarily experience largershifts in phenological synchrony than generalists, their fitness is likely tobe more vulnerable to these shifts. Without direct estimates of the fitnessconsequences of these shifts, it remains uncertain as to whether the shifts inphenological synchrony of the magnitude estimated here would be biologi-cally important for the interacting species. As such, there is a clear needfor studies that consider the ecological and evolutionary consequences ofwarming-driven shifts in phenological synchrony (e.g. Chapter 5).91Chapter 5Direct and indirect effects ofexperimental warming on thephenology and interactionbetween an insect herbivoreand its host plant5.1 SynopsisThe phenology of many species is shifting in response to climatic changes atvarying rates across different species, thereby potentially affecting species?interactions and individual fitness. However, few studies have experimen-tally tested the influence of warming on the timing of species interactions andin particular, the potential for different direct and indirect (via phenologi-cal change) effects of temperature on each species. Here, we experimentallywarmed egg masses and larvae of the western tent caterpillar (Malacosomacalifornicum pluviale) placed on the branches of its host plant, red alder92(Alnus rubra), in the field to test the influence of warming on the phenol-ogy of the two species. Warming advanced the timing of larval but notleaf emergence, leading to a range in phenological mismatch, with larvaeemerging as much as 25 days before to 10 days after the appearance ofleaves. However, even the earliest-emerging larvae had high survival in theabsence of food for up to three weeks prior to leaf emergence. Therefore,early-emerging larvae are surprisingly resistant to starvation. In addition,warming had no net effect on insect fitness, likely driven by opposing directand indirect effects. For example, warming directly accelerated developmentwhereas temperature-driven changes in phenological mismatch decelerateddevelopment. Our results demonstrate that the indirect effects of warmingare as important to consider as the direct effects on insect fitness and thatsome species could be resistant to future climatic warming.5.2 IntroductionClimate change is causing phenological shifts at highly variable rates acrossspecies in different functional groups and trophic levels (e.g. Both et al.,2009; Parmesan, 2007; Thackeray et al., 2010). Such species-specific vari-ation in response to climate has often led to changes in the relative tim-ing of key activities (phenological synchrony) among strongly interactingspecies (e.g. Edwards & Richardson, 2004; Post et al., 2008; Visser & Holle-man, 2001) but not always (Bartomeus et al., 2011; Charmantier et al.,2008). Altered timing of ecological interactions, such as plant-herbivore, orpredator-prey, can influence species? abundances, food-web structure and93ecosystem-level properties such as primary productivity (e.g. Edwards &Richardson, 2004; Liu et al., 2011; Tylianakis et al., 2008).Determining the implications of disruptions in phenological synchronyis critical for predicting how ecological communities will respond to climatechange. In seasonal environments, phenology is thought to be under strongselection, with important consequences for synchrony between species (Fu-tuyma, 1998). Initiating growth or reproduction too early or too late inthe season can reduce fitness via stressful abiotic conditions, reduced re-source availability or quality, or from increases in top-down pressures (e.g.herbivory or predation). For example, premature larval emergence of early-spring feeding insect herbivores might increase the risk of mortality due tosevere weather or from lack of food (Hunter, 1990), whereas late emergencecould reduce fecundity via reduced foliage quality (Feeny, 1970). While thereis some experimental evidence of phenological shifts altering the timing ofinteractions between insect herbivores and their host plants (e.g. Hunter,1990; Martel et al., 2001), such shifts do not always have important fitnessconsequences (e.g. Durant et al., 2005; Kerslake & Hartley, 1997). Our un-derstanding of the issue remains rudimentary, however, given few direct,experimental tests of warming-driven shifts in phenological mismatch on fit-ness or demography (Klapwijk et al., 2010; Liu et al., 2011; Yang & Rudolf,2010).Temperature may influence fitness not only through direct effects (e.g.,temperature-dependent physiological processes), but also indirect effects viaphenological mismatch and altered interactions with other species (e.g. Bar-ton et al., 2009; Harmon et al., 2009; Suttle et al., 2007). Indirect effects94can, in fact, be stronger than the direct effects of temperature and the twomay interact, either amplifying or countering one another (e.g. Barton et al.,2009; Suttle et al., 2007). For example, the net effect of warming on plantproductivity can be due largely to a warming-induced increase in herbivoryrather than direct effects on plant growth (Barton et al., 2009; O?Connor,2009). To our knowledge, no study has compared the direct effects of warm-ing on consumer fitness to indirect effects via warming-driven phenologicalshifts, even though this is essential to predicting the net impact of climatechange on ecological communities.Here, we experimentally manipulated spring temperature around eggmasses and larvae of the western tent caterpillar (Malacosoma californicumpluviale Dyar (Lepidoptera: Lasiocampidae)), an early-spring feeding insectherbivore, placed on the branches of its host plant, red alder (Alnus rubraBongard (Betulaceae)), in the field. Our aim was to test the influence ofwarming on the phenology of the two species. In this system, warmingcould lead to two primary indirect effects: changes in resource availabilityvia phenological mismatch, and changes in resource quality (Figure 5.1).Other indirect effects could include changes in top-down pressures (e.g.predators, pathogens) but we do not measure those effects here. Todetermine the effects of warming-driven changes in leaf quality on insectfitness, we conducted a bioassay in a controlled environment where larvaewere fed leaves from either warmed or control treatments. For proxies ofinsect fitness, we used larval development time, larval weight, and familysurvival (western tent caterpillar females lay 100- 250 eggs in a single eggmass (referred to henceforth as a ?family?)).95Predictive frameworkWe tested a series of predictions concerning the direct and indirect (viaphenological mismatch and leaf quality) effects of warming on insect fitnessproxies (Figure 5.1). First, based on the temperature-body size rule, warm-ing should directly reduce development time and decrease final body size(Atkinson, 1994), while having little or no effect on survival based on theweight of previous evidence, provided the temperature remains below thespecies? upper threshold for viability (Amarasekare & Sifuentes, 2012). Sec-ond, the indirect effect of warming via leaf quality on insect fitness shouldbe negative, given that warming is expected to accelerate leaf maturation(Buse et al., 1998; Dury et al., 1998). Finally, the effect of mismatch shoulddepend on how warming changes the relative timing of larval and leaf emer-gence and on the fitness proxy. In this system, there are likely to be twopossibilities: warming could cause larvae to emerge in advance of leaf emer-gence (negative mismatch) or it could cause larvae to emerge closer to thetime of optimal food quality (reduced positive mismatch, i.e. fewer daysbetween leaf and larval emergence than in control treatments). If warmingleads to a negative mismatch, then development time should increase, aslarvae will not have food available upon emergence, countering the directeffect of warming. Alternatively, if warming reduces positive mismatch, de-velopment time should decrease as food will be available to larvae sooner,thereby amplifying the direct effects of warming by shortening developmenttime. Opposite effects of the two mismatch scenarios are predicted for lar-val weight and family survival (i.e. both should decrease with a negative96mismatch and increase with a reduced positive mismatch).5.3 Methods5.3.1 Study systemWestern tent caterpillars occur in the Pacific Northwest of North America(our study region) where they feed on a variety of plant species including redalder, hawthorn (Crataegus monogyna Jacq), crab apple (Malus diversifoliaL) and wild rose (Rosa nutka) (Myers, 2000). They are univoltine and larvaeemerge from eggs in the early spring at roughly the same time as budburst,generally in late March or early April. Five instars of gregarious larvae feedfor 6-8 weeks and form conspicuous silken tents on which they congregatebetween bouts of feeding. Pupation occurs in June and adults emerge aftertwo weeks and mate. Adults do not feed and only survive for 1-4 days. Oncemated, females deposit a single egg mass. A few weeks later, embryogenesisbegins and larvae develop and remain within the eggs until the followingspring.Red alder is a native deciduous tree and is widely distributed in lowlandsthroughout the Coastal Douglas fir and Coastal Western Hemlock regionsof the Pacific Northwest (Harrington et al., 1994). Budburst occurs in earlyspring and in congenerics, phenophases have been shown to occur earlier inyears with warmer temperatures (Miller-Rushing & Primack, 2008; Sparkset al., 2011).975.3.2 Study design overviewExperiments were carried out at two sites on the University of BritishColumbia (UBC, Vancouver, British Columbia; 49? 15? 42.72N, 123? 14?56.78W) campus: Totem Field (a mowed field research site) and UBC Farm(an experimental farm and forest research site). At both sites, we used thesouth-facing side of fifteen trees adjacent to open fields. In 2010 we con-ducted a pilot experiment at Totem Field to (i) test the influence of ourproposed warming method (polyethylene bags as ?greenhouses?) on branch-level temperature, (ii) provide an initial assessment of warming effects onearly phenology of both species (but not the full life cycle), and (iii) ver-ify the presence of sufficient variation in temperature and mismatch amongbranches within treatments to permit statistical tests of direct vs. indirecteffects of temperature on fitness. On March 10, 2010, single egg masses wereplaced on two branches of each of 15 trees with one branch serving as a con-trol and the other warmed. For the warmed treatment, clear polyethylenebags (35x55cm) with 8-10 holes (5x5cm) for ventilation and water drainagewere attached over the egg masses near the tips of the branches and fastenedwith a tie.A criticism of our experimental warming approach is that the entirephysiological unit (i.e. whole tree) has not been exposed to the warmingtreatment. However, warming only branches has advanced vegetative phe-nology in other tree species (Barton & Jarvis, 1999) and other studies havefound temperature-related shifts in budburst using severed branches in con-trolled environments (Jepsen et al., 2011; Pop et al., 2000). Nevertheless,98to test for potential differences between warming the entire tree and onlythe terminal end (?50 cm) of a branch, in 2011, we compared the date ofbudburst on branches that were cut (i.e. physiologically separated from thetree) to ones that remained attached. Our prediction was that if the controlof budburst on branches were independent from the tree, buds on cut anduncut branches would burst at the same time. For this experiment, we useda 2x2 factorial systematic block design where we had four treatments (cutand attached, warming and control) per tree on 15 trees at Totem Field. Cutand attached branches were paired and placed as close together as possibleon the tree to minimize microclimate variation. Branches were cut afterbags were added and then attached back to the same branch with zip ties.To maintain a water supply to the cut branches, we attached water spikesto the end of the cut branches that were refilled with water semi-daily and athin slice was cut from the branches once a week to optimize water uptake.The experiment was started in conjunction with the main experiment andbudburst was observed daily.The main experiment began on January 27-28 2011, when single eggmasses were deployed to warmed and control branches of each of 15 trees ineach of the two sites (i.e. 30 pairs). At Totem Field, a second warming repli-cate per tree was added. Therefore, in total across the two sites there were30 control and 45 warmed branches. Warming treatments were identical tothose used in the pilot experiment.To determine whether warming indirectly affects insect fitness proxiesthrough changes in leaf quality, we conducted a bioassay in a controlledenvironment in 2011 in which larvae were either fed leaves from warmed99or control treatments from Totem Field. Larvae for this experiment camefrom 20 egg masses placed on trees at Totem Field at the same time themain experiment was started (4 egg masses on each of 5 trees). Ten bagswere added to branches (without egg masses) on each of the five trees toprovide a ?warmed? food supply for the experiment.5.3.3 Warming experimentsBranches varied in the amount of direct sunlight they received given theirposition and height from the ground. While the terminal ends of brancheswere used, some branches were less exposed than others (e.g. found mid-way to the main trunk or under other branches). Branch height ranged fromground level to around 12 feet. Egg masses were randomly assigned to eachbranch. Egg masses in 2010 were taken from an apple orchard on SaturnaIsland (approximately 50 km from UBC, across the Strait of Georgia) onFebruary 17, 2010 and kept outdoors until they were attached to branches.Egg masses used at Totem Field in 2011 were from moths reared from the2010 pilot study (Appendix J). Egg masses used at UBC Farm in 2011were collected along Trueworthy Road on Saturna Island on January 11,2011 and kept outdoors until they were attached to branches. In 2011, eggmasses were split by site to avoid potential viral contamination from theSaturna Island population, where a nucleopolyhedrovirus (Myers, 2000) wasthought to potentially be prevalent that particular year.Small temperature loggers (iButtons R?, Maxim Integrated Products Inc.,USA) were placed next to all egg masses to record the temperature every100hour for the duration of the main experiment. Additionally, two ibuttonsat each location recorded relative humidity (one on a control branch, oneon a warming branch). Only minor differences in relative humidity existedbetween treatments over the span of the whole experiment (mean and maxi-mum daily maximum relative humidity of the warming treatments was 2.34(%RH) and 1.26 (%RH) greater than ambient conditions, respectively).Branches were checked every 1-2 days for budburst, leaf appearance,larval emergence and larval instars. For all phenological events, the firstobservation per branch or family was used for analysis. Budburst was definedwhen the bud scale was broken and leaf tissue visible. Leaf emergencewas defined as a leaf fully separated from the bud (but not yet expanded).Phenological mismatch was the number of days between larval emergenceand first leaf appearance (negative values indicate larvae emerged beforeleaves, and vice versa). Starting at third instar, ten different randomlychosen larvae from each family were weighed every four days (or soonerif they had reached another instar). The field experiment was terminatedbeyond early fourth instar because warmed families had to be moved toother branches to replenish the leaf supply, thereby changing leaf quality.However, until that point (after leaf emergence) larvae were never foliage-limited (i.e. we ensured that bags initially included enough buds). Thecomplete disappearance of larvae precluded some developmental and weightmeasurements on 16 (out of 75) families.We used three proxies for fitness for each family: development time(from larval emergence to fourth instar), final larval weight (based on thefirst weight measurement after molting to fourth instar) and family survival101(number of larvae still living when we terminated the experiment). Fasterdevelopment is beneficial for insects as it reduces the time that larvae arevulnerable to disease and natural enemies (e.g. predators, pathogens) (Ben-rey & Denno, 1997) and larval weight is strongly correlated with fecundityin insects (Hone?k 1993, but see Leather 1988). The latter relationship isexpected to be strong for the western tent caterpillar because they do notconsume food during the adult stage. Given the gregarious behaviour oflarvae and large family size, we could not measure the survival of individuallarvae; however, we were still able to estimate family survival.To better characterize how the growth period differed between treat-ments, we also measured growth rate for each family. Growth rate wasestimated as the slope of the function relating log-transformed weightand time. We included third instar larval weights and the first weightmeasurement after molting to fourth instar.5.3.4 Leaf quality experimentLarvae raised for the leaf quality experiment were kept on the trees undernatural conditions until they reached late second instar, when they werebrought into the lab and maintained in a controlled environmental chamber(day: 22?C, night: 18?C). The first two families (egg masses) to reach latesecond instar on each tree were used. In the lab, families were initiallykept together because young larvae depend on their siblings. They werefed control leaves from their host tree until they reached the late thirdinstar. Twenty-five larvae were then randomly chosen from each family,102individually weighed (day 0) and placed randomly into plastic cups (onelarva per cup) with mature leaves (>5cm) from either warming or controltreatments. Leaves were replaced every 24 hours (i.e. larvae were never food-limited) in approximately the same amount between treatments. Larvaewere weighed every 5 days. For fitness proxies, we used survival (from thirdinstar onwards), female pupal weight and time to reach pupation. We alsoestimated growth rate as the slope of the function relating log-transformedweight (including larval and pupal weight) and time.5.3.5 Statistical analysesThe analysis was divided into four sections (see Appendix K for statisticaldetails). First, to test the effect of warming a branch on budburst, ratherthan a tree, we compared differences in budburst between cut and attachedtreatments. We also tested for the sensitivity of budburst to temperaturewithin treatments. Second, to determine whether warming affected theearly phenology of the two species, we examined treatment differences inphenology from 2010 and 2011. We also tested for treatment differences ininsect fitness proxies from 2011. Third, we tested the effects of temperaturevs. indirect effects via phenological mismatch on insect fitness proxies,using data on within-treatment variation in branch-level temperature andmismatch. However, because the control larvae were not contained in bags,we lost more larvae from the control families than the warmed familiesand only have survival counts for five control families at Totem Field.Therefore, we could not test for treatment differences in family survivalor within-treatment relationships for that treatment-site combination.103Finally, to test the indirect effects of warming via leaf quality on insectfitness proxies, we compared insect fitness proxies between treatments inthe leaf quality experiment. As both environmental conditions and thesource of egg masses differed between Totem Field and UBC Farm, wetreated these sites as two experiments rather than two treatments of anexperimental factor. Consequently, all analyses were conducted separatelyfor Totem Field and UBC Farm and performed using R 2.14.1 (Team, 2012).Within-treatment relationshipsTo test for the direct and indirect effects (via phenological mismatch)of temperature on our proxies of fitness, we explored within-treatment rela-tionships between each proxy of fitness and mismatch and temperature (Ap-pendix K). We used within-treatment variation rather than across-treatmentvariation as preliminary analyses showed different relationships of fitnessproxies with temperature and mismatch in the two treatments. For familysurvival, we accounted for the number of eggs laid the previous year andthe number of larvae that emerged. To assess the relative effects of directand indirect effects of temperature on our proxies of fitness, we comparedstandardized regression coefficients. As the main temperature variable weused mean daily maximum temperature calculated up to May 25 (the aver-age date across all families to reach the fourth instar), which was found tobe the best predictor of larval fitness (Appendix K, Table K.1).1045.4 Results5.4.1 Treatment temperature effectsOver the duration of the pilot warming experiment in 2010, warming treat-ments were 2.93?C (0.18SE) warmer than ambient conditions based on meandaytime temperature. In 2011, mean daytime temperature in warmingtreatments over the duration of the main experiment was 2.97?C (0.22SE)and 1.97?C (0.23SE) warmer than ambient conditions at Totem Field andUBC Farm, respectively. Overall, Totem Field was warmer than UBC Farmwith mean daytime temperatures being 1.45?C (0.36SE) higher for warmedbranches and 0.45?C (0.14SE) higher for controls.5.4.2 Testing the effect of branch-level warming onbudburstControl and warmed branches showed the same patterns in budburst timingregardless of whether they were attached to the tree, showing that budburston branches can respond independently to warming. Buds of branches cutfrom trees did not differ in the date of budburst compared to attachedbranches (?23,2=0.0007, p=0.98), nor was there any difference between onlycontrol (?23,2=0.0008, p=0.98) or warming (?23,2=0.0053, p=0.94) branches(Figure 5.2). Though not significant, budburst was earlier on warmedbranches by 1.3 (0.89SE) days on average (?24,3=2.17, p=0.14) than controlbranches (Figure 5.2).Analysis of the relationship of budburst to the measured within-treatment temperature variation (among both cut and attached branches)105showed that budburst was earlier on warmer replicates for the Totem Fieldwarming treatment (-3.45 days/?C (1.38SE), LRT4,3=5.24, p=0.022) but norelationship existed for any other within treatment-site combination (Totemcontrol: LRT3,2=0.50, p=0.48; Farm control: LRT3,2=1.62, p=0.20; Farmwarming: LRT4,3=1.66, p=0.20). This suggests that short-term changes intemperature at the branch-level can elicit a response in phenology.5.4.3 Treatment differences in phenology and fitness proxiesIn both years, warming treatments led to substantial shifts in larval but notleaf phenology (Table 5.1, Figure 5.3). In 2010, larval emergence was sig-nificantly advanced in warming treatments compared to control treatmentsby an average of 18.1 days (1.19SE) (LRT4,3=65.44, p<0.0001). However,there was no treatment difference in budburst (LRT4,3=0.54, p=0.46). Asa result, a significant difference in the degree of mismatch between treat-ments occurred; on average warmed larvae emerged 17.75 days (1.19SE)after budburst (LRT4,3=67.72, p<0.0001).In 2011, larval emergence was significantly advanced in warming treat-ments compared to controls by an average of 24.0 days (2.39SE) at TotemField (LRT4,3=50.62, p<0.0001) and 22.5 days (1.72SE) at UBC Farm(LRT4,3=58.78, p<0.0001; Figure 5.3). Warming treatments had no sig-nificant effect on the timing of budburst or leaf emergence at either site(Table 5.1). Consequently, there was a significant treatment effect on thedegree of mismatch at both sites (Totem: LRT4,3=44.33, p<0.0001, Farm:LRT4,3=49.03, p<0.0001). At Totem Field, warming treatments caused asubstantial negative mismatch: larval emergence preceded the appearance106of leaves in warmed families by as much as 25 days (Figure 5.4a). Whilewarming treatments also had a substantial effect on mismatch at UBC Farm,almost all larvae emerged after leaf emergence (i.e., there was no negativemismatch, Figure 5.4b).Experimental warming had no net influence on insect development. Nei-ther development time (Table 5.1), final larval weight (Table 5.1) nor growthrate (Totem: LRT13,12=2.10, p=0.15; Farm: LRT14,13=0.03, p=0.86) dif-fered significantly between treatments at either site. In other words, warmedfamilies significantly advanced their entire larval phase (all having gotten tothe fourth instar by May 30 vs June 10; Totem: LRT4,3=34.31, p<0.0001;Farm: LRT4,3=43.16, p<0.0001; Table 5.1) but they required the same num-ber of days as control families to develop (Table 5.1, Figure 5.3).5.4.4 Within-treatment relationshipsMean daily maximum spring temperature and mismatch had varying effectson development time (Table 5.2, Figure 5.4). For both treatments at bothsites, the slope of development time on temperature was negative (Table5.2; Figure 5.4c,d), although only significantly so for warmed families atTotem Field (-2.92 days/?C (0.39SE), LRT5,4=16.37, p<0.0001). This wasthe treatment-site combination with the widest temperature range (8.19?Ccompared to 1.66-5.38?C for the other three treatment-site combinations;Figure 5.4c,d).Temperature-driven shifts in phenological mismatch influenced develop-ment time. The first families to emerge within each treatment took longestto develop, with stronger development-mismatch relationships for warm-107ing than control treatments (Table 5.2; Figure 5.4a,b). Mismatch had alarger effect on development time than temperature for all treatment-sitecombinations (standardized coefficients; Table 5.2). Families in the warm-ing treatment had shorter development times than one would have pre-dicted by extrapolating from the development-mismatch relationship forcontrol families (i.e. warmed families fall below the control-group regres-sion line; Figure 5.4a,b). This suggests that the direct effect of warmingopposes that of mismatch, resulting in no overall treatment differences indevelopment time (Table 5.1). Indeed, development time was significantlyshorter for the warmed families compared to the controls at Totem Field onceleaves were present (8.32 days (2.44SE), LRT4,3=9.21; p=0.0024). The ef-fect of mismatch on development time remained significant after accountingfor temperature (Totem warming: LRT5,4=5.46; p=0.020; Farm warming:LRT4,3=3.54, p=0.06). Therefore, at Totem Field, families took the leastamount of time to develop if they emerged simultaneously with leaves andexperienced relatively warm temperatures over their larval phase (Figure5.4a,c).Temperature and mismatch had weaker effects on fourth instar larvalweight and family survival (Table 5.2). The only clearly significant trendwas that in control treatments at UBC Farm, the later larvae emerged afterleaf emergence, the greater the number of larvae survived to the end of theexperiment (6.13 larvae (2.36SE), LRT5,4=7.24, p=0.0071), suggesting thatlater emergence improved larval survival.1085.4.5 Leaf quality experimentThe sources of leaves (warmed vs. control treatments) had varying effectson insect fitness. Pupal weight of females fed warmed leaves was signifi-cantly lower than those fed control leaves (-0.035g (0.014SE), LRT5,4=4.44,p=0.035), although growth rate did not differ for either sex (females:LRT11,10=0.053, p=0.82; males: LRT11,10=0.24, p=0.62). Finally, neithersurvival (LRT4,3=0.28, p=0.095), nor the amount of time taken to reach pu-pation (LRT5,4=3.44, p=0.064) varied between treatments under these lab-oratory conditions. However, larvae fed warmed leaves took slightly longerto reach pupation than larvae fed control leaves (2.28 days (1.21SE)).109Table 5.1: Differences in phenology and insect fitness proxies between warm-ing and control treatments at Totem Field and UBC Farm. ?Difference? isthe mean difference between warming and control treatments, where a neg-ative difference indicates that warmed families were earlier, faster or smallerthan control families, and vice versa. Model fit was based on likelihood ratiotests (LRT). Significant P values (<0.05) are in bold.Response Site Difference (SE) DF LRT P valueDate of larval emer-gence (days)Totem -24.00 (2.39) 4,3 50.62 <0.0001Farm -22.47 (1.72) 4,3 58.78 <0.0001Date of fourth instar(days)Totem -24.15 (2.88) 4,3 34.31 <0.0001Farm -25.07 (2.33) 4,3 43.16 <0.0001Larval developmenttime (days)Totem -0.68 (2.63) 4,3 0.07 0.79Farm -2.60 (2.91) 4,3 0.84 0.36Larval weight (g) Totem -0.20 (0.11) 5,4 3.17 0.075Farm 0.0680 (0.043) 5,4 2.44 0.12Date of budburst(days)Totem 2.70 (1.45) 4,3 3.4 0.065Farm -1.27 (1.71) 4,3 0.58 0.45Date of leaf emer-gence (days)Totem -0.50 (0.65) 4,3 0.61 0.44Farm -1.93 (1.035) 4,3 3.34 0.068Degree of mismatch(days)Totem -23.50 (2.67) 4,3 44.33 <0.0001Farm -20.53 (1.70) 4,3 49.03 <0.0001110Table 5.2: Within-treatment relationships of insect fitness proxies with mean daily maximum temperature anddegree of mismatch (indirect effect of temperature) at both sites. Larval weight was log-transformed for alltreatments and family survival was square-root transformed for the Totem Field warming treatment. Model fitwas based on likelihood ratio tests (LRT). Significant and marginally significant P values (<0.06) are in bold.Fitness proxy PredictorvariableSite Treatment Coefficient(SE)StandardizedcoefficientDF LRT P valueDevelopmenttime (days)Temperature(?C)Totem Warming -2.92 (0.39) -0.42 5,4 16.37 <0.0001Control -1.14 (2.60) -0.15 3,2 0.22 0.64Farm Warming -1.56 (1.65) -0.27 3,2 1.04 0.31Control -0.66 (3.12) -0.059 3,2 0.052 0.82Mismatch Totem Warming -1.10 (0.21) -0.81 5,4 13.44 <0.0001Control -0.13 (0.14) -0.28 4,3 0.96 0.33Farm Warming -0.86 (0.44) -0.47 3,2 3.76 0.053Control -0.36 (0.22) -0.37 4,3 2.9 0.089Larval weight(g)Temperature(?C)Totem Warming 0.018 (0.028) 0.078 5,4 0.42 0.52Control 0.72 (0.35) 0.48 4,3 3.6 0.058Farm Warming -0.041 (0.025) -0.23 4,3 2.48 0.12Control -0.10 (0.069) -0.26 4,3 2.13 0.14Mismatch Totem Warming 0.0066(0.0058)0.13 5,4 1.26 0.26Control -0.0070 (0.014) -0.14 4,3 0.25 0.62Continued on next page111Fitness proxy PredictorvariableSite Treatment Coefficient(SE)StandardizedcoefficientDF LRT P valueFarm Warming -0.0037(0.0086)-0.065 4,3 0.19 0.6Control 0.0045(0.0089)0.094 4,3 0.26 0.61Familysurvival(number oflarvae)Temperature(?C)Totem Warming 0.80 (0.42) 0.35 6,5 3.66 0.056Control NA NA NA NA NAFarm Warming 9.43 (11.41) 0.33 5,4 0.98 0.32Control -27.56 (23.44) -0.3 5,4 1.81 0.18Mismatch Totem Warming 1.21 (1.39) 0.33 6,5 0.92 0.34Control NA NA NA NA NAFarm Warming 1.30 (3.67) 0.13 5,4 0.19 0.67Control 6.13 (2.36) 0.49 5,4 7.24 0.0071112Figure 5.1: A conceptual diagram of the predicted direct (solid lines) and primary indirect effects (dashed lines)of warming, through shifts in phenological mismatch and changes in leaf quality, on insect fitness proxies. Fitnessproxies are described in the main text.1135.5 DiscussionRecent climatic change has led to shifts in phenological synchrony for manyinteracting species (e.g. Edwards & Richardson, 2004; Post et al., 2008;Visser & Holleman, 2001), but surprisingly few studies have quantified thefitness or demographic consequences of such shifts (e.g. Both et al., 2006;Post et al., 2008). Our study contributes three main findings. First, experi-mental warming advanced the timing of larval emergence and the entire lifecycle of the western tent caterpillar but it did not influence the early vegeta-tive phenophases of red alder. This led to a range in phenological mismatchfrom -20 to +10 days, with larvae emerging before and after the appearanceof leaves, creating a negative mismatch and reducing positive mismatch com-pared to control families (Figure 5.4c,d). Greater temperature sensitivity ofconsumer phenology relative to that of its resource is consistent with somestudies (e.g. Parmesan, 2007; Visser & Holleman, 2001) but not others (e.g.Both et al., 2009; Thackeray et al., 2010). Therefore, it remains difficult topredict trophic differences in phenological sensitivity to temperature.Second, we found that experimental warming had no net effect on theproxies of insect fitness we used. This was likely as a result of opposingdirect and indirect effects of warming. Temperature-driven shifts in pheno-logical mismatch increased insect development time by prolonging the firstinstar (Figure 5.4a,b) while warming shortened development times (Figure5.4c,d), especially once leaves were present (Figure 5.3). Therefore, larvaethat emerged with leaf emergence and higher temperatures were able to de-velop the fastest (Figure 5.4). We also found that warming-driven changes114Figure 5.2: Comparison of the timing of budburst (doy: day of year) onbranches that were cut (i.e. physiologically separated from the tree) to onesthat remained attached for both control and warmed treatments.115Figure 5.3: Development of western tent caterpillar larvae at Totem Field(a,c) and UBC Farm (b,d) represented by the cumulative appearance of eachinstar across trees for families from warmed (a,b) and control treatments(c,d). Each cumulative curve represents a different instar, denoted by adifferent line type, from left to right: larval emergence, second, third, andfourth instars. Vertical dashed lines represent the mean date of budburstand mean date of leaf emergence (from left to right) of each treatment116Figure 5.4: Bivariate relationships of development time (larval emergenceto fourth instar) with mismatch (a,b) and mean daily maximum tempera-ture (c,d) for Totem Field (a,c) and UBC Farm (b,d). Triangles representwarming treatments and circles represent control treatments. Negative (<0)and positive (>0) mismatch indicates that larvae emerged before and afterleaf emergence, respectively. Dashed (warming) and solid (control) linesrepresent a best-fit line from simple linear regression models.117in leaf quality reduced pupal weight of females, and therefore most likelyfecundity, whereas warming increased final larval weight, albeit weakly. To-gether, these results suggest that the indirect effects of warming in thissystem countered the direct effects resulting in no net influence of warmingon insect fitness. This finding reinforces the importance of considering theindirect effects of climate change on species? responses through interspecificinteractions (e.g. Gilman et al., 2010; Harmon et al., 2009; Suttle et al.,2007) and demonstrates that the net effects of climate change may be dif-ficult to predict given the uncertainty of predicting the nature of indirecteffects (McCann, 2007).Third, we found that early-emerging larvae were surprisingly resistant tostarvation. Larvae from some warmed families at Totem Field were withoutfood for three weeks and while we could not estimate individual larval sur-vival, the degree of mismatch did not influence final family size for warmingtreatments. This suggests that emerging early did not lead to mass mortal-ity. The observed starvation tolerance is substantially greater than what hasbeen measured in controlled environments for most other early-spring feed-ing Lepidoptera larvae (forest tent caterpillar (Malacosoma disstria):Smith& Raske 1968; gypsy moth (Lymantria dispar):Hunter 1993; winter moth(Opterophtera brumata):Wint 1983; but see Hunter 1990). This particularlyhigh degree of starvation resistance has been interpreted as an adaptationto variable budburst across years (Parry et al., 1998) where the strategy ofneonate larvae is to manage absent or poor quality foliage by ?persisting?rather than ?escaping? (e.g. wind-assisted dispersal). These results suggestthat the western tent caterpillar is likely to be quite resistant to warming-118driven shifts in phenological mismatch. Furthermore, these findings implythat starvation tolerance is an important trait to consider when predictingthe severity of fitness consequences of warming-driven phenological shifts forspring-feeding insects.The early vegetative phenophases of red alder were insensitive to ourwarming treatments. Temperature is thought to be the primary cue thattriggers budburst in most plant species, especially for early-developingspecies such as red alder (Pau et al., 2011). Given our own qualitativeobservations of inter-annual variation in budburst, other work showing thatphenophases of sister taxa were sensitive to temporal variation in temper-ature (Miller-Rushing & Primack, 2008; Sparks et al., 2011), and that wefound advancement in budburst with warmer temperatures within one ofour treatments, early phenology of red alder is likely to be sensitive to tem-perature. While red alder populations can have high dormancy stability(i.e. prevention of premature budburst) and suffer delayed budburst due toa chilling deficit (having not experienced a threshold period of low tempera-tures; Hamann et al., 2001), chilling requirements for other tree species arenormally met by the end of December in this region (R. Guy pers. comm.).Therefore, it is unlikely that we interfered with its chilling requirements. In-stead, we suspect that alder is simply not as sensitive to spring temperatureas the western tent caterpillar, suggesting that our results are indicative ofthe mismatch - qualitative if not quantitative - likely to be created by cli-mate change in this system. The passive warming of the bags during thefirst weeks of the experiment was rather modest - enough to accelerate lar-val emergence but not enough to significantly advance budburst (mean tem-119perature treatment difference up to the beginning of budburst was 0.92?C(0.023SE) at Totem and 0.43?C (0.18SE) at the Farm). While further workis needed to understand how warming will affect alder phenology, our resultsillustrate a worst-case scenario for the western tent caterpillar, one in whichthere is no shift in early vegetative phenology of its main host plant but amarked shift in the timing of larval emergence.As for virtually all field-warming experiments (e.g. open-top chamberswhich passively warm ground-level temperatures (Marion et al., 1997)), ourwarming treatments likely altered environmental variables besides temper-ature (Wolkovich et al., 2012), including CO2 concentrations and UV-Bradiation, possibly confounding the effects of temperature on insect fitness.However, several lines of reasoning suggest minimal non-temperature re-lated effects of our treatments. Our bags likely reduced CO2 to some extentwhen the leaves were fully developed and air flow was more limited (duringlate third and fourth instars), which potentially affected larval performance(Buse et al., 1998; Lincoln et al., 1993; Robinson et al., 2012) and leaf qual-ity (Robinson et al., 2012). However, given regular wind blowing throughthe bags it is unlikely that CO2 differed significantly among treatments andsince leaves grew faster in the warming treatments than control treatments(Appendix J), it is unlikely the bags limited photosynthetic rates. Later inthe season when overall UV-B exposure was greater, some aspects, but notall, of larval performance and leaf quality may have been affected by reducedUV-B (Buck & Callaghan, 1999; Rousseaux et al., 2004; Zavala et al., 2001).Therefore, it remains unclear what the net effect of reduced UV-B on in-sect fitness might have been. Despite the potential effects of abiotic factors120such as CO2 and UV-B on insect performance, they are thought to be lessimportant overall than temperature (Scriber & Slansky, 1981). Moreover,leaf quality is likely to have the greatest effect on early instars (Scriber &Slansky, 1981) when changes in CO2 and UV-B were likely to be minimal.Finally, before the appearance of leaves, the treatment difference in temper-ature is likely the only factor that could have influenced the emergence oflarvae from their egg mass.5.5.1 ConclusionsTo our knowledge, this is the first study to demonstrate that warming-drivenshifts in phenological mismatch are as important to consider as the directeffects of warming on consumer fitness. Our results revealed that experimen-tal warming, of a magnitude consistent with predicted climate change effectsfor the study region over the next century, led to significant shifts in phe-nological mismatch between the western tent caterpillar and red alder, andadvanced the insect?s life cycle. While warming-driven phenological shiftsmay also influence other factors such as exposure to natural enemies (impor-tant determinants of fitness in natural populations, Fitzgerald, 1995; Myers,2000), we see clear evidence of opposing direct and indirect effects (throughshifts in phenological mismatch and leaf quality) on proxies of insect fitness.By not explicitly accounting for the indirect effects of warming, especiallywhen they counter the direct effects, predictions about the consequences ofclimate change for a particular species could be misleading. Combined withtheir high starvation tolerance, these results suggest that the western tentcaterpillar may be surprisingly resistant to the effects of climate change.121Chapter 6ConclusionsThis thesis tested several hypotheses related to the causes and conse-quences of interspecific and intertaxonomic differences in the sensitivity ofspecies? distributions and phenologies to climate and variation in temper-ature. Chapter 2 focused on characterizing variation among major taxo-nomic groups in the climate-distribution relationship. Chapter 3 focused onunderstanding interspecific variation in the relationship between butterflyphenology and temperature, while Chapter 4 focused on describing potentialimplications of differences between butterflies and plants in the temperaturesensitivity of phenology. Lastly, Chapter 5 experimentally evaluated the fit-ness consequences of differences in phenological sensitivity to temperaturefor one plant-insect interaction. In this concluding chapter, I review thekey findings of each chapter, make some general conclusions, briefly addressthe main limitations of the thesis, and finish with some reflections on futuredirections in this field.In Chapter 2, I tested the hypothesis that species from different taxo-nomic groups (plants, birds, mammals, herptiles and insects) vary in theirclimate-distribution relationships because of differences in life history strate-gies, in particular dispersal ability. Overall, there was a strong climate-122distribution relationship across species; reinforcing the idea that climate isan important determinant of species? distributions at broad scales. How-ever, in support of the hypothesis tested, there was a significant differencein the strength of species climate-distribution relationship across taxonomicgroups. In contrast, I did not find support for the hypothesis that disper-sal ability affects the strength of this relationship. Within and betweentaxonomic differences in the strength of species? climate-distribution rela-tionships suggest a role for ecological traits but it is unclear which trait(s)is important in determining the strength of this relationship.In Chapter 3, I quantified and characterized butterfly phenology-temperature relationships across Canada, testing for earlier flight seasonsunder warmer conditions over space and time. Since flight season timingpredictably advanced in response to warmer temperatures on average, Iconcluded that collection records - seldom used in climate change researchto date - can be used to detect broad- scale relationships between phenologyand climate. I found that the flight season timing in early-season speciesand species with lower dispersal ability was more sensitive to temperaturethan later-season species and those with greater dispersal ability, suggest-ing that ecological traits can explain interspecific variation in phenologicalsensitivity to temperature. Finally, phenological sensitivity to temperaturediffered across space vs. time suggesting that accounting for both of thesedimensions of temperature will be critical for making better predictions ofphenological responses to climate change.In Chapter 4, I compared the sensitivity of flowering and flight seasontiming to temperature across British Columbia over space and time to es-123timate potential shifts in phenological synchrony between butterflies andtheir adult nectar food plants. As predicted, butterfly and plant phenologyadvanced in warmer years and locations; however, the timing of floweringseason was more sensitive and more variable in response to temperaturethan flight season timing. Moreover, the effects of temperature differedacross space vs. time for flowering season timing, suggesting the importanceof other cues for many plant species. Pair-wise interactions revealed sig-nificant butterfly-plant differences in sensitivity suggesting that substantialshifts in phenological synchrony for plant-insect interactions with climatewarming are likely. Finally, the likelihood and magnitude of these shiftsdoes not appear to differ among butterfly specialists and generalists.Lastly, in Chapter 5, I experimentally tested the effects of warming onthe phenology of an insect herbivore and its main host plant and evalu-ated the consequences of warming-driven phenological shifts on the fitnessof the insect herbivore. My results revealed that experimental warming, ofa magnitude consistent with predicted climate change effects for the studyregion over the next century, led to significant shifts in phenological mis-match between the western tent caterpillar and red alder, and advanced theinsect?s life cycle. I also found evidence of opposing direct and indirect ef-fects (through shifts in phenological mismatch and leaf quality) on proxiesof insect fitness, revealing that warming-driven shifts in phenological mis-match (i.e. indirect effects) are as important to consider as the direct effectsof warming on consumer fitness.1246.0.2 General conclusionsFirst, this thesis has shown that collection records can be used to detectbroad-scale phenology-temperature relationships (Chapters 3 and 4). Giventhat collection records are taxonomically, spatially and temporally exten-sive, they represent an invaluable source of data, especially in areas whereother data may be sparse (e.g. tropical ecosystems), for addressing globalchange biology issues (Lavoie, 2013; Vellend et al., 2013). Thus far, they havebeen underutilized in the global change literature but should be particularlyuseful as they can broaden the temporal extent of phenological studies toa century or more and, as I have revealed through this work (Chapter 3and 4), temporal data will be necessary to make accurate predictions aboutphenological shifts over time. Long-term time series are especially usefulgiven inter-annual variations in weather that can make detecting phenology-climate relationships difficult. Finally, this work reinforces the value of mu-seum collections and the efforts to archive and appreciate our natural history(Vellend et al., 2013, http://naturalhistoriesproject.org).Second, this work has illustrated that ecological traits can be used topredict the sensitivity of species? phenology and geographic distributions toclimatic conditions (Chapters 2 and 3). As such, a trait-based approachwill be useful in predicting species responses to future climate change (Al-termatt, 2010; Angert et al., 2011; Diamond et al., 2011). This approachwill be particularly valuable for species where detailed phenology or oc-currence data are not available. Intrinsic differences among species in lifehistory, physiology and other traits already form a central part of the de-125veloping framework for vulnerability assessments (Williams et al., 2008).Therefore, a trait-based approach will improve our understanding of the bi-ological mechanisms underlying vulnerability and will ultimately result in agreater predictive capacity to prioritize vulnerable species to multiple globalchanges.Species? phenological sensitivities to temperature (Chapter 3,4,5) suggestthat phenological shifts in response to future climate change are likely to bewidespread. Changes in species? phenology are already the most observedand documented response to climate change (Parmesan, 2006; Parmesan& Yohe, 2003; Walther et al., 2002) and this work suggests these shifts willcontinue. Since temperature change is predicted to be much greater over thenext century than the last (IPCC, 2007), and that the magnitude of pheno-logical shifts is dependent on the degree of temperature change (Chapter 3),future phenological shifts should be even greater than what we have alreadyobserved. Given the relationship between phenology and other life historytraits (e.g. seed size, Bolmgren & Cowan, 2008), demography (e.g. Miller-Rushing et al., 2010), species? range limits (e.g. Chuine, 2010), ecosystemfunction (e.g. Edwards & Richardson, 2004), and species? interactions (e.g.Yang & Rudolf, 2010), these phenological shifts are likely to have majorecological and evolutionary consequences.Differences across taxa in phenological sensitivities to temperature(Chapters 4 and 5) imply that shifts in phenological synchrony could be sub-stantial for interacting species and lead to important fitness consequences(Chapter 5). Differential changes in phenologies could disrupt the temporalcoordination of longstanding and potentially coevolved species interactions126(Post et al., 2008; Stenseth & Mysterud, 2002; Winder & Schindler, 2004;Yang & Rudolf, 2010). Indeed, climate change-driven shifts in phenologicalsynchrony have already been linked to population declines of a long-distancemigratory bird (Both et al., 2006) and demographic consequences for a largeherbivore (Post et al., 2008). Where timing gaps become too large for in-teracting populations or species, persistence may depend on evolutionaryadaptation, though the rate of adaptation may still not be sufficient (Visser,2008; Visser & Both, 2005). Given that species interactions can stronglydetermine the structure and dynamics of many ecological communities (e.g.Bascompte et al., 2006; Ives & Carpenter, 2007), the effect of these changesare likely to be profound (Tylianakis et al., 2008; Yang & Rudolf, 2010).Finally, the finding that the indirect effects of warming (e.g. phenologi-cal mismatch) can oppose and be greater than the direct effects of warming(Chapter 5) reinforces the idea that the indirect effects of climate change areas important to consider as the direct effects. By not explicitly accountingfor the indirect effects of warming, especially when they counter the directeffects, predictions about the consequences of climate change for a particularspecies could be misleading (Gilman et al., 2010; Tylianakis et al., 2008).More broadly, this finding argues that incorporating the underlying mech-anistic basis of climatic drivers, and thus specific biological knowledge, willbe important to achieving a more complete understanding of the effects ofclimate change (Boggs & Inouye, 2012).1276.0.3 Limitations and future directionsAs with all methods, the ones I used have some limitations. Briefly, speciesdistribution models are only estimates of species? distributions and are notdirect observations of species? presence or population viability. Therefore,my results in Chapter 2 provide only rough estimates of species? climate-distribution relationships. Since collection records are an unconventionaldata source, in that they were not collected for the purposes of long-termecological studies, there are challenges with using them (Vellend et al., 2013).In particular, imprecision in the estimation of a phenological phase andinaccuracies in the location of older and/or non-georeferenced records arelikely to decrease the signal:noise ratio and contribute to greater uncertaintyaround the estimate of the phenology-temperature relationship. Addition-ally, geographic and temporal collection patterns have the potential to biaspatterns of species? range shifts (e.g. Bedford et al., 2012). Finally, fieldexperiments represent imperfect simulations of the real world; however, asfield-warming experiments commonly alter environmental variables besidestemperature (Wolkovich et al., 2012), the warming experiments from Chap-ter 5 are likely to be as realistic as most.The key findings and limitations of this thesis suggest potential avenuesfor future work. First, while collection records have the potential to quan-tify species? phenological responses to recent environmental changes giventheir extensive coverage, they also have their limitations. To maximize thepotential of collection records in global change studies, adopting a hierar-chical Bayesian approach in future studies should allow more explicit ac-128counting for the uncertainties inherent in collection records, thus improvinginferences and ability to make accurate predictions (Clark, 2005). Bayesianapproaches have yet to become widely used in many branches of ecology(Stephens et al., 2007) but have clear potential, especially with data pre-senting multiple sources of uncertainty.Due to the variation in ecological and life history strategies across species,ecological traits have obvious potential in understanding species? recent re-sponses to climate change. However, this potential may be limited in somecases, with relatively small proportions of interspecific variation in climatesensitivity explained by trait differences (Chapter 2 and 3; Angert et al.,2011). It remains unclear whether these limitations are a result of poorknowledge of species? life histories, drawbacks of the traits themselves or alack of understanding of the processes underlying these responses (Angertet al., 2011; Buckley & Kingsolver, 2012). For example, while dispersalability is likely to be a key strategy for dealing with climate change, wecurrently lack reliable data on dispersal rates, either active or passive, ofmany species (Berg et al., 2010). Therefore, our current estimates of disper-sal ability clearly need improvement. Where feasible, estimation of dispersalrates for interacting species would be particularly informative. Additionally,a better understanding of rare long-distance dispersal events (Clark, 1998;Gillespie et al., 2012) or prioritizing trait measurements within leading-edgepopulations, provided populations exhibit local adaptation and genetic dif-ferentiation (Angert et al., 2011), would provide information especially per-tinent to predicting range expansion.While the phenology of most species included in this work were sensi-129tive to temperature, there was substantial interspecific variation across bothinsects and plants (Chapter 3,4,5), even to the same treatment (Chapter5), and the effects of temperature differed across space and time for manyspecies (Chapter 3 and 4). These findings suggest the importance of otherabiotic cues in determining phenology for many species. Most phenologywork thus far has focused on the effects of temperature but it is clear thatphenological cues involve more than temperature. In particular, assessingthe interaction between phenological cues (e.g. temperature vs. day length)may be key to predicting biological responses to climate change (Pau et al.,2011). While single-factor studies are common in the climate change liter-ature (Brown et al., 2011), testing alternative factors concurrently shouldbe a priority in future studies. Most of the phenological studies have beenprimarily from mesic, temperature seasonal systems in the mid-latitudes(Cook et al., 2012; Crimmins et al., 2010); therefore, comparing phenolog-ical sensitivities across diverse regions would increase the generalizabilityof our conclusions about the effects of climate change on phenology (Cooket al., 2012; Diez et al., 2012).While shifts in phenological synchrony are likely to have substantial im-pacts on ecological communities (Yang & Rudolf, 2010), there are relativelyfew examples in which the fitness consequences of these shifts have beenmeasured (Both et al., 2006; Klapwijk et al., 2010; Liu et al., 2011; Postet al., 2008). As such, this represents an area where future research inthis field should be focused. With the goal of being able to make generalpredictions about the sensitivity of species to these shifts, I suggest fourpotential avenues of particular importance. First, most studies have exam-130ined the effects of climate change on the timing of single life-history events,such as flowering time or peak resource demand (Yang & Rudolf, 2010).However, quantifying entire temporal distributions of phenological parame-ters would allow estimation of the phenological overlap between interactingspecies. This in turn would permit meaningful assessments of the likelihoodthat a biologically relevant shift in phenological synchrony will occur for aninteracting pair of species. For example, a species will be more stronglyinfluenced by shifts in synchrony if the resource it interacts with is avail-able only during a narrow temporal window (Durant et al., 2005). Second,since many interactions extend and change over the lifetime of an organismand are strongly affected by the ontogenetic stages of the interacting species(Yang & Rudolf, 2010), considering the implications of shifts in phenologicalsynchrony over multiple life stages is needed. Third, quantifying interactionstrength (the degree to which one species? survival or reproduction is lim-ited by the other species or factors involved in the interaction) for largeassemblages of species is not feasible; therefore, comparing the sensitivityof species that differ in their degree of ecological specialization to shifts inphenological synchrony would be a good start. Lastly, evaluating the po-tential for novel species interactions would help determine the sensitivity ofa species? fitness to climate change-driven shifts in phenological synchrony(Gilman et al., 2010; Miller-Rushing et al., 2010). For example, switching toa new or less-preferred host with different phenology (Alarco?n et al., 2008;Hellmann, 2002; Singer, 1972) or novel interactions among native and ex-otic species could compensate for shifts in synchrony (Gilman et al., 2010;Schweiger et al., 2010).131Finally, phenology has been identified as an especially important traitshaping species? distributions given that it determines the ability to cap-ture variable resources and defines the season and duration of growth andreproduction (Chuine, 2010). For example, beyond northern range limits,growth and reproduction may not be completed because the time period withfavourable conditions is shorter and the conditions are less favourable over-all than at lower latitudes. While there are clear links between phenologyand geographic distributions, there has been little attempt to consider therelationship among recent responses to climate change (Chuine, 2010). Arethose species that are shifting their phenologies also shifting their ranges?Has there been any feedback between these two responses? For example,terrestrial plants might not track climate as fast as their insect herbivores(e.g. Schweiger et al., 2008) given that they are less mobile than insects. Thelower mobility of plants could result in their phenology being more sensitiveto temperature and thus shifting more through time than insects, while atthe same time restricting the rate of their range shift. Therefore, deter-mining the links between species? recent responses to climate change wouldreveal more about the long-term and overall consequences of climate changefor the persistence of species and those with differing life history strategies,than considering each response independently.As a result of substantial effort over the past two decades, we now havea solid understanding of how species? phenologies and distributions are re-sponding to climate change. Moreover, the direct influence of climate, inparticular temperature, on species? distributions and phenologies has beenwell established. However, our understanding of these responses remains132limited given uncertainty about the relative importance of direct and indi-rect effects of climate change, reasons for interspecific variation in species?responses to recent climate change, and the consequences of such variationfor species? interactions and ecological communities.133BibliographyAlarco?n, R., Waser, N. M., & Ollerton, J. (2008). Year-to-year variationin the topology of a plant-pollinator interaction network. Oikos, 117(12),1796?1807.Allen, A., Gillooly, J., & Brown, J. (2005). Linking the global carbon cycleto individual metabolism. Functional Ecology, 19(2), 202?213.Altermatt, F. (2010). Tell me what you eat and I?ll tell you when you fly:diet can predict phenological changes in response to climate change?ll tellyou when you fly: diet can predict phenological changes in response toclimate change. Ecology Letters, 13(12), 1475?1484.Amarasekare, P. & Sifuentes, R. (2012). Elucidating the temperature re-sponse of survivorship in insects. Functional Ecology, 26, 959?968.Angert, A. L., Crozier, L. G., Rissler, L. J., Gilman, S. E., Tewksbury,J. J., & Chunco, A. J. (2011). Do species? traits predict recent shifts atexpanding range edges? Ecology Letters, 14(7), 677?689.Araujo, M. B. & Pearson, R. G. (2005). Equilibrium of species? distributionswith climate. Ecography, 28(5), 693?695.134Araujo, M. B., Pearson, R. G., Thuiller, W., & Erhard, M. (2005a). Val-idation of species-climate impact models under climate change. GlobalChange Biology, 11, 1504?1513.Araujo, M. B., Thuiller, W., Williams, P. H., & Reginster, I. (2005b). Down-scaling European species atlas distributions to a finer resolution: implica-tions for conservation planning. Global Ecology and Biogeography, 14(1),17?30.Atkinson, D. (1994). Temperature and organism size-a biological law forectotherms? Advances in Ecological Research, 25, 1?58.Badeck, F. W., Bondeau, A., Bo?ttcher, K., Doktor, D., Lucht, W., Schaber,J., & Sitch, S. (2004). Responses of spring phenology to climate change.New Phytologist, 162(2), 295?309.Bale, J. S., Masters, G. J., Hodkinson, I. D., Awmack, C., Bezemers, T. M.,Brown, V. K., Butterfield, J., Buse, A., Coulson, J. C., Farrar, J. F.,Good, J. G., Harrington, R., Hartley, S., Jones, T. H., Lindroth, R. L.,Press, M. C., Symrnioudis, I., Watt, A. D., & Whittaker, J. B. (2002).Herbivory in global climate change research: direct effects of rising tem-peratures on insect herbivores. Global Change Biology, 8, 1?16.Bartomeus, I., Ascher, J. S., Wagner, D., Danforth, B. N., Colla, S., Korn-bluth, S., & Winfree, R. (2011). Climate-associated phenological advancesin bee pollinators and bee-pollinated plants. Proceedings of the NationalAcademy of Sciences, 108(51), 20645?20649.135Barton, B. T., Beckerman, A. P., & Schmitz, O. J. (2009). Climate warmingstrengthens indirect interactions in an old-field food web. Ecology, 90(9),2346?2351.Barton, C. & Jarvis, P. (1999). Growth response of branches of Picea sitchen-sis to four years exposure to elevated atmospheric carbon dioxide concen-tration. New Phytologist, 144(2), 233?243.Bascompte, J., Jordano, P., & Olesen, J. M. (2006). Asymmetric coevolu-tionary networks facilitate biodiversity maintenance. Science, 312(5772),431?433.Bates, D., Maechler, M., & Bolker, B. (2011). lme4: Linear mixed-effectsmodels using s4 classes.Beale, C. M., Lennon, J. J., & Gimona, A. (2008). Opening the climate en-velope reveals no macroscale associations with climate in European birds.Proceedings of the National Academy of Sciences, 105(39), 14908?14912.Bedford, F. E., Whittaker, R. J., & Kerr, J. T. (2012). Systemic rangeshift lags among a pollinator species assemblage following rapid climatechange. Botany, 90(7), 587?597.Beerling, D. J., Huntley, B., & Bailey, J. P. (1995). Climate and the dis-tribution of Fallopia japonica: use of an introduced species to test thepredictive capacity of response surfaces. Journal of Vegetation Science,6(2), 269?282.Bellard, C., Bertelsmeier, C., Leadley, P., Thuiller, W., & Courchamp, F.136(2012). Impacts of climate change on the future of biodiversity. EcologyLetters, 15(4), 365?377.Benrey, B. & Denno, R. F. (1997). The slow-growth-high-mortality hypoth-esis: a test using the cabbage butterfly. Ecology, 78(4), 987?999.Bentz, B. J., Re?gnie`re, J., Fettig, C. J., Hansen, E. M., Hayes, J. L., Hicke,J. A., Kelsey, R. G., Negro?n, J. F., & Seybold, S. J. (2010). Climatechange and bark beetles of the western United States and Canada: directand indirect effects. Bioscience, 60(8), 602?613.Berg, M. P., Kiers, E., Driessen, G., Heijden, M. V. D., Kooi, B. W., Kuenen,F., Liefting, M., Verhoef, H. A., & Ellers, J. (2010). Adapt or disperse:understanding species persistence in a changing world. Global ChangeBiology, 16(2), 587?598.Blach-Overgaard, A., Svenning, J., Dransfield, J., Greve, M., & Balslev,H. (2010). Determinants of palm species distributions across Africa: therelative roles of climate, nonclimatic environmental factors, and spatialconstraints. Ecography, 33(2), 380?391.Blomberg, S. P., Jr, T. G., & Ives, A. R. (2003). Testing for phylogeneticsignal in comparative data: behavioral traits are more labile. Evolution,57(4), 717?745.Boggs, C. L. & Inouye, D. W. (2012). A single climate driver has direct andindirect effects on insect population dynamics. Ecology Letters, 15(5),502?508.137Bolmgren, K. & Cowan, P. D. (2008). Time-size tradeoffs: a phylogeneticcomparative study of flowering time, plant height and seed mass in anorth-temperate flora. Oikos, 117(3), 424?429.Bolmgren, K., Eriksson, O., & Linder, H. P. (2003). Contrasting flower-ing phenology and species richness in abiotically and biotically pollinatedangiosperms. Evolution, 57(9), 2001?2011.Both, C., Asch, M. V., Bijlsma, R. G., Burg, A. B. V. D., & Visser, M. E.(2009). Climate change and unequal phenological changes across fourtrophic levels: constraints or adaptations? Journal of Animal Ecology,78(1), 73?83.Both, C., Bouwhuis, S., Lessells, C., & Visser, M. E. (2006). Climatechange and population declines in a long-distance migratory bird. Na-ture, 441(7089), 81?83.Bowman, J., Jaeger, J. A., & Fahrig, L. (2002). Dispersal distance of mam-mals is proportional to home range size. Ecology, 83(7), 2049?2055.Bradley, N. L., Leopold, A. C., Ross, J., & Huffaker, W. (1999). Phenologicalchanges reflect climate change in Wisconsin. Proceedings of the NationalAcademy of Sciences, 96(17), 9701?9704.Brown, C. J., Schoeman, D. S., Sydeman, W. J., Brander, K., Buckley, L. B.,Burrows, M., Duarte, C. M., Moore, P. J., Pandolfi, J. M., & Poloczanska,E. (2011). Quantitative approaches in climate change ecology. GlobalChange Biology, 17, 3697?3713.138Brown, J. H. (1995). Macroecology. Chicago: The University of ChicagoPress.Brown, J. H., Stevens, G. C., & Kaufman, D. M. (1996). The geographicrange: size, shape, boundaries, and internal structure. Annual Review ofEcology and Systematics, 66(3), 597?623.Bryant, S. R., Thomas, C. D., & Bale, J. S. (2002). The influence of thermalecology on the distribution of three nymphalid butterflies. Journal ofApplied Ecology, 39(1), 43?55.Buck, N. & Callaghan, T. V. (1999). The direct and indirect effects ofenhanced UV-B on the moth caterpillar Epirrita autumnata. EcologicalBulletins, 47, 68?76.Buckley, L. B. & Kingsolver, J. G. (2012). Functional and phylogeneticapproaches to forecasting species? responses to climate change. AnnualReview of Ecology, Evolution, and Systematics, 43(1), 205?226.Burke, R. J., Fitzsimmons, J. M., & Kerr, J. T. (2011). A mobility index forCanadian butterfly species based on naturalists? knowledge. Biodiversityand Conservation, 20(10), 2273?2295.Burkle, L. A., Marlin, J. C., & Knight, T. M. (2013). Plant-pollinatorinteractions over 120 years: Loss of species, co-occurrence and function.Science, 339, 1611?1615.Burnham, K. P. & Anderson, D. R. (2002). Model Selection and Multi-Model139Inference. A Practical Information-Theoretic Approach. New York, NY:Springer-Verlag.Buse, A., Good, J., Dury, S., & Perrins, C. (1998). Effects of elevatedtemperature and carbon dioxide on the nutritional quality of leaves of oak(Quercus robur L.) as food for the winter moth (Operophtera brumata L.).Functional Ecology, 12(5), 742?749.Cahill, A. E., Aiello-Lammens, M. E., Fisher-Reid, M. C., Hua, X.,Karanewsky, C. J., Ryu, H. Y., Sbeglia, G. C., Spagnolo, F., Waldron,J. B., & Warsi, O. (2013). How does climate change cause extinction? Pro-ceedings of the Royal Society B: Biological Sciences, 280(1750), 20121890.Cain, M. L., Damman, H., & Muir, A. (1998). Seed dispersal and theHolocene migration of woodland herbs. Ecological Monographs, 68(3),325?347.Carpenter, K. E., Abrar, M., Aeby, G., Aronson, R. B., Banks, S., Bruckner,A., Chiriboga, A., Corte?s, J., Delbeek, J. C., DeVantier, L., et al. (2008).One-third of reef-building corals face elevated extinction risk from climatechange and local impacts. Science, 321(5888), 560?563.Charmantier, A., McCleery, R. H., Cole, L. R., Perrins, C., Kruuk, L. E. B.,& Sheldon, B. C. (2008). Adaptive phenotypic plasticity in response toclimate change in a wild bird population. Science, 320(5877), 800?803.Chen, I., Hill, J. K., Ohlemu?ller, R., Roy, D. B., & Thomas, C. D. (2011).Rapid range shifts of species associated with high levels of climate warm-ing. Science, 333(6045), 1024?1026.140Chuine, I. (2010). Why does phenology drive species distribution?Philosophical Transactions of the Royal Society B: Biological Sciences,365(1555), 3149?3160.Clark, J. S. (1998). Why trees migrate so fast: confronting theory withdispersal biology and the paleorecord. The American Naturalist, 152(2),204?224.Clark, J. S. (2005). Why environmental scientists are becoming bayesians.Ecology Letters, 8(1), 2?14.Cleland, E. E., Allen, J. M., Crimmins, T. M., Dunne, J. A., Pau, S.,Travers, S. E., Zavaleta, E. S., & Wolkovich, E. M. (2012). Phenologicaltracking enables positive species responses to climate change. Ecology,93(8), 1765?1771.Cook, B. I., Cook, E. R., Huth, P. C., Thompson, J. E., Forster, A., &Smiley, D. (2008). A cross-taxa phenological dataset from Mohonk Lake,NY and its relationship to climate. International Journal of Climatology,28(10), 1369?1383.Cook, B. I., Wolkovich, E. M., Davies, T. J., Ault, T. R., Betancourt, J. L.,Allen, J. M., Bolmgren, K., Cleland, E. E., Crimmins, T. M., & Kraft,N. J. (2012). Sensitivity of spring phenology to warming across temporaland spatial climate gradients in two independent databases. Ecosystems,15(8), 1283?1294.Crimmins, T. M., Crimmins, M. A., & Bertelsen, C. D. (2010). Complex141responses to climate drivers in onset of spring flowering across a semiaridelevation gradient. Journal of Ecology, 98(5), 1042?1051.Davis, A. J., Jenkinson, L. S., Lawton, J. H., Shorrocks, B., & Wood, S.(1998). Making mistakes when predicting shifts in species range in re-sponse to global warming. Nature, 391, 783?786.Dell, D., Sparks, T., & Dennis, R. L. H. (2005). Climate change and theeffect of increasing spring temperatures on emergence dates of the but-terfly Apatura iris (Lepidoptera: Nymphalidae). European Journal ofEntomology, 102(2), 161?167.DeLucia, E. H., Nabity, P. D., Zavala, J. A., & Berenbaum, M. R. (2012).Climate change: resetting plant-insect interactions. Plant physiology,160(4), 1677?1685.Dennis, R. L. H. (1993). Butterflies and climate change. Manchester, UK:Manchester University Press.Diamond, S. E., Frame, A. M., Martin, R. A., & Buckley, L. B. (2011).Species? traits predict phenological responses to climate change in butter-flies. Ecology, 92(5), 1005?1012.Diez, J. M., Iba?n?ez, I., MillerRushing, A. J., Mazer, S. J., Crimmins, T. M.,Crimmins, M. A., Bertelsen, C. D., & Inouye, D. W. (2012). Forecast-ing phenology: from species variability to community patterns. EcologyLetters, 15, 545?553.Doi, H., Gordo, O., & Katano, I. (2008). Heterogeneous intra-annual cli-142matic changes drive different phenological responses at two trophic levels.Climate research (Open Access for articles 4 years old and older), 36(3),181.Doi, H. & Takahashi, M. (2008). Latitudinal patterns in the phenologicalresponses of leaf colouring and leaf fall to climate change in Japan. GlobalEcology and Biogeography, 17(4), 556?561.Doligez, B. & Pa?rt, T. (2008). Estimating fitness consequences of dispersal:a road to ?know-where?? Non-random dispersal and the underestimationof dispersers? fitness. Journal of Animal Ecology, 77(6), 1199?1211.Drummond, A., Ashton, B., Buxton, S., Cheung, M., Cooper, A., Duran,C., Field, M., Heled, J., Kearse, M., Markowitz, S., Moir, R., Stones-Havas, S., Sturrock, S., Thierer, T., & Wilson, A. (2012). Geneious v5.6.Available from http://www.geneious.com.Dukes, J. S., Pontius, J., Orwig, D., Garnas, J. R., Rodgers, V. L., Brazee,N., Cooke, B., Theoharides, K. A., Stange, E. E., & Harrington, R. (2009).Responses of insect pests, pathogens, and invasive plant species to climatechange in the forests of northeastern North America: What can we pre-dict? Canadian Journal of Forest Research, 39(2), 231?248.Durant, J. M., Hjermann, D. ?., AnkerNilssen, T., Beaugrand, G., Mys-terud, A., Pettorelli, N., & Stenseth, N. C. (2005). Timing and abundanceas key mechanisms affecting trophic interactions in variable environments.Ecology Letters, 8(9), 952?958.143Dury, S., Good, J., Perrins, C., Buse, A., & Kaye, T. (1998). The effects ofincreasing CO2 and temperature on oak leaf palatability and the implica-tions for herbivorous insects. Global Change Biology, 4(1), 55?61.Edwards, M. & Richardson, A. J. (2004). Impact of climate change onmarine pelagic phenology and trophic mismatch. Nature, 430(7002), 881?884.Elith, J., Graham, C. H., Anderson, R. P., Dudik, M., Ferrier, S., Guisan,A., Hijmans, R. J., Huettmann, F., Leathwick, J. R., Lehmann, A., Li, J.,Lohmann, L. G., Loiselle, B. A., Manion, G., Moritz, C., Nakamura, M.,Nakazawa, Y., Overton, J. M., Peterson, A. T., Phillips, S. J., Richard-son, K. S., Scachetti-Pereira, R., Schapire, R. E., Soberon, J., Williams,S., Wisz, M. S., & Zimmermann, N. E. (2006). Novel methods improveprediction of species? distributions from occurrence data. Ecography, 29,129?151.Ellwood, E. R., Diez, J. M., Iba?n?ez, I., Primack, R. B., Kobori, H., Higuchi,H., & Silander, J. A. (2012). Disentangling the paradox of insect phenol-ogy: are temporal trends reflecting the response to warming? Oecologia,168, 1161?1171.EncinasViso, F., Revilla, T. A., & Etienne, R. S. (2012). Phenology drivesmutualistic network structure and diversity. Ecology Letters, 15(3), 198?208.Feeny, P. (1970). Seasonal changes in oak leaf tannins and nutrients as acause of spring feeding by winter moth caterpillars. Ecology, 51, 565?581.144Felsenstein, J. (1985). Phylogenies and the comparative method. AmericanNaturalist, (pp. 1?15).Fielding, A. H. & Bell, J. F. (1997). A review of methods for the assessmentof prediction errors in conservation presence/absence models. Environ-mental Conservation, 24, 38?49.Fitter, A. & Fitter, R. (2002). Rapid changes in flowering time in britishplants. Science, 296(5573), 1689.Fitzgerald, T. D. (1995). The tent caterpillars. Ithaca, NY: Cornell Univer-sity Press.Forkner, R. E., Marquis, R. J., Lill, J. T., & Corff, J. L. (2008). Tim-ing is everything? Phenological synchrony and population variability inleafchewing herbivores of Quercus. Ecological Entomology, 33(2), 276?285.Forrest, J. R. & Thomson, J. D. (2011). An examination of synchronybetween insect emergence and flowering in Rocky Mountain meadows.Ecological Monographs, 81(3), 469?491.Freckleton, R., Harvey, P., & Pagel, M. (2002). Phylogenetic analysis andcomparative data: a test and review of evidence. The American Naturalist,160(6), 712?726.Freedman, A. H., Buermann, W., Lebreton, M., Chirio, L., & SMITH, T. B.(2009). Modeling the effects of anthropogenic habitat change on savannasnake invasions into African rainforest. Conservation Biology, 23(1), 81?92.145Futuyma, D. J. (1998). Evolutionary Biology. Sunderland, Massachusetts:Sinauer Associates, 3rd edition.Gallagher, R., Hughes, L., & Leishman, M. (2009). Phenological trendsamong Australian alpine species: using herbarium records to identifyclimate-change indicators. Australian Journal of Botany, 57(1), 1?9.Gaston, K. J. (2003). The structure and dynamics of geographic ranges.Oxford, UK: Oxford University Press.Gillespie, R. G., Baldwin, B. G., Waters, J. M., Fraser, C. I., Nikula, R.,& Roderick, G. K. (2012). Long-distance dispersal: a framework for hy-pothesis testing. Trends in Ecology & Evolution, 27(1), 47?56.Gilman, S. E., Urban, M. C., Tewksbury, J., Gilchrist, G. W., & Holt, R. D.(2010). A framework for community interactions under climate change.Trends in Ecology & Evolution, 25(6), 325?331.Good, R. D. O. (1931). A theory of plant geography. New Phytologist, 30,149?149.Gordo, O. & Sanz, J. J. (2005). Phenology and climate change: a long-termstudy in a Mediterranean locality. Oecologia, 146(3), 484?495.Gordo, O. & Sanz, J. J. (2010). Impact of climate change on plant phenologyin Mediterranean ecosystems. Global Change Biology, 16(3), 1082?1106.Graham, C. H., Silva, N., & Vela?squezTibata?, J. (2010). Evaluating thepotential causes of range limits of birds of the Colombian andes. Journalof Biogeography, 37(10), 1863?1875.146Guisan, A., Graham, C. H., Elith, J., Huettmann, F., Dudik, M., Ferrier, S.,Hijmans, R., Lehmann, A., Li, J., Lohmann, L. G., Loiselle, B., Manion,G., Moritz, C., Nakamura, M., Nakazawa, Y., Overton, J. M., Peterson,A. T., Phillips, S. J., Richardson, K., Scachetti-Pereira, R., Schapire,R. E. W., Wisz, M. S., & Zimmermann, N. E. (2007). Sensitivity ofpredictive species distribution models to change in grain size sensitivityof predictive species distribution models to change in grain size. Diversityand Distributions, 13, 332?340.Guisan, A. & Hofer, U. (2003). Predicting reptile distributions at themesoscale: relation to climate and topography. Journal of Biogeography,30(8), 1233?1243.Guisan, A. & Thuiller, W. (2005). Predicting species distribution: offeringmore than simple habitat models. Ecology Letters, 8(9), 993?1009.Guisan, A. & Zimmermann, N. E. (2000). Predictive habitat distributionmodels in ecology. Ecological Modelling, 135(2-3), 147?186.Guppy, C. S. & Shepard, J. H. (2001). Butterflies of British Columbia.Vancouver, BC: University of British Columbia Press.Hamann, A., Namkoong, G., & Koshy, M. (2001). Multiple populationbreeding for uncertain climatic futures with alnus rubra: ecological genet-ics and selection experiments. Forestry Sciences, 70, 321?330.Hampe, A. (2004). Bioclimate envelope models: what they detect and whatthey hide. Global Ecology and Biogeography, 13(5), 469?471.147Hanspach, J., Ku?hn, I., Schweiger, O., Pompe, S., & Klotz, S. (2011). Ge-ographical patterns in prediction errors of species distribution models.Global Ecology and Biogeography, 20(5), 779?788.Harley, C. D. (2011). Climate change, keystone predation, and biodiversityloss. Science, 334(6059), 1124?1127.Harmon, J. P., Moran, N. A., & Ives, A. R. (2009). Species response toenvironmental change: impacts of food web interactions and evolution.Science, 323(5919), 1347?1350.Harrington, C. A., Zasada, J. C., & Allen, E. A. (1994). Biology of red alder(Alnus rubra Bong.), (pp. 3?22). The biology and management of redalder. Oregon State University Press: Corvallis, Oregon.Harrington, R., Woiwod, I., & Sparks, T. (1999). Climate change and trophicinteractions. Trends in Ecology & Evolution, 14(4), 146?150.Harvey, P. H. & Pagel, M. (1991). The Comparative Method in EvolutionaryBiology. Oxford, UK: Oxford University Press.Hassall, C., Thompson, D. J., French, G. C., & Harvey, I. F. (2007). Histor-ical changes in the phenology of British Odonata are related to climate.Global Change Biology, 13(5), 933?941.Hegland, S. J., Nielsen, A., La?zaro, A., Bjerknes, A. L., & Totland, ?.(2009). How does climate warming affect plant-pollinator interactions?Ecology Letters, 12(2), 184?195.148Heikkila?, M., Kaila, L., Mutanen, M., Pen?a, C., & Wahlberg, N. (2012). Cre-taceous origin and repeated tertiary diversification of the redefined butter-flies. Proceedings of the Royal Society B: Biological Sciences, 279(1731),1093?1999.Heikkinen, R. K., Luoto, M., Virkkala, R., Pearson, R. G., & Korber, J. H.(2007). Biotic interactions improve prediction of boreal bird distributionsat macro-scales. Global Ecology and Biogeography, 16, 754?763.Hellmann, J. J. (2002). The effect of an environmental change on mobilebutterfly larvae and the nutritional quality of their hosts. Journal ofAnimal Ecology, 71(6), 925?936.Hellmann, J. J., Prior, K. M., & Pelini, S. L. (2012). The influence of speciesinteractions on geographic range change under climate change. Annals ofthe New York Academy of Sciences, 1249(1), 18?28.Hickling, R., Roy, D. B., Hill, J. K., Fox, R., & Thomas, C. D. (2006). Thedistribution of a wide range of taxonomic groups are expanding polewards.Global Change Biology, 12, 450?455.Hillis, D., Moritz, C., & Mable, B. K. (1996). Molecular Systematics. Sun-derland, MA: Sinauer Associates.Hobbs, R. J., Higgs, E., & Harris, J. A. (2009). Novel ecosystems: impli-cations for conservation and restoration. Trends in Ecology & Evolution,24(11), 599?605.Hodgson, J. A., Thomas, C. D., Oliver, T. H., Anderson, B. J., Brereton, T.,149& Crone, E. (2011). Predicting insect phenology across space and time.Global Change Biology, 17(3), 1289?1300.Hone?k, A. (1993). Intraspecific variation in body size and fecundity in in-sects: a general relationship. Oikos, 66, 483?492.Huey, R. B., Carlson, M., Crozier, L., Frazier, M., Hamilton, H., Harley, C.,Hoang, A., & Kingsolver, J. G. (2002). Plants versus animals: do theydeal with stress in different ways? Integrative and Comparative Biology,42(3), 415?423.Hunter, A. F. (1993). Gypsy moth population sizes and the window ofopportunity in spring. Oikos, 68, 531?538.Hunter, M. D. (1990). Differential susceptibility to variable plant phenol-ogy and its role in competition between two insect herbivores on oak.Ecological Entomology, 15(4), 401?408.Hunter, M. D. & Price, P. W. (1992). Playing chutes and ladders: hetero-geneity and the relative roles of bottom-up and top-down forces in naturalcommunities. Ecology, 73, 724?732.Huntley, B., Collingham, Y. C., Green, R. E., Hilton, G. M., Rahbek, C., &Willis, S. G. (2006). Potential impacts of climatic change upon geograph-ical distributions of birds. Ibis, 148(s1), 8?28.Huntley, B., Collingham, Y. C., Willis, S. G., & Green, R. E. (2008). Poten-tial impacts of climatic change on european breeding birds. PLoS One,3(1), e1439.150Huntley, B., Green, R. E., Hill, J. K., Willis, S. G., Bartlein, P. J., Cramer,W., Hagemaijer, W. J. M., & Thomas, C. D. (2004). The performance ofmodels relating species geographical distributions to climate is indepen-dent of trophic level. Ecology Letters, 7, 417?426.Hurlbert, A. H. & Liang, Z. (2012). Spatiotemporal variation in avian mi-gration phenology: citizen science reveals effects of climate change. PloSone, 7(2), e31662.Hutchinson, G. E. (1957). Concluding remarks. Cold Spring Habor Symposiaon Quantitative Biology, 22, 415?427.Iler, A. M., H?ye, T. T., Inouye, D. W., & Schmidt, N. M. (2013). Long-term trends mask variation in the direction and magnitude of short-termphenological shifts. American Journal of Botany, 100, 1398?1406.Illan, J. G., Gutierrez, D., Diez, S. B., & Wilson, R. J. (2012). Elevationaltrends in butterfly phenology: implications for species responses to climatechange. Ecological Entomology, 37(2), 134?144.Inouye, D. (2000). The ecological and evolutionary significance of frost inthe context of climate change. Ecology Letters, 3(5), 457?463.Inouye, D. W. (2008). Effects of climate change on phenology, frost damage,and floral abundance of montane wildflowers. Ecology, 89(2), 353?362.IPCC (2007). Climate Change 2007: The Physical Science Basis. Contri-bution of Working Group I to the Fourth Assessment Report of the In-151tergovernmental Panel on Climate Change. Technical report, CambridgeUniversity Press.Ives, A. R. & Carpenter, S. R. (2007). Stability and diversity of ecosystems.Science, 317(5834), 58?62.Jepsen, J. U., Kapari, L., Hagen, S. B., Schott, T., Vinstad, O. P. L., Nilssen,A. C., & Ims, R. A. (2011). Rapid northwards expansion of a forest insectpest attributed to spring phenology matching with subarctic birch. GlobalChange Biology, 17(6), 2071?2083.Katoh, K. & Toh, H. (2008). Recent developments in the mafft multiplesequence alignment program. Briefings in bioinformatics, 9(4), 286?298.Keane, T. M., Creevey, C. J., Pentony, M. M., Naughton, T. J., & Mclnerney,J. O. (2006). Assessment of methods for amino acid matrix selection andtheir use on empirical data shows that ad hoc assumptions for choice ofmatrix are not justified. BMC evolutionary biology, 6(1), 29.Kerslake, J. & Hartley, S. (1997). Phenology of winter moth feeding on com-mon heather: effects of source population and experimental manipulationof hatch dates. Journal of Animal Ecology, 66, 375?385.Klapwijk, M. J., Grobler, B. C., Ward, K., Wheeler, D., & Lewis, O. T.(2010). Influence of experimental warming and shading on host?parasitoidsynchrony. Global Change Biology, 16(1), 102?112.Klinkenberg, B. (2012). E-fauna BC: Electronic Atlas of the Fauna of BritishColumbia.152Lane, J. E., Kruuk, L. E., Charmantier, A., Murie, J. O., & Dobson, F. S.(2012). Delayed phenology and reduced fitness associated with climatechange in a wild hibernator. Nature, 489, 554?557.Lavergne, S., Mouquet, N., Thuiller, W., & Ronce, O. (2010). Biodiversityand climate change: integrating evolutionary and ecological responses ofspecies and communities. Annual Review of Ecology, Evolution, and Sys-tematics, 41, 321?350.Lavoie, C. (2013). Biological collections in an ever changing world: Herbariaas tools for biogeographical and environmental studies. Perspectives inPlant Ecology, Evolution and Systematics, 15, 68?76.Lavoie, C. & Lachance, D. (2006). A new herbarium-based method forreconstructing the phenology of plant species across large areas. AmericanJournal of Botany, 93(4), 512.Layberry, R. A., Hall, P. W., & Lafontaine, J. D. (1998). The Butterflies ofCanada. Toronto: NRC Research Press, Canada Institute for Scientificand Technical Information, University of Toronto Press.Leather, S. R. (1988). Size, reproductive potential and fecundity in insects:things aren?t as simple as they seem. Oikos, 51(3), 386?389.Lester, S. E., Ruttenberg, B. I., Gaines, S. D., & Kinlan, B. P. (2007). Therelationship between dispersal ability and geographic range size. EcologyLetters, 10(8), 745?758.Lincoln, D. E., Fajer, E. D., & Johnson, R. H. (1993). Plant-insect her-153bivore interactions in elevated CO2 environments. Trends in Ecology &Evolution, 8(2), 64?68.Liu, Y., Reich, P. B., Li, G., & Sun, S. (2011). Shifting phenology andabundance under experimental warming alters trophic relationships andplant reproductive capacity. Ecology, 92(6), 1201?1207.Lobo, J. M., Jime?nezValverde, A., & Real, R. (2008). AUC: a mislead-ing measure of the performance of predictive distribution models. GlobalEcology and Biogeography, 17(2), 145?151.Luoto, M., PoI`yry, J., Heikkinen, R. K., & Saarinen, K. (2005). Uncertaintyof bioclimate envelope models based on the geographical distribution ofspecies. Global Ecology and Biogeography, 14, 575?584.Luoto, M., Virkkala, R., & Heikkinen, R. K. (2007). The role of land coverin bioclimatic models depends on spatial resolution. Global Ecology andBiogeography, 16(1), 34?42.Maclean, I. M. & Wilson, R. J. (2011). Recent ecological responses to climatechange support predictions of high extinction risk. Proceedings of theNational Academy of Sciences, 108(30), 12337?12342.Marion, G., Henry, G., Freckman, D., Johnstone, J., Jones, G., Jones, M.,Levesque, E., Molau, U., M?lgaard, P., & Parsons, A. (1997). Open-topdesigns for manipulating field temperature in highlatitude ecosystems.Global Change Biology, 3(S1), 20?32.Marsico, T. D. & Hellmann, J. J. (2009). Dispersal limitation inferred from154an experimental translocation of Lomatium (Apiaceae) species outsidetheir geographic ranges. Oikos, 118(12), 1783?1792.Martel, J., Hanhima?ki, S., Kause, A., & Haukioja, E. (2001). Diversity ofbirch sawfly responses to seasonally atypical diets. Entomologia Experi-mentalis et Applicata, 100(3), 301?309.McCann, K. (2007). Protecting biostructure. Nature, 446(7131), 29?29.McLaughlin, J. F., Hellmann, J. J., Boggs, C. L., & Ehrlich, P. R. (2002).Climate change hastens population extinctions. Proceedings of the Na-tional Academy of Sciences, 99(9), 6070?6074.McPherson, J., Jetz, W., & Rogers, D. J. (2004). The effects of species? rangesizes on the accuracy of distribution models: ecological phenomenon orstatistical artefact? Journal of Applied Ecology, 41(5), 811?823.McPherson, J. M. & Jetz, W. (2007). Effects of species? ecology on theaccuracy of distribution models. Ecography, 30(1), 135?151.Memmott, J., Craze, P. G., Waser, N. M., & Price, M. V. (2007). Globalwarming and the disruption of plant-pollinator interactions. Ecology Let-ters, 10, 710?717.Menzel, A., Sparks, T. H., Estrella, N., Koch, E., Aasa, A., Ahas, R.,ALMKU?BLER, K., Bissolli, P., Braslavska?, O., & Briede, A. (2006).European phenological response to climate change matches the warmingpattern. Global Change Biology, 12(10), 1969?1976.155Merriam, C. H. (1894). Laws of temperature control of the geographicdistribution of terrestrial animals and plants. National Geographic, 6,229?238.Miller-Rushing, A. J., H?ye, T. T., Inouye, D. W., & Post, E. (2010). Theeffects of phenological mismatches on demography. Philosophical Trans-actions of the Royal Society B: Biological Sciences, 365(1555), 3177?3186.Miller-Rushing, A. J., Inouye, D. W., & Primack, R. B. (2008). How welldo first flowering dates measure plant responses to climate change? theeffects of population size and sampling frequency. Journal of Ecology,96(6), 1289?1296.Miller-Rushing, A. J. & Primack, R. B. (2008). Global warming and flow-ering times in Thoreau?s Concord: a community perspective. Ecology,89(2), 332?341.Moritz, C., Patton, J. L., Conroy, C. J., Parra, J. L., White, G. C., &Beissinger, S. R. (2008). Impact of a century of climate change onsmall-mammal communities in Yosemite National Park, USA. Science,322(5899), 261?264.Moussus, J. P., Clavel, J., Jiguet, F., & Julliard, R. (2011). Which arethe phenologically flexible species? A case study with common passerinebirds. Oikos, 120(7), 991?998.Myers, J. H. (2000). Population fluctuations of the western tent caterpillarin southwestern British Columbia. Population Ecology, 42(3), 231?241.156O?Connor, M. (2009). Warming strengthens an herbivore-plant interaction.Ecology, 90(2), 388?398.O?Connor, M. I., Gilbert, B., & Brown, C. J. (2011). Theoretical predictionsfor how temperature affects the dynamics of interacting herbivores andplants. The American Naturalist, 178(5), 626?638.Opler, P. A., Lotts, K., & Naberhaus, T. (2012). Butterflies and Moths ofNorth America.Paradis, E., Baillie, S. R., Sutherland, W. J., & Gregory, R. D. (1998).Patterns of natal and breeding dispersal in birds. Journal of AnimalEcology, 67(4), 518?536.Paradis, E. & Claude, J. (2002). Analysis of comparative data using general-ized estimating equations. Journal of theoretical biology, 218(2), 175?185.Paradis, E., Claude, J., & Strimmer, K. (2004). APE: analyses of phyloge-netics and evolution in R language. Bioinformatics, 20(2), 289?290.Parmesan, C. (2006). Ecological and evolutionary responses to recent cli-mate change. Annual Review of Ecology, Evolution, and Systematics, 37,637?669.Parmesan, C. (2007). Influences of species, latitudes and methodologieson estimates of phenological response to global warming. Global ChangeBiology, 13(9), 1860?1872.Parmesan, C., Duarte, C., Poloczanska, E., Richardson, A. J., & Singer,157M. C. (2011). Overstretching attribution. Nature Climate Change, 1(1),2?4.Parmesan, C. & Yohe, G. (2003). A globally coherent fingerprint of climatechange impacts across natural systems. Nature, 421, 37?43.Parra, J. L. & Monahan, W. B. (2008). Variability in 20th century climatechange reconstructions and its consequences for predicting geographic re-sponses of California mammals. Global Change Biology, 14(10), 2215?2231.Parry, D., Spence, J. R., & Volney, W. J. A. (1998). Budbreak phenologyand natural enemies mediate survival of first-instar forest tent caterpillar(Lepidoptera: Lasiocampidae). Environmental Entomology, 27(6), 1368?1374.Pateman, R. M., Hill, J. K., Roy, D. B., Fox, R., & Thomas, C. D. (2012).Temperature-dependent alterations in host use drive rapid range expan-sion in a butterfly. Science, 336(6084), 1028?1030.Pau, S., Wolkovich, E. M., Cook, B. I., Davies, T. J., Kraft, N. J. B.,Bolmgren, K., Betancourt, J. L., & Cleland, E. E. (2011). Predictingphenology by integrating ecology, evolution and climate science. GlobalChange Biology, 17, 3633?3643.Paul, J. R., Morton, C., Taylor, C. M., & Tonsor, S. J. (2009). Evolutionarytime for dispersal limits the extent but not the occupancy of species?potential ranges in the tropical plant genus Psychotria (Rubiaceae). TheAmerican Naturalist, 173(2), 188?199.158Pearce, J. & Ferrier, S. (2000). An evaluation of alternative algorithmsfor fitting species distribution models using logistic regression. EcologicalModelling, 128(2), 127?147.Pearson, R. G. & Dawson, T. P. (2003). Predicting the impacts of climatechange on the distribution of species: Are bioclimate envelope modelsuseful? Global Ecology and Biogeography, 12(5), 361?371.Pearson, R. G., Thuiller, W., ArauI`jo, M. B., Martinez-Meyer, E., Brotons,L., McClean, C., Miles, L., Segurado, P., Dawson, T. P., & Lees, D. C.(2006). Model-based uncertainty in species range prediction. Journal ofBiogeogrpahy, 33, 1704?1711.Pelham, J. P. (2011). A catalogue of the butterflies of the United States andCanada. Technical report, University of Washington.Pelini, S. L., Dzurisin, J. D., Prior, K. M., Williams, C. M., Marsico, T. D.,Sinclair, B. J., & Hellmann, J. J. (2009). Translocation experiments withbutterflies reveal limits to enhancement of poleward populations underclimate change. Proceedings of the National Academy of Sciences, 106(27),11160?11165.Phillips, S. J., Anderson, R. P., & Schapire, R. E. (2006). Maximum entropymodeling of species geographic distributions. Ecological Modelling, 190,231?259.Pinheiro, J. C., Bates, D., DebRoy, S., & Sarkar, D. (2012). nlme: Linearand nonlinear mixed effects models.159Pinheiro, J. C. & Bates, D. M. (2000). Mixed-effects models in S and S-PLUS. New York, NY: Springer.Polgar, C. A., Primack, R. B., Williams, E. H., Stichter, S., & Hitchcock,C. (2013). Climate effects on the flight period of Lycaenid butterflies inMassachusetts. Biological Conservation, 160, 25?31.Pollard, E. (1988). Temperature, rainfall and butterfly numbers. Journal ofApplied Ecology, (pp. 819?828).Pop, E. W., Oberbauer, S. F., & Starr, G. (2000). Predicting vegetativebud break in two arctic deciduous shrub species, Salix pulchra and Betulanana. Oecologia, 124(2), 176?184.Post, E., Pedersen, C., Wilmers, C. C., & Forchhammer, M. C. (2008).Warming, plant phenology and the spatial dimension of trophic mismatchfor large herbivores. Proceedings of the Royal Society B: Biological Sci-ences, 275(1646), 2005?2013.Po?yry, J., Luoto, M., Heikkinen, R. K., & Saarinen, K. (2008). Species traitsare associated with the quality of bioclimatic models. Global Ecology andBiogeography, 17(3), 403?414.Primack, D., Imbres, C., Primack, R. B., Miller-Rushing, A. J., & Tredici,P. D. (2004). Herbarium specimens demonstrate earlier flowering timesin response to warming in Boston. American Journal of Botany, 91(8),1260?1264.Primack, R. B., Iba?n?ez, I., Higuchi, H., Lee, S. D., Miller-Rushing, A. J.,160Wilson, A. M., & Silander, J. A. (2009). Spatial and interspecific vari-ability in phenological responses to warming temperatures. BiologicalConservation, 142(11), 2569?2577.Pulliam, H. R. (2000). On the relationship between niche and distribution.Ecology Letters, 3(4), 349?361.Pyle, R. M. (2002). The Butterflies of Cascadia. Seattle, USA: SeattleAudubon Society.Revell, L. J. (2010). Phylogenetic signal and linear regression on speciesdata. Methods in Ecology and Evolution, 1(4), 319?329.Robbirt, K. M., Davy, A. J., Hutchings, M. J., & Roberts, D. L. (2011).Validation of biological collections as a source of phenological data for usein climate change studies: a case study with the orchid Ophrys sphegodes.Journal of Ecology, 99(1), 235?241.Robinson, E. A., Ryan, G. D., & Newman, J. A. (2012). A metaanalyti-cal review of the effects of elevated CO2 on plant-arthropod interactionshighlights the importance of interacting environmental and biological vari-ables. New Phytologist, 194, 321?336.Rock, G., Stinner, R., Bacheler, J., Hull, L., & Hogmire, H. (1993). Predict-ing geographical and within-season variation in male flights of four fruitpests. Environmental Entomology, 22(4), 716?725.Romo, J. & Eddleman, L. (1995). Use of degree-days in multiple-temperature experiments. Journal of Range Management, (pp. 410?416).161Ronquist, F., Maxim, T., van der Mark, P., Ayres, D., Darling, A., Ho?hna,S., Larget, B., Liu, L., Suchard, M., & Huelsenbeck, J. (2012). Mrbayes3.2: efficient Bayesian phylogenetic inference and model choice across alarge model space. Systematic Biology, 61(3), 539?542.Rousseaux, M. C., Julkunen-Tiitto, R., Searles, P. S., Scopel, A. L., Aphalo,P. J., & Ballare?, C. L. (2004). Solar UV-B radiation affects leaf qualityand insect herbivory in the southern beech tree Nothofagus antarctica.Oecologia, 138(4), 505?512.Roux, P. C. L. & McGeoch, M. A. (2008). Rapid range expansion andcommunity reorganization in response to warming. Global Change Biology,14(12), 2950?2962.Roy, D. & Sparks, T. (2000). Phenology of British butterflies and climatechange. Global Change Biology, 6(4), 407?416.Sala, O. E., Chapin, F. S., Armesto, J. J., Berlow, E., Bloomfield, J., Dirzo,R., Huber-Sanwald, E., Huenneke, L. F., Jackson, R. B., Kinzig, A., Lee-mans, R., Lodge, D. M., Mooney, H. A., Oesterheld, M., Poff, N. L.,Sykes, M. T., Walker, B. H., Walker, M., & Wall, D. H. (2000). Globalbiodiversity scenarios for the year 2100. Science, 287, 1770?1774.Sanders, N. J. & Gordon, D. M. (2003). Resource-dependent interactionsand the organization of desert ant communities. Ecology, 84, 1024?1031.Santika, T. (2011). Assessing the effect of prevalence on the predictiveperformance of species distribution models using simulated data. GlobalEcology and Biogeography, 20(1), 181?192.162Schliep, K. (2011). PHANGORN: phylogenetic analysis in R. Bioinformat-ics, 27(4), 592?593.Schneider, C. (2003). The influence of spatial scale on quantifying insectdispersal: an analysis of butterfly data. Ecological Entomology, 28(2),252?256.Schweiger, O., Biesmeijer, J. C., Bommarco, R., Hickler, T., Hulme, P. E.,Klotz, S., Ku?hn, I., Moora, M., Nielsen, A., & Ohlemu?ller, R. (2010).Multiple stressors on biotic interactions: how climate change and alienspecies interact to affect pollination. Biological Reviews, 85(4), 777?795.Schweiger, O., Settele, J., Kudrna, O., Klotz, S., & Ku?hn, I. (2008). Cli-mate change can cause spatial mismatch of trophically interacting species.Ecology, 89(12), 3472?3479.Scriber, J. & Slansky, F. (1981). The nutritional ecology of immature insects.Annual Review of Entomology, 26(1), 183?211.Segurado, P., Arau?jo, M. B., & Kunin, W. (2006). Consequences of spatialautocorrelation for niche-based models. Journal of Applied Ecology, 43(3),433?444.Sexton, J. P., McIntyre, P. J., Angert, A. L., & Rice, K. J. (2009). Evolutionand ecology of species range limits. Annual Review of Ecology, Evolution,and Systematics, 40, 415?436.Sinervo, B., Mendez-De-La-Cruz, F., Miles, D. B., Heulin, B., Bastiaans,E., Cruz, M. V.-S., Lara-Resendiz, R., Mart??nez-Me?ndez, N., Caldero?n-163Espinosa, M. L., & Meza-La?zaro, R. N. (2010). Erosion of lizard diversityby climate change and altered thermal niches. Science, 328(5980), 894?899.Singer, M. C. (1972). Complex components of habitat suitability within abutterfly colony. Science, 176(4030), 75?77.Skaug, H., Fournier, D., Nielsen, A., Magnusson, A., & Bolker, B. (2012).glmmADMB: Generalized Linear Mixed Models Using AD Model Builder.Smith, G. & Raske, A. (1968). Starvation experiments with first instar foresttent caterpillar larvae. Canadian Department of Fisheries and ForestryBimonthly Research Note, 24, 39.Smith, M. A. & Green, D. M. (2005). Dispersal and the metapopulationparadigm in amphibian ecology and conservation: are all amphibian pop-ulations metapopulations? Ecography, 28(1), 110?128.Sparks, T., Huber, K., & Dennis, R. L. H. (2006). Complex phenologicalresponses to climate warming trends? Lessons from history. EuropeanJournal of Entomology, 103(2), 379?386.Sparks, T. & Yates, T. (1997). The effect of spring temperature on theappearance dates of British butterflies 1883?1993. Ecography, 20(4), 368?374.Sparks, T. H., Gorska-Zajaczkowska, M., Wojtowicz, W., & Tryjanowski, P.(2011). Phenological changes and reduced seasonal synchrony in westernPoland. International Journal of Biometeorology, 55(3), 447?453.164Stefanescu, C., Penuelas, J., & Filella, I. (2003). Effects of climatic changeon the phenology of butterflies in the northwest Mediterranean Basin.Global Change Biology, 9, 1494?1506.Stenseth, N. C. & Mysterud, A. (2002). Climate, changing phenology, andother life history traits: nonlinearity and match?mismatch to the environ-ment. Proceedings of the National Academy of Sciences, 99(21), 13379?13381.Stephens, P. A., Buskirk, S. W., & del Rio, C. M. (2007). Inference inecology and evolution. Trends in Ecology & Evolution, 22(4), 192?197.Suding, K. N., Lavorel, S., Chapin, F. S., Cornelissen, J. H., DIAz, S.,Garnier, E., Goldberg, D., Hooper, D. U., Jackson, S. T., & NAVAS,M. (2008). Scaling environmental change through the community-level:a trait-based response-and-effect framework for plants. Global ChangeBiology, 14(5), 1125?1140.Sutherland, G. D., Harestad, A. S., Price, K., & Lertzman, K. P. (2000).Scaling of natal dispersal distances in terrestrial birds and mammals. Con-servation Ecology, 4(1), 16.Suttle, K. B., Thomsen, M. A., & Power, M. E. (2007). Species interactionsreverse grassland responses to changing climate. Science, 315, 640?642.Swets, J. A. (1988). Measuring the accuracy of diagnostic systems. Science,240, 1285?1293.Syphard, A. D. & Franklin, J. (2010). Species traits affect the performance165of species distribution models for plants in southern California. Journalof Vegetation Science, 21(1), 177?189.Team, R. D. C. (2012). R: A language and environment for statistical com-puting.Thackeray, S. J., Sparks, T. H., Frederiksen, M., Burthe, S., Bacon, P. J.,Bell, J. R., Botham, M. S., Brereton, T. M., Bright, P. W., & Carvalho,L. (2010). Trophic level asynchrony in rates of phenological change formarine, freshwater and terrestrial environments. Global Change Biology,16(12), 3304?3313.Thomas, C. D., Cameron, A., Green, R. E., Bakkenes, M., Beaumont, L. J.,Collingham, Y. C., Erasmus, B. F., Siqueira, M. F. D., Grainger, A.,& Hannah, L. (2004). Extinction risk from climate change. Nature,427(6970), 145?148.Thuiller, W. (2003). BIOMOD - Optimizing predictions of species distribu-tions and projecting potential future shifts under global change. GlobalChange Biology, 9, 1353?1362.Thuiller, W., Araujo, M. B., & Lavorel, S. (2003). Generalized models vs.classification tree analysis: predicting spatial distributions of plant speciesat different scales. Journal of Vegetation Science, 14(5), 669?680.Thuiller, W., Broennimann, O., Hughes, G., ALKEMADE, J. R. M., MIDG-LEY, G. F., & Corsi, F. (2006). Vulnerability of African mammals toanthropogenic climate change under conservative land transformation as-sumptions. Global Change Biology, 12(3), 424?440.166Thuiller, W., Brotons, L., AraA?jo, M. B., & Lavorel, S. (2004a). Effectsof restricting environmental range of data to project current and futurespecies distributions. Ecography, 27(2), 165?172.Thuiller, W., Lavorel, S., Araujo, M. B., Sykes, M. T., & Prentice, I. C.(2005). Climate change threats to plant diversity in Europe. Proceedingsof the National Academy of Sciences, 102(23), 8245?8250.Thuiller, W., Lavorel, S., Midgley, G., Lavergne, S., & Rebelo, T. (2004b).Relating plant traits and species distributions along bioclimatic gradientsfor 88 Leucadendron taxa. Ecology, 85(6), 1688?1699.Tobin, P. C., Nagarkatti, S., Loeb, G., & Saunders, M. C. (2008). Historicaland projected interactions between climate change and insect voltinismin a multivoltine species. Global Change Biology, 14(5), 951?957.Tsoar, A., Allouche, O., Steinitz, O., Rotem, D., & Kadmon, R. (2007).A comparative evaluation of presence-only methods for modelling speciesdistribution. Diversity and Distributions, 13(4), 397?405.Tylianakis, J. M., Didham, R. K., Bascompte, J., & Wardle, D. A. (2008).Global change and species interactions in terrestrial ecosystems. EcologyLetters, 11(12), 1351?1363.Valtonen, A., Ayres, M. P., Roininen, H., Po?yry, J., & Leinonen, R. (2011).Environmental controls on the phenology of moths: predicting plasticityand constraint under climate change. Oecologia, 165(1), 237?248.van Riper III, C., van Riper, S. G., Goff, M. L., & Laird, M. (1986). The167epizootiology and ecological significance of malaria in hawaiian land birds.Ecological Monographs, 56(4), 327?344.Vellend, M., Brown, C., Kharouba, H. M., McCune, J., & Myers-Smith,I. (2013). Historical ecology: using unconventional data sources to studytemporal change in plant populations and communities. American Journalof Botany, 100, 1294?1305.Venier, L., Pearce, J., McKee, J., McKenney, D., & Niemi, G. (2004). Cli-mate and satellite-derived land cover for predicting breeding bird distri-bution in the great lakes basin. Journal of Biogeography, 31(2), 315?331.Visser, M. E. (2008). Keeping up with a warming world; assessing therate of adaptation to climate change. Proceedings of the Royal Society B:Biological Sciences, 275(1635), 649?659.Visser, M. E. & Both, C. (2005). Shifts in phenology due to global climatechange: the need for a yardstick. Proceedings of the Royal Society B:Biological Sciences, 272(1581), 2561?2569.Visser, M. E. & Holleman, L. J. M. (2001). Warmer springs disrupt thesynchrony of oak and winter moth phenology. Proceedings of the RoyalSociety of London.Series B: Biological Sciences, 268, 289?294.Voigt, W., Perner, J., Davis, A. J., Eggers, T., Schumacher, J., Ba?hrmann,R., Fabian, B., Heinrich, W., Ko?hler, G., & Lichter, D. (2003). Trophiclevels are differentially sensitive to climate. Ecology, 84(9), 2444?2453.Walther, G., Post, E., Convey, P., Menzel, A., Parmesan, C., Beebee, T.,168Fromentin, J., Hoegh-Guldberg, O., & Bairlein, F. (2002). Ecologicalresponses to recent climate change. Nature, 413, 389?396.Walther, G.-R. (2010). Community and ecosystem responses to recent cli-mate change. Philosophical Transactions of the Royal Society B: BiologicalSciences, 365(1549), 2019?2024.Warren, A., Ogawa, J., & Brower, A. (2008). Phylogenetic relationships ofsubfamilies and circumscription of tribes in the family Hesperiidae (Lep-idoptera: Hesperioidea). Cladistics, 24(5), 642?676.Warren, A., Ogawa, J., & Brower, A. (2009). Revised classification of thefamily Hesperiidae (Lepidoptera: Hesperioidea) based on combined molec-ular and morphological data. Systematic Entomology, 34(3), 467?523.Watkinson, A. R. & Gill, J. A. (2002). Climate change and dispersal, (pp.45?56). Blackwell Science Publishing: Oxford, UK.Webb, C. T., Hoeting, J. A., Ames, G. M., Pyne, M. I., & Poff, N. L. (2010).A structured and dynamic framework to advance traits-based theory andprediction in ecology. Ecology Letters, 13(3), 267?283.Webb III, T. (1986). Is vegetation in equilibrium with climate? How tointerpret late-Quaternary pollen data. Vegetatio, 67(75-91).Wiens, J. J. (2011). The niche, biogeography and species interactions.Philosophical Transactions of the Royal Society B: Biological Sciences,366(1576), 2336?2350.169Williams, J. W. & Jackson, S. T. (2007). Novel climates, no-analog commu-nities, and ecological surprises. Frontiers in Ecology and the Environment,5(9), 475?482.Williams, S. E., Shoo, L. P., Isaac, J. L., Hoffmann, A. A., & Langham, G.(2008). Towards an integrated framework for assessing the vulnerabilityof species to climate change. PLoS biology, 6(12), 2621?2626.Willis, C. G., Ruhfel, B., Primack, R. B., Miller-Rushing, A. J., & Davis,C. C. (2008). Phylogenetic patterns of species loss in Thoreau?s woodsare driven by climate change. Proceedings of the National Academy ofSciences, 105(44), 17029.Wilson, J. (2011). Assessing the value of DNA barcodes for molecular phy-logenetics: Effect of increased taxon sampling in Lepidoptera. PloS one,6(9), e24769.Winder, M. & Schindler, D. E. (2004). Climatic effects on the phenology oflake processes. Global Change Biology, 10(11), 1844?1856.Wint, W. (1983). The role of alternative host-plant species in the life ofa polyphagous moth, Operophtera brumata (Lepidoptera: Geometridae).The Journal of Animal Ecology, 52, 439?450.Wisz, M. S., Hijmans, R., Li, J., Peterson, A. T., Graham, C., & Guisan, A.(2008). Effects of sample size on the performance of species distributionmodels. Diversity and Distributions, 14(5), 763?773.Wolkovich, E. M. & Cleland, E. E. (2010). The phenology of plant inva-170sions: a community ecology perspective. Frontiers in Ecology and theEnvironment, 9(5), 287?294.Wolkovich, E. M., Cook, B., Allen, J., Crimmins, T., Betancourt, J.,Travers, S., Pau, S., Regetz, J., Davies, T., & Kraft, N. (2012). Warmingexperiments underpredict plant phenological responses to climate change.Nature, 485, 494?497.Wood, S. N. (2011). Fast stable restricted maximum likelihood and marginallikelihood estimation of semiparametric generalized linear models. Journalof the Royal Statistical Society (B), 73, 3?36.Yang, L. H. & Rudolf, V. (2010). Phenology, ontogeny and the effects ofclimate change on the timing of species interactions. Ecology Letters,13(1), 1?10.Zavala, J., Scopel, A., & Ballare?, C. (2001). Effects of ambient UV-B radia-tion on soybean crops: Impact on leaf herbivory by Anticarsia gemmatalis.Plant Ecology, 156(2), 121?130.Zuur, A. F., Ieno, E. N., Walker, N. J., Saveliev, A. A., & Smith, G. M.(2009). Mixed effects models and extensions in ecology with R. New York,NY: Springer.171Appendix AMethodological attributes ofstudies172Table A.1: Methodological attributes used to build SDMs for each study used in the analysis. Presented arethe model types, number of climatic variables used in the model, resolution (km2) of the model, total spatialextent over which the model was built (km2) and average absolute latitude of the region for which the modelwas built.Study Model type (s) Variables Resolution(km2)Spatial extent(km2)Latitude(o)Araujo et al. 2005b GAM 7 2500 1.105x107 47Huntley et al. 2006 GAM, locally weighted regres-sion4 12227 2.40x107 0Huntley et al. 2008 locally weighted regression 3 2500 1.105x107 47Huntley et al. 2004 locally weighted regression 3 2500 1.105x107 47Beale et al. 2008 ANN 3 2500 6.04x106 47Thuiller et al. 2006 GAM 6 256 3.02x107 0McPherson & Jetz2007autologistic regression 1-28 2975 8.27x106 47Elith et al. 2006 Mars, gdm, maxent, brt, domain,bruito, GAM, GARP, GLM, bio-clim, lives11 1 1.465 x107 14Heikkinen et al. 2007 GAM 3 100,1600 3.381 x105 64Luoto et al. 2005 GAM 3 100 3.381 x105 47Parra & Monahan2008maxent 19 16 4.240 x105 15Phillips et al. 2006 GARP, maxent 13 30.25 19621904 15Thuiller 2003 GLM, CART 7 2500 5222500 15Continued on next page173Study Model type (s) Variables Resolution(km2)Spatial extent(km2)Latitude(o)Freedman et al. 2009 maxent 7 1 475442 6Guisan & Hofer 2003 GLM 12 1 4.10 x104 47Venier et al. 2004 logistic regression 10 25 8.0 x105 15Pearson et al. 2006 ANN, GARP, GAM, CGM 5 2.56 1.22 x106 15McPherson et al.2004logistic regression, discriminant 61 648 2.77 x106 15Thuiller et al. 2004a GAM 4 2500 6.525 x106 15Thuiller et al. 2003 GLM 7 2500 1.105x107 15174Appendix BTaxonomic attributes ofstudies175Table B.1: Taxonomic attributes of studies used in the analysis. Presented are the total number of unique speciesused, the number of birds, herptiles, butterflies, mammals, and plants.Study Total number of species Birds Herptiles Butterflies Mammals PlantsAraujo et al. 2005b 1778 157 103 0 152 1366Huntley et al. 2006 1457 1457 0 0 0 0Huntley et al. 2008 214 214 0 0 0 0Huntley et al. 2004 173 36 0 37 0 100Beale et al. 2008 42 42 0 0 0 0Thuiller et al. 2006 272 0 0 0 272 0McPherson & Jetz 2007 176 176 0 0 0 0Elith et al. 2006 30 0 0 0 0 30Heikkinen et al. 2007 2 2 0 0 0 0Luoto et al. 2005 79 0 0 79 0 0Parra & Monahan 2008 57 0 0 0 57 0Phillips et al. 2006 2 0 0 0 2 0Thuiller 2003 2 0 0 0 0 2Freedman et al. 2009 3 0 3 0 0 0Guisan & Hofer 2003 8 0 8 0 0 0Venier et al. 2004 10 10 0 0 0 0Pearson et al. 2006 4 0 0 0 0 4McPherson et al. 2004 5 5 0 0 0 0Thuiller et al. 2004a 1 0 0 0 0 1Thuiller et al. 2003 2 0 0 0 0 2176Appendix CAttributes of studies withdispersal distances177Table C.1: Attributes of studies that contained dispersal distances for species in our datasetStudy Taxonomic Group DetailsBowman et al. 2002 mammals Conducted a literature review to find maximum distancemoved by adult mammals after translocation.Cain et al. 1998 plants Measured dispersal distance for Asarum canadense via di-rect observations of seed movement by ants; searched theliterature for measured dispersal distances for other wood-land herbs. Some of these were directly observed and otherswere based on measured fall rates of seeds combined withtypical wind speeds.Paradis et al. 1998 birds Used survey data from the ringing scheme of the BritishTrust for Ornithology 1909-1994. Included only birds ringedand recovered during the breeding season (i.e. excluded mi-gration distances). Estimated both natal and breeding dis-persal distances.Schneider 2003 butterflies Compiled mean distances reported in mark-release-recapture studies.Smith & Green 2005 amphibians Compiled a list of the longest distances moved in both mark-recapture and displacement studies.Sutherland et al. 2000 mammals and birds Compiled data on natal dispersal distances from a literaturesearch. Most data were based on incidental observations.Did not accept data from ?likely migrants?.178Appendix DCollinearity correlationsTable D.1: Collinearity (Spearman r coefficients) between all continuouscovariates (n=4317).latitude area resolutionarea -0.903resolution -0.589 0.503variables 0.277 -0.29 -0.479179Appendix EAnalysis of deviance resultsfor ?large studies?180Table E.1: Analysis of deviance table for the relationship between model accuracy, covariates and taxonomic groupwhen studies that contributed more than half of the total number of species in one taxonomic group were removed.Presented are the differences in AIC and deviance between full and reduced models as well as the associated pvalue. The difference in degrees of freedom between full and reduced models was four for all comparisons andsubsets. The full model includes number of variables, log(spatial extent), resolution and model type. WhenHuntley et al. 2006 is removed there are 2860 species from nineteen published studies, without Araujo et al. 2005there are 2539 species from nineteen published studies and without Luoto et al. 2005 there are 4238 species fromnineteen published studies.Studies subset Model for comparison Model terms ?AIC ?Deviance pWithout Huntley et al. 2006 Just intercept +Taxonomic group 33.46 41.46 <0.0001Full model +Taxonomic group 31.88 39.88 <0.0001Without Araujo et al. 2005 Just intercept +Taxonomic group 14.96 22.96 0.000129Full model* +Taxonomic group 15.18 23.18 0.000117Without Luoto et al. 2005 Just intercept +Taxonomic group 41.6 49.6 <0.0001Full model +Taxonomic group 41.76 49.76 <0.0001* Only includes spatial extent and resolution, none of the other covariates in the model led to estimation181Appendix FDescription of datapreparationButterfly datasetOur main data source was the National Collection of Butterflies database(updated as of January 2011, Layberry et al., 1998). Each collection recordincludes a specimen preserved in one of 40 Canadian natural history col-lections: specimens were collected and identified initially by lepidopteristsand re-verified recently by Lepidopteran systematists (see Layberry et al.,1998). We supplemented the database with additional data for BritishColumbia, Alberta, Yukon and Ontario from the Spencer Entomological Col-lection (UBC; Papilionidae) and the personal and professional collectionsof Canadian butterfly experts. However, not all observations in the non-institutional databases are associated with a specimen i.e. not all individualswere collected after identification. Nevertheless, expert Lepidopterists iden-tified these individuals so data quality should not be greatly affected. Datafrom the main database can be found online (http://www.cbif.gc.ca/spp -pages/butterflies/index e.php) and effort is underway to make the supple-mental database available as well, conditional on the approval of the collec-tors. Before data management, the total database included 250,950 georef-erenced records.182To avoid pseudo-replication, duplicate records of a species on a particu-lar day at a specific location (within ?2km) were considered a single record.We worked at the species level so records identified only to genus level or hy-brids were not used and sub-species were grouped to the species level. Basedon known geographic ranges (Bedford et al., 2012), we removed species sus-pected of having large geographical collection biases (large areas of severeundersampling): Plebejus glandon, Boloria chariclae, Boloria frigga, Ly-caena phlaeas, and Oeneis melissa. Finally, to avoid including erroneousrecords, single records (per location and year) that occurred in December,January, or February on days that had a daily maximum temperature lessthan 10?C were eliminated. We assumed that these dates were incorrect asnone of the species in our database are known to be in flight during thosemonths (Layberry et al., 1998).To isolate the effect of spring temperature on the same year?s adultflight season, we excluded all non-resident species in Canada (migratory,rare strays etc.). In summary, the analysis included ?48,000 georeferencedrecords for 204 species from 1873 to 2010. There were 1792 weather stationsand 6172 unique butterfly sites included that spanned the entire country(9,984,670 km2; Figure F.1). Species varied in the total number of records(11-1475), total number of years with data (4-113) and range of years(30-135) with data.Overall statistical approachUnless otherwise specified, we used linear mixed-effects models (?nlme?package in R (Pinheiro et al., 2012) and lmer; ?lme4? package in R (Bates183et al. 2011). All models were fitted using restricted maximum likelihood(REML), which deals appropriately with correlated structures, unbalanceddesigns and gives less biased variance estimates than the method ofmaximum likelihood (Pinheiro & Bates, 2000). Likelihood ratio testswere used for model comparison and AIC was used for model selectionfollowing Burnham & Anderson (2002). If sample size was less than 40,AICc (AIC corrected for small sample size) was used. We judged significantimprovement in model fit when AIC or AICc improved by greater than two(Burnham & Anderson, 2002).Estimating temperatureBased on studies conducted in similar climates, the best predictors ofbutterfly phenology are temperatures from late winter to the end of sum-mer, depending on the region (Dell et al., 2005; Pollard, 1988; Roy & Sparks,2000; Sparks & Yates, 1997). Therefore, we first determined which threemonth block from January to July best predicted the timing of flight seasonon average for species by comparing the relationship between temperatureand timing of flight season across all 205 non-migratory species for eachseasonal block (Table F.1). Weather station and year were included as ran-dom effects. To estimate temperature for each block, we used mean dailymaximum temperature (for all 24 hours). Daily maximum temperatureswere used rather than mean temperatures because they provide a more rel-evant characterization of conditions when butterflies are developing duringthe day (Dell et al., 2005). If more than 10% of the days were missing datafor a given three-month block, the temperature value for that year was con-184sidered missing. Of the four seasonal blocks we tested, no block was theclear best model (i.e. all ?AIC<2). However, late spring (April 1 to June30) was ranked the best model most often across species (Table F.1), so weproceeded with late spring temperature.Cumulative degree-days have also been shown to be an important pre-dictor of insect phenology (Bryant et al., 2002; Hodgson et al., 2011; Romo& Eddleman, 1995). However, since calculating cumulative degree-days forall records across Canada is prohibitively computationally intensive, wefirst wanted to determine whether degree-day was a significantly betterpredictor of the timing of flight season than late spring temperature acrossa subset of sites. Therefore, we calculated cumulative degree-days above5?C from April 1 to June 30 for the province of British Columbia, whichis especially climatically heterogeneous. For each species we compareda model with day of year as a function of temperature calculated basedon mean daily maximum spring temperature to one based on degree-daytemperature (number of species = 130). Based on the average median AICdifference across all species, the model with late spring temperature was abetter fit than the model with degree-day temperature (mean ?AIC= 7.63(0.24SE)). Therefore, we proceeded with mean daily maximum temperaturefor April 1 to June 30.185Table F.1: Relationship between seasonal temperature and timing of flight season across species (n=205). Shownis the mean temperature sensitivity of flight season timing (days/?C), the number of species with that seasonalblock as the best fit, median ?AIC between models across species.Season Months Mean temperature sensitivity(SE) (days/?C)Number of species withmodel ranked as bestMedian ?AICWinter January-March -0.46 (0.08) 22 1.28Early spring March- May -1.49 (0.11) 43 0.54Late spring April-June -1.85 (0.13) 80 0.98Summer May- July -1.56 (0.13) 60 0186Table F.2: Predicted relationships between relevant ecological traits and the sensitivity of flight season timing totemperature.Trait Prediction RationaleOverwintering stage(egg, larvae, chrysalis,adult)More advanced stages will bemore sensitiveAdults are more mobile than other developmentalstages and can readily respond without further de-velopment (Dennis, 1993; Diamond et al., 2011)Flight season length Species with shorter flight sea-sons will be more sensitiveSpecies with shorter flight seasons should haveheavier fitness costs of mis-timing their seasonVoltinism Species with fewer generationswill be more sensitiveSpecies with more generations are more dependenton photoperiod as a cue for diapause (Tobin et al.,2008)Timing of flight season Early season species will bemore sensitiveEarly season species will have heavier fitness costsof mis-timing because of increased likelihood offrost and more variable weather and they are likelyto be more attuned to abiotic cues (Pau et al., 2011)Dispersal ability Species with lower dispersalability will be more sensitiveSpecies with greater dispersal ability are better ableto track suitable climatic conditions and are lessdependent on tracking the phenology of their hostplants (Diamond et al., 2011)Continued on next page187Trait Prediction RationaleLarval host breadth(monophagous,oligophagous orpolyphagous)Species with a broader dietwill be more sensitiveSpecies with a broader diet will be less dependenton their host plants and better able to track suit-able temperature (Altermatt, 2010; Diamond et al.,2011)Range size Species with a smaller rangewill be more sensitivePopulations of species that occur over a wide geo-graphic range are adapted to a wide range of cli-mates and are better able to track suitable climaticconditions (Altermatt, 2010; Diamond et al., 2011)188a)# b)#Figure F.1: Pattern of butterfly collecting efforts in the Canadian National Collection of Butterflies database(Layberry et al., 1998) and supplemental database used in the analysis, from each decade over the past 130 years(a) and across Canada (b).189Appendix GDescription of mainstatistical analysesOverall statistical approachSee Appendix F for details on overall statistical approach.Testing for interactions, non-linearity and temporal autocorrelationTo test for a potential interaction between spatial and temporal dimen-sions of temperature in predicting the timing of flight season, we evaluatedthe contribution of an interaction term to each species? model. However,on average, the interaction term did not improve model fit (?AIC= 0.265;n=204) so no interaction term was included.To test for potential nonlinear relationships in each species? model (i.e.day of year as a function of temporal and spatial dimensions of temperature,temperature as a function of year and day of year as a function of year), weconducted generalized additive mixed-effects models (gamm: ?mgcv? pack-age in R (Wood, 2011)) with a smooth function. We used the estimateddegrees of freedom of the smooth term (which allows for the best fit) todetermine the degree of nonlinearity. This dimension equaled 1, which indi-190cates linearity, for the majority of species and models (68-80%, number ofspecies=207, number of models=3), so we continued with lme models.To test for the presence of temporal autocorrelation in each species?model (i.e. day of year as a function of temporal and spatial dimensions oftemperature, temperature as a function of year, and day of year as a functionof year), we specified the correlated within-group error structure with thelme function in R using an autoregressive order of one, which allows thecorrelation function to decrease in absolute value exponentially with lag. Wealso experimented with different correlation structures, (a moving averageparameter, increasing the order of the autoregressive parameter) but theyresulted in problems with model convergence. For each species, a modelwith an autoregressive correlated error structure was compared to a modelwith no specified correlation term.The model fit for spring temperature as a function of year for themajority of species was significantly improved by accounting for temporalautocorrelation (median ?AIC= 2.88 (1.72SE), 111/207 of species had?AIC>2). Therefore, we specified a serial correlation structure for springtemperature as a function of year for all species. Since models witha temporal error structure did not improve model fit for the majorityof species for day of year as a function of spring temperature nor as afunction of year (median ?AIC=-0.40 (1.03SE), 78/207 of species had?AIC>2; median ?AIC=-0.35 (0.66SE), 76/207 of species had ?AIC>2,respectively), we concluded that temporal autocorrelation was not presentin the majority of cases and proceeded with lme models for these two models.191Evaluating trends in phenology over timeUsing the regression coefficients from each analysis, we evaluatedwhether a species? temporal phenological shift could be predicted by i) thedegree of temperature change for the set of locations and times at whichthere were collection records for that species, ii) its phenological sensitivityto spatial temperature and iii) its phenological sensitivity to temporaltemperature. We used generalized least squares with maximum likelihoodto test the significance of these relationships. For all models, an outlierwas removed to meet the assumption of a normal error distribution. Theremoval of the outlier did not affect the results so the results are presentedwithout the outlier.192Appendix HPhylogenetic data andanalysesPhylogenetic treeIn order to account for potential phylogenetic non-independence in ouranalyses, a phylogenetic tree of all species included in our study was re-constructed (Figure H.1). Sequences from the mitochondrial cytochromeoxydase I gene (COI; the ?barcoding? gene) were obtained from Genbank.Complete sequences (?1475 base pairs) were selected when available, butpartial sequences were also used. For each taxon, sequences from two differ-ent sources were obtained when possible. All genera were represented in thefinal dataset, which included 284 sequences. Sequences were retrieved fromGenbank and aligned in GENEIOUS v 5.6 (Drummond et al., 2012). Align-ment was performed with the MAFFT (Katoh & Toh, 2008) plugin usingdefault settings and was visually inspected. The final trimmed alignmentincluded 1473 sites.The program MODELGENERATOR (Keane et al., 2006) was employedto determine the best substitution model to use for tree reconstruction.The Akaike information criterion (AIC) and the hierarchical likelihood ratio193tests implemented in the software both selected the GTR + I + ? modelas the best one (general time reversible model with a proportion of invari-able sites and a gamma-distributed rate variation).This model was used forphylogenetic inference using the software MRBAYES v 3.2 (Ronquist et al.,2012). The dataset was partitioned in three subsets according to codonposition. All model parameters were unlinked for analysis, i.e. each parti-tion had distinct parameters (nucleotide frequencies, rates of transitions andtransversions, shape of gamma distribution, proportion of invariable sites).Preliminary runs without topological constraints were performed to evaluateif constraints were required, as homoplasy and saturation in COI may hin-der the reconstruction of very deep nodes in phylogenies (Hillis et al., 1996),which has been found to occur in Lepidoptera (Wilson, 2011). These prelim-inary results revealed that our COI dataset failed to recover the monophylyof Papilionoidea, placing Papilionoidae as a sister taxon to a clade compris-ing Hedylidae + Hesperiidae + other papilionoid families (Lycaenidae +Pieridae + Nymphalidae). However, this result was replicated recently in athorough investigation of butterfly phylogenetics that included 191 morpho-logical characters and 7 nuclear genes in addition to COI (Heikkila? et al.,2012). We also obtained paraphyly of the hesperiid subfamily Pyrginae,which is consistent with current knowledge of butterfly phylogenetics (War-ren et al., 2008, 2009). No additional discrepancy from recognized highertaxa in butterfly systematics was noted. Since our results were congruentwith other studies, we decided not to enforce any topological constraint. Fi-nal computations were performed with four simultaneous independent runs,each comprising four chains, for 2 million generations, discarding the first19425% as burn-in. The average standard deviation of split frequencies acrossruns was 0.017 and the effective sample size (ESS) of all model parameterswas above 100 (range: 132-2772). Two genera (Cupido and Erora), whichare included here, to the best of our knowledge, for the first time as in-group taxa in a phylogenetic analysis, did not cluster with their respectivesubfamily.A matrix of cophenetic distances (distances between pairs of tips ona tree) was calculated from the final tree (with branch lengths) with theR package APE v 3.0 (Paradis et al., 2004). Species included in the phy-logenetic analysis that were missing from our trait dataset were removed.For species in the trait dataset but missing from the final tree, speciesfrom the same genus were used instead to obtain inter-generic distances.Intra-generic distances for missing species were considered to be the meanpairwise distances among species from the same genus. If only a singlespecies from that genus was included (i.e. it was impossible to calculatea mean intra-generic distance for that genus), intra-generic distances wereinstead approximated as the mean intra-generic distances from all genera inthe tree (excluding the few genera that were revealed to be paraphyletic).These operations were performed with custom R scripts written by theauthors and applied on the cophenetic matrix. The final matrix of pairwisecophenetic distances was used to calculate a new tree with the UPGMAfunction implemented in the R package PHANGORN (Schliep, 2011). Thistree is shown in Figure H.1 and was used for subsequent statistical analyses.Phylogenetic analyses195Following Revell (2010), we first assessed whether a phylogenetic analysiswas necessary by comparing residuals from linear models to phylogenetically-adjusted linear models. To do so, we used a phylogenetic generalized least-squares model framework (pgls, (Freckleton et al., 2002); gls function inthe ?nlme? package (Pinheiro et al., 2012) in R) to account for potentialphylogenetic non-independence across species. This method uses the phylo-genetic variance/covariance matrix estimated from the phylogeny to adjustfor correlated error structure. For each trait, we conducted a standard glsmodel (not accounting for phylogenetic relationships) and then applied sev-eral common branch length transformations including a Brownian motionmodel (no transformation), Ornstein-Uhlenbeck, and Pagel?s lambda. Allgls and pgls models were fit with maximum likelihood. To confirm that apgls model was the most appropriate (rather than a standard gls model), wecompared standard gls and pgls models with AIC (Table H.2, H.3, H.4) andthen reported the results from the model that best fits the data (?AIC>2).If there was no clear best model (i.e. ?AIC<2), the results from the glsmodel were reported. Using pgls models instead of standard regressionswhen the residual error is non-phylogenetic can increase type I error (Rev-ell, 2010). In any case, results were similar between gls and pgls modelsfor these traits. If the pgls model fitted the data better than the standardgls model, then the residual variation in the dependent variable exhibitedphylogenetic signal (Blomberg et al., 2003; Freckleton et al., 2002).For each trait, the absolute value of species? phenological sensitivity totemperature was the response variable. To meet assumptions of normalityand homoscedasticity, we square-root transformed species? phenological sen-196sitivity to temperature and range size, as well as log-transformed averageflight season length and wingspan. To improve homoscedasticity betweenspecies? phenological sensitivity to temperature and timing of flight seasonand mobility, we used an exponential variance structure and for larval hostbreadth, we used a power variance structure. Larval host plant breadth andoverwintering strategy were defined as integers.Since the variance increased with the mean for the trait ?average numberof generations? and standard variance transformations (Pinheiro & Bates,2000) did not improve heteroscedasticity, generalized estimating equations(GEE; Paradis & Claude, 2002) with a Gamma error structure were used(compar.gee function in the ?ape? package (Paradis et al., 2004) in R). Thismethod incorporates species relatedness as a correlation matrix and uses ageneralised linear model approach, allowing data to be analysed using non-normal response variables. GEE models assume a Brownian motion modelof evolution so to incorporate other models of evolution, we transformedthe tree branch-length using the estimated from the pgls model (fitted bymaximum likelihood) and used the transformed tree in a GEE model. Wecompared the GEE models to one without phylogeny by setting ? to 0 (nophylogenetic signal) and re-running the GEE model. We compared modelsusing QIC, a quasi-information criterion. Finally, an F test was used todetermine the significance of overall model fit (between a full and reducedmodel).197Table H.1: Model comparison between phylogenetic and non-phylogenetic gls models for the relationship betweenspecies? traits and phenological sensitivity to overall, spatial and temporal temperature for only those species withnegative sensitivities (warmer temperatures lead to earlier flight season timing). Only traits found to be significantin the main analysis (Table 3.1) are shown here. Models in bold represent the best model. The response variableis the square root of the absolute value of phenological sensitivity to temperature as defined in the main text.Trait Temperaturesensitivity(residualdf)Model ?AIC TransformationparameterCoefficient(SE)LRT df pvalueTimingof flightseasonOverall(181)No phylogeneticcorrection21.18 NA -0.0029(0.0013)Brownian 179.35 NA -0.014 (0.0026)Ornstein-Uhlenbeck5.93 1562.03 -0.0030(0.0013)Pagel 0 0.16 -0.0038(0.0015)6.49 5,4 0.011Spatial(166)No phylogeneticcorrection16.18 NA -0.0051(0.0018)Brownian 217.36 NA -0.019 (0.0055)Ornstein-UhlenbeckNA NA NAPagel 0 0.15 -0.0054(0.0019)7.56 5,4 0.006Continued on next page198Trait Temperaturesensitivity(residualdf)Model ?AIC TransformationparameterCoefficient(SE)LRT df pvalueTemporal(180)No phylogeneticcorrection15.26 NA -0.00048(0.0013)Brownian 177.46 NA -0.0037(0.0012)Ornstein-Uhlenbeck0 277.46 -0.00082(0.0013)0.38 5,4 0.54Pagel 0.68 0.026 -0.00038(0.0014)Mobility Overall(178)No phy-logeneticcorrection0 NA -0.073(0.025)8.33 4,3 0.0039Brownian 247.83 NA -0.37 (0.062)Ornstein-Uhlenbeck1.38 1566.95 -0.073 (0.025)Pagel 1.3 0.043 -0.068 (0.027)Spatial(163)No phy-logeneticcorrection0.42 NA -0.068(0.035)3.75 4,3 0.053Brownian 204.87 NA -0.52 (0.12)Ornstein-UhlenbeckNA NA NAContinued on next page199Trait Temperaturesensitivity(residualdf)Model ?AIC TransformationparameterCoefficient(SE)LRT df pvaluePagel 0 0.14 -0.062 (0.041)Temporal(177)No phy-logeneticcorrection0 NA -0.033(0.026)1.64 4,3 0.2Brownian 198.78 NA -0.035 (0.048)Ornstein-Uhlenbeck1.73 703.05 -0.033 (0.026)Pagel 1.99 0.0072 -0.032 (0.026)200Table H.2: Model comparison between phylogenetic and non-phylogenetic gls models for phenological sensitivityto overall temperature with each species? traits. Models in bold represent the best model that is reported inthe main text. Range size was square-root transformed, and average flight season length and wingspan werelog-transformed. ?*? means QIC rather than AIC.Trait Model ?AIC TransformationparameterCoefficient (SE)Average number ofgenerations (residualdf=203)No phylogenetic correction 0 NA -0.028 (0.044)Brownian 1.18* NA 0.12 (0.020)Pagel 1.93* ?=0.063 -0.028 (0.044)Average length offlight season (residualdf=203)No phylogenetic correction 0.18 NA -0.054 (0.066)Brownian 358.07 NA -0.098 (0.066)Ornstein-Uhlenbeck 1.79 ?=2055.95 -0.054 (0.066)Pagel 0 ?=0.059 -0.036 (0.068)Timing of flight season(residual df=203)No phylogenetic correction 19.41 NA -0.0031 (0.0013)Brownian 221.82 NA -0.0019 (0.0030)Ornstein-Uhlenbeck 4.44 ?=1875.14 -0.0032 (0.0013)Pagel 0 ?=0.102 -0.0037(0.0014)Continued on next page201Trait Model ?AIC TransformationparameterCoefficient (SE)Mobility (residualdf=199)No phylogenetic correction 0 NA -0.064 (0.025)Brownian 348.27 NA 0.041 (0.093)Ornstein-Uhlenbeck 1.69 ?=2466.46 -0.063 (0.025)Pagel 0.97 ?=0.046 -0.057 (0.027)Larval host breadth(residual df=189)No phylogenetic correction 0.12 NA 0.040 (0.063)Brownian 342.45 NA 0.032 (0.13)Ornstein-Uhlenbeck 1.49 ?=1844.67 0.043 (0.063)Pagel 0 ?=0.061 0.042 (0.064)Range size (residualdf=167)No phylogenetic correction 0 NA -3.04e-5 (7.4e-5)Brownian 191.51 NA 3.72e-4 (5.81e-5)Ornstein-Uhlenbeck NA NA NAPagel 0.079 ?=0.074 -6.5e-6 (7.36e-5)Wingspan (residualdf=196)No phylogenetic correction 0 NA -0.16 (0.101)Brownian 319.95 NA 1.065 (0.31)Ornstein-Uhlenbeck 1.42 ?=1781.32 -0.16 (0.10)Pagel 1.13 ?=0.048 -0.096 (0.12)Overwintering (residualdf=188)No phylogenetic correction 0 NA -0.0019 (0.049)Brownian 342.89 NA 0.15 (0.14)Continued on next page202Trait Model ?AIC TransformationparameterCoefficient (SE)Ornstein-Uhlenbeck 1.64 2053.41 0.0020 (0.049)Pagel 0.83 ?=0.041 0.0043 (0.052)203Table H.3: Model comparison between phylogenetic and non-phylogenetic gls models for phenological sensitivityto temporal temperature with each species? traits. Models in bold represent the best model that is reportedin the main text. Range size was square-root transformed, and average flight season length and wingspan werelog-transformed. ?*? means QIC rather than AIC.Trait Model ?AIC TransformationparameterCoefficient (SE)Average number ofgenerations (residualdf=202)No phylogenetic correction 14.27 NA 0.047 (0.035)Brownian 0 NA 0.061 (0.0097)Pagel 14.05* 0.0063 0.047 (0.35)Average length offlight season (residualdf=202)No phylogenetic correction 0 NA -0.021 (0.070)Brownian 410.31 NA 0.31 (0.079)Ornstein-Uhlenbeck NA NA NAPagel 2 0.0022 -0.021 (0.070)Timing of flight season(residual df=202)No phylogenetic correction 14.18 NA -0.0030 (0.0014)Brownian 321.24 NA 0.0086 (0.0048)Ornstein-Uhlenbeck NA NA NAPagel 0 -0.00024 -0.0020(0.0014)Continued on next page204Trait Model ?AIC TransformationparameterCoefficient (SE)Mobility (residualdf=198)No phylogenetic correction 0 NA -0.039 (0.028)Brownian 396.8 NA -0.57 (0.13)Ornstein-Uhlenbeck 1.99 7229.49 -0.039 (0.028)Pagel 1.77 -0.024 -0.043 (0.026)Larval host breadth(residual df=188)No phylogenetic correction 0 NA 0.023 (0.062)Brownian 412.07 NA 0.17 (0.078)Ornstein-Uhlenbeck NA NA NAPagel 1.86 -0.016 0.018 (0.061)Range size (residualdf=166)No phylogenetic correction 0 NA -6.7e10-5(7.89e10-5)Brownian 312.59 NA 8.61e10-4(8.84e10-5)Ornstein-Uhlenbeck NA NA NAPagel 1.98 -0.0094 6.94e10-5(7.88e10-5)Wingspan (residualdf=195)No phylogenetic correction 0 NA -0.089 (0.11)Brownian 409.39 NA 1.603 (0.41)Ornstein-Uhlenbeck NA NA NAPagel 1.72 -0.023 -0.109 (0.094)Overwintering (residualdf=187)No phylogenetic correction 0 NA 0.039 (0.052)Continued on next page205Trait Model ?AIC TransformationparameterCoefficient (SE)Brownian 400.85 NA 0.032 (0.18)Ornstein-Uhlenbeck NA NA NAPagel 1.99 0.0012 0.039 (0.052)206Table H.4: Model comparison between phylogenetic and non-phylogenetic gls models for phenological sensitivityto spatial temperature with each species? traits. Models in bold represent the best model that is reported inthe main text. Range size was square-root transformed, and average flight season length and wingspan werelog-transformed. ?*? means QIC rather than AICTrait Model ?AIC TransformationparameterCoefficient (SE)Average number ofgenerations (residualdf=202)No phylogenetic correction 3.30* NA -0.095 (0.055)Brownian 11.74* NA -0.18 (0.040)Pagel 0 0.148 -0.087 (0.053)Average length offlight season (residualdf=202)No phylogenetic correction 2.45 NA -0.061 (0.079)Brownian 329.78 NA 0.33 (0.074)Ornstein-Uhlenbeck NA NA NAPagel 0 0.143 -0.034 (0.082)Timing of flight season(residual df=202)No phylogenetic correction 18.29 NA -0.0045 (0.0015)Brownian 306.24 NA -0.019 (0.0043)Ornstein-Uhlenbeck NA NA NAPagel 0 0.17 -0.0039(0.0018)Continued on next page207Trait Model ?AIC TransformationparameterCoefficient (SE)Mobility (residualdf=198)No phylogenetic correction 0.33 NA -0.073 (0.020)Brownian 340.19 NA 0.097 (0.073)Ornstein-Uhlenbeck NA NA NAPagel 0 0.11 -0.062 (0.035)Larval host breadth(residual df=188)No phylogenetic correction 0.36 NA -0.079 (0.076)Brownian 333.71 NA 0.15 (0.091)Ornstein-Uhlenbeck NA NA NAPagel 0 0.096 -0.066 (0.078)Range size (residualdf=166)No phylogenetic correction 1.64 NA -1.33e10-4(9.0e10-5)Brownian 156.55 NA -1.06e10-4(6.34e10-5)Ornstein-Uhlenbeck 3.64 4237.45 -1.33e10-4(9.1e10-5)Pagel 0 0.197 -9.16e10-5(8.8e10-5)Wingspan (residualdf=195)No phylogenetic correction 1.47 NA -0.14 (0.12)Brownian 270.58 NA 2.76 (0.33)Ornstein-Uhlenbeck NA NA NAPagel 0 0.147 -3.98e10-5 (0.17)Continued on next page208Trait Model ?AIC TransformationparameterCoefficient (SE)Overwintering (residualdf=187)No phylogenetic correction 1.53 NA -0.030 (0.060)Brownian 307.28 NA 0.71 (0.16)Ornstein-Uhlenbeck NA NA NAPagel 0 0.121 -0.033 (0.070)209Figure H.1: Phylogenetic tree of all Lepidoptera species included in thisstudy. This tree was obtained by pruning the tree obtained from theBayesian analysis of a more extensive mtDNA COI sequence dataset andadding species for which no sequence data were available (see Appendix Hfor details). Branch colors indicate different subfamilies. Outgroups wereremoved for visual representation, but the tree is shown as rooted, withoutgroups placed at the most basal node.210Appendix IAdditional analysesBest temperature predictor of phenology analysis and resultsTo determine the best temperature predictor of phenology, we used a3-month moving window approach (January-March, February- April etc.)from January to August covering the time period of 95% of the butterflycollection records and herbarium specimens. Following results from similarstudies, we focused on three aspects of temperature (Dell et al., 2005; Gordo& Sanz, 2010; Hodgson et al., 2011). For each of these time blocks, we tookthe mean temperature, mean daily maximum temperature and cumulativedegree-days above a base threshold above 0?C (for plants only) and 5?C. Ifmore than 10% of the days were missing information for that three-monthperiod, the temperature value for that year was considered missing.The metric of temperature that best predicted overall phenological sen-sitivity differed between groups. For plants, there was no model that was asignificant improvement over other models (Table I.1); therefore, we evalu-ated the sensitivity of our results using two different metrics of temperature(Table I.2). We took the seasonal block-temperature combination with thelowest median AIC (mean maximum for February to April) and the one withthe greatest mean phenological sensitivity (mean for March to May) (Table211I.1,I.2). Since both dimensions of temperature (spatial and temporal) pre-dicted the timing of plant flowering season using both temperature metrics(Table I.2), we only report a single metric of temperature for plants in themain analyses - mean temperature for March 1 to May 31 - the metric withgreater phenological sensitivity and more species with data.For butterflies, there was no single model that was a significantimprovement over other models; however, the best model (lowest median?AIC) also had the greatest mean phenological sensitivity, so we proceededwith that model - mean daily temperature for May to July (Table I.3).212Table I.1: Relationship between the timing of flowering season and temperature for 56 plant species. Mean dailyand mean maximum daily temperature, and cumulative degree-days above 0?C and 5?C for six 3-month blocks isshown. Shown is the mean sensitivity of flowering season timing to temperature (days/?C), the number of specieswith that seasonal block-temperature combination as the best fit and the median AIC between models. The bestmodels (lowest ?AIC) and greatest mean phenological sensitivity are in bold.TemperaturevariableSeason Months Mean sensitivity (SE)(days/?C)Number of specieswhere model isbest?AICMean Winter January-March -0.99 (0.55) 7 1.62Late winter February- April -1.90 (0.77) 7 0.18Early spring March- May -2.27 (0.77) 5 0.32Late spring April- June -2.10 (0.68) 4 0.87Early summer May- July -2.19 (0.55) 6 1.47Late summer June- August -2.01 (0.56) 7 1.66Mean dailymaxWinter January-March -1.26 (0.55) 3 0.22Late winter February- April -1.94 (0.76) 5 0Early spring March- May -1.91 (0.51) 2 0.16Late spring April- June -1.36 (0.44) 3 0.67Early summer May- July -1.11 (0.42) 5 1.36Late summer June- August -0.94 (0.46) 2 1.24Degreeday model(0?C)Winter January-March -0.033 (0.0096) 0 8.34Continued on next page213TemperaturevariableSeason Months Mean sensitivity (SE)(days/?C)Number of specieswhere model isbest?AICLate winter February- April -0.037 (0.010) 0 8.71Early spring March- May -0.033 (0.010) 0 9.4Late spring April- June -0.028 (0.0068) 0 9.73Early summer May- July -0.024 (0.0062) 0 10.37Late summer June- August -0.023 (0.0066) 0 10.54Degreeday model(5?C)Winter January-March -0.095 (0.025) 0 7.14Late winter February- April -0.065 (0.018) 0 8.45Early spring March- May -0.038 (0.011) 0 9.09Late spring April- June -0.030 (0.0072) 0 9.62Early summer May- July -0.025 (0.0063) 0 10.32Late summer June- August -0.023 (0.0065) 0 10.59214Table I.2: Differences in plant phenological sensitivity to temporal and spatial temperature based on two differentmetrics of temperature. Where a sign test instead of a one-sample t-test was used, the s test statistic and medianslope value is shown. Significant p values are in bold.Temperature variable orgroup of species (number ofspecies)Model (units of coefficient) Mean coefficient(SE)T teststatisticPvalueMean max for February- April(n=59)Sensitivity to temporal temperature(days/?C)-2.40 (0.82) -2.92 0.0025Sensitivity to spatial temperature(days/?C)-1.57 (1.83) S=24 0.096Mean for March- May (n=60) Sensitivity to temporal temperature(days/?C)-3.69 (1.57) -2.35 0.011Sensitivity to spatial temperature(days/?C)-1.82 (1.31) S=22 0.026215Table I.3: Relationship between the timing of flight season and temperature for 71 non-migratory butterfly species.Mean daily and mean maximum daily temperature, and cumulative degree-days above 5?C for six 3-month blocksis shown. Shown is the mean sensitivity of flowering season timing to temperature (days/?C), the number ofspecies with that seasonal block-temperature combination as the best fit and the median ?AIC between models.The best model (lowest ?AIC) is in bold.TemperaturevariableSeason Months Mean sensitivity (SE)(days/?C)Number of specieswhere model isbest?AICMean Winter January-March -1.02 (0.33) 1 2.84Late winter February- April -1.29 (0.50) 9 3.21Early spring March- May -2.14 (0.54) 9 2.38Late spring April- June -2.10 (0.51) 6 1.78Early summer May- July -2.33 (0.64) 9 0Late summer June- August -2.16 (0.63) 10 2.4Mean dailymaxWinter January-March -1.16 (0.36) 3 3.57Late winter February- April -1.50 (0.44) 9 3.82Early spring March- May -1.81 (0.39) 1 5.45Late spring April- June -1.49 (0.38) 4 1.58Early summer May- July -1.54 (0.48) 4 1.98Late summer June- August -1.54 (0.42) 2 6Degree daymodelWinter January-March -0.12 (0.0041) 0 6.78Late winter February- April -0.057 (0.017) 1 9.22Continued on next page216TemperaturevariableSeason Months Mean sensitivity (SE)(days/?C)Number of specieswhere model isbest?AICEarly spring March- May -0.015 (0.053) 0 16.05Late spring April- June -0.02 (0.0077) 0 8.14Early summer May- July -0.027 (0.0069) 0 11.57Late summer June- August -0.024 (0.0066) 0 12.6217Testing for non-linearity and temporal autocorrelationTo test for potential nonlinear relationships, we conducted generalizedadditive mixed-effects models (gamm function in the ?mgcv? package in R(Wood, 2011)) with a smooth function. We used the estimated degrees offreedom of the smooth term (which allows for the best fit) to determine thedegree of nonlinearity; a dimension of 1 indicates linearity. For each species,a model with a smoothing term was compared to a model with no smoothingterm. However, to the best of our knowledge, there is currently no mixed-effects model function in R that can deal with two non-nested random effectsand a smoothing term. Given this, we assessed the difference in model fitfor day of year as a function of temperature based on a model that includedone random effect at a time (site and year).For the model of day of year as a function of temperature, models werelinear for the majority of species and average model fit did not improve withthe smoothing term (Table I.4). This suggests that non-linearity is unlikelyto be an issue in the full model (with both site and year random effects). Forthe model of spring temperature as a function of year, the dimension equaled1 for the majority of species (plants: n=61, 75.41%; butterflies: n=120,50.83%) and on average, model fit did not improve with the smoothingterm (plants: median ?AIC= 2.0 (0.45SE); butterflies: median ?AIC= 2.0(1.04SE)). For the model of timing of flowering or flight season as a functionof year, models were linear for the majority of species (plants: n=66, 83.33%;butterflies: n=119, 81.51%) and average model fit did not improve with thesmoothing term (plants: median ?AIC= 2.0 (0.099SE); butterflies: median?AIC= 2.0 (0.26SE)). Therefore, we continued with lme models.218To test for the presence of temporal autocorrelation, we specified thecorrelated within-group error structure with the lme function in R usingan autoregressive order of one, which allows the correlation function to de-crease in absolute value exponentially with lag. We also experimented withdifferent correlation structures, (a moving average parameter, increasing theorder of the autoregressive parameter) but they resulted in problems withmodel convergence. For each species, a model with an autoregressive corre-lated error structure was compared to a model with no specified correlationterm. To the best of our knowledge, there is currently no mixed-effectsmodel function in R that can deal with two non-nested random effects anda correlation structure. Given this, we assessed difference in model fit forday of year as a function of temperature based on a model that includedone random effect at a time (site and year).For both butterflies and plants, model fit either did not significantly im-prove with a temporal error structure or was a worse fit (plants: year model:mean ?AIC=5.72 (0.39SE), n=70; station model: median ?AIC=5.49(0.54SE), n=71; butterflies: year model: median ?AIC=0.16 (1.06SE),n=117; station model: median ?AIC=1.14 (0.31SE), n=119), suggestingthat temporal autocorrelation would not be an issue in the full model (withboth site and year random effects). The model fit for spring temperature as afunction of year for the majority of species was a better fit without account-ing for temporal autocorrelation (plants: n=64, mean ?AIC= 4.09 (0.23SE);butterflies: n=121, median ?AIC=1.40 (0.52SE)). Similarly, models with atemporal error structure did not improve model fit for the majority of speciesfor timing of flowering or flight season as a function of year (plants: n=69,219median ?AIC=4.11 (0.32SE); butterflies: n=120, ?AIC=1.076 (0.29SE)).Therefore, we concluded that temporal autocorrelation was not present inthe majority of cases and proceeded with lme models for all models.220Table I.4: A comparison of mixed-effects models that evaluate the phenological sensitivity of butterflies and plantsto temperature based on the inclusion of a smoothing term to assess the potential for non-linearity. Shown is thecomparison when a single random effect is included in the model, the percentage of species with a linear smoothingterm for both temporal and spatial terms, and median difference in AIC between models. Mean temperature forMarch to May and mean temperature for May to July is used for plant and butterflies, respectively.Trophic group Number of species Random effect Linearity (%) (temporal, spatial) ?AIC (SE)Plant 34 Station 94.11, 97.059 -4 (0.050)37 Year 89.19, 100 -4 (0.056)Butterfly 115 Station 86.96, 79.13 -4 (0.26)115 Year 89.38, 59.29 -4 (1.24)221Appendix JAdditional methods,statistical analyses andresultsWe used linear mixed-effects models (?nlme? package in R (Pinheiro et al.,2012)) for analyses unless specified otherwise. Overall model residuals andwithin-group residuals were checked for normality and homoscedasticityfollowing Pinheiro & Bates (2000). To improve homoscedasticity, varianceswere allowed to vary by a fixed or random effect if the model fit wassignificantly improved (see below for details). To evaluate model fit, thefull model was compared to a reduced model using a likelihood ratio test.Restricted maximum likelihood was used to compare nested models inwhich only the random effects differed, but maximum likelihood was usedwhen only the fixed effects differed (Pinheiro & Bates, 2000). Restrictedmaximum likelihood was used for parameter estimation once the best modelwas selected. The model selection process followed a top-down strategyfollowing Zuur et al. (2009).222i) RearingTo rear adults from the pilot experiment in 2010, five females and fivemales were placed together as one group in paper containers with fresh twigsthat were changed daily. We randomized the family and assigned-treatmentof individuals within mating groups whenever possible but early and latein the season we were limited by the number of individuals that pupated atthe same time (males eclose first and generally only survive for 2-4 days)and by the sex ratio of each family increasing the potential for inbreeding.ii) Leaf growth rateLeaf length (from base of leaf to tip) on three leaves from every treatmentwas measured weekly from March 28 to June 19, 2011. To determine whethergrowth rate of leaves differed between treatments, we used a repeated mea-sures approach with length as our response variable and time (day of year)and treatment as fixed effects. Tree, branch and day (to account for mul-tiple leaves measured) were random effects. To improve homoscedasticity,variance was allowed to vary exponentially with time. We used a compoundsymmetry correlation structure to account for temporal non-independencefor both sites.Growth rate was faster in warming treatments at Totem Field (0.0228cm/day (0.0048SE), LRT10,9=21.83, p<0.0001) and UBC Farm (0.0121cm/day (0.0058SE), LRT10,9=4.34, p=0.0373). Leaves in warming treat-ments at Totem Field were smaller by 1.86 cm (0.52SE) across the experi-ment (t29=3.60, p=0.0012), whereas leaves did not differ between treatmentsat UBC Farm (t14=-1.28, p=0.2199).223Appendix KDescription of mainstatistical analysesWe used linear mixed-effects models (?nlme? package in R (Pinheiro et al.,2012)) for analyses unless specified otherwise (see Appendix J for generalstatistical approach).i) Testing the effect of branch-level warming on budburstTo test for differences in budburst between cut and attached treatments,we used a generalized mixed-effects model (lmer function, lme4 package in R(Bates et al., 2011)) with tree as a random effect and a Poisson distribution,and we combined warming and control treatments.To test for within-treatment relationships between temperature andbudburst, we used mixed-effects models for warming treatments at bothsites and generalized least squares for both control treatments. We usedmean daily maximum temperatures up to mean larval emergence across alltrees.ii) Treatment differences in phenology and fitness proxies224To determine whether warming treatments influenced phenology (forboth 2010 and 2011) and development time, we used tree as a random factor.For final larval weight, tree and family (to account for multiple individualsweighed per day) were random effects. Variances were allowed to vary bytreatment for development time at Totem Field.To determine whether growth rate differed between treatments, weused tree, family, instar, and day of measurement (to account for multipleindividuals weighed per day) as random effects. Time was a fixed factor andweight was the response variable. To improve homoscedasticity, variancewas allowed to vary for each instar-treatment combination at both sites.To allow for temporal non-independence, the best model also includedan autoregressive correlation structure of order 1 at Totem Field and anautoregressive-moving average correlation structure of order 1 at UBCFarm. Weight was log transformed for growth rate and final larval weightanalyses to improve normality.iii) Within-treatment relationshipsTo test for the direct and indirect effects of temperature on our proxiesof fitness, we explored within-treatment relationships between each proxyof fitness and mismatch and temperature. For development time and fam-ily survival, mixed-effects models were used for the warming treatment atTotem Field (tree as a random effect) and generalized least squares modelswith restricted maximum likelihood were used for all other treatments. Forlarval weight, mixed-effects models were used for all treatments. We usedtree (only for the warming treatment at Totem Field) and family (to account225for multiple individuals weighed per day) as random effects.To improve homoscedasticity between development time and mismatch,we used an exponential variance structure for all treatments except warmingtreatments at UBC Farm. To improve homoscedasticity between develop-ment time and temperature, a power variance structure was used for thewarming treatment at Totem Field. For warmed families at Totem Field,we improved homoscedasticity between number of larvae and mismatchby allowing variance to vary with degree of mismatch. We square roottransformed number of larvae to improve normality between number oflarvae and temperature. Larval weight was log transformed to improvenormality.iv) Leaf quality experimentWe used tree and family as random effects to determine whether femalepupal weight and time to pupation differed between treatments. Only treewas used as a random effect to determine whether there were significantdifferences in survival between treatments. To determine whether individualgrowth rate differed between treatments, we used tree, family, individuallarva and instar as random effects. To improve homoscedasticity, we logtransformed weight. We also included an autoregressive-moving averagecorrelation structure of order 1 to allow for temporal non-independence (forboth sexes).v) Temperature analysisDaytime temperatures (from 6:00 to 18:00) were calculated up to May 25226(the average date across all families to reach the fourth instar). Maximumdaily, mean daily maximum, mean daily temperature and cumulative degreedays above 10?C (above which larvae are active) were calculated.We then compared models to determine the best temperature predictorof development time and larval weight (see main text for description of thefitness proxies) at each site. Models were fit as described above. Each modelincluded treatment and was fit individually to avoid multicollinearity (i.e.there were four models for the four temperature variables). We could notcompare models for family survival (number of larvae) because we could notcompare treatments (see main text for details). At Totem Field, variancewas allowed to vary with treatment for development time to improve ho-moscedasticity. There was no significant interaction between treatment andtemperature in any model.For each site, the best model was based on the one with the smallest AIC(or AICc for sample sizes less than 40) and we judged significant improve-ment in model fit when AIC improved by greater than two (i.e. betweenthe best model and the model with the next smallest AIC; Burnham &Anderson, 2002). Since mean daily maximum temperature correlated withdevelopment time at Totem Field was the only case where there was a sig-nificantly improvement in model fit (Table K.1), we proceeded using meandaily maximum temperature in the within-treatment analysis.227Table K.1: Model comparison of different temperature variables and fitnessproxies (development time and larval weight) at each site. Shown is thedifference in AIC between models at each site.Fitness proxy Site Temperature variable ?AICDevelopment time Totem maximum daily 3.07mean daily maximum 0mean daily 12.39degree day 2.23Farm maximum daily 1.45mean daily maximum 0.75mean daily 8.54degree day 0Larval weight Totem maximum daily 0mean daily maximum 0.78mean daily 1.44degree day 1.3Farm maximum daily 0.21mean daily maximum 0mean daily 0.17degree day 0.19228

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

Usage Statistics

Share

Embed

Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                        
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.24.1-0058479/manifest

Comment

Related Items