Marine Phytoplankton Temperature versus Growth Responses from Polar to Tropical Waters – Outcome of a Scientific Community-Wide Study

“It takes a village to finish (marine) science these days” Paraphrased from Curtis Huttenhower (the Human Microbiome project) The rapidity and complexity of climate change and its potential effects on ocean biota are challenging how ocean scientists conduct research. One way in which we can begin to better tackle these challenges is to conduct community-wide scientific studies. This study provides physiological datasets fundamental to understanding functional responses of phytoplankton growth rates to temperature. While physiological experiments are not new, our experiments were conducted in many laboratories using agreed upon protocols and 25 strains of eukaryotic and prokaryotic phytoplankton isolated across a wide range of marine environments from polar to tropical, and from nearshore waters to the open ocean. This community-wide approach provides both comprehensive and internally consistent datasets produced over considerably shorter time scales than conventional individual and often uncoordinated lab efforts. Such datasets can be used to parameterise global ocean model projections of environmental change and to provide initial insights into the magnitude of regional biogeographic change in ocean biota in the coming decades. Here, we compare our datasets with a compilation of literature data on phytoplankton growth responses to temperature. A comparison with prior published data suggests that the optimal temperatures of individual species and, to a lesser degree, thermal niches were similar across studies. However, a comparison of the maximum growth rate across studies revealed significant departures between this and previously collected datasets, which may be due to differences in the cultured isolates, temporal changes in the clonal isolates in cultures, and/or differences in culture conditions. Such methodological differences mean that using particular trait measurements from the prior literature might introduce unknown errors and bias into modelling projections. Using our community-wide approach we can reduce such protocol-driven variability in culture studies, and can begin to address more complex issues such as the effect of multiple environmental drivers on ocean biota.


