Phytoplankton life strategies, phenological shifts and climate change in the North Atlantic Ocean from 1850 to 2100

Significant phenological shifts induced by climate change are projected within the phytoplankton community. However, projections from current Earth System Models (ESMs) understandably rely on simplified community responses that do not consider evolutionary strategies manifested as various phenotypes and trait groups. Here, we use a species‐based modelling approach, combined with large‐scale plankton observations, to investigate past, contemporary and future phenological shifts in diatoms (grouped by their morphological traits) and dinoflagellates in three key areas of the North Atlantic Ocean (North Sea, North‐East Atlantic and Labrador Sea) from 1850 to 2100. Our study reveals that the three phytoplanktonic groups exhibit coherent and different shifts in phenology and abundance throughout the North Atlantic Ocean. The seasonal duration of large flattened (i.e. oblate) diatoms is predicted to shrink and their abundance to decline, whereas the phenology of slow‐sinking elongated (i.e. prolate) diatoms and of dinoflagellates is expected to expand and their abundance to rise, which may alter carbon export in this important sink region. The increase in prolates and dinoflagellates, two groups currently not considered in ESMs, may alleviate the negative influence of global climate change on oblates, which are responsible of massive peaks of biomass and carbon export in spring. We suggest that including prolates and dinoflagellates in models may improve our understanding of the influence of global climate change on the biological carbon cycle in the oceans.

