Low complexity model to study scale dependence of phytoplankton dynamics in the tropical Pacific
KKeep it simple: A low-complexity model to study scale dependence of phytoplanktondynamics in the tropical Pacific
Jozef Sk´akala , and Paolo Lazzari Plymouth Marine Laboratory, Prospect Place, PL1 3DH Plymouth, United Kingdom, National Centre for Earth Observation, PL1 3DH Plymouth, United Kingdom, and Instituto Nazionale di Oceanografia e di Geofisica Sperimentale, Trieste, 34151, Italy. (Dated: September 9, 2020)We demonstrate that a simple model based on Fisher-Kolmogorov-Petrovski-Piskunov (FKPP)equation forced by realistic surface velocities and nutrients is skilled in reproducing the distributionsof the surface phytoplankton chlorophyll in the tropical Pacific. We use the low-complexity FKPPmodel to investigate the scale-relationships in the impact of different drivers (diffusion, mean andeddy advection, primary productivity) on the phytoplankton chlorophyll concentrations. We findthat in the 1 / ◦ model, advection has a substantial impact on the rate of primary productivity,whilst the diffusion term has a fairly negligible impact. Diffusion has an impact on the phytoplank-ton variability, with the impact being scale-propagated and amplified by the larger scale surfacecurrents. We investigate the impact of a surface nutrient decline and some changes to mesoscaleeddy kinetic energy (climate change projections) on the surface phytoplankton concentrations. TheFKPP model suggests that unless mesoscale eddies radically change, phytoplankton chlorophyllscales sub-linearly with the nutrients, and it is relatively stable with respect to the nutrient con-centrations. Furthermore we explore how a white multiplicative Gaussian noise introduced into theFKPP model on its resolution scale propagates across spatial scales through the non-linear modeldynamics under different sets of phytoplankton drivers. The unifying message of this work is thatthe low complexity (e.g. FKPP) models can be successfully used to realistically model some specificaspects of marine ecosystem dynamics and by using those models one can explore many questionsthat would be beyond computational affordability of the higher-complexity ecosystem models. I. INTRODUCTION
There is no effective scale in ecology [1]. New struc-tures and processes appear with every new scale down tothe fundamental spatial scale of molecular biology, whichis far beyond the reach of our ecosystem models. Ourmodels provide merely an effective description for theecosystem dynamics (e.g.[2, 3]), so that any impact of themodel sub-grid processes is either parametrized, or can berepresented by a stochastic noise. To be able to correctlydescribe the impact of sub-grid processes on the modelgrid scale, it is beneficial to have some understandingof how the ecosystem equations, or ecosystem variablesevolve with the spatial, or temporal scale ( (cid:96) ). We will callan approach that provides such understanding a ”scalinganalysis”. Scaling analysis has been largely advancedwithin the framework of renormalization group (e.g.[4–6]) with many fundamental applications across particlephysics, solid state physics and complex dynamical sys-tems (e.g.[7–10]). Interestingly, renormalization groupmethods were also applied to Navier-Stokes equations(e.g.[11]) and reaction-diffusion models (for a recent com-prehensive review see [12]). The renormalization groupturned out to be particularly well suited to describe scale-invariant properties of the examined system and has beenwidely applied to study critical phenomena and univer-sality (e.g.[13–15]).With the increased ecosystem model resolution, aswell as the increased model complexity, more phenom-ena are included into the ecosystem model. However,a model does not necessarily provide good understand- ing for all the phenomena it represents. Indeed, un-derstanding phenomena often requires a specific scale:for example to understand oceanic gyres it is desirableto look at a long-time, spatially large-scale oceanic andatmospheric behavior. Although model that capturesocean mesoscale, or sub-mesoscale dynamics representsalso ocean gyres, their behavior remains hidden behindthe dominant short-time small-spatial scale eddy signal.Similarly if we managed to run a model on a molecu-lar scale, the eddy behavior would remain hidden behindthe thermal fluctuations of the molecules and atoms (andthe same type of situation happens if we switched fur-ther from the atomic scales to the scales of the currentelementary particle theory). Here lies another benefit ofthe scaling analysis: it provides us with a natural toolto understand diverse phenomena with a wide range ofcharacteristic spatio-temporal scales (e.g. turbulence, ge-ologic processes, climate, financial markets, [16–20]), as itsimultaneously compares processes across different scales.Apart of improving our ecosystem models and under-standing processes, there is also a third potential bene-fit of scaling analysis: Fine resolution models representbroad range of ecosystem phenomena, but they are com-putationally expensive. Understanding dynamics acrossrange of scales might also optimize the performance ofhigh resolution models by converting them into multi-scale models (e.g.[21]). This means each “separate” partof the model dynamics could be represented at the max-imum scale where it occurs, eventually leading to sub-stantial reduction in the model computational cost.The main point of this work is to develop a schematicmulti-scale understanding for some essential aspects a r X i v : . [ q - b i o . P E ] S e p of ecosystem dynamics. This provides different view-point from the standard ecosystem modelling, where theecosystem model is understood at some specific (fine-resolution) scale, whilst the larger-scale phenomena al-ways “emerge” from the model small-scale complex dy-namics (e.g.[22]). We will show that to get a sufficientlyrealistic representation of the primary productivity in thetropical Pacific, the high complexity model can be forspecific purposes bypassed by a simplified “toy” reaction-diffusion model based on Fisher-Kolmogorov-Petrovski-Piskunov (FKPP) equation [23, 24] (we will further callsuch model “the FKPP model”).We forced the FKPP model by the realistic Coperni-cus Marine Environment Monitoring Service (CMEMS)reanalyses for the surface currents and nutrients, andthe model remarkably successfully captures the dynam-ics of chlorophyll also provided by a CMEMS reanaly-sis. This is a surprising result: Although simple (oftenone-dimensional) models based on the FKPP equationwere often used to address conceptual problems, such ashow species survival depends on diffusion rate, advec-tion rate, or on the characteristic patch size occupied bythe specie (e.g.[25–30] and for an overview see [31]), onewould assume that a sufficiently realistic marine modelmust be much more complex than the FKPP-like mod-els. The relationship between the ecosystem model skilland the model complexity is non-trivial ([32]), howevercertain minimum amount of model complexity is alwaysassumed; the real world marine biogeochemistry is ad-dressed either by the medium-complexity models ([33–35]), or by the high complexity models ([36, 37]) that have often tens of state variables and more than hundredparameters. Such assumptions are without any doubtfounded, but this paper shows that for some suitablychosen problems of high scientific interest, even the sim-plest model based on FKPP equation is capable to pro-duce surprisingly good approximation to the selected realworld ecosystem data. One then gets the best of bothworlds: the advantage of the FKPP model is that itis cheap to run, and it depends only on three free pa-rameters, whose impact on chlorophyll distributions canbe easily understood, modified and studied across widerange of spatial and temporal scales. In the same timethe FKPP model appears (within its constrained frameof reference) to be sufficiently realistic for the results ofsuch analyses to be taken seriously. II. METHODSA. The FKPP model
The growth of biomass starts with the photosynthe-sis in the autotrophic species, and for marine ecosystemsthese are the diverse species of phytoplankton. The fre-quently used proxy quantity for phytoplankton biomassis chlorophyll a , with a clear advantage of large volume ofocean-color derived observations available for the oceansurface concentrations of chlorophyll (e.g.[38]). In thiswork we focus on the chlorophyll dynamics, modelled bya FKPP equation expressed as ∂ρ ( t, (cid:126)x ) ∂t = − (cid:126)u ( t, (cid:126)x ) · ∇ ρ ( t, (cid:126)x ) + κ · ∇ ρ ( t, (cid:126)x ) + P · N ( t, (cid:126)x ) ρ ( t, (cid:126)x ) − D · ρ ( t, (cid:126)x ) , (1)where ρ ( t, (cid:126)x ) represents chlorophyll concentrations, N ( t, (cid:126)x ) nutrients, (cid:126)u ( t, (cid:126)x ) is the current velocity, κ isthe diffusion parameter, P the net primary productiv-ity (growth) rate and D is the damping (mortality) rate.The diffusion parameter ( κ ) describes the integrated ef-fect of sub-grid eddy mixing and determines the rate ofsmall-scale chlorophyll smoothing. The damping rate D integrates phytoplankton loss due to the limitation inresources, mortality, respiration and grazing by highertrophic-level species. D also impacts the degree to whichchlorophyll and nutrients are correlated: if substantialphytoplankton concentrations get advected into the low-nutrient areas they die off quickly if the damping rate D is high. Conversely, in the nutrient-rich areas the highrate of damping D will not allow phytoplankton to growabove certain threshold in concentrations, constrainingthe correlation between ρ and N . Finally, the growthparameter ( P ) describes the rate of photosynthesis. The P parameter determines (for a fixed D ) the average levelsof chlorophyll ( (cid:104) ρ (cid:105) ) on the domain. For the purpose of this study the FKPP modelis constrained to a two-dimensional horizontal plane,representing the ocean surface in the Pacific cen-tral tropical region (155E - 110W, 30S - 30N, seeFig.1). The selected region spans most of tropicalPacific with meridional dimension ∼ ∼ / ◦ . The oceansurface current velocity ( (cid:126)u ( t, (cid:126)x ), see Eq.1) and nutri-ents ( N ( t, (cid:126)x )) were provided for the FKPP model ex-ternally; the ocean surface current velocity was takenfrom the 2017-2018 daily resolution CMEMS reanal-ysis (GLOBAL ANALYSIS FORECAST PHY 001 024, http://marine. - copernicus.eu ), which is based on as-similation of satellite sea surface temperature, sealevel anomaly, as well as in situ temperature andsalinity into 1 / ◦ ORCA012 model configuration ofthe Nucleus for European Modelling of the Ocean(NEMO, v3.1, [39], for details on the reanalysis,see http://marine.copernicus.eu/ documents/QUID/ - FIG. 1. The CMEMS horizontal 2017-2018 mean surface current velocity. Left hand panel shows the zonal mean currentvelocity component ( U ), the middle panel the meridional mean current velocity component ( V ) and the right hand panel showsthe mean surface eddy speed. CMEMS-GLO-QUID-001-024.pdf ). To represent thesurface currents (cid:126)u on the 1 / ◦ FKPP model gridwe upscaled the CMEMS data from their original(1 / ◦ ) scale of resolution. In the Fig.1 we showthe 2017-2018 mean values of the surface current ve-locity vector components and also the mean surfaceeddy speed. The nutrients N ( t, (cid:126)x ) have been es-timated as a sum of nitrate and phosphate usingthe outputs of 2017-2018 CMEMS hindcast based on1 / ◦ resolution NEMO coupled with the biogeochemi-cal model Pelagic Interactions Scheme for Carbon andEcosystem Studies (PISCES, GLOBAL REANALYSIS -BIO 001 029, http://marine.copernicus.eu/ ). No assim-ilation was used in the biogeochemical run. Phosphateand nitrate were the only nutrient data available withthe desired resolution, however taking the sum of nitrateand phosphate is only one of multiple seemingly equiv-alent choices of how to represent the nutrients. Sincephosphate and nitrate concentrations are shaped by sim-ilar drivers, the two nutrients have been found to bereasonably highly Pearson correlated ( R = 0 . ±
30% vari-ance) around the 0 . mg/m mean. The tests (not shownhere) have demonstrated that on the timescale of (cid:38) emissions to the atmosphere [40, 41]. Largeparts of the region, such as the eastern equatorial Pa-cific, are characteristic of high-nutrient (nitrate and phos-phate) concentrations due to the equatorial upwelling,but comparably low chlorophyll concentrations (the so-called “high chlorophyll low-nutrient regions” [42–44]).The comparably low primary productivity around theequator is often understood to be caused by the limitedresources of iron [41, 45, 46], although the elevated levelsof grazing also may play a role [43, 47]. In the olig-otrophic regions further away from the equator (higherlatitude than 10 ◦ ) the conditions are very different andphytoplankton is mostly nutrient-limited [41].There are several reasons why tropical Pacific is an optimal choice for our experiment:a) It is an open ocean region with little impact ofbathymetry on the ecosystem dynamics.b) The 1-st baroclinic Rossby radius is in the tropicalPacific on the scale of 100-s of km [48] and the eddyscales can get close to ∼
500 km [49] (see also Fig.2), sothe 1 / ◦ model resolution allows us to see a wide rangeof interesting scales for the chlorophyll dynamics [50].c) Phytoplankton dynamics can be fairly complex andhave strong seasonal signatures (e.g. spring blooms) dueto seasonal variability in the upper ocean mixing and sun-light that drives photosynthesis. However, the seasonalcycles in the tropical Pacific are weak, and phytoplank-ton production is primarily regulated by the availablenutrients (Fig.3), with additional impact of advection bythe surface currents. Fig.3 shows that in the tropicalPacific nutrients and chlorophyll are strongly correlated(Pearson correlation, R=0.77), with nutrient spatial ge-ography playing an essential role for the phytoplanktondistributions. In such case one can represent the bio-logical productivity as a simple function of the nutrientconcentrations, as is done in the FKPP model (Eq.1).d) The FKPP model used in this study is a single-equation model with externally supplied nutrients. Onecould argue that the FKPP model needs adding a similardynamical equation for the nutrients, as has been donemany times in the literature (e.g.[31]). In the tropicalPacific the nutrient sinks and sources depend largely onthe vertical mixing (e.g. equatorial upwelling) and suffi-ciently near the coastline could reflect other forcing fields,such as the river discharge. In such case it becomes dif-ficult to implement a two-equation nutrient-chlorophyllmodel without substantially increasing the model com-plexity. However, we argue that for the purpose of thisstudy the single equation model (Eq.1) is in the tropicalPacific a reasonable approximation to the phytoplank-ton chlorophyll dynamics. There are two issues here thatneed to be raised: Firstly, within this study we will ex-plore the impact of the modified CMEMS data for thesurface currents ( (cid:126)u ) and the diffusion parameter ( κ ) onthe phytoplankton chlorophyll concentrations ( ρ ). Thechanged surface advection and diffusion can potentiallychange the nutrient concentrations ( N ) relative to theirexternally supplied CMEMS values. Secondly, the phy-toplankton concentrations change as a function of themodified advection and the changes to the nutrient up-take by the changed phytoplankton ( ρ ) could be anothersource that modifies the nutrients relative to their sup-plied CMEMS values. There are, however, two argu-ments why we could reasonably neglect those changes tothe supplied nutrients and still use the CMEMS product:i) The nutrient distributions are much more geographi-cally stable than the chlorophyll (Fig.3 and Fig.4), bywhich we mean that the nutrient anomalies are relativelysmall when compared to the nutrient spatial geographyestimated from the 2017-2018 mean values (Fig.3). Thenutrient geographic sinks and sources, which largely cor-respond to the upwelling and downwelling zones, thenconsequently play a key role in the representation ofthe nutrient distributions, with other drivers (such aseddy mixing, or time-fluctuations in the uptake by phy-toplankton) playing mostly a secondary role. Moreover,this study will explicitly demonstrate that it makes lit-tle difference to the simulated chlorophyll, whether weforce the FKPP model with a time-changing, or 2017-2018 time-averaged nutrient distributions. ii) A substan-tial change to the CMEMS phytoplankton chlorophyllconcentrations might indeed introduce some changes tothe CMEMS nutrients through the uptake. However, theconcern of this study are not the changes to the nutrients,but the impact of those nutrient changes on the phy-toplankton distributions. Although the single-equation FKPP model does not represent the changes to the nu-trients, the quadratic damping term in the FKPP model(the D -term in Eq.1) effectively integrates into the phyto-plankton dynamics the impact of the resource limitationdue to nutrient uptake.Although iron is an important limiting factor in someareas of tropical Pacific [43], the daily products for irondistributions were unavailable and could not be used aspart of the FKPP model forcing. The limitation by ironwas, similarly to the nutrient uptake, included into theFKPP model only implicitly as part of the quadraticdamping term. The FKPP model assumes that anydamping effect included in the quadratic term is pro-portional to the chlorophyll concentration. This can beeasily justified for the rate of phytoplankton mortality,nutrient limitation, or for phytoplankton grazers (theirdensity is expected to be proportional to phytoplanktondensity), and to some degree it can be justified also forthe iron limitation, as the chlorophyll concentrations arehighest in the iron-limiting equatorial upwelling region(Fig.3). However, we acknowledge that representing theiron limitation only implicitly is definitely a shortcomingof the FKPP model. B. Some analytical results about the FKPP modelsolutions
In this section we briefly outline some analytical prop-erties of the FKPP model, which will be later used tobetter understand the results of the study. Since advec-tion and diffusion do not change the spatially averagedchlorophyll concentration (cid:104) ρ (cid:105) , the FKPP model (Eq.1)has a simple stochastic steady state ( ∂ (cid:104) ρ (cid:105) /∂t = 0) solu-tion: (cid:104) N (cid:48) ( (cid:126)x ) ρ (cid:48) ( (cid:126)x ) (cid:105) = (cid:104) ρ (cid:48) ( (cid:126)x ) (cid:105) , (2)where N (cid:48) = P.N and ρ (cid:48) = D.ρ . (For a region withboundaries, we assume in Eq.2 also constant Dirichletboundary conditions.) Applying (cid:104) N (cid:48) .ρ (cid:48) (cid:105) = (cid:104) N (cid:48) (cid:105) . (cid:104) ρ (cid:48) (cid:105) + Cov ( N (cid:48) , ρ (cid:48) ) , where “Cov” is covariance, or consequently (cid:104) ρ (cid:48) (cid:105) = (cid:104) ρ (cid:48) (cid:105) + V ar ( ρ (cid:48) ) , where “Var” stands for variance, one can transform thestochastic steady state solution (Eq.2) into a quadradicpolynomial equation for (cid:104) ρ (cid:105) : (cid:104) ρ (cid:48) (cid:105) − (cid:104) N (cid:48) (cid:105)(cid:104) ρ (cid:48) (cid:105) + V ar ( ρ (cid:48) ) − Cov ( N (cid:48) , ρ (cid:48) ) = 0 . (3)By solving Eq.3 we obtain a relationship between theaverage chlorophyll and the nutrient concentrations as: (cid:104) ρ (cid:48) (cid:105) = (cid:104) N (cid:48) (cid:105) ± (cid:118)(cid:117)(cid:117)(cid:116)(cid:18) (cid:104) N (cid:48) (cid:105) (cid:19) + V ar ( ρ (cid:48) ) (cid:32)(cid:115) V ar ( N (cid:48) ) V ar ( ρ (cid:48) ) · R ( ρ (cid:48) , N (cid:48) ) − (cid:33) , (4)with R being the Pearson correlation coefficient. A sim-pler relationship between (cid:104) ρ (cid:48) (cid:105) and (cid:104) N (cid:48) (cid:105) can be derived,if we assume that the standard deviation of both ρ (cid:48) and N (cid:48) is directly proportional to their mean values: (cid:112) V ar ( ρ (cid:48) ) = c ρ (cid:48) . (cid:104) ρ (cid:48) (cid:105) and (cid:112) V ar ( N (cid:48) ) = c N (cid:48) . (cid:104) N (cid:48) (cid:105) . ThenEq.3 leads directly to a linear relationship: (cid:104) ρ (cid:48) (cid:105) = 1 + R ( ρ (cid:48) , N (cid:48) ) .c ρ (cid:48) .c N (cid:48) c ρ (cid:48) · (cid:104) N (cid:48) (cid:105) . (5)If we lower advection, chlorophyll becomes highly cor-related with nutrients and the Pearson correlation R ( ρ (cid:48) , N (cid:48) ) in Eq.5 approaches R ( ρ (cid:48) , N (cid:48) ) = 1, whereaswith the high levels of mixing ρ (cid:48) and N (cid:48) decorrelate( R ( ρ (cid:48) , N (cid:48) ) → (cid:104) N (cid:48) (cid:105) , lowers the meanchlorophyll concentrations.If there is neither advection, nor diffusion ( (cid:126)u = κ = 0),and N (cid:48) does not depend on time, Eq.1 has the followingexact solutions: ρ (cid:48) ( t, (cid:126)x ) = N (cid:48) ( (cid:126)x )1 + ρ o · exp {− N (cid:48) ( (cid:126)x ) · t } , (6)which converge for ρ (cid:48) > ρ (cid:48) ( (cid:126)x ) = N (cid:48) ( (cid:126)x ) , (7)whilst for ρ (cid:48) < −∞ . Thesolutions from Eq.7 approach the steady state attractor(Eq.6) as: ∆( t, (cid:126)x ) (cid:39) exp {− N (cid:48) ( (cid:126)x ) · t } , (8)where ∆ is the distance measured on the real line be-tween the approaching solution and the attractor. Eq.8means that the higher nutrient concentration, the fasterthe chlorophyll distributions converge to the steady statesolution from Eq.7.For the exact steady state solution (Eq.7) chlorophyllis maximally correlated with nutrients, R ( ρ (cid:48) , N (cid:48) ) = 1. Asimple consistency check shows that for R ( ρ (cid:48) , N (cid:48) ) = 1,Eq.4 is solved by the averaged form of the linear rela-tionship in Eq.7: (cid:104) ρ (cid:48) (cid:105) = (cid:104) N (cid:48) (cid:105) . (9)together with V ar ( ρ (cid:48) ) = V ar ( N (cid:48) ) . (10)Eq.7 and Eq.9-10 imply that if the first two statisticalmoments of ρ (cid:48) and N (cid:48) are equal, then the steady state ρ (cid:48) and N (cid:48) are maximally correlated ( R ( ρ (cid:48) , N (cid:48) ) = 1). Withthe increased advection N (cid:48) and ρ (cid:48) decorrelate, and in thelimit of R ( ρ (cid:48) , N (cid:48) ) = 0, one obtains (cid:104) ρ (cid:48) (cid:105) = (cid:104) N (cid:48) (cid:105) ± (cid:115)(cid:18) (cid:104) N (cid:48) (cid:105) (cid:19) − V ar ( ρ (cid:48) ) . (11) C. The scaling analysis
In this work, we borrow insights from the long his-tory of the studies on turbulence and multifractals([16, 20, 51–57]), and use a simple measure for the scaledependence of the system variables (∆ (cid:96) ρ ) as:∆ (cid:96) ρ = (cid:104)| ρ ( x + (cid:96) ) − ρ ( x ) |(cid:105) . (12)Here ∆ (cid:96) ρ represents a (scale-dependent) magnitude ofspatial and temporal variability of ρ , x is the spatial, ortemporal variable (spatial vector for spatial variability,or time for temporal variability), (cid:96) is the scale of interestand the averaging in Eq.12 runs through the relevantspatial domain, or the time interval. ∆ (cid:96) ρ corresponds tothe first statistical moment of what is in the multifractalliterature often called “increments” (e.g.[57]).∆ (cid:96) ρ has the advantage of being methodologically sim-ple and has been many times proven fruitful in the lit-erature (e.g.[20, 56, 58–64]): For the scale-invariant sys-tems the scaling of ∆ (cid:96) ρ follows a power law and it hasbeen found that its power law exponent is often an im-portant indicator of the system dynamics, e.g. specificinterval for the values of the power law exponent hasbeen found for tracers passively advected by a turbulentflow [59, 64, 65]. In the recent work of [66] it has beenshown that the scaling described by Eq.12 is frequentlya piece-wise power law with the scaling transition be-tween different power laws corresponding to a transitionbetween different dynamical regimes (see also [60, 62]).The power law exponents correspond to the scaling slopeof ∆ (cid:96) ρ (we use ˜∆ (cid:96) ρ notation), which can be analysed bynormalizing the ∆ (cid:96) ρ value as˜∆ (cid:96) ρ = ∆ (cid:96) ρ/ ∆ L ρ, (13)where L is some maximum spatial, or temporal scale ofinterest [66]. ˜∆ (cid:96) ρ can be then used as a simple “probe”to test the impact of dynamical drivers (e.g. eddy andmean advection, diffusion, biological productivity) on thevariable of interest (e.g. chlorophyll) across a wide rangeof spatio-temporal scales. FIG. 2. The spatial and temporal eddy scales. The left hand side panels (A,C) show the Pearson correlation ( R , y-axis) in theeddy surface velocity as a function of spatial (A) and temporal (C) scale (x-axis), the right hand side panels (B,D) show thechlorophyll magnitude of spatial (B) and temporal (D) variability (∆ (cid:96) ρ , Eq.12, y-axis) as a function of temporal scale (D), orspatial log-scale (B, x-axis). The chlorophyll from the panels B and D was a FKPP model output with κ = N =0 and with (cid:126)u represented only by the eddy field (the eddy (cid:126)u was estimated by subtracting the 2017-2018 mean CMEMS currents from theCMEMS daily output). Since the panels B and D focus only on the scaling slope of ∆ (cid:96) ρ , the values of ∆ (cid:96) ρ are not shown.Both analyses point consistently to the maximum eddy spatial scale of 500 km and the maximum time scale of 50 days (thiscan however be much shorter than eddy life-time, as eddies move). The spatial large-scale correlation ( R ∼ .
1) that can beseen in the panel A has been found (not shown here) to be caused by a meridional cross-correlation across the equator due toseasonal variations in the currents. Similarly, the ∆ (cid:96) ρ scaling slope within the intermediate time scale between 50-150 days inthe panel D has been found (not shown here) to correspond to the seasonal variability in the currents. III. VALIDATION OF THE FKPP MODEL
An ensemble of FKPP model simulations was run un-til the optimal set of P , D , κ (Eq.1) values was deter-mined to be: κ = 300 m s − , P = 7 . − m mmol − s − , D = 1 . . − m mg − s − . The set of optimal parametervalues was chosen based on the match-ups between theFKPP model and the CMEMS data using three metricsshown in Fig.5: a) the 2-year mean spatial distributionof chlorophyll, b) the magnitude of spatial and temporalvariability ∆ (cid:96) ρ (Eq.12) across 25-2500 km and 1 day -1 year range of scales. The first metric (a, Fig.5:A,C)measures the FKPP model skill to estimate the averagechlorophyll concentrations and to represent the dominantchlorophyll patterns. Since the spatial chlorophyll pat-terns dominate over the temporal chlorophyll patterns(Fig.3:A-B) the metric entirely focuses on the chloro-phyll spatial distributions. The two remaining metrics(b, Fig.5:B,D) measure how well the FKPP model re-produces the CMEMS magnitude of chlorophyll spatialand temporal variability. The magnitude of chlorophyllspatial and temporal variability will be used to identifythe impact of drivers on the chlorophyll concentrationsacross a wide range of spatial and temporal scales. Since the impact analysis for the chlorophyll drivers relies fullyon the FKPP model, it is essential that the FKPP modelreproduces realistically the scaling of the magnitude ofchlorophyll spatial and temporal variability.All the three metrics in Fig.5 show that the FKPPmodel is skilled in representing the CMEMS chlorophylldata, i.e. the magnitude of spatial and temporal vari-ablity match on most scales within 10% and on all scaleswithin 20%, with the exception of the magnitude of tem-poral variability on the annual scale. The sudden dropin CMEMS data temporal variability on the annual scaleis due to the bi-annual periodicity in the chlorophyll dis-tributions (see Fig.3) driven by the bi-annual seasonalitypattern in the solar radiation (at the equator the seasonalpattern has bi-annual periodicity, because the seasons inthe Southern and Northern Hemisphere have identicalimpact on the equator). Since the FKPP model does notrepresent the solar cycle, it is understandable that it failsto capture the bi-annual, or annual periodicity in the ∆ (cid:96) ρ of the CMEMS data.The FKPP model parameters can be characterized bythe relative magnitude of three types of drivers: diffu-sion, advection and biological activity. The Damk¨ohlernumber Da (see [29]) gives the scale ( (cid:96) ) dependent ratio FIG. 3. The CMEMS 2017-2018 mean surface concentrations for the chlorophyll a (panel A) and nutrients (panel C). Thepanels B and D show the 2017-2018 time series for the spatial mean of surface chlorophyll (B) and nutrients (D). It is shownthat chlorophyll has a modest bi-annual periodicity (panel C), which is driven by the seasonal solar cycle (since the region ismeridionally symmetric across the equator, the solar seasonal cycle here is bi-annual).FIG. 4. The distributions expressed by the Probability Density Function (PDF) and the Cummulative Density Function (CDF)for the chlorophyll (panel A) and the nutrient (panel C) anomalies calculated relative to the 2017-2018 mean concentrations.The values on the x-axis are scaled (in %) relative to the 2017-2018 spatio-temporal mean. The plots show that the relativespread of chlorophyll distribution is substantially larger than the relative spread of nutrients. FIG. 5. The skill of the FKPP model to reproduce the CMEMS 2018 mean chlorophyll, as well as the chlorophyll a magnitudeof spatial (panel C) and temporal (panel D) variability across a range of spatial and temporal scales (∆ (cid:96) ρ calculated as anappropriate average of the 2018 daily data). The panels show that the model is skilled in reproducing both CMEMS chlorophyll a spatial distributions (panels A and C) and the magnitude of variability (panels B and D), except for the magnitude oftemporal variability around the half year-to-annual scale. This can be easily explained: the FKPP model does not include thetime variability in the solar input and hence does not reproduce adequately the bi-annual periodicity of the CMEMS data.The chlorophyll magnitude of spatial variability over 2500 km starts decreasing (panel B), since the chlorophyll distributionshave a meridional symmetry across the equator. Similarly, as mentioned before, the local minimum of the CMEMS chlorophyllmagnitude of temporal variability at the annual scale (panel D) is due to the annual cycle (annual cycle seems more pronouncedthan the bi-annual cycle). between the biological rate of the process and the ad-vection rate: Da = biological rate / advection rate = (cid:96).P. (cid:104) N (cid:105) / (cid:104)| (cid:126)u |(cid:105) . We can then easily calculate the scale (cid:96) da where biological rate ≈ advection rate as (cid:96) da ≈ km .At the scales (cid:96) << (cid:96) da advection dominates biologicalprocesses and vice versa. If we interpret “the muchsmaller” as a separation by two orders of magnitude, weconclude that advection is expected to dominate biolog-ical processes at the O(10) km scales. Similarly to theDamk¨ohler number, we can introduce P´eclet number [29] as P e = advection rate / diffusion rate = (cid:96) (cid:104)| (cid:126)u |(cid:105) /κ . Thenfor the scale (cid:96) pe where advection rate ≈ diffusion rate,we obtain (cid:96) pe ≈ m . At the scales (cid:96) >> (cid:96) pe advectiondominates over diffusion and vice versa. The (cid:96) pe scalesuggests that advection should be dominant over diffu-sion on the scales of O(100) km. The estimates usingDamk¨ohler and Peclet numbers are broadly consistentwith the results of this study, however we will show thatadvection can propagate the impact of diffusion to re-markably large scales. IV. IMPACT OF PRIMARY PRODUCTIVITYDRIVERS ACROSS DIFFERENT SCALESA. Spatial analysis
What will be the impact on the chlorophyll concentra-tions if we switch off horizontal advection or diffusion in the FKPP model? We have done multiple experimentswith: i) switched off mesoscale eddies, in which case (cid:126)u (Eq.1) was taken as mean currents only, estimated froma 2017-2018 average of the CMEMS data (see Fig.1), ii)switched off mean currents, in which case the 2017-2018means were subtracted from the CMEMS data for (cid:126)u toestimate the eddy field, iii) no advection at all ( (cid:126)u = 0).In each of these cases (i-iii) and also in the case forced
FIG. 6. The impact of nutrient concentrations on the mean chlorophyll. The values shown in the Figure are the averagesthrough the FKPP model spatial domain and the year 2018. It is shown that eddies stabilize the chlorophyll concentration:without eddies the 50% decrease in nutrients leads to almost 50% decrease in chlorophyll, whilst in the presence of CMEMSeddies the 50% decrease in nutrients lowers chlorophyll only by ∼ by CMEMS data for (cid:126)u we ran two separate simulations,with and without diffusion (diffusion was removed bysetting κ = 0). For the simulation with switched offmesoscale eddies (i), it is desirable to remove the eddysignatures also from the nutrient ( N ) data. We havecompared two simulations with the eddy advection (cid:126)u : a)one that used the CMEMS product for nutrients ( N ) andb) another simulation, which used for N the 2017-2018mean CMEMS nutrient concentrations. The two simu-lations produced very similar results for the chlorophyll(not shown here), e.g. the differences in the magnitude ofspatial variability were on all scales < n L , n S and ch N , ch S where n L , ch L are “large” nutrient and chlorophyll concentra-tions, whilst n S , ch S are “small” nutrient and chlorophyllconcentrations. Since “large” is larger than “small” wehave: ( n L − n S ) · ( ch L − ch S ) > n L ch L + n S ch S > n L ch S + n S ch L . (15)Advection (e.g. eddy mixing) brings large chlorophyllconcentrations ch L to areas with worse growth condi-tions (small nutrient concentrations n S ) and vice versa,the growth term then corresponds to the right side ofEq.15, whereas if there was no advection the growth termis described by the left side of Eq.15. This means whenthere is advection (eddy or mean) the growth term issmaller than if there is no advection. However, focus-ing purely on eddies, their size matters: the eddies thatimpact primary productivity have to act on a scale withsubstantial nutrient variability. Otherwise the inequalitybetween the two growth terms in Eq.15 has small impactsince ch L and ch S ( n L and n S ) are of comparable size.For example the diffusion term representing eddy mixingbeneath the 25 km scale has been found to have very littleimpact ( ∼ (cid:96) ρ ) increases onall scales roughly fourfold (Fig.7:B and Tab.I). This is nothugely surprising, since the removal of mesoscale eddiesincreased primary productivity and doubled the meanchlorophyll concentrations. The increased chlorophyllconcentrations then usually imply a higher chlorophyllvariability. However, the different scales of eddy impacton the chlorophyll distributions can be estimated fromthe chlorophyll scaling slopes ˜∆ (cid:96) ρ (Eq.13 and see also[66]), rather than directly from the magnitude of chloro-phyll spatial variability (∆ (cid:96) ρ , Eq.12). Eddies shouldlower variability (steepen the scaling slope) above thecharacteristic eddy scale (they mix, therefore smooth)and increase variability (flatten the scaling slope) at therange of scales with eddies (due to characteristic eddypatchiness). Assuming that any smoothing effect abovethe eddy scale goes away at a sufficiently large scale, onecan determine the range of scales where the ˜∆ (cid:96) ρ differsbetween the case with and without eddies. Given thatabove (cid:38)
500 km ˜∆ (cid:96) ρ scales with similar slope in bothcases (the case with eddies vs the case without eddies),it is natural to assume that the only important impactof eddy patchiness, or eddy mixing, on the chlorophyllvariability happens at (cid:46)
500 km where the removal ofeddies steepens the ∆ (cid:96) ρ scaling slope (10% increase invariability under ∼
250 km, due to eddy patchiness, seeFig.7 and also Tab.I). It is interesting to analyze the interaction betweenexplicit advection terms and the sub-grid eddy mixingcaptured by the diffusion term (Fig.8, Tab.II). Due tomesoscale eddies and large scale currents (“mean” flows)the smoothing impact of diffusion spreads to the largespatial scales, i.e. at the resolution ∼
25 km scale re-moving diffusion more than triples the chlorophyll vari-ability (Tab.II) and it increases variability by at least10% up to 2000 km scale (see Fig.8:D). We can thenseparate out the relative impact of the mesoscale ed-dies and the mean currents on the large scale smooth-ing (see Fig.8). With the model advection completelyturned off, removing diffusion increased the magnitudeof chlorophyll spatial variability by a maximum 10% atthe resolution scale (Fig.8:A, Tab.II), with a detectableimpact on the chlorophyll variability constrained to the (cid:46)
70 km scales. By switching on mean currents, but nomesoscale eddies, removing the diffusion increased thechlorophyll variability by about 100% at the resolutionscale, and the impact of diffusion on chlorophyll variabil-ity lasted up to ∼ TABLE I. We show the impact of different drivers on the chlorophyll magnitude of spatial variability (∆ (cid:96) ρ ). The Table showsthe values displayed in Fig.7 and Fig.8. The first column shows the percentage change in the magnitude of spatial variabilityat 2500 km after we removed a specific driver (diffusion, eddy and mean advection) from the fully forced FKPP run. Thenumbers in the first column amount to the comparison of the different curves from Fig.7:B at 2500 km, and the purpose of thosenumbers is to show the overall change to the spatial variability at the regional scale. The second-to-fourth column display thepercentage change to the spatial scaling slopes ˜∆ (cid:96) ρ (the scaling slopes are understood as a ratio ∆ (cid:96) ρ/ ∆ L ρ with L = 2500 km,see Fig.7:A) in the situation without a specific driver when compared to the fully forced FKPP model. The percentage changeis shown for a range of values within three intervals of spatial scales: 25-100 km, 100-500 km and 500-2500 km. The ↑↓ symbolsbefore the numbers indicate whether the FKPP value increases ( ↑ ), or decreases ( ↓ ) when the specific driver is removed. removed driver ∆ (cid:96) ρ at 2500 km 25-100 km 100-500 km 500-2500 kmdiffusion ↑ ↑ ↑ ↑ ↑ ↓ ↓ ↓ ↑ ↓ ↓ ↓ ↑ ↓ ↓ ↓ B. Temporal analysis
Fig.9 shows an exact analogue of Fig.8 with the tempo-ral scaling replacing the spatial scaling. It is shown thatwith no advection, the diffusion term has a negligible ef-fect on the magnitude of chlorophyll temporal variabilityabove the daily scale (Fig.9:A). The impact of advectionon the chlorophyll diffusive smoothing (Fig.9:C-D) ap-pears highly non-linear: The largest effect is observed due to the mean currents and this effect is perhaps surpris-ingly reduced when also eddies are removed (Fig.9:D).However, more broadly the conclusions based on the tem-poral analysis (Fig.9) are consistent with the spatial anal-ysis (Fig.8). Fig.9 confirms that advection substantiallyincreases the impact of diffusion on the chlorophyll vari-ability on a large range of scales ( >
10% for up to the180 day time-scale).1
FIG. 7. The panel A shows the percentage reduction in the magnitude of chlorophyll spatial variability (∆ (cid:96) ρ ) when comparedto the magnitude of its spatial variability at the scale L=2500km (∆ L ρ ), or equivalently it compares the spatial scaling slopes( ˜∆ (cid:96) ρ ) for the different simulations. The panel B shows the absolute values for the magnitude of spatial variability ∆ (cid:96) ρ . The∆ (cid:96) ρ (panel B) and ˜∆ (cid:96) ρ (panel A) curves represent the 2018 annual averages of the spatial scaling of the daily data. Both x andy axes are on a log-scale. We show the relative (panel A) and absolute (panel B) chlorophyll magnitude of spatial variability forthe different dynamical scenarios of the FKPP model: a) model forced by both mean and eddy surface currents (“FKPP, adv”),b) model forced only by the mean currents (“FKPP, only mean adv”) and c) model with all the (eddy and mean) advectionremoved (“FKPP, no adv”). In addition to the chlorophyll variability, the cyan line in the panel A shows the magnitude ofspatial variability for the 2017-2018 averaged nutrient concentrations (“CMEMS, nutrient clim”). The dashed lines parallel tothe variability curves mark a 100% and 300% increase in the magnitude of spatial variability with respect to the FKPP modelforced by both eddy and mean advection. The vertical lines show the scales from which the relative scaling remains within10% from the fully (eddy & mean advection) forced FKPP model.TABLE II. We show how the different drivers (eddy and mean advection) propagate the impact of diffusion on the chlorophyllmagnitude of spatial variability (∆ (cid:96) ρ ), as displayed in the Fig.8. The first column shows the percentage change in the magnitudeof spatial variability at 2500 km between the runs with and without diffusion, after we removed a specific driver (diffusion, eddyand mean advection) from the fully forced FKPP run. The numbers in the first column amount to the comparison of the pairsof curves from Fig.8:A-D at 2500 km, and the purpose of those numbers is to show the overall change to the spatial variabilityat the regional scale. The second-to-fourth column display the percentage change to the spatial scaling slopes ˜∆ (cid:96) ρ (the scalingslopes are understood as a ratio ∆ (cid:96) ρ/ ∆ L ρ with L = 2500 km, see Fig.7:A) in the situation with and without diffusion after aspecific driver was removed from the fully forced FKPP model. The percentage change is shown for a range of values withinthree intervals of spatial scales: 25-100 km, 100-500 km and 500-2500 km. The ↑↓ symbols before the numbers indicate whetherthe FKPP value increases ( ↑ ), or decreases ( ↓ ) when the specific driver is removed. removed driver ∆ (cid:96) ρ at 2500 km 25-100 km 100-500 km 500-2500 kmnone ↑ ↑ ↑ ↑ ↑ ↑ ↑ ↑ ↑ ↑ ↑ ↑ ↑ ↑ FIG. 8. The impact of the diffusion term in the FKPP model on the chlorophyll magnitude of spatial variability (2018 averagesfrom daily ∆ (cid:96) ρ ) depending on the advection input: no advection (panel A), only mesoscale eddy advection (panel B), onlymean advection (panel C), both mean and mesoscale eddy advection (panel D). The dashed lines parallel to the ∆ (cid:96) ρ curvesmark a 10% and 100% increase in the magnitude of spatial variability with respect to the FKPP model with the diffusion.The vertical lines show the scale from which ∆ (cid:96) ρ matches the fully forced model within 10%. When there is no advection thediffusion term has a spatially limited impact up to ∼
70 km scale. With mesoscale eddies and/or mean currents the impact ofdiffusion on the chlorophyll magnitude of spatial variability increases 2-4 times at the resolution scale and becomes substantialup to 600-700 km scale with mean currents having non-negligible impact up to the largest scale ( ∼ C. A relationship between chlorophyll spatial andtemporal scales
The chlorophyll distributions are influenced by thecomplex dynamics occurring at wide ranges of spatialand temporal scales. To have a simultaneous under-standing of the ecosystem processes across a range ofspatio-temporal scales, it is of general interest to find arelationship between the characteristic spatial and tem-poral scales for the processes driving surface chlorophyll.In this short section we will not distinguish between thespecific processes driving chlorophyll, but we will demon-strate (Fig.10) a methodology (developed in [66]) on howto find a relationship between the spatial and the tem-poral scales for the magnitude of chlorophyll variability.In essence, the relationship is defined by computing the magnitude of temporal variability for a sequence of lowpass filtered CMEMS chlorophyll spatial distributions ata range of spatial scales (125 km, 500 km, 2000 km).Spatial filtering removes processes that occur on the sub-filter spatial scales and those processes typically influencethe chlorophyll dynamics on some specific range of tem-poral scales. For example, the processes removed by thespatial filtering may lead to a substantial decrease in theCMEMS chlorophyll daily variability (Fig.10). As oneincreases the temporal scale, the spatial high-resolutionscale processes that were removed by the low pass filterplay lesser role in the magnitude of chlorophyll temporalvariability (∆ (cid:96) ρ ) and the ∆ (cid:96) ρ curves of the spatially fil-tered and the unfiltered chlorophyll start converging toeach other. This means that the difference in the dailyvariability between the filtered and the unfiltered chloro-3 FIG. 9. The impact of the diffusion term in the FKPP model on the chlorophyll magnitude of temporal variability dependingon the advection input: no advection (panel A), only mesoscale eddy advection (panel B), only mean advection (panel C),both mean and mesoscale eddy advection (panel D). The dashed lines parallel to the variability curves mark a 10% and 100%increase in the magnitude of temporal variability with respect to the FKPP model with the diffusion. The vertical lines showthe scale from which the variability matches the fully forced model within 10%. The combined impact of mesoscale eddiesand mean currents on how diffusion spreads across time-scales is highly non-linear: removing mean currents close to the dailytime-scale impacts chlorophyll temporal variability more (in both absolute numbers and proportionally) than removing bothmean currents and mesoscale eddies. If there is no advection, diffusion has no impact on the chlorophyll temporal variabilityabove the daily scale (upper left panel). By including mean currents, but no mesoscale eddies, diffusion increases the magnitudeof chlorophyll temporal variability by >
10% on the full range of scales (1 day - 1 year). By including eddy advection, but nomean currents, diffusion has <
10% impact on the magnitude of chlorophyll temporal variability above the ∼
80 day scale, whichis broadly consistent with Fig.6. Removing both mesoscale eddies and mean currents has <
10% impact on the magnitude ofchlorophyll temporal variability above the scale of a half year. phyll (we will call it “Missing Daily Variability of theFiltered Data” and abbreviate it with MDVFD) is re-duced when we increase the temporal scales. The connec-tion between spatial and temporal variability is providedas follows: For each spatial filter (at spatial scale (cid:96) ) wesubdivide the temporal scales into different ranges ( < − > N % at a certainrange of temporal scales (e.g 1-6 months) then we saythat this specific range of temporal scales contains N % of MDVFD. It is then clear that as one increases thespatial scale of the low pass filter one removes processeswith longer temporal scales and larger fraction of MD-VFD will be concentrated at larger temporal scales (e.gabove the scale of 6 months). This provides a connectionbetween the spatial scale of the low pass filter and theranges of temporal scales of MDVFD. The Fig.10 showsthis spatio-temporal relationship: while the 125 km spa-tial filter has 50% of MDVFD on sub-monthly scales andonly 2% of MDVFD on scales larger than half year, the2000 km spatial filter corresponds to 17% of MDVFD on4the sub-monthly scale and almost 50% of MDVFD on the scales larger than half year. D. Impact of nutrients and eddies on chlorophyll
In the climate change scenarios the upper ocean warmsup, leading to increased ocean stratification. The in-creasingly stratified ocean acts as a barrier to verticalnutrient mixing and lowers the surface nutrient concen-trations [74, 75]. Besides nutrients, the increased verticalstratification influences the first baroclinic Rossby radiusimpacting on the mesoscale eddy kinetic energy (EKE,[76, 77]). In this section we will use the FKPP modelto explore the impact of the changed nutrients and EKEon the surface chlorophyll. Although the FKPP modelis a major simplification, we believe it might offer atleast some qualitative insights into how phytoplanktonmight respond to some of the environmental changes. Aform of analytical relationship between the mean chloro-phyll and the mean nutrients has been derived for thestochastic steady state of the FKPP model in Eq.4. How-ever, in reality chlorophyll might be far from a stochasticsteady state prediction described by Eq.4 and we havefound (not shown here) that the stochastic steady statemodel does not approximate well the simulations fromthis study.In Fig.6 we show the spatio-temporal means (for 2018and the FKPP spatial domain) of chlorophyll and nu-trients plotted against each other in a series of experi-ments, where the CMEMS nutrients and EKE (forcingthe FKPP model) were re-scaled by constant factors, i.e.as k.N ( t, (cid:126)x ), where N are the nutrients from the CMEMSmodel. The constant ( k ) factors were for nutrients takenfrom the k ∈ (0 . , .
4) interval and for EKE from the k ∈ (0 , .
7) interval (in case of EKE we rescaled eacheddy velocity component with the same factor). Fig.6demonstrates that changing EKE by ± ∼
85% and more) can have a large impact on the meanchlorophyll concentrations. Fig.6 also demonstrates thatunder the increased EKE, phytoplankton becomes in-creasingly insensitive to the changing nutrients: withzero EKE chlorophyll scales almost linearly with nutri-ents (e.g. 50% decrease in nutrients amounts to 50% decrease in chlorophyll concentrations, similar to Eq.5),while with increased EKE the scaling becomes increas-ingly sub-linear (e.g. for CMEMS EKE 50% decreasein nutrients there is about 20% decrease in chlorophyllconcentrations). This is an interesting result implyingthat in the increased EKE scenario the phytoplanktonconcentrations become more stable. In particular, theFKPP model suggests (Fig.6) that within ±
70% of thecurrent EKE levels, a dramatic decline of nutrients hascomparably small impact on chlorophyll. It is not en-tirely clear how to interpret this result, neither how se-riously it should be taken: we would recommend to takeit with a lot of caution, unless it is reconfirmed in morerealistic simulations using higher complexity models.It is also interesting to explore how the chlorophyll sur-face distributions respond to the changes imposed on thenutrients, or EKE. Since in tropical Pacific the chloro-phyll spatial variability dominates over temporal vari-ability (Fig.3) it is useful to understand how the spa-tial regional patterns of chlorophyll change under thechanged chlorophyll mean. In Fig.11 we show the im-pact of halved nutrient concentrations (Fig.11:A-B) anddecreased EKE by 88% (Fig.11:C-D) on the chlorophyllannual mean spatial distributions. The changed nutrientconcentrations and the eddy velocities were re-scaled ver-sions of the original CMEMS data, where by “re-scaled”we mean the original CMEMS distributions multiplied bya spatio-temporally constant factor. The Fig.11 showsthat the resulting chlorophyll 2018 mean spatial distri-butions are far from being the re-scaled versions of the2018 mean chlorophyll forced by the CMEMS data. Inparticular Fig.11:A shows that under the nutrient de-cline chlorophyll changes by substantially larger propor-tion in the areas with higher chlorophyll concentrations(eastern tropical Pacific). This indicates that areas withthe highest biological activity are also most vulnerableto change. It is perhaps surprising that reducing nutri-ents (Fig.11:A-B) has proportionally largest impact onchlorophyll in the chlorophyll-rich areas, since the sameareas have correspondingly highest eddy activity (Fig.1)and chlorophyll is less sensitive to the nutrient concen-trations in the presence of eddies (Fig.6).
E. Scale-propagation of a multiplicative stochasticnoise
In the last part of our analysis we investigate the im-pact of a stochastic Gaussian white noise on the chloro-phyll dynamics across range of spatial scales. Such whitenoise usually represents a number of higher-complexity, scale-constrained processes that were not explicitly in-cluded into the dynamical model. If such processeshave linear relationship to the dynamical model vari-ables, their impact on the model variables will remainconstrained to the (spatio-temporal) scales of those pro-cesses. However if the relationship between those pro-cesses and dynamical model variables is highly non-5
FIG. 10. The upper panel shows the magnitude of temporal variability for the spatially filtered CMEMS chlorophyll data(moving median filter) at a range of scales: 25 km (original resolution), 125 km, 500 km and 2000 km, with the purpleline being time variability of the regional mean value. The two local minima in the purple line correspond to the bi-annualperiodicity of the CMEMS chlorophyll. The scale where the temporal variability of the spatially filtered data meets with thetemporal variability of the original CMEMS (25 km) data is the scale where the processes removed by the spatial filtering haveno longer impact on the magnitude of chlorophyll temporal variability. The upper panel then provides connection between thespatial and the temporal scales shown in the bottom panel. The bottom panel demonstrates how the temporal variability of thespatially filtered data (at 125 km, 500 km and 2000 km) splits (in %) into three categories: < > > FIG. 11. Panels A-B show the impact of 50% nutrient decrease on the annual 2018 mean chlorophyll concentrations. Similarlyto A-B, the panels C-D show the impact of 88% decrease in the EKE on the mean annual 2018 chlorophyll concentrations.The panels B and D show the absolute change (in mg/m ) in chlorophyll concentrations when compared to the simulationforced by the CMEMS nutrients and EKE. The panels A and C show the same change, but relative (in %) to the values of thesimulation using the CMEMS data. linear, the impact of those processes on the model vari-ables may propagate beyond the original scale of the pro-cess. A simple example is the impact of wind stress onthe vertical mixing and primary productivity in the wa-ter column: the phenomena observable on weekly time-scales, such as phytoplankton blooms (e.g. see the criticalturbulence hypothesis in [78]), may be sensitive to suchdetails, as to whether we capture wind stress with anhourly, or 3-hourly resolution [79].We have run the 2017-2018 model simulation with amultiplicative white noise [80, 81] to account for a ran-dom variability in the growth rate parameter P (Eq.1).The multiplicative Gaussian noise has already proven tobe both realistic and useful in the population dynam-ics models [82, 83]. The Fig.12 compares simulationsin which the growth parameter ( P ) was perturbed bythe Gaussian noise with 20% standard deviation (cor-responding to ∆ P = ± . . − m mmol − s − ). The random perturbations were applied as a white noise onthe FKPP model-grid spatio-temporal scale (25 km and1 day). The Fig.12 shows the magnitude and scale-propagation of the stochastic noise impact on chlorophyllin simulations using different sets of dynamical drivers(the FKPP model, the FKPP model without mean cur-rents, the FKPP model without eddies, the FKPP modelwithout any advection). The outputs for the stochasticsimulations were low-pass filtered at different scales (25km, 100 km, 400 km, 1600 km and at the “regional” scale,6400 km, where only total spatial averages were calcu-lated) and compared with the corresponding low-pass fil-tered deterministic simulations (with the fixed P value).The chosen metric for the comparison was the Root MeanSquare Difference (RMSD). The Fig.12:A shows the per-centage of the 25 km scale RMSD that remains on scales >
25 km. The larger is the percentage, the more is the25 km white noise propagated to the larger scales by the7model dynamics. The reduction of chlorophyll RMSDas a function of scale is compared to the scaling of themean absolute value of the white noise originally appliedto the FKPP model (shown by the black dashed curvein Fig.12:A labeled as “Noise”). The white noise is bydefinition uncorrelated on the scales above 25 km, butit remains visible also on the 100-1000 km spatial scales(on the level of < ≥
100 km scales effectively averages outthe white noise over a finite number of samples, so thelow-pass filtered mean will differ from the theoretical zeromean of the sampling Gaussian distribution. The numberof samples increases with the spatial scale of the low-passfilter and in the limit of infinite scale the mean absolutevalue of the noise is precisely zero. Since the non-lineardynamics of the FKPP model is expected to propagatethe white noise to larger scales, the mean absolute valueof the white noise applied to the FKPP model is ex-pected to reduce faster than the RMSD of chlorophyll.The black dashed curve in Fig.12:A can be then inter-preted as a “theoretical maximum” for the chlorophyllRMSD reduction as a function of scale, such theoreti-cal maximum being reached when the FKPP does notscale-propagate the stochastic noise.For a non-advective FKPP model ( (cid:126)u = κ = 0) the mul-tiplicative noise generates at each spatial point a typeof random-walk solution which is constrained to someneighborhood of the steady state solution (Eq.7). The steps of the random walk are larger in nutrient-rich areas,however this might be compensated by the fact that theconvergence of a perturbed solution to the unperturbedsolution might be faster in the areas with larger nutrientconcentrations (Eq.8). The Fig.12 shows that the multi-plicative noise with 20% standard deviation leads to 4%RMSD in chlorophyll when the model has no advection,or runs with only mean advection (Fig.12:B). The FKPPmodel without mesoscale eddies (and diffusion) does notpropagate the noise to the ≥
100 km scales, as the noisereduction in those simulations is close to its “theoreticalmaximum” (Fig.12:B). When mesoscale eddies (and dif-fusion) are included, the fluctuations in chlorophyll intro-duced by the stochastic noise on the 25 km scale, decreaseto 2%, or 1% depending on whether we include also themean currents (Fig.12:B). However, mesoscale eddies andthe diffusion term introduce scale-propagation into thechlorophyll noise, with 10-30% of the 25-km fluctuationsvisible on the 100-500 km scales (Fig.12:A). The reasonfor this scale-propagation of the chlorophyll noise is theeddy and diffusive mixing, which smooths the chlorophyllnoise, lowering the size of the chlorophyll random fluctua-tions (see the lower RMSD at the 25 km scale, Fig.12:B),but introducing larger-scale correlations to the randomfluctuations. These larger scale correlations explain whythe RMSD reduces comparably slowly as a function ofscale (Fig.12:A).
F. Summary
Low complexity models based on the Fischer-Kolmogorov-Petrovski-Piskunov (FKPP) equation havebeen often used to study conceptual questions in popu-lation biology, such as the critical patch size for popu-lation survival [31]. However, a realistic simulation ofphytoplankton dynamics in a specific global region istypically assumed to require a medium, or high com-plexity models. Here we demonstrate that for very spe-cific purposes in a suitably tailored choice of region (e.g.tropical Pacific), the FKPP model forced by a higher-complexity model outputs for nutrients and surface cur-rents, provides a sufficiently realistic simulation for thephytoplankton chlorophyll concentrations (a proxy forprimary productivity and phytoplankton biomass). Theadvantage of the FKPP model is that the model dependsonly on few external inputs and model parameters, all ofwhich are straightforward to interpret and modify. Sincethe model is computationally cheap to run and can beeasily perturbed with a stochastic noise, one can pro-duce almost arbitrary number of both deterministic andstochastic simulations.We use the FKPP model to develop a multi-scale viewof a driver (eddy and mean advection, diffusion) impacton the chlorophyll distributions. The impact of different drivers on chlorophyll is explored in a series of simula-tions, where we remove specific set of drivers and anal-yse the changes to the chlorophyll variability on a rangeof spatial (25-2500 km) and temporal (1 day - 1 year)scales. We show that for the 1 / ◦ model, advection hasa major impact on the mean chlorophyll concentrations.Diffusion has a negligible impact on the mean chloro-phyll concentrations, but it is propagated by the largerscale currents and influences chlorophyll variability on awide range of spatial and temporal scales. We analysethe impact of surface nutrient decline and changes to themesoscale eddy kinetic energy (EKE) on the mean sur-face chlorophyll concentrations (some changes to nutri-ents and EKE are projected in the future climate scenar-ios). The FKPP model indicates that unless EKE radi-cally changes from its current levels, chlorophyll tends toscale sub-linearly with nutrients, which implies that thechlorophyll concentrations are relatively stable with re-spect to the nutrient decline. However, the FKPP modelalso shows that the chlorophyll sensitivity to nutrientsgoes through a sudden transition and becomes substan-tially larger if we minimise the EKE to 0-15% from itscurrent value. In the limit of vanishing EKE, chloro-phyll scales with nutrients approximately linearly. Wealso investigate the spatial scale-propagation of a whitemultiplicative Gaussian noise, introduced into the FKPP8 FIG. 12. The two panels show the impact of different drivers (e.g. mean and eddy advection, biological activity) on thepropagation of a white stochastic noise in the FKPP model. The panel B shows the Root Mean Square Difference (RMSD)in chlorophyll between the stochastic run and the corresponding deterministic run (y-axis) vs spatial log-scale (x-axis). TheRMSD values are divided by the mean 2018 chlorophyll of the deterministic run and shown in %. The panel A shows thesame quantity, only compared (in %) to its own value at the lowest, 25 km scale. The purpose of the panel A is to show howthe impact of the stochastic noise propagates through the spatial scales under different dynamical scenarios. The differentscenarios are: a) the FKPP model configuration with the mean and eddy currents (”FKPP, adv”), the FKPP model withthe mean currents removed (”FKPP, only eddy adv”), the same FKPP configuration with the mesoscale eddies and diffusionremoved (”FKPP, only mean adv”), and the FKPP model without any advection (”FKPP, no adv”). The panel A comparesthe chlorophyll RMSD to the scaling of a white noise (”Noise”) applied at the model resolution scale. model on the model resolution scale. We demonstratethat the impact of the stochastic noise on the chloro-phyll concentrations propagates to 100-500 km spatialscales through the mixing by eddy advection and diffu-sion term.This study aims to provide an inspiration for re-searchers to further explore specific contexts in whichlow-complexity models could serve as a sufficiently real-istic tool to address questions that would often be beyondthe computational affordability of the higher-complexitymodels. The limitations of the low-complexity modelneed to be always recognized, but this should not meanthat low-complexity models have to be always discardedas a tool of realistic modelling. Eventually the fu- ture modelling could become a multi-complexity effort,where high and medium complexity models become inte-grated with low-complexity models, each serving its op-timal purpose while mutually achieving the desired goalwith a reduced computational cost. Moreover, the low-complexity models such as the FKPP model used in thisstudy, could provide a priceless public educational toolto enhance the understanding of marine biogeochemistryin different realistic situations.
Acknowledgments:
The data used to force theFKPP model were obtained from Copernicus Marine En-vironment Service (CMEMS) and can be freely down-loaded from https://resources.marine.copernicus.eu/. [1] S. A. Levin, The problem of pattern and scale in ecology:the robert h. macarthur award lecture, Ecology , 1943(1992).[2] C. Heinze and M. Gehlen, Modeling ocean biogeochem-ical processes and the resulting tracer distributions, in International Geophysics , Vol. 103 (Elsevier, 2013) pp.667–694.[3] D. Ford, S. Key, R. McEwan, I. Totterdell, andM. Gehlen, Marine biogeochemical modelling and dataassimilation for operational forecasting, reanalysis, andclimate research, New Frontiers in Operational Oceanog-raphy , 625 (2018). [4] K. G. Wilson, Renormalization group and critical phe-nomena. i. renormalization group and the kadanoff scal-ing picture, Physical Review B , 3174 (1971).[5] K. G. Wilson and J. Kogut, The renormalization groupand the expansion, Physics reports , 75 (1974).[6] K. G. Wilson, The renormalization group: Critical phe-nomena and the kondo problem, Reviews of modernphysics , 773 (1975).[7] M. Kardar, G. Parisi, and Y.-C. Zhang, Dynamic scalingof growing interfaces, Physical Review Letters , 889(1986).[8] Y. G. Sinai and K. Khanin, Renormalization group method in the theory of dynamical systems, InternationalJournal of Modern Physics B , 147 (1988).[9] J. Cardy, Scaling and renormalization in statisticalphysics , Vol. 5 (Cambridge university press, 1996).[10] H. D. Rozenfeld, C. Song, and H. A. Makse, Small-worldto fractal transition in complex networks: a renormal-ization group approach, Physical Review Letters ,025701 (2010).[11] D. Forster, D. R. Nelson, and M. J. Stephen, Large-distance and long-time properties of a randomly stirredfluid, Physical Review A , 732 (1977).[12] N. Dupuis, L. Canet, A. Eichhorn, W. Metzner,J. Pawlowski, M. Tissier, and N. Wschebor, The nonper-turbative functional renormalization group and its appli-cations, arXiv preprint arXiv:2006.04853 (2020).[13] S.-k. Ma and G. F. Mazenko, Critical dynamics of ferro-magnets in 6- ε dimensions: general discussion and de-tailed calculation, Physical Review B , 4077 (1975).[14] M. E. Fisher, Scaling, universality and renormalizationgroup theory, in Critical phenomena (Springer, 1983) pp.1–139.[15] L. P. Kadanoff, Scaling and universality in statisticalphysics, Physica A: Statistical Mechanics and its Appli-cations , 1 (1990).[16] A. N. Kolmogorov, The local structure of turbulence inincompressible viscous fluid for very large reynolds num-bers, Cr Acad. Sci. URSS , 301 (1941).[17] B. B. Mandelbrot, The fractal geometry of nature , Vol.173 (WH freeman New York, 1983).[18] D. L. Turcotte,
Fractals and chaos in geology and geo-physics (Cambridge university press, 1997).[19] B. B. Mandelbrot,
Fractals and scaling in finance:Discontinuity, concentration, risk. Selecta volume E (Springer Science & Business Media, 2013).[20] S. Lovejoy and D. Schertzer,
The weather and climate:emergent laws and multifractal cascades (Cambridge Uni-versity Press, 2013).[21] I. Grooms and K. Julien, Multiscale models in geophys-ical fluid dynamics, Earth and Space Science , 668(2018).[22] S. A. Levin, Self-organization and the emergence of com-plexity in ecological systems, Bioscience , 1075 (2005).[23] R. A. Fisher, The wave of advance of advantageous genes,Annals of eugenics , 355 (1937).[24] A. Kolmogorov, I. Petrovskii, and N. Piskunov, Studyof a diffusion equation that is related to the growth ofa quality of matter, and its application to a biologicalproblem, Byul. Mosk. Gos. Univ. Ser. A Mat. Mekh ,26 (1937).[25] J. G. Skellam, Random dispersal in theoretical popula-tions, Biometrika , 196 (1951).[26] S. Petrovskii, On the plankton front waves acceleratedby marine turbulence, Journal of marine systems , 179(1999).[27] R. Cantrell and C. Cosner, Spatial heterogeneity and crit-ical patch size: area effects via diffusion in closed environ-ments, Journal of Theoretical Biology , 161 (2001).[28] J. Joo and J. L. Lebowitz, Population dynamics in spa-tially heterogeneous systems with drift: The generalizedcontact process, Physical Review E , 036112 (2005).[29] D. A. Birch, Y.-K. Tsang, and W. R. Young, Boundingbiomass in the fisher equation, Physical Review E ,066304 (2007).[30] Z. Neufeld, Stirring effects in models of oceanic plank- ton populations, Chaos: An Interdisciplinary Journal ofNonlinear Science , 037102 (2012).[31] A. B. Ryabov and B. Blasius, Population growth andpersistence in a heterogeneous environment: the role ofdiffusion and advection, Mathematical Modelling of Nat-ural Phenomena , 42 (2008).[32] E. A. Fulton, A. D. Smith, and C. R. Johnson, Effect ofcomplexity on marine ecosystem models, Marine EcologyProgress Series , 1 (2003).[33] C. Schrum, I. Alekseeva, and M. S. John, Developmentof a coupled physical–biological ecosystem model ecosmo:part i: model description and validation for the north sea,Journal of Marine Systems , 79 (2006).[34] A. Yool, E. Popova, and T. Anderson, Medusa-2.0: an in-termediate complexity biogeochemical model of the ma-rine carbon cycle for climate change and ocean acidifi-cation studies, Geoscientific Model Development , 1767(2013).[35] O. Aumont, C. ´Eth´e, A. Tagliabue, L. Bopp, andM. Gehlen, Pisces-v2: An ocean biogeochemical modelfor carbon and ecosystem studies., Geoscientific ModelDevelopment Discussions (2015).[36] M. Vichi, N. Pinardi, and S. Masina, A generalized modelof pelagic biogeochemistry for the global ocean ecosys-tem. part i: Theory, Journal of Marine Systems , 89(2007).[37] M. Butensch¨on, J. Clark, J. N. Aldridge, J. I. Allen,Y. Artioli, J. Blackford, J. Bruggeman, P. Cazenave,S. Ciavatta, S. Kay, et al. , Ersem 15.06: a generic modelfor marine biogeochemistry and the ecosystem dynamicsof the lower trophic levels, Geoscientific Model Develop-ment , 1293 (2016).[38] S. B. Groom, S. Sathyendranath, Y. Ban, S. Bernard,B. Brewin, V. Brotas, C. Brockmann, P. Chauhan, J.-k.Choi, A. Chuprin, et al. , Satellite ocean colour: currentstatus and future perspective, Frontiers in Marine Science , 485 (2019).[39] G. Madec et al. , Nemo ocean engine, (2015).[40] R. H. Stanley, J. B. Kirkpatrick, N. Cassar, B. A. Bar-nett, and M. L. Bender, Net community production andgross primary production rates in the western equatorialpacific, Global Biogeochemical Cycles (2010).[41] M. J. Behrenfeld, K. Worthington, R. M. Sherrell, F. P.Chavez, P. Strutton, M. McPhaden, and D. M. Shea,Controls on tropical pacific ocean productivity revealedthrough nutrient stress diagnostics, Nature , 1025(2006).[42] S. W. Chisholm, F. M. Morel, et al. , What controls phy-toplankton production in nutrient-rich areas of the opensea?, (1991).[43] R. Dugdale and F. Wilkerson, Low specific nitrate uptakerate: A common feature of high-nutrient, low-chlorophyllmarine ecosystems, Limnology and Oceanography ,1678 (1991).[44] R. Eppley and E. Renger, Nitrate utilization by planktonin the equatorial pacific march 1988 along 150 w, Journalof Geophysical Research: Oceans , 663 (1992).[45] Z. I. Johnson, R. Shyam, A. E. Ritchie, C. Mioni, V. P.Lance, J. W. Murray, and E. R. Zinser, The effect ofiron-and light-limitation on phytoplankton communitiesof deep chlorophyll maxima of the western pacific ocean,Journal of Marine Research , 283 (2010).[46] J.-Y. Park, J. P. Dunne, and C. A. Stock, Ocean chloro- phyll as a precursor of enso: An earth system modelingstudy, Geophysical Research Letters , 1939 (2018).[47] B. W. Frost, The role of grazing in nutrient-rich areasof the open sea, Limnology and Oceanography , 1616(1991).[48] D. B. Chelton, R. A. DeSzoeke, M. G. Schlax, K. El Nag-gar, and N. Siwertz, Geographical variability of the firstbaroclinic rossby radius of deformation, Journal of Phys-ical Oceanography , 433 (1998).[49] D. B. Chelton, M. G. Schlax, and R. M. Samelson, Globalobservations of nonlinear mesoscale eddies, Progress inoceanography , 167 (2011).[50] R. Hallberg, Using a resolution function to regulate pa-rameterizations of oceanic mesoscale eddy effects, OceanModelling , 92 (2013).[51] E. Novikov and R. Stewart, Intermittency of turbulenceand spectrum of fluctuations in energy-disspation, Izv.Akad. Nauk. SSSR. Ser. Geofiz , 408 (1964).[52] A. Yaglom, The influence of fluctuations in energy dissi-pation on the shape of turbulence characteristics in theinertial interval, in Soviet Physics Doklady , Vol. 11 (1966)p. 26.[53] A. Obukhov,
Structure of the temperature field in tur-bulent flow , Tech. Rep. (ARMY BIOLOGICAL LABSFREDERICK MD, 1968).[54] B. B. Mandelbrot, Intermittent turbulence in self-similarcascades: divergence of high moments and dimension ofthe carrier, Journal of Fluid Mechanics , 331 (1974).[55] U. Frisch, P.-L. Sulem, and M. Nelkin, A simple dy-namical model of intermittent fully developed turbulence,Journal of Fluid Mechanics , 719 (1978).[56] D. Schertzer and S. Lovejoy, Physical modeling and anal-ysis of rain and clouds by anisotropic scaling multiplica-tive processes, Journal of Geophysical Research: Atmo-spheres , 9693 (1987).[57] D. Schertzer and S. Lovejoy, Multifractals, generalizedscale invariance and complexity in geophysics, Interna-tional Journal of Bifurcation and Chaos , 3417 (2011).[58] D. Schertzer and S. Lovejoy, Multifractal simulationsand analysis of clouds by multiplicative processes, At-mospheric research , 337 (1988).[59] F. Schmitt, D. Schertzer, S. Lovejoy, and Y. Brunet, Mul-tifractal temperature and flux of temperature variance infully developed turbulence, EPL (Europhysics Letters) , 195 (1996).[60] L. Seuront, F. Schmitt, Y. Lagadeuc, D. Schertzer,S. Lovejoy, and S. Frontier, Multifractal analysis of phy-toplankton biomass and temperature in the ocean, Geo-physical Research Letters , 3591 (1996).[61] S. Lovejoy, D. Schertzer, and J. Stanway, Direct evi-dence of multifractal atmospheric cascades from plane-tary scales down to 1 km, Physical Review Letters ,5200 (2001).[62] S. Lovejoy, W. Currie, Y. Tessier, M. Claereboudt,E. Bourget, J. Roff, and D. Schertzer, Universal mul-tifractals and ocean patchiness: phytoplankton, physicalfields and coastal heterogeneity, Journal of Plankton Re-search , 117 (2001).[63] J.-S. Gagnon, S. Lovejoy, and D. Schertzer, Multifractalearth topography, Nonlinear Processes in Geophysics ,541 (2006).[64] L. De Montera, M. Jouini, S. Verrier, S. Thiria, andM. Cr´epon, Multifractal analysis of oceanic chlorophyllmaps remotely sensed from space, Ocean Science , 219 (2011).[65] L. De Montera, S. Verrier, C. Mallet, and L. Barth`es, Apassive scalar-like model for rain applicable up to stormscale, Atmospheric Research , 140 (2010).[66] J. Sk´akala, T. Smyth, R. Torres, C. Buckingham,A. Brearley, P. Hyder, and A. Coward, SST dynamicsat different scales: Evaluating the oceanographic modelresolution skill to represent SST processes in the southernocean, Journal of Geophysical Research: Oceans (2019).[67] H. Kierstead and L. B. Slobodkin, The size of watermasses containing plankton blooms, J. mar. Res , 141(1953).[68] A. Okubo, Advection - diffusion in the presence of sur-face convergence, in Oceanic fronts in coastal processes (Springer, 1978) pp. 23–28.[69] A. Okubo, Horizontal dispersion and critical scales forphytoplankton patches, in
Spatial pattern in planktoncommunities (Springer, 1978) pp. 21–42.[70] J. D. Murray and R. P. Sperb, Minimum domains forspatial patterns in a class of reaction diffusion equations,Journal of Mathematical Biology , 169 (1983).[71] E. E. Holmes, M. A. Lewis, J. Banks, and R. Veit, Partialdifferential equations in ecology: spatial interactions andpopulation dynamics, Ecology , 17 (1994).[72] D. C. Speirs and W. S. Gurney, Population persistencein rivers and estuaries, Ecology , 1219 (2001).[73] A. Okubo and S. A. Levin, Diffusion and ecological prob-lems: modern perspectives , Vol. 14 (Springer Science &Business Media, 2013).[74] M. S. Lozier, A. C. Dave, J. B. Palter, L. M. Gerber,and R. T. Barber, On the relationship between stratifi-cation and primary productivity in the North Atlantic,Geophysical Research Letters (2011).[75] A. Capotondi, M. A. Alexander, N. A. Bond, E. N. Cur-chitser, and J. D. Scott, Enhanced upper ocean stratifi-cation with climate change in the cmip3 models, Journalof Geophysical Research: Oceans (2012).[76] O. A. Saenko, Influence of global warming on baroclinicRossby radius in the ocean: A model intercomparison,Journal of Climate , 1354 (2006).[77] M. Sueyoshi and T. Yasuda, Reproducibility and fu-ture projection of the ocean first baroclinic rossby radiusbased on the cmip3 multi-model dataset, Journal of theMeteorological Society of Japan. Ser. II , 821 (2009).[78] J. Huisman, P. van Oostveen, and F. J. Weissing, Criticaldepth and critical turbulence: two different mechanismsfor the development of phytoplankton blooms, Limnologyand Oceanography , 1781 (1999).[79] H. R. Powley, J. Bruggeman, J. Hopkins, T. Smyth, andJ. Blackford, Sensitivity of shelf sea marine ecosystems totemporal resolution of meteorological forcing, Journal ofGeophysical Research: Oceans , e2019JC015922 (2020).[80] A. Schenzle and H. Brand, Multiplicative stochastic pro-cesses in statistical physics, Physical Review A , 1628(1979).[81] R. Graham and A. Schenzle, Carleman imbedding of mul-tiplicative stochastic processes, Physical Review A ,1731 (1982).[82] G. Denaro, D. Valenti, B. Spagnolo, G. Basilone, S. Maz-zola, S. W. Zgozi, S. Aronica, and A. Bonanno, Dynam-ics of two picophytoplankton groups in mediterraneansea: Analysis of the deep chlorophyll maximum by astochastic advection-reaction-diffusion model, PLoS One , e66765 (2013). [83] D. Valenti, G. Denaro, B. Spagnolo, S. Mazzola,G. Basilone, F. Conversano, C. Brunet, and A. Bo- nanno, Stochastic models for phytoplankton dynamics inmediterranean sea, Ecological Complexity27