Introduction
To date, much of the progress in understanding how climate change will manifest itself in the ocean has come from projections obtained from global modelling experiments using general circulation or coupled ocean atmosphere models [1,2]. These types of models provide predictions of how bulk properties such as phytoplankton stocks (based on the proxy chlorophyll a) or ecosystem-level properties, such as downward export flux, will be altered by climate in the coming decades. However, at present many environmental projections from models -for example global maps of altered upper ocean temperature or nutrient concentrations -cannot be put into the urgently-needed wider context of the biological or ecological implications resulting from climate change. Such under-utilisation of model outputs is due to the current dearth of information on the physiological performance (often expressed as fitness versus environment, [3]) of many phytoplankton groups, species or ecotypes that are key players in the biogeochemical cycling of major (C, N, P) and minor (Fe, Zn, Co) elements in the ocean [4]. Given that phytoplankton photosynthesis and nitrogen fixation make major contributions to global C [5,6] and N [7,8] inventories, respectively, they can potentially drive significant feedbacks on climate change. However, the sign and magnitude of most biologically-mediated feedbacks are not yet known [9]. Both global modelling experiments [10] and timeseries data [11] reveal that the upper ocean is already warming in many regions (excluding upwelling regions, [12] and is highly likely to continue warming as a result of a changing climate. The temperature of the upper ocean is a fundamental control on phytoplankton metabolic processes [13,14] and sets the biogeographical boundaries or biomes of major phytoplankton groups [15,16]. As such, temperature response functions of phytoplankton are included in several widely used productivity models [17], leading to sensitivity in global primary production models to projected warmer surface ocean temperatures [6]. Furthermore, experiments in which multiple environmental properties are manipulated reveal that temperature has significant interactive (i.e. synergistic or antagonistic) effects with other properties such as carbon dioxide and/or iron concentrations, and on phytoplankton processes such as growth, photo-physiology, and calcification [18,19,20,21]. There have been several syntheses of temperature versus growth rate relationships for a range of laboratory-cultured phytoplankton [22] to look at generic relationships across a wide range of ocean temperatures. However, most comparative studies examining the influence of temperature on phytoplankton physiology were conducted several decades ago [23,24] and so such syntheses have had to rely upon available datasets that were collected using a range of protocols, both for culturing and estimation of growth rates. Moreover, many of the lab-cultured isolates used in these studies were isolated from the coastal ocean and/or were long-institutionalised (i.e. decades), and easily cultured ''weed'' species [22].
At present, there is a growing disjoint between the proliferation, improved accuracy and resolution of model projections of how oceanic conditions will change in the next 4-5 decades [25], and the availability of physiological datasets needed to contextualise such environmental projections. For example, the availability of datasets describing both the temperature optima for phytoplankton growth and the thermal limits (termed thermal niche width) at which pronounced decreases in physiological performance occurs, is limited, particularly for open ocean phytoplankton groups that drive major biogeochemical cycles [26]. Moreover, because different experimental approaches have been used to conduct laboratory culture experiments and calculate growth rates, the validity of comparing many of the datasets is questionable. Obtaining datasets that reveal the fundamental responses of organisms to altered ocean conditions, such as temperature, that are well-replicated across relevant physiological or environmental ranges for specific organisms is time-consuming, and funding may be difficult to obtain due to perceptions that such research has already been conducted in previous decades [22] or is unnecessarily simplistic.
However, these datasets are fundamental for parameterising models to make informed projections of how oceanic biota and ecosystems will respond to change in the many biomes that make up the global ocean [15]. Given that global change is occurring, and is highly likely to continue do so in the coming decades, we cannot afford to delay obtaining and employing such datasets to advise our models [27].
A recent development in other major disciplines, faced with similarly complex systems, has been the adoption of a communitybased approach to tackle the issues associated with a daunting number of permutations. In the disciplines of astronomy [28] and biochemistry/protein-folding [29] unprecedented rapid progress has been made through implementation and fostering of such community-wide initiatives. We maintain that in order to address issues of complex system science -e.g., multiple oceanic biomes with many different phytoplankton species, and a wide range of potential algal responses to environmental change -such community-wide efforts are necessary. However, because sufficient data relating to phytoplankton responses to climate change variables do not exist we carried out a pilot community-wide laboratory study. As part of this study we established a common laboratory approach for conducting experiments and then using this approach, we measured growth responses of cultured representatives from key phytoplankton groups to temperature.
Cultured isolates were examined originating from polar to tropical oceanic regions, and from coastal to remote offshore waters ( Table 1). The phytoplankton included eukaryotes such as a Southern Ocean diatom, isolated in waters of ,3uC and used for iron and photo-physiology studies [30,31] to prokaryotic nitrogen fixers (diazotrophs), isolated from tropical oligotrophic waters, being investigated for their response to ''greenhouse'' ocean conditions (i.e. a higher CO 2 warmer ocean, Fu et al., unpublished data, Hutchins et al. unpublished data) (Figure 1a). Other nearshore non-diazotrophic cyanobacterial species also being examined for their response to future ocean conditions (higher temperature and pCO 2 ; Ozmon et al., unpublished data), provided a contrast with offshore species. In addition, eukaryotic and mixotrophic dinoflagellate species involved in estuarine and coastal harmful algal blooms (HABs) and nearshore eutrophication were examined [32] as blooms of these groups of phytoplankton are thought to be favored under future climate scenarios [33]. A further contrast was provided between species in the diatom genus Thalassiosira that have been comprehensively studied physiologically and genetically, such as Thalassiosira pseudonana [34,35,36], Thalassiosira weissflogii [23,37,38] and Thalassiosira rotula [39,40] to those about which relatively little is known -such as the small diazotroph Crocosphaera watsonii [41] and the polar diatom Proboscia enermis [42] that have only recently been isolated.
In addition to testing the efficacy of this community-wide pilot study, this research provides insights as to how regional projections of warming might alter the physiological performance of phytoplankton groups/species that reside in distinctly different biomes across the world ocean ( Figure 1). Our study also enables a comparison with a recent collation [43] of published data on phytoplankton temperature versus growth rate relationships that were made with a wide range of experimental protocols. It is critical to establish the value of such syntheses of the earlier literature, whether we can use their parameterizations describing physiological growth response of different phytoplankton to temperature, and assess if we can relate these parameterizations to modelled projections of primary production and biogeographic distributions in response to regional warming. Our dataset provides an example of what is required to better parameterise models and predict how phytoplankton (as well as other microbial and planktonic groups) will be affected by changing temperatures [44,45,46] and hence the degree to which pelagic ecosystems will be potentially restructured in the future.

Temperature Growth Curves
The growth curves (hereafter termed reaction norms) for each species or strain are presented for polar and temperate species in Figures 2 and 3, and for tropical species in Figure 4 (see also Figures S1 and S2). The maximum growth rates ranged from 0.3 to 1.4 d 21 for the polar and temperate species, and were ,0.3 d 21 for all of the tropical nitrogen fixers (Figures 2, 3, and 4). The shapes of the reaction norms varied considerably from strongly asymmetric for the polar diatom to a more symmetric response for warmer water species and strains. Interspecific differences in the reaction norms were determined using species for which multiple strains were examined and a bootstrapping approach. The reaction norms of remaining species, where just a single strain was examined, were compared qualitatively. We first describe the main features of the reaction norms for all species, and then provide a statistical comparison of species with multiple strains.
The reaction norms of the diatoms species examined differed considerably. For the polar diatom ( Figure 3D), a temperature increase of 3uC (relative to the isolation temperature ,3uC, c.f. Figure 1B) resulted in a 25% increase in growth rate (m d 21 ), but then a further 1uC of warming caused a rapid decrease in m followed by mortality. All three temperate diatom species (Figures 2 and 3 C) had broadly similar shaped reaction norms, exhibiting up to a four-fold increase in m as temperatures increased but no further increases in growth rates at temperatures .20uC, or 25uC in the case of Thalassiosira rotula (Figure 2).
For the tropical species, the unicellular nitrogen fixing isolates (Crocosphaera watsonii) had a similar reaction norm curve to the ubiquitous cyanobacterium Synechococcus ( Figure 3A), a nondiazotrophic unicellular cyanobacteria, as did the three diazotrophic Trichodesmium isolates (Figure 4). Each diazotroph had a plateau in growth rate at 24 to 28uC. Growth rates dropped off at the upper end of this temperature range in a similar manner for both the unicellular and colonial nitrogen fixers (i.e. to zero by 35uC). However, compared to Crocosphaera the biggest difference between these reaction norms is that Trichodesmium has lower minimum temperature limits and lower maximum growth rates ( Figure 4). The dinoflagellate, A. sanguinea, displayed maximum growth rates at 25uC and an upper temperature limit of 30 to 33uC ( Figure 2). Attempts to grow the strains at 35uC failed repeatedly with some strains remaining vegetative (but with no discernible growth) followed by mortality. All strains had robust growth rates at 15uC. Attempts to grow the strains at 10uC resulted in little/no growth but with all strains maintaining vegetative populations, indicating a lower temperature tolerance (for growth) of between 10-15uC and an ability to survive at relatively low temperatures compared to the thermal optimum. The dinoflagellate P. donghaiense exhibited a ,3-fold increase in growth rate between 15 and 20uC, then exhibited ,15% decreases in growth rate at 30uC ( Figure 3B).

Intraspecific Variations in Growth Rate
Multiple strains were examined for the dinoflagellate A. sanguinea and the diatoms Thalassiosira rotula and Thalassiosira pseudonana at temperatures ranging from 4uC to 35uC. The coefficient of variation (CV) was examined to assess variability within species at each temperature. For A. sanguinea, the CV ranged from 8-82% over a temperature range of 15-33uC ( Table 2). The CV for Thalassiosira pseudonana ranged from 8% to 133% across a range of 10 to 32.5uC. For Thalassiosira rotula, the CV had a range of 9-53% across a temperature range of 4-25uC (no growth at 30uC). In general, these species exhibited more variation in growth rates among strains at the low and high temperature extremes ( Table 2). For example, in A. sanguinea the upper temperature boundary (33uC) had a CV that was an order of magnitude higher than the  Ocean Survey of Microalgal Thermal Reaction Norms PLOS ONE | www.plosone.org mid-range of the reaction norm (20uC). As noted above, strains also exhibited varying survival at the lowest and highest temperature ranges. Contrary to what one might expect, the least and greatest reduction in growth at elevated temperatures (33uC) occurred in the most northerly strains (BWA and GBB) ( Figure 2 middle panel), both isolated from the state of Washington, while the most southerly strain (YRB) showed the least sensitivity to low temperature (15uC), suggesting that geographical location does not necessarily correlate with optimal conditions for these coastal strains. This highlights the importance of evaluating multiple strains for a given region, since substantially different results would be obtained for the two isolates from Washington if one of those strains were used as representative of the species.
The N 2 fixers, C. watsonii and T. erythraeum generally had lower ranges of CV's across their temperature range (17-24% and 8-33%, respectively) than for the temperate diatoms. With few exceptions, there were significant differences in growth rates among strains within a species (p,0.05, Table 2). At every temperature tested, there were significant intraspecific differences in m amongst Thalassiosira rotula strains and C. watsonii strains. For Thalassiosira pseudonana, A. sanguinea and T. erythraeum, significant differences in growth rate occurred at every temperature except 15, 20 and 24uC, respectively ( Table 2).
The reaction norm curves for each species provided the opportunity to examine growth response over a range of temperatures. For example, each of the strains of the cosmopolitan diatom species Thalassiosira rotula (Table 1) was characterised by  Table 1) exhibited strikingly different reaction norms with respect to temperature. For example strain YRB continued to grow at rates close to m max at temperatures .25uC, whereas strains BWA and GBB had a pronounced decrease in m above 26uC. In the case of Thalassiosira pseudonana, the six strains appear to split into two groups in their temperature responses. Strains CCMP 1011 and 1012 have lower optima and narrower niches than CCMP 1013, 1014, 1015 and 1335 ( Figure 2, right panel). However, these groupings appear to be unrelated to the different locales that they were isolated from.

