Plankton biogeography in the North Atlantic Ocean and its adjacent seas: Species assemblages and environmental signatures

Abstract Plankton biodiversity is a key component of marine pelagic ecosystems. They are at the base of the food web, control the productivity of marine ecosystems, and provide many provisioning and regulating ecological services. It is therefore important to understand how plankton are organized in both space and time. Here, we use data of varying taxonomic resolution, collected by the Continuous Plankton Recorder (CPR) survey, to map phytoplankton and zooplankton biodiversity in the North Atlantic and its adjacent seas. We then decompose biodiversity into 24 species assemblages and investigate their spatial distribution using ecological units and ecoregions recently proposed. Finally, we propose a descriptive method, which we call the environmental chromatogram, to characterize the environmental signature of each plankton assemblage. The method is based on a graphic that identifies where species of an assemblage aggregate along an environmental gradient composed of multiple ecological dimensions. The decomposition of the biodiversity into species assemblages allows us to show (a) that most marine regions of the North Atlantic are composed of coenoclines (i.e., gradients of biocoenoses or communities) and (b) that the overlapping spatial distribution of assemblages is the result of their environmental signatures. It follows that neither the ecoregions nor the ecological units identified in the North Atlantic are characterized by a unique assemblage but instead by a mosaic of assemblages that overlap in many places.