Phenological shifts in phytoplankton are expected because of ocean warming (and also freshening in the poles) and its strengthening influence on water column stability, which is expected to diminish nutrients supply in the euphotic zone, favouring small phytoplankton (e.g. nanoplankton and picoplankton) at the expense of diatoms (Bopp et al., 2005;Marinov et al., 2010), a major group thought to be responsible for 40% of total marine primary production (Field et al., 1998;Tréguer et al., 2018).
Although shifts in phytoplankton phenology are now widely observed among marine and freshwater ecosystems (Friedland et al., 2018;Poloczanska et al., 2013;Thackeray et al., 2016), the examination of diatom seasonality in some regions of the North Atlantic suggests a relative stability (Chivers et al., 2020;Edwards & Richardson, 2004). These unexpected results may originate from the range of strategies that have been developed by diatoms and enable them to occur in diverse environments (Kemp & Villareal, 2013, 2018. Among those strategies, diatom morphological traits might have an important influence on species phenology . Indeed, it has recently been shown that diatom cell shape has evolved as a key adaptation that confers to a species a specific phenology, oblate (i.e. flattened) diatoms being dominant in well-mixed, nutrients-rich, low-stratified waters, whereas prolate (i.e. elongated) diatoms dominate in the stratified low-nutrient waters . Cell elongation of prolate diatoms enhances their buoyancy without altering their capacity to absorb nutrients, conferring them a selective advantage in stratified lownutrient waters .
Most projected phenological shifts (Asch et al., 2019;Henson et al., 2013Henson et al., , 2018Yamaguchi et al., 2022), which rely on the analyses of variables such as chlorophyll concentration or net primary production, originate from ESMs that only consider a few phytoplankton types (only one for diatoms and none for dinoflagellates). Therefore, those models cannot anticipate subtle changes in phytoplankton community (Kemp & Villareal, 2018;Séférian et al., 2020;Tréguer et al., 2018) and some discrepancies have been reported between satellite observations and model projections (Cabré et al., 2016).
Here we applied a species-based modelling approach, combined with large-scale plankton observations, to investigate past, contemporary and future long-term changes  in the phenology of oblates, prolates and dinoflagellates using six ESMs and two climate warming scenarios (a medium and a high emission scenario, i.e. ).

| Biological data
Phytoplankton data came from the continuous plankton recorder (CPR) survey. The CPR is a long-term plankton monitoring programme that has collected plankton on a monthly basis in the North Atlantic Ocean and its adjacent seas since 1946. The sampling instrument is a high-speed plankton recorder towed behind voluntary merchant ships (called 'ships of opportunity') at a depth of approximately 7 m (Beaugrand et al., 2003;Reid et al., 2003). As phytoplankton sampling remained unchanged since 1958 (Warner & Hays, 1994) and historical climatic simulations ended in 2014 (see Section 2.2 below), we used the abundance data collected for the diatoms and the dinoflagellates between 1958 and 2014.
Each value of abundance corresponds to a number of cells per CPR sample, which corresponds to ~3 m 3 of seawater filtered (Jonas et al., 2004). Species that were first identified after 1958 were discarded from the analyses. To account for the climatic variability that was observed in the North Atlantic sector during the period (Beaugrand et al., 2019;Edwards et al., 2002) and minimise the noise associated with the CPR sampling (e.g. exceptional high abundance caused by local hydro-meteorological events, patchiness or abundance misestimation), we calculated a daily mean abundance climatology (based on 365 days) for each species for three 18-year periods: 1958-1976, 1977-1995 and 1996-2014 (a total of 365 days × 3 = 1095 days). The three periods were chosen to be of even length. Therefore, we were able to consider the seasonal signal and the long-term trend, which would have been impossible by using a single climatology based on the entire available time period (Mannocci et al., 2017). Climatologies were estimated in three distinct oceanic regions of the North Atlantic: (i) the North Sea (51°N to 60°N, −3°E to 9°E), (ii) the North-East Atlantic (61°N to 63.5°N, −21°E to −8°E) and (iii) the Labrador Sea (48°N to 60°N, −55°E to −40°E). Climatologies were calculated in a given region for the species/taxa that were present in more than 100 CPR samples. A total of 44 species/taxa were therefore selected (Table S1).
Climatologies (for the three periods) were smoothed in each region by means of a double 6-order simple moving average (i.e. 13-day window) and subsequently standardised between 0 and 1, applying the approach of Caracciolo et al. (2021). Standardisation was performed as follows: with A * (i,j) the standardised abundance of species j on day i , A ij the abundance of species j on day i and max A j the maximal abundance of species j in the three regions (i.e. North Sea, North-East Atlantic and Labrador Sea) between 1958 and 2014. By doing so, the variations in A * (i,j) reflect both the differences within and between the three regions. We then estimated the 90th percentile of the abundance (P90) of each species using the non-standardised daily climatologies in the three oceanic regions (Table S1).
Species were divided into three groups: oblate and prolate diatoms and dinoflagellates, which are known to have distinct environmental requirements (Irwin et al., 2012;Kléparski et al., 2022; Table S1). The oblate group gathered diatoms (species or taxa) that have a mean cell diameter greater than their mean cell height and the prolate group the diatoms (species or taxa) that have a mean cell diameter smaller than their mean cell height. Information on diatom cell shapes were retrieved from Kléparski et al. (2022;their figure 3).
Mean monthly abundances between 1958 and 2014 were estimated for each group and in each oceanic region (Figure 1a-i). Data were subsequently smoothed by means of a first-order simple moving average (i.e. 3-month smoothing window) and then standardised between 0 and 1 using the 90th percentile of the monthly abundance of each taxonomic group (P90 m ): With N * m the standardised monthly abundance of a taxonomic group for month m, N m the monthly abundance for month m and P90 m the 90th percentile of the monthly abundance of the taxonomic group. In rare cases, standardised monthly abundance above 1 were fixed to 1. We chose to use the 90th percentile here instead of the maximum abundance value (as in Equation 1) because the latter would have been too sensitive to outliers that may originate from multiple causes (e.g. exceptional abundance, patchiness or abundance misestimated from only 3 m 3 of seawater filtered; see Figure S1).

| CMIP6 climate simulations 1850-2100
Species responses to environmental variability were modelled using three environmental variables known to influence phytoplankton phenology at high latitudes (Beaugrand & Kirby, 2018;Boyce et al., 2017;Caracciolo et al., 2021;Lewandowska & Sommer, 2010;Miller, 2004): that is sea surface temperature (SST; °C), surface downwelling shortwave radiation (SDSR; W m −2 ) and dissolved nitrate concentration (mol m −3 ). Climate projections for the three variables originated from the Coupled Model Intercomparison Project Phase 6 (CMIP6; Eyring et al., 2016) and were obtained from the Earth System Grid Federation. We used the shared socioeconomic pathways (SSP) 2-4.5 and 5-8.5 corresponding, respectively, to a medium and a high radiative forcing by 2100 (4.5 and 8.5 W m −2 ; O'Neill et al., 2016O'Neill et al., , 2017. The simulations of six different ESMs (i.e. CNRM-ESM2-1, GFDL-ESM4, IPSL-CM6A-LR, MPI-ESM1-2-LR, NorESM2-LM and UKESM1-0-LL) covering the time period 1850-2014 (historical simulation) and 2015-2100 (future projections for the two SSP scenarios) were used. An historical perspective is important to better understand the magnitude of the present and future shifts . The six ESMs were chosen on the basis of data availability for 1850-2100 and the two warming scenarios. All the data were interpolated into daily on a 0.5 × 0.5° regular grid. Key references (i.e. DOI and dataset version) are provided in Text S1. Long-term projected changes of the three variables, for the six models, the two scenarios and the three regions are dis- Comparisons between observed and modelled long-term monthly abundance of the three taxonomic groups (oblates, prolates and dinoflagellates). (a-i) Long-term changes in the monthly abundance of (a-c) oblates, (d-f) prolates and (g-i) dinoflagellates in the North Sea (a, d and g), the North-East Atlantic (b, e and h) and the Labrador Sea (c, f and i). (j-r) Modelled long-term changes in the monthly abundance of (j-l) oblates, (m-o) prolates and (p-r) dinoflagellates in the North Sea (j, m and p), the North-East Atlantic (k, n and q) and the Labrador Sea (l, o and r). Differences between observed and modelled mean monthly abundance of (s-u) oblates, (v-x) prolates and (yα) dinoflagellates in the North Sea (s, v and y), the North-East Atlantic (t, w and z) and the Labrador Sea (u, x and α). In a-r, colours denote the mean monthly observed or modelled abundance of a taxonomic group. The abundance in each panel was standardised between 0 and 1 for comparison purpose (Section 2). Spearman correlation coefficients between observed and modelled abundance are displayed in Table 1. In (sα), colours denote the difference between observed and modelled abundance, red and blue colours indicating an underestimation and an overestimation of the modelled abundance, respectively. Modelled abundances are the average based on the six Earth System Models (ESMs).

| The MacroEcological Theory on the Arrangement of Life
A framework from the MacroEcological theory on the arrangement of life (METAL) theory was used to investigate past, contemporary and future phenological changes in the three taxonomic groups.
METAL is a theory that attempts to explain how biodiversity, from the individual to the community level, is organised in space and time and how it responds to environmental changes (Beaugrand, 2015).
Based on the concept of the ecological niche sensu Hutchinson (i.e. the set of conditions enabling a species to growth and reproduce; Hutchinson, 1957), one fundamental assumption of METAL is that the niche-environment interaction is a fundamental interaction that enables one to unify and predict a large number of ecological and biogeographical phenomena, as well as the effect of climateinduced environmental changes on individuals, species and communities (Beaugrand, 2015;Beaugrand & Kirby, 2016, 2018. More information on METAL can be found in https://biodi versi te.macro ecolo gie.climat.cnrs.fr/. Recently, it has been shown that a framework originating from METAL and using the data collected by the CPR survey could also explain phytoplankton phenology and the resulting annual plankton succession in the North Sea, based on SST, SDSR and nitrate concentration (Caracciolo et al., 2021). By creating a local pseudo-community, composed of fictive species (i.e. pseudospecies), it was possible to reconstruct at a species-level phytoplankton phenology and at a community level the annual plankton succession (Caracciolo et al., 2021). We therefore applied the same framework here to explore how diatom (oblate and prolate) and dinoflagellate phenology might be modified in the coming decades with global climate change.

| Niche characterisation
We first estimated daily climatologies (based on 365 days) of the three environmental variables described above (i.e. SST, SDSR and nitrate concentration) for three time periods (1958-1976, 1977-1995 and 1996-2014 i.e. 365 days × 3 = 1095 days), three oceanic re-  SST and SDSR, because of data availability, daily means were estimated in the three regions between 1959 and 2014 and then converted into daily climatologies for the three time periods (1959-1976, 1977-1995 and 1996-2014 Then, we generated a set of 1,755,000 pseudo-species (i.e. fictive species) and we calculated their abundance along the daily climatologies using a three-dimensional Gaussian niche (Caracciolo et al., 2021): with B the abundance of a pseudo-species along a given environmental gradient x, c the maximum abundance of a pseudo-species (here c was fixed to 1), x 1 to x n the environmental gradients, with n = 3 (SST, SDSR and nitrate concentration). x opt1 to x optn are the ecological niche optima along x 1 to x n and t 1 to t n the niche amplitude along x 1 to x n . Niche op-  (Gause, 1934;Hutchinson, 1978). Therefore, the total number of pseudo-species (1,755,000) corresponds to the total number of unique combinations. Niche intervals were chosen to optimise computational cost while covering the largest set of combinations.
We compared the modelled and observed (CPR) abundance at a species and a group (oblates, prolates and dinoflagellates) level.
At a species level, we calculated the mean absolute error (MAE) between the daily abundance climatology (based on 365 days × 3 temporal periods = 1095 days) of each species and each of the 1,755,000 pseudo-species we created, for each oceanic region and ESM. Some species exhibited erratic short peaks that were difficult to explain in some regions. To prevent our model to be influenced by these events that might be related to misidentification, CPR silk contamination or species expatriation (i.e. species occurring in unsuitable environ- ( To assess whether the use of climatologies based on three temporal periods (1958-1976, 1977-1995 and 1996-2014) might influence model skill and performance, we also estimated the MAEs between the observed and modelled daily abundance climatologies based on a single  and two temporal periods (1958-1985and 1986-2014Figures S8 and S9).

| Long-term changes in the abundance of the three taxonomic groups
Daily abundances of the 44 species were estimated from Equation (3) between 1850 and 2100, using environmental variables assessed from the six models and for the two scenarios (Section 2.4 above).
To assess the mean abundance of each group, we weighted species abundance using their P90 (Section 2.1 and Table S1) to account for species contemporary difference in abundance within the three groups (oblates, prolates and dinoflagellates). We used the 90th percentile instead of the maximum abundance value because the latter would have been too sensitive to outliers that may originate from multiple causes (e.g. exceptional abundance, patchiness or abundance misestimated from only 3 m 3 of seawater filtered). Predicted abundances were finally standardised between 0 (the lowest abundance) and 1 (the highest) for each taxonomic group, all ESMs and scenarios as follows: with D * m the standardised predicted abundance of a group for a given ESM, scenario and month m, D m the predicted abundance for month m and min(D) and max(D) the minimal and maximal predicted abundance between 1850 and 2100 in a given oceanic region, respectively. Finally, the relationships between modelled and observed standardised monthly abundances were quantified for 1958-2014 (i.e. overlapping period of observed and modelled data) using the Spearman correlation coefficient (Figure 1a-r and Table 1). Differences between observed and modelled standardised monthly abundances were also calculated to compare the seasonal and long-term trends (Figure 1sα).
To test whether the long-term changes in phenology and abundance of the three groups were well modelled by our approach, we also estimated the annual (long-term changes) and monthly (d-f) Relationships between the Spearman rank correlation (r Spearman ) calculated between the observed and the modelled abundances, and the corresponding 90th percentile of the abundance (P90; see Table S1) values for (d) oblates (blue dots), (e) prolates (red dots) and (f) dinoflagellates (black dots) in the three oceanic regions. Displayed MAEs and Spearman rank correlations were calculated for each Earth System Model (ESM) (see Tables S2-S5). Comparison between the modelled abundance forced by the six ESMs and the observed abundance by the Continuous Plankton Recorder (CPR) survey are shown in Animations S1-S3. Mean monthly observed and predicted standardised abundances of the three taxonomic groups in the three oceanic regions are shown in Figure 1.
subsequently compared by means of a Spearman rank correlation coefficient (Tables S6 and S7).

| Long-term phenological changes in the three taxonomic groups
We examined phenological shifts of the three taxonomic groups be- Tables S8-S16; Kwiatkowski et al., 2020).

| RE SULTS
Visual comparisons between observed and reconstructed SST, SDSR and nitrate concentration showed that the seasonal cycle and the long-term trend of each variable was in general well reconstructed by the ESMs in the three oceanic regions and for the three time periods (except for nitrate long-term trend, see Section 2 and Figures S5-S7).

| Model performance
Model performance was first assessed by visually comparing the observed and modelled daily abundance of each species for the three periods and North Atlantic regions (Animations S1-S3). The phenology of most species was well modelled by our approach and the MAE between the observed and modelled abundance was generally low  Table S2; see Section 2). MAEs calculated between observed and modelled abundance based on a single  or two temporal periods (1958-1985 and 1986-2014) climatologies were higher, indicating that our model performed better when based on the three temporal periods climatologies ( Figures S8 and S9). We then estimated the Spearman rank correlation coefficient between daily observed and modelled patterns, which revealed that for most species, correlations were highly positive (Figure 2d-f and Tables S3-S5). A few low correlations were observed, however, indicating that for some species our approach did not work well in some North Atlantic regions and for some ESMs. Low correlations between observations and model reconstructions may be caused by some regional hydro-meteorological events that are not well reproduced by an ESM or by some methodological limitations of our approach (see Section 4). However, their weight in the subsequent analyses was low, as shown by the value of their P90 (i.e. the 90th percentile of their abundance, Section 2, Figure 2d-f and Table S1); the higher the P90, the higher the abundance of a species in a North Atlantic region and therefore the greater its weight in subsequent analyses at the group level (Figures 4-6).
Modelled species abundances were finally aggregated at a monthly scale for each group (oblates, prolates and dinoflagellates) by weighting the abundance of each species according to the value of their P90 (Table S1). A monthly scale was first chosen for comparison of the modelling approach with the in situ monthly observations collected by the CPR survey, whom sampling is carried out at a monthly interval (Richardson et al., 2006). Differences between observed and modelled long-term changes in monthly abundance showed that model perfor- Atlantic, and 0.5 and 0.58 in the Labrador Sea; Figure 1 and Table 1).
The comparison between long-term changes in the annual abundance of prolates in the three regions, of oblates in the North-East Atlantic and of oblates and dinoflagellates in the Labrador Sea were globally TA B L E 1 Spearman rank correlations between long-term changes in observed and modelled monthly abundance of the three taxonomic groups (oblates, prolates and dinoflagellates) in the three studied regions. Note: All correlations were significant (p < .01). The degree of freedom associated with each oceanic region is also indicated. Long-term monthly changes in observed and modelled abundance are displayed in Figure 1.
well reconstructed ( Figure S10), but the correlations were not always significant, suggesting that some periods and years did not match well (Table S6). On the contrary, our model did not fully reconstruct the long-term trends in oblates and dinoflagellates in the North Sea and of dinoflagellates in the North-East Atlantic ( Figure S10), which was confirmed by Spearman rank correlations (Table S6). Such difficulties to explain observed long-term changes in abundance were expected since ESMs are known to have some difficulties to reconstruct highfrequency natural variability in climate well (Stock et al., 2011). Finally, the comparison between observed and modelled monthly changes in abundance revealed that the seasonal cycles of the three groups in the three regions were correctly reconstructed (significant Spearman rank correlations above 0.85; Figure S11 and Table S7).

| Phenological shifts in the North Atlantic
In  (Figures 4-6). Expectedly, the magnitude of the projected changes in phenology and abundance appeared to be more important under the high warming scenario (Figures 4-6).

| Phenological shifts in the North Sea
In the North Sea (Figure 4a), modelled long-term daily phenological changes exhibited similar patterns for both SSP scenarios, with a shift in the diatom community composition, although changes were highest under SSP5-8.5 (Figure 4b- ,e and Table S8). In contrast, prolate MAA rose  Table S10). The duration of their SRP is expected to rise throughout the 21st cen-  Table S10).

F I G U R E 3
Theoretical diagram explaining graphically the six phenological indices used in this study. Six phenological indices were defined: the Maximum Annual Abundance (MAA; Index 1), the day where MAA is reached (Index 2), the initiation (Index 3) and the termination (Index 4) of the Seasonal Reproductive Period (SRP), that is the first and the last days where the abundance is ≥50% of MAA, the seasonal duration (Index 5), that is the number of consecutive days where abundance is ≥50% of MAA, and the Integrated Mean Annual Abundance (IMAA; Index 6).

| Phenological shifts in the North-East Atlantic
In the North-East Atlantic region (Figure 5a), phenological shifts were also predicted by our approach, the magnitude of which depending on warming intensity (Figure 5b  in 1850-1859 and 204 ± 16 or 197 ± 19 for SSP2-4.5 and 5-8.5, respectively) and a later termination (around day number 254 ± 9 in 1850-1859 and 262 ± 5 or 268 ± 9 in 2090-2099 for scenarios SSP2-4.5 and SSP5-8.5, respectively, Table S13) of their SRP. Their IMAA rose above those of oblates during the first half of the 21st century, although there was a high inter-model variance (Figure 5p-q and Table S13).

| Phenological shifts in the Labrador Sea
In the Labrador Sea (Figure 6a), distinct phenological shifts were found, also depending on the SSP scenarios (Figure 6b Figure 6h and Table S15) accompanied by an earlier phenology ( Figure 6k and  Figure 6d and Table S16) without a significant increase in the duration of their SRP (Figure 6h and Table S16) or a change in initiation or termination (Figure 6l and Table S16). Their IMAA slightly increased but remained above that of oblates (Figure 6p and Tables S14 and S16).
Under SSP5-8.5, our approach predicted a shift from a system dominated by oblates to a system dominated by dinoflagellates (Figure 6b-q). The oblates exhibited a decrease in their MAA from   Figure 6g and Table S16). These changes were more prominent at the end of the 21st century, from 66 ± 12 days in 2000-2009 to 82 ± 18 in 2090-2099 ( Figure 6i and Table S16), resulting in a phenological dilatation (Figure 6o and

| DISCUSS ION
Some discrepancies between modelled and observed abundance indicate that our model may have some limitations (Figure 1, Figures S10 and S11 and Animations S1-S3). For example, the decrease in the abundance of dinoflagellates in the North Sea and the North-East Atlantic, and the increase of oblates in the North Sea, was not captured well by our approach (Figure 1 and Figure S10).
In the same way, the fit between observed and reconstructed daily abundances varied among the three different temporal periods (Animations S1-S3). For example, the abundance of Leptocylindrus danicus in the North Sea was correctly reconstructed for the temporal periods 1977-1995 and 1996-2014 but not for 1958-1977 (Animation S1). These inconsistencies may originate from (i) our framework, (ii) the CPR sampling and (iii) the ESMs.

| Limitations of our approach
In this study, we made a series of choices and assumptions that might have affected model performance. First, we used an approach based on METAL rather than a species distribution model (SDM).
This choice was made because we showed previously that METAL correctly reconstructs phytoplankton phenology and annual succession when applied to the CPR data (Caracciolo et al., 2021) while it has been shown that most SDM poorly performed with those data (Brun et al., 2016). In addition, using METAL enables us to convert the monthly CPR observations at a daily scale, which is more adapted (i) to reconstruct/project the long-term phenological changes from 1850 to 2100, (ii) to examine their trends  and more importantly (iii) to demonstrate that the niche-environment interaction governs long-term phenological phytoplankton shifts (Beaugrand & Kirby, 2018).
Second, we used a symmetric Gaussian function to model the thermal niche of the 44 phytoplanktonic species. Phytoplankton are ectotherms and may therefore have asymmetric thermal niches Martin & Huey, 2008;Thomas et al., 2012); in particular, they might be able to tolerate a larger range of cold temperatures (i.e. below their thermal optimum) and a narrower North Sea using symmetrical niches. Furthermore, the use of asymmetrical Gaussian functions would have also increased computational cost. Nevertheless, we acknowledge that the consideration of asymmetrical niches might perhaps improve our results. For example, the existence of an upper thermal limit above which the abundance rapidly declined substantially (e.g. the events observed with the dinoflagellates in the North Sea and the North-East Atlantic; Figure 1g,h and Figure S10c,f) was not reproduced by our model; indeed, the modelled shifts were almost always gradual (Figure 1p-q and Figure S10c,f).
Last, we supposed that species niches will be conserved over the time period 1850-2100 (i.e. the assumption of niche conservatism), an assumption that is generally formulated by modellers that project species response to climate change (Soberon & Nakamura, 2009).
The large population size and fast generation time of phytoplankton are supposed to allow a rapid alteration of their niches  but some studies (on zooplankton however) have not seen any niche alteration at a decadal scale (Helaouët & Beaugrand, 2007) and there must exist fundamental constraints that limit adaptation (e.g. phytoplankton morphological traits, which affect nutrient uptake, are supposed to be highly conserved; Beaugrand & Kirby, 2018;Kléparski et al., 2022;Litchman et al., 2007).