Interspecific Variation in Thermal Traits
The five species for which multiple strains were tested (A. sanguinea, C. watsonii, Thalassiosira pseudonana, Thalassiosira rotula, and Trichodesmium erythraeum) differed significantly in their trait distributions ( Table 3). The distribution of delta AICc values enabled us to evaluate whether interspecific differences in thermal traits existed, and permitted us to incorporate our uncertainty in the model fits. For both optimal temperature and maximum growth rate, the 95% confidence intervals on the delta AICc were .0 indicating significant differences among species for these traits ( Table 4). The confidence intervals for temperature niche width overlapped with zero, indicating that distributions of niche width cannot be easily distinguished based on species identity. However, this may be because our uncertainty in this trait is large (as estimated from the distribution of this trait in the bootstraps, Table 3), particularly in the case of T. pseudonana, T. rotula and A. sanguinea.

Intraspecific Variations in Thermal Traits -Comparison with the Literature
Of the eight species investigated in this study, we found previously published growth-temperature data for five common species (Table 5). For thermal niche width, there was less agreement between our study and prior published data ( Figure 5). For example see C. watsonii, where the niche width is twofold larger from the prior literature relative to our study. In contrast, the thermal trait T opt ( Figure 6) was similar to previously published values for T. erythraeum and Crochosphaera, lower than previously reported for Thalassiosira rotula, and higher for Thalassiosira pseudonana and A. sanguinea. A comparison of maximum growth rates between the present and prior studies revealed significant differences for four of the five common species (Figure 7). In most cases, we investigated more strains than investigators did in the prior literature and hence we see greater variability in each of these thermal traits in the present study (Figures 5, 6, and 7; Figures S3, S4, S5, and S6).

Discussion
Although we were able to characterise the thermal reaction norms of only a small subset of the resident phytoplankton in the global ocean in this illustrative community-wide experiment, they provide valuable insights into how a warming ocean could influence marine floristics. They also enabled us to conduct the first (as far as we are aware) evaluation of how much confidence we should have in exploiting the rich datasets of the prior physiological literature, and raise some issues about how best to relate such physiological data to the future temperature projections from climate change models. Hence they offer some lessons on where we should focus our efforts in future community-wide experiments.