the biological carbon pump (Brun et al., 2019). Plankton are also good indicators of climate change impacts because they are in general not commercially exploited and have a relatively short life cycle.
Understanding the influence of climate change on plankton necessitates having a good understanding of its spatial distribution and its natural annual variations (Reid & Edwards, 2001). However, large-scale plankton-monitoring programs are rare (Richardson & Poloczanska, 2008). In the North Atlantic and some of its adjacent seas, plankton have been investigated by the Continuous Plankton Recorder (CPR) survey on a monthly basis since 1946 (Reid et al., 2003).
Recently, the survey has been implicated in the European research project ATLANTOS, gathering 18 countries with the aims of bringing together all the existing Atlantic Ocean observing activities into a more integrated wide observation system (AtlantOS, 2019).
As part of the ATLANTOS project, a new partition of the North Atlantic Ocean and adjacent seas has been proposed . This new partition divides the North Atlantic into 13 ecological units and 40 ecoregions (see Figures S1 and S2).
An ecological unit is a group of observations with a homogeneous environmental regime, similar biodiversity and seasonal variability.
It can be subsequently divided into ecoregions that are a group of interconnected observations .
The partition was primarily based on spatiotemporal fluctuations in plankton biodiversity but also physical data (e.g., bathymetry, SST, mean surface current). This new partition was built at a relatively fine spatial resolution (i.e., 0.5° latitude × 0.5° longitude) in contrast to previous global studies (Longhurst, 1998) and was the result of the analysis of six key plankton groups (dinoflagellates, diatoms, small and large copepods, small and large zooplankton other than copepods). However, the study did not decompose the biodiversity into species groups, being based exclusively on spatial and temporal (i.e., 2-month time periods) changes in taxonomic richness.
The main objective of this study was to decompose the biodiversity into species assemblages and to relate them with the partition proposed by . Such a work is needed to provide an information on the biotic composition of the ecological units and ecoregions of the North Atlantic and its adjacent seas. Although a study has identified assemblages of calanoid copepods (108 species or taxa), no joint analyses have been performed on both phytoplankton and zooplankton collected from the CPR survey (Beaugrand et al., 2002). We identified taxonomic assemblages using all taxa (species, genera, or higher taxonomic resolution) recorded by the CPR survey. All phytoplankton and zooplankton were considered and combined together to be associated into groups of taxa. First, we gridded spatially the abundance of all taxa for each 2-month period . We assumed that spatial variance was much more pronounced than temporal (i.e., year-to-year to decadal) variance. To minimize the potential influence of year-to-year and decadal variability, we performed data smoothing. Then, we mapped the biodiversity of phytoplankton and zooplankton. We used a cluster analysis to examine the relationships among ~300 plankton taxa that had a level of abundance sufficiently high to be analyzed. Then, we examined the spatial distribution and annual changes in each assemblage. We used nine environmental parameters to characterize the environmental signature of each assemblage (i.e., a combination of environmental variable that characterizes a species group) using a procedure that leads to the display of what we propose to call an environmental chromatogram (i.e., a graphic that identifies where species of an assemblage aggregate along an environmental gradient composed of multiple ecological dimensions). We prefer the term environmental signature (i.e., the environmental regime where species aggregation within an assemblage is highest) instead of the term ecological niche because the latter is usually restricted to the species level (Chase & Leibold, 2003). Finally, we examined the composition of each ecoregion and ecological unit (sensu  in terms of assemblage.

| Physical data
We used a set of nine physical variables: sea surface temperatures (SST), bathymetry, monthly mean downward solar radiation flux at surface (DSRF), macronutrients (nitrate, phosphate, and silicate), distance to nearest coastlines, density mixed layer depth (DMLD), and sea surface salinity (SSS). Information on data sources can be found in Appendix 1.

| Biological data, the CPR survey
Biological data originated from the Continuous Plankton Recorder (CPR) survey. It is a long-term plankton-monitoring program currently operated by the Marine Biological association of Plymouth.
The CPR is the longest and most extensive program of that kind in the world. The machine is a high-speed plankton recorder towed behind voluntary merchant ships, called "ships of opportunity," that filters phytoplankton and zooplankton at a depth of ~7 m (Hays & Warner, 1993). The taxonomic resolution of the data used in this study is shown in Table S1. More information about the CPR survey (advantages and limitations) can be found in Text S1 and in Section 4.1 of the Discussion. Information on data sources can be found in Appendix 1.

| Mathematical analyses
Six main analyses, all coded in MATLAB, were performed in this study ( Figure S3).

| Analysis 1. Spatial regularization
Spatial sampling by the CPR survey is heterogeneous. Therefore, we first created a spatial grid for every plankton species/taxa sampled by the CPR survey for every two-month period using data collected between 1948 and 2016. The spatial grids were identical to the ones used in . For each geographical cell of 0.5° of latitude and 0.5° of longitude from 80.5°W to 9.5°E and from 40.5°N to 65.5°N, we calculated average abundance values of each plankton taxa for every two-month period based on the period 1948-2016. The procedure led to a three-dimensional matrix of 304 taxa × 9,231 geographical cells (latitudes × longitudes) × 6 two-month periods.
2.3.2 | Analysis 2. Spatial smoothing of the gridded data Abundance data can be highly variable from a geographical cell to another which can be in part attributed to the CPR sampling (e.g., variable seawater filtered by the machine, ship speed) and year-toyear variability (Jonas et al., 2004). We therefore smoothed spatially all two-month abundance grids for each species/taxa. Smoothing was performed by applying a c-order spatial simple moving average: y i,j is the smoothed abundance value of a species/taxa at the geographical cell corresponding to latitude i and longitude j and x is the original abundance value for a species/taxa. Threshold selection c depends upon the size of geographical cells (here a cell is 0.5° latitude × 0.5° longitude), the noise inherent to the data, and the type and the location of spatial structures. By trial and error, we fixed c to 2 as a compromise between noise reduction and potential numerical artifacts. Many values can be missing in some areas, and we fixed to 20 the maximum number of missing values allowed to have an estimation.

| Analysis 3. Spatial distribution of the taxonomic richness
The spatial distribution of taxonomic richness was mapped using the 3-dimensional matrix resulting from Analyses 1 and 2. For this analysis, abundance was converted into presence-absence using a threshold of 0 (i.e., abundance > 0 means an occurrence). The use of other thresholds did not affect our conclusions. We then summed all presence for phytoplankton and/or zooplankton species or taxa and smoothed the resulting matrix by applying a first-order spatial triple moving average (c = 1) to obtain a map of the spatial distribution of taxonomic richness.
The analysis was performed for total (i.e., phytoplankton and zooplankton) taxonomic richness, phytoplankton taxonomic richness, and zooplankton taxonomic richness. Among the 304 taxa we used, 149 were phytoplankton, 155 were zooplankton; ~60% of plankton taxa were identified at the species level ( Figure 1 and Table S1).

| Analysis 4. Identification of taxonomic assemblages
We then decomposed the biodiversity into taxonomic assemblages.
We first calculated a squared distance matrix (taxa × taxa) using the Hellinger distance coefficient, which is robust to a high number of double zeros (Legendre & Legendre, 1998): with y 1j the abundance of the species/taxa 1 in geographical cell j, y 1+ , the total abundance of species/taxa 1 across all geographical cells, y 2j the abundance of species/taxa 2 in geographical cell j, y 2+ the total abundance of species/taxa 2 across all geographical cells.
The distance coefficient was calculated between species/taxa on the basis of their patterns of abundance in space (i.e., all geographical cells covering the Atlantic Ocean) and time (i.e., two-month period).
Prior to the use of Hellinger's distance coefficient, abundance data were transformed using the function log 10 ( x + 1 ), a procedure frequently applied to the CPR data that also clearly limit the Euclidean's distance paradox.
A cluster analysis was subsequently applied using Ward's minimum variance method (Legendre & Legendre, 1998). The resulting dendrogram is in Figure S4. Spatial distribution of all assemblages was subsequently mapped by applying the same procedure to represent the spatial distribution of taxonomic richness (Analysis 3) but at the assemblage level. For each assemblage, we mapped the percentage of species/taxa aggregation (i.e., percentage of co-occurring taxa of a given taxonomic assemblage in a cell) by averaging the maps built for each 2-month period (1948-2016; Figure 2). The procedure was also applied for each two-month period to examine seasonal changes in each taxonomic assemblage ( Figure 3).
Taxonomic composition of each assemblage was characterized by the use of pie charts in two ways: • First, we determined the number of species/taxa belonging to six categories: (i) diatoms, (ii) dinoflagellates, (iii) other phytoplankton, (iv) large copepods, (v) small copepods, and (vi) other zooplankton ( Figure S5).
• Second, we determined the number of taxa identified at (i) a species level and (ii) a higher taxonomic resolution ( Figure S6).
2.3.5 | Analysis 5. Relationships between taxonomic assemblages and both ecoregions and ecological units identified by  From Analysis 4, we calculated the percentage of geographical cells with a percentage of taxonomic aggregation (i.e., percentage of co-occurring taxa of a given taxonomic assemblage in a cell) higher than 10% (Figure 4 and Figure S7) and 50% ( Figures S8 and S9) in  coast, all data were linearly interpolated for each two-month period. We therefore obtained a matrix 9,231 geographical cells × 6 by two-month period for each environmental variable; note that for bathymetry and distance to coast, the same values were repeated for each two-month period. By this way, it was possible to relate the abundance of each species/taxa with any environmental variable in space and time. All values of each environmental variable were then standardized between 0 (i.e., lowest value) and 1 (highest value).
We then divided all environmental values between 0 (lowest values) and 1 (highest values) into 100 categories and calculated the abundance of each species/taxa that corresponded to each environmental category between 0 and 1. The choice of the categories (100) resulted from a compromise between the resolution of the F I G U R E 1 Spatial distribution of total taxonomic richness (a), phytoplankton taxonomic richness (b), and zooplankton taxonomic richness (c) sampled by the CPR survey in the North Atlantic Ocean. The magenta line corresponds to isotherm 10°C. The blue lines correspond to current velocities from 0.5 to 2 m/s. The black lines denote isobath 200 m. We caution that phytoplankton and zooplankton biodiversity maps were based on taxa identified at the species (72.5% for phytoplankton and 47.1% for zooplankton), genus (26.2% and 23.9%), and higher taxonomic ( We also characterized the environmental signature of all taxonomic assemblages by considering phytoplankton and zooplankton separately ( Figure S10). This comparison was made to check whether signatures were identical. By this way, we checked the homogeneity of phytoplankton and zooplankton signatures among each species assemblage. We then compared patterns of environmental signature for phytoplankton and zooplankton by means of a correlation analysis. We used a Spearman correlation coefficient tested by a Monte Carlo test using 10,000 permutations to estimate the probability (Jackson & Somers, 1989). Figure 1a shows the spatial distribution of taxonomic richness in the North Atlantic Ocean and its adjacent seas based on 304 phytoplankton and zooplankton species/taxa. A strong gradient was observed between the northern and the southern part of the oceanic basin, with an increasing taxonomic richness toward the equator (we used data with different degrees of taxonomic resolution, see Table S1). A second gradient was also observed between the western and the eastern sides of the Atlantic Ocean, with an increasing taxonomic richness eastward. Biodiversity was greatest south of the 10°C isotherm (magenta line in Figure 1) with taxonomic richness varying between 100 and 140 taxa in the south and between 60 and 80 in the north. Higher taxonomic richness occurred in areas F I G U R E 2 Mean spatial distribution of taxonomic assemblages identified by the cluster analysis (see Figure S4 and Table S1) and based on an average of six 2-month period . On the top right of each panel, the number of species/taxa (n) in the assemblage is indicated. The assemblage number, corresponding to the numbers in Figure S4, is displayed in bold on the top left. The green line corresponds to isotherm 10°C. The blue lines correspond to current velocities from 0.5 to 2 m/s. The black lines denote the isobath 200 m. Panels are classified from a to x by decreasing taxonomic richness characterized by warmer temperatures (e.g., West European Basin) and high current velocity (e.g., the Gulf-Stream). Similar patterns in taxonomic richness were observed for phytoplankton ( Figure 1b) and zooplankton ( Figure 1c) although maximal taxonomic richness was observed in the North Sea for phytoplankton and in the West European Basin for zooplankton. Zooplanktonic taxonomic richness was higher than phytoplanktonic taxonomic richness in the Gulf-Stream Extension.

| RE SULTS
Plankton biodiversity was then decomposed into 24 plankton assemblages ( Figure 2) by means of a cluster analysis (see Figure S4).
We chose a cutoff level of 1.7 in the dendrogram to select a maximum of groups. When thresholds were too high (i.e., >1.7), groups remained too spatially and/or temporally heterogeneous, and when thresholds were too small, we had too many isolated species. Five assemblages (Assemblages 16, 23, 12, 19, 22, and 21; Figure 2a-d,f-g) occurred nearly everywhere in the studied area; they were therefore eurygraphic. All the remaining assemblages were located south of the  Figure S1. For this analysis, we used a threshold of taxonomic aggregation of 10% (i.e., species or taxa that had a percentage of cooccurrence higher than 10% in a given geographical cell for a given assemblage). Other thresholds were tried and did not alter substantially our conclusions when they were fixed below 50% (see Figures S8-S9).  Figure 4i,j,l,s,u) occurred in EUs covering the continental shelf (e.g., the Cold-Temperate Shallow Neritic, the F I G U R E 4 Percentage of taxonomic aggregation (i.e., number of species or taxa of the same assemblage) greater than 10% in each ecological unit as defined by ).
Here, the percentage of species/ taxa aggregation is used to identify an assemblage characteristic of an ecological unit. On the top right of each panel, the number of species/taxa (n) in the assemblage is indicated. The assemblage number, corresponding to the numbers in Figure S4 and in Figure 2, Figure 2 and Figure S4), and the number at the top right (n) indicates the taxonomic richness of each assemblage Cold-Temperate Neritic, and the Ocean-Influenced Cold-Temperate Neritic EUs) as well as the shelf-edge for Assemblage 5 (Figure 4s) 1,18,4,9,11,and 7 (Figure 4k,o,p,q,t,v,w,x), occurred mostly in EUs close to or over the European shelf-edge, for example, the Diverse and productive Oceanic Temperate EU, the Pseudo-Oceanic Warm-Temperate EU, and the Mixed Coastal-Oceanic Highly-Seasonally dynamical EU. Therefore, all assemblages were linked to specific EUs of the North Atlantic Ocean and its adjacent seas.
Similar conclusions were reached for ecoregions, although results were more difficult to interpret as ecological units originate from the aggregation of ecoregions (see Figures S8-S9); therefore, the number of ecoregions (40) is higher than the number of ecological units (13) and it is more difficult to find a pattern on a higher number of units.
We Then, we compared the environmental signature of phytoplankton and zooplankton species/taxa composing each assemblage in order to check whether signatures were identical. By this way, we verified the homogeneity of phytoplankton and zooplankton signatures among each species assemblage. Expectedly, we found the same environmental signature for phytoplankton and zooplankton of the same assemblage ( Figure S10). We also compared the environmental signature of phytoplankton and zooplankton for each assemblage by means of a Spearman correlation coefficient (Table S2 and Figure S10). All correlations were significant (p < .01; Table S2). With the exception of four assemblages (1, 4, 9, and 11) that had a correlation lower than 0.5, all assemblages had correlations above 0.5. This result confirmed that the environmental signature of phytoplankton and zooplankton was similar for all assemblages, although the strength of these similarities varied from one assemblage to another.