| Limitations related to the CPR sampling
Another source of discrepancy between modelled and observed abundance may originate from the sampling by the CPR. First, sampling is carried out at a monthly scale over large areas (Richardson et al., 2006) and it is likely that such a sampling frequency is not al-  (Tables S3-S5 and Animations S1-S3) suggesting that the problem might not be related with a sampling bias affecting model performance. Third, the CPR sampling is heterogeneous in space and time and some regions were not continuously sampled between 1958 and 2014 (Richardson et al., 2006). For example, no CPR sample was available in the Labrador Sea at the end of the 1980s (white bands in Figure 1c,f,i,u,x,α). However, it is thought that the CPR collects an important fraction of the abundance of each taxon and is therefore though to be robust to examine seasonal and interannual patterns (Richardson et al., 2006).

| Limitations from ESMs
The reconstructed regional climatic variability between 1958 and 2014 by ESMs is not supposed to exactly match the variability that was observed during that time (Stock et al., 2011), a bias that might have also affected niche estimates and that might therefore explain the discrepancies between observed and modelled mean long-term changes in abundance (Figure 1a-r and Figure S10). For example, the abundance of C. closterium was well modelled in the North Sea with CNRM-ESM2-1, IPSL-CM6A-LR and NorESM2-LM but not with MPI-ESM1-2-LR, GFDL-ESM4 and UKESM1-0-LL (Table S3)