Reaction Norm Shape -Implications for Biogeographical Change
The two open ocean end-members are the polar diatoms and tropical diazotrophs, which span a temperature range of 2.4 to 27.5uC that is comparable to that presented in well-cited collations [22,23]. They represent two groups that play a key role in the biogeochemical cycles of N, C, Fe and Si [47,48,49]. The first comparison of reaction norms of unicellular and colonial diazotrophs reveals marked differences between these groups in maximum growth rates and thermal optima that will have implications for future floristic shifts within the diazotroph assemblage (Fu et al., unpublished data). Such a range of temperature reaction norms (Figure 4) supports indirect evidence from oceanographic surveys of the role of environmental conditions in setting the relative distributions of different diazotroph groups [50]. Such floristic shifts (Fig. 3) will certainly have ecological consequences (different pathways for N fixed by uni-cells versus colonial diazotrophs, [8]. Moreover, as thermallytolerant non-diazotrophic groups may replace N 2 fixers, as annual temperature ranges begin to exceed the maximum limits of the diazotrophs, there may also be biogeochemical ramifications of such floristic shifts 51,52] (Figure 4). These datasets provide some of the first direct physiological evidence of the potential to alter biome boundaries, as predicted by models based on previously available relatively poor physiological details for different phytoplankton functional groups [53].
Although we have data for only one polar diatom (P. inermis, Table 1), the reaction norm reveals that warming could have a detrimental effect on Southern Ocean diatoms by the end of the century if this species is representative (Table 6) as this species would be at or beyond its thermal limit with a 3.5uC rise in temperature ( Table 6, Figure 1B). Although there are few prior studies of the temperature ranges of Southern Ocean diatoms [53], at least one study [54] reports a similar growth temperature range of 2 to 9.2uC for a Chaetoceros sp. isolated from the Southern Ocean. The reaction norm for the Southern Ocean diatom highlights the lack of oceanic refugia (i.e. colder waters) for these polar species, and raises uncertainties as to whether such a geographically-isolated phytoplankton group (i.e. South of a major physical barrier, the Polar Front) which has a very different photo-physiology [31] can acclimate or adapt on a timescale of decades to warming temperatures. By extension, it is unclear which phytoplankton group(s) might emerge or replace diatoms if they become extinct and how such a floristic shift will affect trophodynamics. Potentially coccolithophores could play expand their range and there is some recent evidence that in the Subantarctic this group is extending their southerly extent [55]. Shifts between diatoms and calcifying phytoplankton could also have major biogeochemical implications for the Southern Ocean Si and C cycles. There are major uncertainties in global environmental change research that are difficult to tackle in the confines of laboratory research. For example, far less than 1% of phytoplankton taxa are in culture or culturable and many of those in culture have been cultivated for thousands of generations (decades). Choosing representative species appears daunting as over 40,000 phytoplankton species are described, with thousands more to be discovered [56]. However, by focusing on the major contributors to primary production i.e. phytoplanktonic functional groups in the ocean (cyanobacteria, diatoms, coccolithophorids and dinoflagellates), the list can be narrowed [57]. Within each of those groups, representative species that are important contributors to bloom formation, downward carbon flux, and enhancement of upper ocean N inventories are logical targets for laboratory research. For example, we focused on diatoms that comprise the Thalassiosira spp. such as Thalassiosira rotula, because they are an important contributor to bloom formation [39,58] and Thalassiosira pseudonana, because it is one of the beststudied physiological model species [37,59] with a fully sequenced genome [36]. Understanding the growth responses of additional genera that are central to biogeochemical cycling such as the coccolithophores in the genus Emiliania and diatoms in the genus Chaetoceros will provide important data for understanding phytoplankton response to climate change.
It is also likely that species dominance will alter with changing environmental conditions, leading to the appearance of new and likely unstudied organisms that will need careful physiological characterization. For example, the diatom species Neodenticula seminae was sampled in large numbers from North Atlantic waters in 1999 and has become established there over the last decade [60]. This is the first recorded presence of N. seminae in the Atlantic in over 800,000 yrs. The ''invasion'' of this diatom is thought to reflect increasing transport of Pacific waters into North Atlantic waters via the Arctic [60].
One confounding issue of using representative species for phytoplankton functional types is the major differences in reaction norms (Figure 2), between strains for each species we considered. These trends of different thermal reaction norms potentially have parallels in recent studies that report differential responses to CO 2 enrichment among strains of coccolithophores [61,62], and diazotrophic cyanobacteria (Hutchins et al. unpublished data). In diatoms, significant inter-and intra-specific variation exists in response to many factors, including light intensity [63,38,64]. Our datasets on temperature and diazotrophs suggest that within each genus (i.e. unicellular or colonial), reactions norms differed among strains, albeit to a small degree (Table 3). There were also marked differences in thermal optima between genera, illustrating the dangers of using a sole strain of a cultured phytoplankton species to represent an entire phytoplankton functional group in ecosystem or biogeochemical models. However, we may be more optimistic about parameterizing models at the species level: there is broad similarity in the shape of strains' reaction norms within a species (Figures 2 and 4, Table 3). This may reflect selection on different taxa to maintain particular 'shapes' or underlying constraints in their ability to adapt to different temperature environments. More importantly, it suggests that we may justifiably parameterise biogeochemical models with traits of A. Growth rates were determined at a minimum of six temperature conditions. B. A minimum of three replicate growth rates were determined. C. All other environmental variables were held constant within each individual experiment, other than temperature. These include day length, culture medium, and culture protocols. Saturating nutrients were used to avoid nutrient-induced growth limitation. D. Isolates grown at saturating light intensity (mmol quanta m 22 s 21 ). E. Semicontinuous cultures were diluted using media that was previously adjusted to the appropriate temperature. Dilution frequencies were set so that cells were maintained in constant exponential growth phase and growth rates were reported when cultures were fully acclimated to the experimental conditions, after statistically invariant growth rates were recorded for at least 3-5 generations [98] F. Upper and lower thermal limits were tested repeatedly (at least 3 times) G. Multiple biomass parameter/ proxies were used to determine daily abundance and included cell counts, extracted chlorophyll a, and in vivo chlorophyll a fluorescence. Each method could be used reliably to determine steady-state acclimation. H. At each temperature, the maximum acclimated specific growth rate (d 21 ) for each isolate was determined by regressing the change in the log of fluorescence, cell count or chlorophyll a over time and testing the equality of slopes from at least three serial cultures (a = 0.05) [99]. If slopes of serial growth curves were not significantly different, the average regression coefficient was used to estimate the common slope, which represented the average acclimated growth rate and the standard error. Footnotes: 1 Growth rates were determined at five temperatures. 2 Growth rates were determined at 4 temperatures. Cells failed to grow at 35uC and reliable growth estimates could not be obtained at 10uC. 3 Growth rates were determined for 4 temperatures, cells exhibited no growth at 35uC. 4 Isolates also grown in 4 different nitrogen species (nitrate, ammonium, urea, and glutamate). 5 Carbonate system also held within a specific range at ambient conditions. 6 Recorded growth for 5 days after acclimation (not necessarily 3-5 generations). 7 After .8 generations. 8 This was performed for upper limit only. 9 Upper limited tested repeatedly. Lower limit was below 0uC, the lowest temperature tested. 10 Upper and lower limits were tested for P. donghaiense, however only lower limit was tested for Synechococcus. Upper and lower limit tests were performed twice, not three times. 11 This was performed for upper limit in all isolates, and in one isolate for the lower limit. 12 Fluorescence alone was used. 13 In vivo fluorescence reported, but Chl measured at the first and last culture days, as well as PN and PC. 14 Growth rates were determined using in vivo fluorescence, but were not significantly different from growth rates determined using cell counts. 15 Significant differences between slopes of replicate cultures were not tested. Instead, the mean slope was used. Variation within a temperature treatment was much smaller than variation between temperatures. doi:10.1371/journal.pone.0063091.t005 Figure 5. A comparison of the thermal trait, niche width (6C) using box and whisker plots, between previously published studies (using a wide range of experimental protocols, see [43]) and the species/strains used in the present study. The black bands denote the median value, the bottom and top of the red/blue boxes represent the 1st and 3rd quartile of the data respectively. The 'whiskers' extending from the boxes indicate the positions of the lowest & highest values in the data. If the sample size is small enough, the whiskers may not appear (e.g. if there are only 3 equally spaced points, the value represented as the 1st quartile is the lowest value). doi:10.1371/journal.pone.0063091.g005 Figure 6. A comparison of the thermal trait, T opt (6C) (box and whisker plots), between previously published studies (using a wide range of experimental protocols, see [43]) and the species/strains used in the present study. For details see Figure 5