| Potential limitations
Our study has potential limitations related to the nature of the CPR data and to our methodology. First, the CPR is not a perfect sampling mechanism. It underestimates some components of the plankton, for example, large plankton such as fish larvae and delicate gelatinous plankton. Due to the mesh size (~270 µm) of the CPR silk (Jonas et al., 2004), some organisms are only semiquantitatively recorded and abundance of small species is probably underestimated when compared to other water sampling methods. Despite this bias, the proportion of the population that is retained by the CPR silk reflects the major changes in abundance, distribution, and specific composition; that is, the percentage retention is roughly constant within each species even with very small-celled species (Edwards et al., 2006).
Second, taxonomic identification of plankton has evolved since 1948, which might have affected our results. Nevertheless, on the 304 species/taxa considered in our study, only 27 underwent a change in taxonomic resolution (Batten et al., 2003), which represented less than 9%. Moreover, the taxonomic names recorded in the database have mostly remained unchanged throughout years . Some zooplankton species/taxa can also be recorded twice by the CPR taxonomists during the identification laboratory process. The only species recorded in eyecount and traverse analyses in our study was Centropages chierchiae, and the two categories were grouped in the same assemblage (i.e., Assemblage 4; Figure 2t). When the same taxon is recorded twice (e.g., Chaetognatha), it is often because it includes different species or developmental stages (e.g., Euphausiacea nauplii and calyptopis) . Because this information was also important for taxonomists and scientists working on the CPR survey, we chose to keep them in our analyses. The choice to keep all categories had no effect on the species assemblages because of the use of the Ward algorithm that minimizes the intragroup variance (Legendre & Legendre, 1998).
Third, sampling by the CPR survey is also restricted to the surface water (~6.5 m in depth; Hays & Warner, 1993), which might affect our perception of how ecosystems are organized. However, a recent study showed that the seasonal and diel patterns in the abundance of Calanus finmarchicus at surface were positively correlated with patterns of abundance observed at 100 m (Hélaouët et al., 2016). Despite these limitations, it has been shown that the CPR gives a correct picture of both temporal (i.e., seasonal and diel scale) and spatial (i.e., regional to basin scale) changes in plankton (Batten et al., 2003;Richardson et al., 2006).
It is widely recognized that all plankton sampling systems have their own limitations and nuances and that all underestimate abundance to some degree. However, it is important to note that the CPR survey is the only scientific monitoring program of that kind in the world, with no equivalent existing program. It covers an important timescale from 1948 to present, still active, with a large spatial scale covering the whole North Atlantic Ocean and its adjacent seas such as the Channel, the Celtic Sea, and the North Sea. The CPR is now "the most extensive long-term survey of marine organisms in the world" (Reid et al., 2003).
Potential limitations also arose because of the methodological choices we made. First, we assumed that the spatial variance was more pronounced than temporal (i.e., year-to-year and decadal) variance; we did not consider year-to-year to decadal variability.
This assumption was needed to cover as fully as possible the North Atlantic Ocean and its adjacent seas. The same assumption was made in  to propose a new partition of the North Atlantic Ocean. The effect of temporal variability was to inflate local spatial variance. To reduce this effect, we spatially smoothed the data (see Analysis 2 in Section 2) prior to conducting other analyses such as the cluster analysis and the identification of environmental signatures. This assumption had no effect on the identification of species assemblages because the environmental signature of assemblages was based upon their spatial and 2-monthly aggregation for the period 1948-2016.
Environmental signatures should be stable at a decadal scale because of niche conservatism (Crisp et al., 2009 Figure 2l). It was also important to consider different developmental stages (e.g., Euphausiacea nauplii and calyptopis) because they are key for ecosystem trophodynamics (Kirby et al., 2008). For example, most eggs and larvae were clustered in Assemblage 19, which was composed of species or taxa mainly abundant south of the Polar Front (Figure 2d). Although the analysis that leads to the identification of species assemblages was not affected by merging different types of entities, this was not so for the maps of biodiversity we present. Our maps therefore display taxonomic richness and not species richness. However, merging different taxonomic entities inside the same analysis did not alter our perception of the spatial difference in biodiversity among regions because patterns observed in this study are close to those observed in other studies that focussed only on calanoid copepods (Beaugrand et al., 2000), (Beaugrand et al., 2010), and phytoplankton (Righetti et al., 2019).
Third to produce abundance maps, we spatially regularized and smoothed the CPR data for each two-month period. Our procedure gave estimations of gridded abundance (e.g., Calanus finmarchicus) similar to that obtained from kriging (Planque et al., 1997), the inverse distance method (Beaugrand et al., 2002), or spatial regularization (Helaouët & Beaugrand, 2007). In addition, our biodiversity maps gave similar biodiversity patterns to those originated from kriging, the inverse square distance method, and spatial regularization using a first-order jackknife procedure (Beaugrand, 1999;Beaugrand et al., 2000Beaugrand et al., , 2010.

| Factors contributing to the large-scale pelagic biodiversity patterns
We provide a map of plankton taxonomic richness based on both phytoplankton and zooplankton at the North Atlantic basin scale.
We caution that phytoplankton and zooplankton biodiversity maps were based on taxa identified at the species (72.5% for phytoplankton and 47.1% for zooplankton), genus (26.2% and 23.9%), and higher taxonomic (1.3% and 29%) levels. Similar maps have been shown in , but the maps we provide here have an improved spatial resolution (2° × 2° vs. 0.5° × 0.5°) due to the implementation of the smoothing algorithm into our procedure of spatial regularization (see Analysis 2 in Section 2).
The examination of the biodiversity maps ( Figure 1) has revealed two gradients in taxonomic richness in the North Atlantic: (a) a meridional gradient (south to north) corresponding to the latitudinal biodiversity gradient and (b) a zonal (west to east) gradient. Thus, we observed the greatest biodiversity in the southeastern part of our studied zone and the lowest biodiversity in the Subarctic Gyre.
Higher taxonomic richness observed south of our studied area coincided with warmer sea surface temperatures and to a lesser extent with oceanic circulation, for example, the Gulf-Stream and its extension the North Atlantic Current. Some studies have investigated the relationships between temperature and plankton biodiversity using three taxonomic groups and found clear nonlinear relationships between mean and annual variability in temperature and biodiversity (Beaugrand et al., 2010). We have also seen a strong influence of the OPF (identified here by the 10°C isotherm) on biodiversity. These results tend to confirm the influence of the 9-10°C isotherm on plankton as revealed by Beaugrand et al. (2008) using three trophic levels: phytoplankton, zooplankton, and fish. These results also reflect the strong biodiversity difference between the Polar and Westerlies-Wind biomes sensu Longhurst (1998).
Because biodiversity is low in the subarctic gyre and higher in the south of the OPF, our study suggests that temperature may be an important factor controlling plankton biodiversity in the North Atlantic. During the eighties, Colebrook (Colebrook, 1982(Colebrook, , 1985 already suggested a relationship between plankton abundance and temperature. More recently, temperature has been identified as a key driver for both phytoplankton biodiversity and zooplankton biodiversity (Beaugrand et al., 2010). Phytoplankton biodiversity is three times higher in the tropics than in the higher latitudes (Righetti et al., 2019). Such a pattern is related to the Latitudinal Biodiversity Gradient (LBG). The LBG has also been observed for zooplankton (Rombouts et al., 2009).
Other secondary factors may also influence biodiversity patterns at a regional scale. Investigating foraminifera biodiversity in the North Atlantic, Ruddiman (1969) stressed that the biodiversity gradient was virtually erased by the strength of the diverse subtropical North Atlantic gyre. To explain phytoplankton biodiversity patterns in the North Atlantic, Righetti et al. (2019) proposed that phytoplankton biodiversity was influenced by a great species turnover resulting from high seasonal variability in wind stress, turbulence, and light limitation. Margalef (1978) highlighted that there are two important parameters in phytoplankton biology: turbulence that controls sedimentation rates and variance in current velocity that affects β diversity, that is, the differences between local community.
An increase in β diversity may also explain why there is a very high taxonomic richness over the Celtic Sea and the Bay of Biscay at the boundary between the continental shelf and the open ocean. The strong spatial variability in the bathymetry in these areas enables the coexistence of oceanic and neritic (meroplankton and holoplankton) species but also pseudo-oceanic species (i.e., species that occur in the ocean and over the continental shelf but are mainly abundant along the shelf-edge).
Oceanic circulation has a strong regional or local influence on biodiversity. Many authors have provided evidence that regional biodiversity can be highly influenced by surface currents (Beaugrand et al., 2001;Longhurst, 1998;Ruddiman, 1969;Van der Spoel, 1994).
The influence of oceanic circulation is crucial at the Atlantic Basin scale, and the zonal difference in North Atlantic biodiversity is clearly explained by the warmer North Atlantic Current that flows northwards in the eastern part of the North Atlantic (Beaugrand et al., 2002). More locally, the Gulf-Stream and its northward extension the North Atlantic Current bring more warm-water species polewards. The OPF (Dietrich, 1964) has a major influence by separating the low biodiversity of the Polar biome (and the subarctic gyre) from regions of higher biodiversity in the Westerlies-Wind biome (sensu Longhurst, 1998). This influence is stronger on warm-water species than on cold-water species. The European shelf-edge current has also a major influence on the biodiversity of the Bay of Biscay and the western regions of the British Isles in modulating local upwelling and warm-water advection northwards . The Labrador Current is characterized by poor biodiversity and the occurrence of a few species such as Calanus hyperboreus and C. glacialis. In the North Sea, the Flamborough Front has a strong influence on biodiversity by separating stratified waters in the north (i.e., lower taxonomic richness) from mixed waters to the south (i.e., higher taxonomic richness). In the northeastern part of Georges Bank (shelf-edge at the southeastern part of Newfoundland), Flemish Cape is characterized by a lower taxonomic richness (i.e., ~70 species/taxa; Figure 1a) and the Northwest Corner (Worthington, 1976), which is located to the northeastern part of Newfoundland, has a higher taxonomic richness (i.e. ~100 species/taxa). North of the British Isles, the Faroe Current, and the associated Iceland Faroe Front (Read & Pollard, 1992) limit the spatial location of annual maximum taxonomic richness ( Figure 1a). Southwest of the British Isles, maximum taxonomic richness recorded in the Southwest European Basin (west of the Bay of Biscay) may also originate from hydrographic processes. The area closely corresponds with the northward spreading of Intermediate Mediterranean Water (Käse & Zenk, 1996). This water mass has a maximum influence at ~1,000 m depth, but its range extends from 600 to 2,500 m (Käse & Zenk, 1996). Therefore, the influence of hydrographical features resulting from oceanic circulation and in some cases topographic features exerts a strong secondary influence on plankton biodiversity in the North Atlantic by controlling temperature or by controlling more directly plankton advection. Therefore, plankton biodiversity in the North Atlantic Ocean is mainly driven by temperature, oceanic circulation, and bathymetry, which have a more local/regional influence.

| Decomposition of plankton biodiversity into assemblages: relationships with North Atlantic ecological units
In a similar study on calanoid copepods, Beaugrand et al. (2002) applied the IndVal method (Dufrêne & Legendre, 1997)  the Polar, the Westerlies, and the Coastal Boundary Zone Biomes (Longhurst, 1998  Although many assemblages were found in a specific biome (sensu Longhurst, 1998), none were characteristic of a specific ecological unit sensu . In contrast to the terrestrial realm, lack of strong physical barriers associated with changing hydroclimatic conditions may partly explain this observation ( Van der Spoel, 1994). Spatial coenoclines (i.e., gradients of biocoenoses or communities) we observed in our study (Figures 2-4) can be compared with annual succession observed in ecosystems during the year (Romagnan et al., 2015). These biological gradients in space and time are a perfect illustration of the ecological principle of impermanence resulting from the constant biological adjustments that originate from niche-environment interaction (Beaugrand, Conversi, et al., 2019). Sharp or gradual environmental gradients interact with the niche of each species within a multidimensional space to generate a variety of biogeographic patterns.
We do not think that the overlapping spatial distribution of some assemblages with the partition proposed by  is due to differences in spatial resolution. Partitioning the pelagic ocean is difficult because of the absence of geographical barriers and because any partition changes over time from daily to decadal scales (Reygondeau et al., 2013). In addition, because of the principle of competitive exclusion of Gause (1934), we should expect all spatial ranges and phenologies to be unique (Caracciolo et al., 2021). The necessary corollary of this principle is that any synthetic partition should not be expected to work for all species or species assemblages. Temporal dynamics and plankton dispersal or expatriation are also mechanisms adding further complexity.
We found no obvious relationships between the spatial distribution of the species assemblages, seasonal patterns, and taxonomic composition in most cases (Figure 2 vs. Figure S5). We also found no obvious relationships between assemblages and taxonomic resolution ( Figure S6). However, assemblages occurring in the same eco- The environmental signature of a group was characterized by means of an environmental chromatogram. This procedure has enabled a rapid display of the 9-dimensional environmental signature of the 24 assemblages into a two-dimensional space (Figures 5 and S10).
The environmental chromatogram separated plankton assemblages according to their degree of environmental tolerance (i.e., degree of euryoecy) and optima (e.g., degree of thermophily). We suggest that the method could also characterize the (multidimensional) ecological niche (sensu Hutchinson, 1957) of a species where both the amplitude and the optimum of an environmental factor would be represented as a function of the number of ecological dimensions, with the abundance (instead of a percentage of aggregation) as a third dimension (i.e., the color of the contour plot). The environmental chromatogram is a new method that enables one to display all multidimensional niches composing an assemblage (there are as many niches as species). Such chromatograms allow a rapid display of the environmental signature of an assemblage. Furthermore, the visual comparison of the chromatograms can immediately reveal the degree of overlapping of two environmental signatures and explain why the spatial distribution of two assemblages overlaps. For example, Assemblage 8 ( Figure 2o) has a spatial distribution that is included inside Assemblage 5 (Figure 2s), a result that is expected from the environmental chromatogram of the two assemblages (Figure 5o,s).
To our knowledge, there is no graphical method enabling a clear representation of a multidimensional ecological niche into more than two dimensions. A radar plot can be used but results are less easy to interpret, especially when a large number of species are considered (Reygondeau et al., 2012). Multivariate analyses have been sometimes used (Helaouët & Beaugrand, 2007), but because each component is a linear combination of different ecological factors, the resulting multidimensional niche is more difficult to interpret.
It has been shown that cold-water species (e.g., Calanus finmarchicus, Assemblage 12) have a larger niche breadth than their warmer-water counterparts (e.g., Calanus helgolandicus, Assemblage 19) (Sunday et al., 2011). Indeed, in our 24 assemblages, we noticed that few were specific of northern ecological units (e.g., Assemblages 10 and 6 in Figure 4n,r). Among them, many had a large spatial range covering sometimes the whole North Atlantic Ocean and its adjacent seas. We suggest that the spatial range of these assemblages illustrates well Rapoport's effect, which postulates that high-latitude species have a larger geographical range than mid-latitude species (Stevens, 1989).

| CON CLUS IONS
Spatial patterns in plankton biodiversity are the result of environmental factors acting at both global and local scales (Beaugrand et al., 2000). We have decomposed CPR-sampled plankton biodiversity into 24 taxonomic assemblages using both phytoplankton and zooplankton. These assemblages are interesting because they characterize specific hydroclimatic conditions and so they can be used as biological indicators, either of a substrate-biotope (water mass) component or of a stable-biotope (key area) component (Beaugrand et al., 2002). Their associated specific environmental signatures have allowed us to better characterize the ecosystems of the North Atlantic Ocean and its adjacent seas and to improve our understanding of the arrangement of plankton biodiversity.
Although some assemblages were characteristic of the three biomes observed in our studied areas, none of them were observed specifically in the ecological units identified by . This important result can in part be explained by the regional complexity of the northern part of the North Atlantic Ocean, which exhibits a pronounced hydrodynamic variability with interwoven substrate and stable-biotope com- ponents. An alternative hypothesis is that it could be a universal feature. All species having a unique ecological niche after the principle of competitive exclusion of Gause (1934), the corollary is that they are expected to exhibit distinct spatial distributions. We have clearly shown that each taxonomic assemblage had a unique environmental signature and that therefore Gause's corollary can be extended at the community level. In the same way that coenoclines are observed during annual succession and that community reorganization takes place all the time from year-to-year to multidecadal scales, coenoclines may occur everywhere in the ocean.
Pelagic ecosystems are therefore likely to be more complex than previously envisioned.

S U PP O RTI N G I N FO R M ATI O N
Additional supporting information may be found online in the Supporting Information section.
How to cite this article: Kléparski L, Beaugrand G, Edwards M.