| Long-term changes in the phenology of diatoms and dinoflagellates
Modelled shifts are in agreement with theoretical projections, that is that spring species (i.e. oblates) experience a phenological contraction and an earlier phenology, whereas late summer species (i.e. prolates and dinoflagellates) exhibit a phenological expansion with both an earlier initiation and a later termination (Beaugrand & Kirby, 2018). The phenological contraction of oblates, associated with a decline in abundance, is indicative of species/taxa that reach the limits of their phenological plasticity (Beaugrand & Kirby, 2018), probably because an earlier phenology in spring is impossible as winter photoperiod and light levels limits primary production at high latitudes (Boyce et al., 2017;Caracciolo et al., 2021) and/or because the decrease in nitrate concentration ( Figures S2-S4), which are an essential constituent of proteins and nucleic acids that are, in turn, critical for the synthesis of organic matter (Miller, 2004), limits the size of the bloom. Therefore, these diatoms might encounter more favourable environmental conditions at high latitudes, causing a biogeographical shift, a postulate that would be in agreement with the increased primary production projected by ESMs in the Arctic region  and the observed increase and decrease in diatom populations in the North-Sea and the bay of Biscay, respectively . On the contrary, the phenological expansion of prolates and dinoflagellates is indicative of species/ taxa that are adapted to the reinforced stratified and nutrientdepleted conditions in the future North Atlantic Ocean (Beaugrand & Kirby, 2018;Kwiatkowski et al., 2020). The elongated cell shape of prolates reduces their sinking speed while not affecting their capacity to uptake nutrients, enabling them to remain in the well-lit upper part of the water column Padisak et al., 2003).
Furthermore, it has been shown that some taxa belonging to this taxonomic group (e.g. Rhizosolenia) were able to passively migrate between the surface layers and the nutricline, enabling them to harvest nutrients more efficiently (Kemp & Villareal, 2013, 2018. Similarly, the dinoflagellates have flagella which enable them to move and actively exploit nutrients in the water column. They are also capable of switching nutritional strategies to mixotrophy when less nutrients are available (Miller, 2004).
Because of the shifts in phytoplankton phenology, changes in species succession are expected to desynchronise species interaction and trigger trophic mismatch (Beaugrand et al., 2003;Cushing, 1990;Edwards & Richardson, 2004). For a given season, mismatch may occur when some species exhibit different range of tolerance for a given environmental variable. Such alterations of the trophodynamics have already been observed in freshwater lakes (Winder & Schindler, 2004) but remain difficult to anticipate in the marine environment because each species has a specific response to environmental variability, which is driven by its life strategy and its niche (Beaugrand & Kirby, 2018). Furthermore, it has also been shown that only seasonally heterogeneous environmental changes can lead to a desynchronisation in species interaction (e.g. warmer temperature in spring but not in summer; Straile et al., 2015). Results from the METAL theory suggest that trophic mismatch might be limited however, because interacting species have at least a part of their ecological niche in common, providing that no other important ecological dimensions or habitat requirements differ among them (Beaugrand & Kirby, 2018). Indeed, species are more expected to track changes along the environmental variables for which they have the narrowest tolerance (Ackerly, 2003).
Therefore, in the pelagic environment, one could expect a mismatch among species that are not linked to the same habitat components (i.e. the stable-and substrate-biotope components; van der Spoel, 1994).
For example, it has been shown that a shift in phytoplankton phenology is more likely to affect the fish that spawn in geographically fixed areas rather than those that have their spawning grounds moving according to the environmental changes (Asch et al., 2019).
Diatoms and dinoflagellates have been historically separated into two distinct functional types (Kemp & Villareal, 2018;Margalef, 1978). According to Margalef's mandala (Margalef, 1978), diatoms thrive when the water column is poorly stratified and turbulent and when nutrients concentration is relatively high. In contrast, dinoflagellates dominate in stratified conditions and lower nutrients concentration. By simplifying Margalef's framework (Margalef, 1978), a decline in primary production and carbon export (therefore the biological pump; Passow & Carlson, 2012) is anticipated in the North Atlantic because ocean warming will enhance water column stratification and reduce nutrients, which are believed to negatively affect the diatoms as a whole (Bopp et al., 2005;Kemp & Villareal, 2018;Kléparski et al., 2022;Marinov et al., 2010;Tréguer et al., 2018). Our results nuance these projections. Although they show that diatom abundance as a whole will decrease, they also suggest that prolates may indeed persist and increase in a more stratified ocean, albeit probably not to an abundance level that might compensate the large diminution in the abundance experienced by oblates. Indeed, massive blooms of those last diatoms in spring generate high carbon export through the generation of marine snow Raven & Waite, 2004;Smetacek, 1985) and it is therefore unlikely that the slight increase in prolates will overcome the decrease in carbon exportation induced by the reduction in oblate abundance. Furthermore, the transfer efficiency of carbon to the deep ocean is reduced with smaller sized slow-sinking cells, and therefore the shift from large round oblates to smaller elongated prolates will also diminish the rate of carbon export (Henson et al., 2022). Whether this shift might also affect species aggregation remains an open question as the causes of this phenomenon are not clear (e.g. defence against predation, response to changes in turbulence or buoyancy regulation) (Lürling, 2021;Pančić & Kiørboe, 2018;Smayda, 1970;Sournia, 1982). Our results also suggest that dinoflagellate abundance will rise in the three North Atlantic regions.
Although the contribution of this group to the biological pump is poorly understood, a study has shown that dinoflagellate blooms sometimes coincide with a rise in carbon export at an annual scale in the North-East Atlantic at 3000 m (Henson et al., 2012). It has also been shown that the consideration of mixotrophic organisms (e.g. dinoflagellates) within biogeochemical models can increase the flux of carbon export (Ward & Follows, 2016). Therefore, this group might alleviate to some extent the consequences of the reduction in oblates on the biological pump in the North Atlantic, albeit export from the surface ocean is also strongly influenced by synergistic effects such as nutrient delivery through remineralisation or direct inputs, processes that are difficult to predict (Passow & Carlson, 2012).
ESMs used as part of CMIP5 included a few functional types to represent the phytoplanktonic community, often with a single grouping for diatoms Tréguer et al., 2018).
The last generations of ESMs in CMIP6 have improved the representation of phytoplankton ecology by including more functional types (e.g. diazotrophs for some models), which has enhanced model performance. However, diatoms still remain represented by a single group . Results from model inter-comparisons highlight the great sensitivity of models to the types of functional groups they include, as well as the ecological traits they represent (Fu et al., 2016;Taucher & Oschlies, 2011). Dutkiewicz et al. (2020) have recently used morphological traits (i.e. cell size), biogeochemical functions and thermal tolerances to define more than 350 phytoplankton types. Their results showed that species richness projections were altered by the number as well as the nature of phytoplankton traits, each controlling a distinct aspect of the species' responses (Dutkiewicz et al., 2020).
They concluded that the addition of traits, such as cell shape, might improve model projections (Dutkiewicz et al., 2020). In situ studies have also highlighted the importance of phytoplankton composition for carbon exportation (Dybwad et al., 2021;Henson et al., 2012). However, such models with hundreds of phytoplankton types cannot realistically be used to perform centennial climatic simulations because of computational cost. Hence, although our results with oblates give some support to the ESMs that only considered a single diatom group, they also emphasise that prolates and dinoflagellates exhibit distinct alterations of their phenology and abundance. Therefore, the distinction of oblate and prolate diatoms as well as dinoflagellates in ESM models may be a simple way to improve model projections because these three groups have different phenologies and responses to climate change, which may affect our current understanding of the consequences of phytoplankton shifts on the biological carbon pump (Karp-Boss & Boss, 2016;Kléparski et al., 2022;Litchman & Klausmeier, 2008;Naselli-Flores et al., 2021;Ryabov et al., 2021).

AUTH O R CO NTR I B UTI O N S
L.K. and G.B. designed the study, performed the analyses and wrote the first draft. L.K., G.B., M.E. and C.O. discussed the results and reviewed the manuscript.

ACK N OWLED G M ENTS
We acknowledge the World Climate Research Programme, which, through its Working Group on Coupled Modelling, coordinated and promoted CMIP6. We thank the climate modelling groups for producing and making available their model output, the Earth System Grid Federation (ESGF) for archiving the data and providing access, and the multiple funding agencies who support CMIP6 and ESGF.
We also thank the owners and crews that have towed the CPRs on a

CO N FLI C T O F I NTE R E S T S TATE M E NT
The authors declare no conflict of interest.

DATA AVA I L A B I L I T Y S TAT E M E N T
All the data used in this article are already freely available (see Section 2 and Text S1).