Differential Physiological Responses to Warming between Strains
A further major challenge in global environmental change research is to understand why the reaction norms of multiple strains of the diatoms Thalassiosira rotula and Thalassiosira pseudonana and the dinoflagellate A. sanguinea differed. For example, one of six Thalassiosira rotula strains was unable to grow at 4uC and had a significantly higher optimal growth temperature than other strains of Thalassiosira rotula. This strain was isolated (,20 years ago) from the Mediterranean and it is tempting to link its growth response to its provenance. However, the reaction norms of the other two Mediterranean strains tested were more similar to strains sampled from northern temperate regions. This highlights the importance of examining multiple strains from a single region, particularly when testing for region-specific responses to environmental change.
Within species, the variation of growth rates amongst strains changed considerably with temperature. For example, the coefficient of variation (CV) amongst Thalassiosira rotula strains was fivefold greater at the low end of its temperature tolerance than at the upper end. In A. sanguinea, higher CV's were observed at the high end of its temperature tolerance than the middle. In general, variation among strains in both diatoms and dinoflagellates was lowest (7-17%) near the optimal growth temperature. This is comparable to previously observed intra-specific CVs of 5-15% in diatoms [63] and 13-39% in dinoflagellates [65] at nearoptimal growth temperatures. In contrast, variation amongst strains was often highest at the extremes of temperature tolerance suggesting that genotypic selection pressures would have the largest influence at these points.
In the case of A. sanguinea, the temperature-growth response is similar to the trends presented in [66] who reported optimal growth rates at 25uC and a salinity of 20, but with positive growth from 10-30uC and salinities from 10-40. Dinoflagellates are classified as ''modified latitudinal cosmopolitan'' [67] and true endemism is rare [68]. A. sanguinea follows this pattern and is widespread, observed along the west coast of the United States [69,70,71,72,73], the Gulf of Mexico [74], Brazil [75], Peru [76], Hong Kong [77], Japan [66], Korea [78], and the Black Sea [79]. While several hypotheses have been put forward regarding their range expansion including dispersion of vegetative cells either naturally or by ballast-water transport, there have been no genetic studies to date that provide population structure or gene flow, although microsatellite markers have been developed [80]. While the relatively uniform temperature-growth responses [66,69] (this study) are consistent with the classification of this species as eurythermal, we note again the strain variability at the low and high end of the temperature range. This strain variability suggests Figure 7. A comparison of the maximum specific growth rate (day 21 ), using box and whisker plots, between previously published studies (using a wide range of experimental protocols, see [43]) and the species/strains used in the present study. doi:10.1371/journal.pone.0063091.g007 Table 6. A summary of projected increases in global sea surface temperature for 2020-2029 (relative to 1980 to 1999) and for 2090-2099 (relative to 1980 to 1999) from three IPCC scenarios [104]). The A2 and A1B scenarios for CO 2 emissions are very similar to that for observed global emissions [106] and hence were used here for the comparison in Figure 7. doi:10.1371/journal.pone.0063091.t006 that the widespread success of this organism in neritic waters is due at least in part to substantial strain (genetic) variability within the species. Successful expansion of this genera due to either transport or changes in ocean climate may be largely dependent on the strain(s) present within a given region.

Comparison of Our Findings to the Temperature Versus Growth Literature
There is a large body of literature collected between 1962-2010 [43], but its value may be limited if there is too much confounding interference due to the differing protocols employed across individual studies. Our datasets provide some preliminary checks and balances to appraise the worth of such prior datasets. The thermal trait analysis reported in [43] differed from that carried out on the species in our study in several ways.
Departures between the literature collation relative to the present study may be related to the quality of the data used in this thermal trait analysis (the literature provides a large repository of data but of more uneven quality relative to the present study) and allows us to characterize areas where efforts such as this study would be most fruitful. For example, we found that differences in optimum temperature for growth between our study and literature estimates were small (Figure 6), suggesting that this trait is robust across a range of experimental methodologies. In contrast, a comparison of the maximum growth rates showed that this trait varied more than for the others examined (Figure 7), likely due to its higher sensitivity to experimental conditions, such as light intensity [69] or composition of lab culture media [81]. Temperature niche width is a more difficult trait to characterize, requiring measurements across a larger range of temperatures than is typically feasible. Therefore, although we see large differences in the distributions of this trait ( Figure 5), this is likely due to the high uncertainty in our estimates, rather than absolute differences.
Despite some of the above methodological limitations, these prior published estimates provide a 'parameter envelope' reflecting uncertainty introduced through experimental and statistical methods. This uncertainty provides valuable information that can be incorporated when modelling phytoplankton populations [82] and communities [83]. Models focusing on the community level pose the greater challenge, as the uncertainty introduced by measurements made under different environmental conditions must be carefully accounted for (though this may require the collection of additional data). We hope that in future, interactions between temperature and other important parameters (e.g. light, nutrients) will be better characterized, thereby constraining this uncertainty further and improving our ability to model community interactions.

Trait-environment Relationships
A recent synthesis study [42] reported found that optimum temperature was strongly related to mean environmental temperature, indicative of past adaptation. This pattern was not evident in the present study, probably because of insufficient data from high latitude isolates. In contrast, Thomas et al. did not observe any clear relationship between temperature traits and environmental variability. In our study we do see a positive saturating relationship between temperature niche width and annual temperature range, as predicted by the climate variability hypothesis [84] (see Figure S7). However, intraspecific variation in temperature traits is unrelated to environmental temperature variation and we do not have evidence for local adaptation within a species (see Table 3). This trait-environment relationship is consistent with a simple ecological interpretation: that changes in the temperature variability gradient, rather than local adaptation, drives species turnover. Nevertheless, our study may also lack the power to detect any effects of local adaptation.
Another obvious uncertainty in global environmental change research is the capacity of phytoplankton to evolve new thermal windows, or the lack of such evolutionary resilience [85]. In a thermal experimental evolution study [26], used ''ratchet'' experiments (incremental increases in growth temperature) to demonstrate differential abilities of phytoplankton groups to adapt to increasing temperature. In general, they found that species from continental lakes were able to adapt to large increases in temperature compared to open ocean phytoplankton groups; for instance, the globally distributed coccolithophore species Emiliania huxleyi was one of the species that displayed little or no ability to adapt to warming. Further evolutionary studies in regard to thermal tolerance are called for, including studies focused on the diversity of intraspecific responses. Long term evolutionary inferences can also be made from biogeographical observations, such as the lack of cyanobacteria or coccolithophores in polar marine waters [86]. Evolution of higher temperature tolerance may also be constrained by the existence of correlations and tradeoffs among traits and possible opposing selective pressures on correlated traits [87].
Inferences about species range limits can be drawn using their thermal reaction norms. These can be estimated under current and future climate predictions, and aggregated to make predictions about changes to diversity patterns and biogeochemical processes. A recent analysis [43] found that ocean warming this century is likely to lead to a decline in tropical phytoplankton diversity, as many tropical strains, in the absence of evolution, will be unable to survive even small increases in temperature. Merging this approach with ocean biogeochemical models and data on species' growth vs. nutrient curves will allow us to make even more fine-grained predictions of growth and possibly productivity of different groups in the future. However, evolutionary responses to temperature change need to be accounted for; experiments subjecting different phytoplankton species to different thermal regimes will prove valuable in estimating the rapidity of this process in single species. Future experiments that subject mixed phytoplankton communities to elevated temperature conditions will also be needed, as species interactions can negatively influence both rates of adaptation and productivity [88,89].

Merging Physiological Data with Climate Change Model Projections
One of our aims was to exploit the growing datasets of climate change projections from models which provide regional or global maps of how ocean conditions, such as temperature, irradiance, nutrient and trace metal supply, pH and carbon dioxide will be altered in the coming decades [2,90]). In Figure 8 we attempt to merge some of our physiological data for the two thermal endmembers -polar diatoms and tropical diazotrophs -with those of model projections for warming from the most realistic IPCC scenario -3.5uC warming by 2100 (Table 6; Figure 1B). This comparison is instructive as it reveals that additional information is required to merge the physiological and the model projections with confidence. The red bars in Figure 8 denote the projected increase in temperature based on the annual minimum and maximum temperature for the regions in which the phytoplankton were isolated [91]. Although information is available regarding the in situ temperature at the locales where cells were isolated -for example 3uC for the polar diatom, we do not know whether this species successfully subsists from austral winter to austral summer, or whether it is succeeded by other species over the course of the annual growth season.
Of course, beyond warming effects there are many other uncertainties about how the coastal ocean will respond to climate change. Additional anthropogenic pressures such as eutrophication, acidification, and alterations in salinity and irradiance regimes will all impact HABs and other coastal phytoplankton groups [92]. Consideration of these multivariate influences greatly complicates predictions of the responses of complex natural phytoplankton assemblages to environmental change [86]. In addition, in the laboratory we are limited to examining individual cultured isolates so we have little understanding of how the diverse range of strains/ecotypes present in natural systems might confer some stability or buffer ecosystems against such environmental change [85].

Next Steps and Lessons Learnt
Although we devoted considerable effort to developing a consensus-based experimental protocol, there were minor departures from it by most participants, that this is one area that can be improved upon. Such improvements, and the development of a more educated community with respect to protocol developments, are essential before we can tackle methodologically more challenging manipulations -trace metals, pH, or manipulation of multiple drivers concurrently, or to increase the number of participating laboratories to .10. This approach builds upon major efforts by the Ocean Acidification community to advocate best methodological practices [93] by being more proactive in developing community-wide protocols. Such community-wide initiatives would benefit from 'hands-on'' workshops similar to that conducted regularly by the Group for Aquatic Primary Productivity (GAP) workshops [94].
In climate change research there is an urgent need to alter what sort of experiments we do -more emphasis should be focused on community-wide experiments (in this short study we conducted measurements of 675 growth rates from 9 species). If we are to develop such community-wide initiatives we must also deal with some of the intrinsic limitations of our existing science culture, in that while collaborative efforts are often warranted and even encouraged, all investigators also have a need to demonstrate their individual productivity and creativity to their peers, their employers, and their funding agencies. Thus, it may be difficult to motivate researchers sufficiently in the long run (i.e. years) to devote considerable effort and scarce funding to obtaining collective datasets. The problem of inadequate motivation to participate in community-based efforts could be partially alleviated if there was a demonstrated interest by funding agencies in funding for and recognition of the needs and benefits of such community-wide research. These types of funding issues are especially problematic for international collaborative groups, due to the lack of coordination and communication between national funding agencies. These pressing issues and impasses have also been recognised recently for the field of macro-ecology [95].

Rationale for Selection of Temperature Versus Growth Rate
We designed a community-wide experiment examining the effect of temperature on phytoplankton growth for two general reasons. First, there is a high degree of consensus that there will be significant changes in temperature based on projections from climate change models (relative to other environmental properties such as incident irradiance or trace metal supply where changes are more uncertain, [86]. Second, it is relatively straightforward to manipulate temperature in laboratory experiments using cultured phytoplankton isolates (versus pH or CO 2 , [96,93] and hence this variable is more amenable to developing and implementing a common protocol for community-wide experiments. In many Figure 8. Thermal reaction norms for the two end-members species from our study compared with predicted ocean warming trends. Projected warming by 2020-2029 and 2090-2099 red bars (see Figure 1B) and temperature range (arrow) (from [91]), over the annual cycle is overlaid on the two reaction-norms (0.9 to 4.3uC for the polar diatom and 23.7 to 28.3uC for the Crocosphaera strains). doi:10.1371/journal.pone.0063091.g008 cases, temperature is considered to be a master variable or ultimately limiting property for phytoplankton growth [23,26].
The choice of growth rate as a physiological metric to measure the effect of temperature was selected because it is fundamental to the success of phytoplankton species and it ultimately integrates phytoplankton physiology [37], since for an asexually reproducing unicellular organism it represents the most direct measure of reproductive success and therefore of fitness in the environment. By directly measuring fitness (e.g. growth rate) over a range of temperatures, a reaction norm can be described for each genotype (culture isolate) with respect to temperature [97]. Ultimately, one would describe reaction norms for each physiological variable under a range of environmental conditions to understand both the direct and interactive effects of factors that affect phytoplankton growth.

Phytoplankton Species
The species used for this community-wide experiment in some cases were those already being maintained in culture for other experimental studies by each participating lab. They included both eukaryotes (diatoms and dinoflagellates) and prokaryotes (cyanobacteria) collected from a range of environments including polar, temperate and tropical waters ( Figure 1A). A summary of the provenance, distribution and environmental relevance of each species/strain is provided in Table 1. For five species, up to six strains were tested to get an estimate of intraspecific variation in thermal traits (Table 2).

Lab Culturing Protocol
The following protocol was agreed upon by the participating labs and generally implemented here. Inevitably there were minor departures from this consensus-based protocol which are detailed in Table 5. These departures were largely due to laboratoryspecific logistics and provide insights as to how future communitywide manipulation experiments can be improved upon.

1)
For each isolate, growth rates were determined at a minimum of six temperatures for each growth curve.

2)
A minimum of three replicate growth rates were determined at each temperature.

3)
All other environmental variables were held constant within each individual experiment, other than temperature. These include light intensity, day length, culture medium and culture protocols. Saturating nutrient concentrations were used in all experiments to avoid nutrient-induced growth limitation. Light intensity varied among experiments but was most often at saturating levels (determined based on a range of approaches including PE curves and prior investigation of irradiance versus growth responses) for the particular culture isolate (Table 5). Each laboratory used protocols that were most appropriate for the cultured isolates they were working with. A description of culturing conditions for each species is provided in Table 5.

4)
Semi-continuous cultures were diluted using media that had been pre-conditioned to the appropriate temperature. Dilution frequencies were set so that cells were maintained in exponential growth phase and the growth rates reported are for cultures fully acclimated to the experimental conditions, (e.g. after statistically invariant growth rates were recorded for at least 4 generations [98]).

5)
The upper and lower thermal limits for each species/strain were determined by repeatedly (i.e. up to three times) incubating cultures at temperatures at which they did not grow (see Discussion).

6)
Multiple biomass parameters/proxies were used to determine growth rates. These included: cell counts, extracted chlorophyll a, and in vivo chlorophyll a fluorescence. Each method reliably estimated biomass and growth rates, and was used to ascertain when steady-state growth rates were achieved, and to monitor the acclimation of the cultures as outlined above.

7)
At each temperature, the mean steady-state specific growth rate (d 21 ) for each isolate acclimated for at least four generations, was determined by regressing and testing the equality of slopes (i.e. temporal change in the log of fluorescence, cell count or chlorophyll a) from at least three serial cultures (a = 0.05) [99]. If the slopes of serial growth curves were not significantly different, the average regression coefficient was used to estimate the common slope, which represented the average steady-state growth rate in cultures acclimated to a specific temperature and the standard error. We did not impose a standard method to generate the temperature growth rate response curves. The individual methods used are presented in Table 5.

8)
For each species where multiple strains were examined, analysis of variance was used to test for intraspecific differences in growth rates at each temperature examined (a = 0.05) and the coefficient of variation (CV) was calculated to determine the extent of growth rate variation among strains. Statistical analyses were performed in SPSS V19 (IBM, Inc).

Estimation of the Thermal Reaction Norm Parameters
We used the following equation (modified from [100]) to describe the thermal reaction norms of each strain: where specific growth rate f depends on temperature, T, as well as parameters z, w, a, and b. w is the temperature niche width, while the other three possess no explicit biological meaning. We fit (1) to the growth data for each strain using maximum likelihood to obtain estimates for parameters z, w, a and b. In addition, we estimated the optimum temperature for growth and maximum growth rate by numerically maximizing the equation after estimating the parameter values. Note, estimates of temperature niche width for three T. pseudonana strains (CCMP 1014, 1015 & 1335) and one strain of A. sanguinea (YRB) are inflated because of poor resolution of the lower limit for growth (see Results).

Thermal Trait Analyses
These analyses focused on three traits that describe the thermal tolerance curve: optimum temperature for growth, temperature niche width and maximum growth rate.
Comparison with literature traits. We tested for heterogeneity in intraspecific variance of all three traits using Levene's test centered by the median [101]. We also compared the estimated species traits to those estimated from previously published thermal reaction norms of the same species. These data were collated from ten publications (Appendix S1). In cases where different prior published studies measured growth rates on the same strain, we pooled the measurements before estimating the thermal reaction norm parameters. This has the drawback of ignoring interactive effects caused by differences in experimental conditions between previous studies, but was necessary for two reasons. Firstly, in several prior studies only growth rates below the optimum were measured and these data could not be used without pooling, and secondly we have a poor understanding of how the differences in experimental conditions (irradiance, day length, nutrient concentrations, etc.) interact with temperature to determine growth rate, preventing us from incorporating covariates in any meaningful way.
Testing for interspecific differences in thermal traits. We also tested whether there were significant differences among species in thermal trait distributions. In order to account for uncertainty in our estimates of these traits, we used a parametric bootstrapping approach coupled with information criterion-based model selection [102]. This comparison was performed on the five species for which we measured more than a single strain (Thalassiosira pseudonana, Thalassiosira rotula, T. erythraeum, A. sanguinea and C. watsonii). For each strain, we fitted the function described earlier to the growth rate measurements and extracted the residuals from this fit. We then performed 1000 residual bootstraps, a procedure in which the residuals are randomly 'reassigned' to predicted values (each of which corresponds to a growth rate measurement) and added to them, thereby generating a slightly different thermal reaction norm [103]. For each iteration, we refitted the function and estimated the parameters (z, w, a, b) and also the two derived traits, maximum growth rate and optimum temperature for growth. Examining the distribution of these parameters and traits over the 1000 bootstraps allows us to quantify the uncertainty in our estimates, which we can then incorporate in models seeking to explain variation in them.
For each set of bootstraps of all 21 strains (i.e. a relatively small dataset for interspecific comparisons), we fitted two linear models to each trait (optimum temperature, temperature niche width and maximum growth rate) and compared their explanatory ability using AICc values (Akaike Information Criterion, corrected for small sample size) models rather than AIC. One of these models included just an intercept term while the other had both an intercept and species identity as an explanatory variable. We then examined the distribution of delta AICc (AICc [model with species identity] minus AICc [model without]) values over the 1000 bootstraps to determine whether species identity was a useful predictor of each of the three thermal traits. If the entire 95% confidence interval of AICc values exceeded zero, we concluded that species identity was a useful predictor and that species differed in the distribution of the trait.
Comparison of traits with environmental data. We obtained estimates for mean annual and monthly mean temperatures at locations close to the isolation location of each strain using the World Ocean Atlas [91]. We approximated the annual temperature ranges for each species using the difference between the maximum monthly mean temperature and the minimum monthly mean temperature at each location, i.e. thermal variability sensu [46]. Thereafter, we attempted to explain variation in three temperature traits (optimum temperature for growth, temperature niche width, and maximum growth rate) with linear models containing these environmental parameters for each species. The best model was chosen using AIC values. Statistical analyses were performed using a combination of least squares and maximum likelihood estimation techniques using the R statistical computing environment (version 2.15.0). Figure S1 The individual reactions norms (specific growth rates per day) of all measured strains and cultures obtained by fitting a thermal tolerance function (see Methods) to these data. (TIF) Figure S2 Intraspecific variation in thermal reaction norms for species in which only one strain was available.

Supporting Information
(TIF) Figure S3 A summary of temperature optima, maximum growth rates and niche width -expressed as box and whiskers plots -for each of the species used in our study. The black bands denote the median value, the bottom and top of the box represent the 1st and 3rd quartile of the data, respectively. The 'whiskers' extending from the boxes indicate the positions of the lowest & highest values in the data. If the sample size is small enough, the whiskers may not appear (e.g. if there are only 3 equally spaced points, the value represented as the 1st quartile is the lowest value). (TIF) Figure S4 A summary of temperature optima (6C) obtained here (red boxes) and from the literature (blue boxes), expressed as box and whiskers plots, for each of the species used in our study (red). The thick black line in each box represents the median temperature. (TIF) Figure S5 A summary of niche width (6C), expressed as box and whiskers plots -for each of the species used in our study. (TIF) Figure S6 A summary of maximum growth rate expressed as box and whiskers plots -for each of the species used in our study. (TIF) Figure S7 A plot of niche widths versus annual temperature range. Data points are coloured by species. This plot omits the niche widths that are poorly resolved (i.e. the 6 T. pseudonana +1 A. sanguinea strain). Niche widths increase as the annual temperature range increases, in accordance with the climate variability hypothesis. (TIF) Appendix S1 (DOCX)