Dissecting mesopelagic particulate organic carbon budgets in the North Atlantic: A mechanistic diagnosis and evaluation of PISCESv2_RC

Biogeochemical and physical processes in the mesopelagic layer regulate long-term carbon storage in the ocean interior. However, uncertainties in particulate organic carbon (POC) budgets limit quantitative understanding of the biological carbon pump and its representation in models. Here we analyse POC budgets simulated by the NEMO4–PISCESv2_RC model in the upper 1000 m of the North Atlantic over a climatological seasonal cycle, evaluating model performance and diagnosing the mechanisms that regulate POC export and transfer efficiency. Comparison with satellite, Argo-float, shipboard, and sediment-trap data indicates generally realistic POC stocks and fluxes. However, the model exhibits systematic biases, including (i) excessive diatom dominance at mid–high latitudes, (ii) compensation between too-low primary production and excessive epipelagic export in the subtropics, and (iii) underestimated mesopelagic POC stocks at mid and low latitudes. In the model, the mesopelagic POC supply driven by gravitational export and large detritus increases poleward, while winter–spring mixing supplies an additional 37 % in the subpolar region. Up to 60 % of mesopelagic POC supply is intercepted by zooplankton, yet most is recycled to detritus through fragmentation and trophic processing. This detrital loop modulates particle size, sinking speed, and degradation pathways, and is seasonally reinforced at mid–high latitudes. The lability of exported POC increases with latitude, whereas temperature-driven decay rates decrease with latitude. These opposing gradients result in maximal mesopelagic POC flux attenuation at midlatitudes, coincident with peak productivity. Small detritus plays a central role in the model, accounting for ∼ 55 % of mesopelagic POC decay and driving 33 %–50 % of vertical flux at 1000 m. These large modelled contributions reflect substantial biological supply and a pronounced vertical decline in lability, together favouring mesopelagic transfer of small POC. We propose POC budget analysis as a mechanistic framework for identifying structural biases and constraining inter-model spread in projections of the biological particulate carbon pump.

1 Introduction

Through photosynthesis, phytoplankton fix inorganic carbon (C) into organic matter that fuels marine food webs (Lindeman, 1942). This organic matter flows in particulate form through trophic interactions between microbial primary producers (phytoplankton), their protist and metazoan predators (zooplankton), and detrital particles – collectively termed particulate organic carbon (POC). A fraction of primary production is released as dissolved organic carbon (DOC); in addition to being remineralised to dissolved inorganic carbon (DIC) and nutrients, DOC can re-enter the particulate pool via the microbial loop (Pomeroy, 1974).

A portion of surface POC escapes rapid recycling in the upper ocean and is transferred to depth via gravitational sinking and other transport mechanisms (Eppley and Peterson, 1979; Volk and Hoffert, 1985; Kwon et al., 2009; Boyd et al., 2019). Additional vertical fluxes include biogenic particulate inorganic carbon sinking (Honjo, 1980; Neukermans et al., 2023) and DOC transport (Jiao et al., 2010). The vertical transfer of biogenic carbon pools, and their gradual transformation to DIC and subsequent storage in the ocean interior, constitute the biological carbon pump (Volk and Hoffert, 1985; Marinov et al., 2008b; Legendre, 2024; Frenger et al., 2024). Here, we focus on its particulate organic component.

POC export to depth occurs through multiple pathways, including “passive” gravitational sinking (Volk and Hoffert, 1985), turbulent mixing and diffusion (Gardner et al., 1995; Bol et al., 2018), advection (Levy et al., 2013), and “active” transport by diel and seasonal zooplankton migrations (Jónasdóttir et al., 2015; Brun et al., 2019; Gorgues et al., 2019). Once in the mesopelagic, POC is transformed by zooplankton detritivory (Steinberg and Landry, 2016; Giering et al., 2014) and particle aggregation-disaggregation processes (Mayor et al., 2014; Takeuchi et al., 2019; Briggs et al., 2020), and ultimately removed through bacterial degradation (Belcher et al., 2016). As a result, the vertical POC flux is predominantly attenuated in the mesopelagic zone (∼ 100–1000 m; Martin et al., 1987). The balance between vertical transfer and remineralisation regulates carbon storage, nutrient regeneration (Weber et al., 2016; Rodgers et al., 2024), and sea-air CO₂ exchange (Ricour et al., 2023), thereby influencing climate feedbacks (Kwon et al., 2009; Le Quéré et al., 2016).

POC spans a broad size spectrum, from hundred nm to mm (Stemmann and Boss, 2012). Large particles (> 100 µm), including phytodetrital aggregates and faecal pellets with their microbial colonisers (Kiørboe et al., 2003) dominate fast gravitational export due to their fast settling velocities in the order of tens to hundreds of m d⁻¹ (Stemmann and Boss, 2012). In contrast, small particles (< 100 µm), which comprise > 85 % of the POC stock (Baker et al., 2017; Galí et al., 2022) are mostly suspended or sink slowly (< 10 m d⁻¹). However, owing to their large standing stocks, small particles contribute substantially to non-gravitational export pathways, such as turbulent diffusion – the so-called “mixed layer pump (Giering et al., 2016; Lacour et al., 2019) – and advection (Alonso-González et al., 2009).

Particle size is linked to reactivity, reflecting variations in composition, age and bioavailability (Walker et al., 2016; Kharbush et al., 2020; Johnson et al., 2020). In surface waters, living organisms may dominate the POC pool, especially in productive settings (Gasol et al., 1997), whereas the proportion of detrital POC increases with depth. However, uncertainties in POC apportionment persist, especially in the mesopelagic zone (Galí et al., 2022; Koestner et al., 2024). As particles undergo biotic and abiotic transformations, their size, sinking velocity and reactivity may vary in concert along the water column. Slow-sinking detritus experiences longer residence times in the upper ocean and is more exposed to biological transformation. Therefore, sinking-reactivity models (Aumont et al., 2017) predict slow-sinking detritus to become refractory at shallower depths than fast-sinking detritus.

Observationally, POC export and mesopelagic transformation are difficult to quantify due to the intermittency (Berger and Wefer, 1990), seasonality (de Melo Viríssimo et al., 2024), and fine spatial scales (Bol et al., 2018; Briggs et al., 2020; Lacour et al., 2023) of POC export events. For instance, annual mean export estimates can carry uncertainties of up to ±60 % in strongly seasonal regions (Henson et al., 2015). Additional challenges arise from limited quantitative knowledge about mesopelagic bacteria and zooplankton activities (Arístegui et al., 2009; Hernández-León et al., 2019), missing or poorly quantified POC inputs (Boyd et al., 2019), biases in metabolic measurements (Burd et al., 2010), and temporal mismatches between POC inputs and losses (Uchimiya et al., 2018). Together, these limitations hinder closure of mesopelagic POC budgets using observations alone (Giering et al., 2014).

Despite recent advances in observing POC stocks (Claustre et al., 2021), export fluxes (Henson et al., 2024), and biological turnover rates (Bressac et al., 2024), comprehensive quantification of POC fluxes across regions and seasons remains elusive. Therefore, ocean biogeochemistry models are essential tools to mechanistically integrate physical transport, ecosystem structure and biogeochemical transformations (Fennel et al., 2022; Henson et al., 2022), enabling exploration of spatiotemporal scales and processes that are inaccessible to observations alone. Here, we investigate mesopelagic POC dynamics and budgets using the NEMO 4.0.4–PISCESv2_RC coupled model, focusing on the North Atlantic – a data-rich region with sharp biogeochemical gradients.

Compared to other models participating in the Coupled Model Intercomparison Project (CMIP) (Henson et al., 2022), PISCESv2_RC incorporates several distinctive features relevant to mesopelagic POC cycling: (i) a variable reactivity or reactivity continuum (RC) scheme (Boudreau and Ruddick, 1991) that simulates vertical changes in small and large POC reactivity (Aumont et al., 2017); (ii) explicit zooplankton detritivory through both phagotrophic particle ingestion and flux feeding; (iii) parameterisations for fragmentation of large detritus by mesozooplankton and bacterial processes; and (iv) abiotic aggregation-disaggregation dynamics that interconvert between DOC and detrital POC.

Future projections of productivity, export efficiency (EE), and transfer efficiency (TE) under climate change diverge widely in both magnitude and spatial patterns (Wilson et al., 2022; Doney et al., 2024; Wang and Fennel, 2024; Walker and Palevsky, 2025; Doléac et al., 2025). These discrepancies reflect differences in the representation of physical circulation and particle export, remineralisation, and food-web structure (e.g., Walker and Palevsky, 2025; Brabson et al., 2025). As PISCES underpins several CMIP model configurations, a detailed diagnosis of its detrital POC budgets in a relatively well-observed region like the North Atlantic is essential for interpreting inter-model spread and identifying structural and parametric biases relative to observations.

The objectives of this study are to:

• Analyse the variability of modelled POC stocks and fluxes in the upper 1000 m of the North Atlantic, from subtropical to subpolar latitudes, and evaluate them against available observations from satellites, autonomous floats, shipboard surveys and sediment traps.

• Quantify the budgets of small and large detrital POC in the dynamically consistent framework provided by NEMO–PISCESv2_RC, accounting for the full suite of simulated food-web and physicochemical process rates (sources and sinks) and various physical transport modes.

• Dissect the mechanistic drivers of the biological particulate carbon pump across contrasting North Atlantic biogeochemical regimes and reveal their influence on the epipelagic export efficiency (EE) and mesopelagic transfer efficiency (TE) summary metrics.

2 Methods

2.1 Model configuration and simulation setup

We use the dynamical ocean model NEMO 4.0.4 (Nucleus for European Modelling of the Ocean; Madec, 2008), coupled with the sea-ice model SI³ (Vancoppenolle et al., 2009) and the biogeochemical model PISCESv2_RC (Pelagic-Interactions Scheme for Carbon and Ecosystem Studies; Aumont et al., 2015, 2017). Hereinafter, we will use the terms NEMO and PISCES to refer to NEMO 4.0.4 and the version of PISCESv2_RC released alongside it (Table 1). NEMO and PISCES are coupled through the Tracers in Ocean Paradigm (TOP) tracer manager (NEMO TOP Working Group, 2018), which performs advective and diffusive transports of both physical and biogeochemical tracers (see details in Supplement Sect. S1).

Table 1 Abbreviations used in the text and their definitions.

^* Biogeochemical model tracer and fraction of POC.

Global simulations were run on the ORCA1L75 grid (∼ 1° horizontal resolution, 75 non-linear vertical levels) with a 45-minute time step for ocean dynamics and biogeochemistry. The simulation followed the Ocean Model Intercomparison Project protocol (OMIP-2; Tsujino et al., 2018) and was forced by the JRA-55 atmospheric reanalysis (Kobayashi et al., 2015). After initialisation, biogeochemistry was spun up for 1302 years under recurring circulation and fixed preindustrial atmospheric CO₂. The subsequent transient simulation included rising atmospheric CO₂ concentrations over 1850–2019. Only results from the final 1958–2019 forcing period are analysed here, focusing on the North Atlantic. For comparison with satellite products, a monthly climatology over 1998–2019 was used. Further details are given in Supplement Sect. S1.

2.2 Biogeochemical model description

PISCES simulates marine planktonic food webs and the biogeochemical cycling of key elements (carbon, nitrogen, phosphorus, silicon, and iron) using 24 tracers that represent living, detrital and inorganic compartments. Here we briefly describe the representation of POC cycling in the model, introducing tracer names in italics.

POC is represented by six tracers: two phytoplankton types – miscellaneous small phytoplankton (phymisc) and diatoms (phydiat); two zooplankton types – microzooplankton (zmicro) and mesozooplankton (zmeso); and two detrital pools – small (sdetoc) and large (ldetoc) particles. Of these, only detrital tracers sink, with fixed sinking velocities of 2 and 50 m d⁻¹ for sdetoc and ldetoc, respectively. Phytoplankton growth depends on light availability and uptake of nitrate, ammonium, phosphate, dissolved iron and, for diatoms, silicate. Phytoplankton is lost through zooplankton grazing and linear and quadratic mortality terms, the latter representing aggregation upon bloom demise. Both zooplankton types prey on phytoplankton and sdetoc with different preferences. Mesozooplankton additionally predates on microzooplankton and consumes sdetoc and ldetoc via flux feeding.

Detritus is produced mainly through plankton mortality and zooplankton inefficient (“sloppy”) feeding, with partitioning between sdetoc and ldetoc depending on the organism and process involved. Detrital particles are degraded to DOC using the variable reactivity scheme of Aumont et al. (2017), which allows decay rate constants (k, d⁻¹) to vary in space and time. This scheme prescribes a maximum k of 0.035 d⁻¹ at 0 °C for freshly produced detritus in the upper mixed layer. With increasing depth, k decreases as sinking particles lose their most reactive components and become increasingly refractory. This decrease is not prescribed explicitly but emerges from the balance between particle sinking, the local supply of fresh detritus, and the progressive degradation of 15 lability classes implicitly represented by a gamma distribution. As a result, the vertical profile of k varies among regions, seasons, and detrital size classes. Detritus decay rates also increase as a power function of temperature, approximately doubling every 10 °C (Q10 of 1.9). Following PISCES terminology, degradation of sdetoc and ldetoc is referred to as remineralisation, although detrital POC is first converted to DOC before reaching dissolved inorganic carbon (DIC).

Mesozooplankton flux feeding and fragmentation are key features in PISCESv2 (Aumont et al., 2017). Given that the current parameterisation differs markedly from previous model versions (Gehlen et al., 2006), here we provide an updated description. At each time step, a variable fraction of mesozooplankton biomass engages in flux feeding. This fraction is diagnosed as the ratio of flux-feeding ingestion to total mesozooplankton ingestion. Flux feeding rates are calculated as the product of flux-feeding mesozooplankton biomass, a flux-feeding cross-section parameter, and the sinking flux of each detrital tracer (Supplement Sect. S1.2, Eq. S4). The flux-feeding parameter corresponds to a biomass-normalised area of 3 m² (mol C)⁻¹, roughly representative of temperate and subpolar flux feeders (Jackson, 1993; Stemmann et al., 2004b; Stukel et al., 2019).

Fragmentation of large detritus by mesozooplankton is expressed as a variable proportion of flux feeding rates, resulting in additional attenuation of ldetoc fluxes and transfer to sdetoc. The ratio between fragmentation and flux feeding increases, from 0.2 to around 2, as a function of the ratio of biogenic silica to ldetoc concentrations. This dependence is intended to represent the greater fragility of large detritus composed of diatom aggregates compared to more compact forms such as faecal pellets.

2.3 Detrital POC budgets

We diagnosed the full budgets of sdetoc and ldetoc by outputting all individual terms from their prognostic equations, including sources (production processes), sinks (consumption processes) and transports (Fig. 1). Transports include gravitational sinking and tridimensional advection and diffusion. Diffusion – equivalent to mixing – represents transport by subgrid scale ocean dynamics. Budgets were computed separately for the epipelagic and mesopelagic layer in three regions outlined in Fig. 2: the Subpolar North Atlantic (SPNA), the Transition Area (Trans_Area) and the Subtropical North Atlantic (STNA). The boundary between the two layers, Z_prod, is defined as the depth at which the annual mean vertical profile of primary production (PP) falls below 1 % of its maximum. For each region, Z_prod was fixed to ensure that budgets were computed over a consistent vertical domain over the seasonal cycle. The choice of a biome-specific Z_prod, rather than a fixed reference depth, avoids conflating productive-layer structure with mesopelagic transfer when comparing regions. For instance, a uniform reference depth of 100 m would underestimate export in the SPNA, where substantial net POC removal occurs above 100 m, and overestimate export in the STNA, where photosynthetic production extends below 100 m (Palevsky and Doney, 2018, 2021; Wilson et al., 2022; Walker and Palevsky, 2025).

For a given region and layer, the budget equation takes the generic form

$$\begin{array}{}\text{(1)}& {\displaystyle \frac{\mathit{\delta }\left[X\right]}{\mathit{\delta }t}}=\mathrm{Sources}-\mathrm{Sinks}-\mathrm{GS}-\mathrm{ADV}-\mathrm{DIFF}\end{array}$$

where Sources and Sinks represent vertically integrated rates averaged over each region. Transport terms, namely gravitational sinking (GS), advection (ADV) and diffusion (DIFF), are expressed as the divergence of the fluxes along the lateral and vertical directions. Divergence corresponds to the net difference between input and output fluxes through the domain boundaries, divided by the domain length in the corresponding direction. Hereafter, fluxes across a bounding surface are denoted by F followed by a subscript (e.g., F_GS, F_ADV, F_DIFF), whereas GS, DIFF and ADV denote the associated flux divergences. Subscripts Z and L indicate the vertical and lateral components, respectively (e.g., F_ADVZ, F_ADVL, DIFF_Z, DIFF_L).

We verified budget closure by comparing the simulated tracer tendency $\frac{\mathit{\delta }\left[X\right]}{\mathit{\delta }t}$ with the sum of all diagnosed source, sink and transport terms over each domain. Budgets were closed both over successive model time steps (within $\pm \mathrm{6}\times {\mathrm{10}}^{-\mathrm{4}}$ % error) and over monthly periods.

Below, we provide a qualitative representation of the detrital POC budgets (Eqs. 2 and 3). The full prognostic equations for sdetoc and ldetoc, explicitly showing their controlling parameters, are provided in Supplement Sect. S1.2 and Table S2.1, as these equations have been described previously by Aumont et al. (2015, 2017). Superscript numbers are used as labels to link terms to Table S2.1 and Fig. 1:

$$\begin{array}{}\text{(2)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{\mathit{\delta }\left[\mathrm{sdetoc}\right]}{\mathit{\delta }t}}& =({\mathrm{mortalities}}^{\mathrm{1}}+{\text{zmicro non assimilation}}^{\mathrm{2}}\\ & +{\text{bacterial breakdown of ldetoc}}^{\mathrm{3}}\\ & +{\text{DOC aggregation}}^{\mathrm{4}}\\ & +{\text{zmeso fragmentation of ldetoc}}^{\mathrm{5}})\\ & -({\text{zmicro ingestion}}^{\mathrm{6}}\\ & +{\text{zmeso ingestion and flux feeding}}^{\mathrm{7}}\\ & +{\text{remineralization}}^{\mathrm{8}}+{\text{sdetoc aggregation}}^{\mathrm{9}})\\ & -({\text{ADV}}_L^{\mathrm{10}}+{\text{ADV}}_Z^{\mathrm{11}}+{\text{DIFF}}_L^{\mathrm{12}}\\ & +{\text{DIFF}}_Z^{\mathrm{13}}+{\text{GS}}^{\mathrm{14}})\end{array}\text{(3)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{\mathit{\delta }\left[\text{ldetoc}\right]}{\mathit{\delta }t}}& =({\text{mortalities}}^{\mathrm{15}}+{\text{meso non assimilation}}^{\mathrm{16}}\\ & +{\text{sdetoc aggregation}}^{\mathrm{9}}\\ & +{\text{DOC aggregation}}^{\mathrm{17}})\\ & -({\text{zmeso flux feeding}}^{\mathrm{18}}\\ & +{\text{zmeso fragmentation of ldetoc}}^{\mathrm{5}}\\ & +{\text{remineralization}}^{\mathrm{19}}\\ & +{\text{bacterial breakdown of ldetoc}}^{\mathrm{3}})\\ & -({\text{ADV}}_L^{\mathrm{20}}+{\text{ADV}}_Z^{\mathrm{21}}+{\text{DIFF}}_L^{\mathrm{22}}\\ & +{\text{DIFF}}_Z^{\mathrm{23}}+{\text{GS}}^{\mathrm{24}})\end{array}\end{array}$$

Figure 1 Schematic of detrital POC dynamics and budgets in the PISCES model. Red arrows represent sources (production processes); blue arrows represent sinks (consumption processes). Purple arrows and tags represent the transport processes, which can supply or remove detritus from the budget. The top panel illustrates the six tracers that form POC in PISCES. Full budgets are analysed only for the detrital tracers (sdetoc and ldetoc), but all POC tracers are included in calculations of vertical non-gravitational fluxes and comparison to observational estimates.

For clarity, we lumped together some budget terms that are small and/or follow similar dynamics: (i) sdetoc consumption by phagotrophic ingestion and flux feeding; (ii) physical aggregation processes, and (iii) mortality terms. Finally, note that processes that interconvert between sdetoc and ldetoc appear in both budgets with opposite signs.

2.4 Export Efficiency (EE) and Transfer Efficiency (TE) metrics

We summarised spatial variations in vertical POC fluxes using three common metrics. Export efficiency (EE; Buesseler, 1998), also referred to as export ratio, is defined as the ratio between POC export flux through the bottom of the productive (epipelagic) layer and the net primary production rate integrated over that layer (NPP):

$$\begin{array}{}\text{(4)}& \mathrm{EE}={\displaystyle \frac{F_{Z_{\mathrm{prod}}}}{\mathrm{NPP}}}\end{array}$$

Transfer efficiency (TE; François et al., 2002) quantifies the fractional transfer of exported POC from the base of the productive layer to a deeper reference depth – in our case, the base of the mesopelagic layer (1000 m):

$$\begin{array}{}\text{(5)}& \mathrm{TE}={\displaystyle \frac{F_{Z\mathrm{1000}}}{F_{Z_{\mathrm{prod}}}}}\end{array}$$

The product of EE and TE is the particulate carbon pump efficiency (Buesseler and Boyd, 2009) at 1000 m:

$$\begin{array}{}\text{(6)}& \mathrm{EE}\times \mathrm{TE}={\displaystyle \frac{F_{Z_{\mathrm{prod}}}}{\mathrm{NPP}}}\times {\displaystyle \frac{F_{\mathrm{1000}}}{F_{Z_{\mathrm{prod}}}}}={\displaystyle \frac{F_{\mathrm{1000}}}{\mathrm{NPP}}}\end{array}$$

Vertical POC fluxes at the same two depths bounding the mesopelagic domain were also used to calculate Martin's b exponent (Martin et al., 1987), often referred to as the flux attenuation coefficient:

$$\begin{array}{}\text{(7)}& b=-{\displaystyle \frac{\ln (F_{Z\mathrm{1000}}/F_{Z_{\mathrm{prod}}})}{\ln (\mathrm{1000}/Z_{\mathrm{prod}})}}\end{array}$$

Stronger vertical attenuation of POC fluxes corresponds to lower TE and larger positive values of b. All metrics were calculated using annual-mean climatological fluxes for each region to avoid issues associated with non-steady-state seasonal dynamics (Giering et al., 2017) (see Sects. 3.4 and 4.3).

To assess differences in attenuation among POC pools, Eqs. (4)–(6) were applied both to detrital POC fluxes alone (gravitational, diffusive and advective fluxes of sdetoc and ldetoc) and to total POC fluxes, which additionally include the diffusive and advective vertical transports of plankton tracers.

2.5 Model evaluation

We evaluated model performance using gridded observational products for mixed layer depth (MLD), vertically integrated net primary production (NPP), phytoplankton carbon biomass, mesozooplankton carbon and small POC (sPOC) stocks, the latter assessed in both the epi- and mesopelagic layers. After regridding the observational products onto the model grid (ORCA1_L75), model skill was assessed using standard metrics applied to (i) spatial patterns of annual-mean fields (Fig. 2) and (ii) monthly climatological seasonal cycles averaged over each region (Fig. 3). Key characteristics and references for each observational dataset are summarised in Table 2. Further evaluation results, including skill metrics (Table S2.3) and seasonal maps comparing simulated and observed fields (Figs. S3.2–3.3), are presented in the Supplement Sect. S3.

While some simulated variables have a good correspondence in the observations – most notably MLD and NPP – other comparisons require simplifying assumptions, which are described here. First, modelled phytoplankton is evaluated using a satellite product that resolves three phytoplankton size classes. For this comparison, we assume that satellite microphytoplankton corresponds to model diatoms, and that the sum of pico- and nanophytoplankton corresponds to miscellaneous phytoplankton in the model. This group includes mostly haptophytes and cyanobacteria and, to some extent, dinoflagellates. Although this assumption introduces some uncertainty, it is broadly consistent with wide evidence from size-fractionated biomass measurements, ecophysiological studies, pigment signatures, and bio-optical data, all of which underpin satellite algorithms (Uitz et al., 2006, 2008) and model parameterisations. Second, we assume that the sPOC estimated from bio-optical data corresponds to the sum of four POC tracers with nominal sizes smaller than 100 µm: the two phytoplankton classes, microzooplankton, and small detritus, as previously discussed (Galí et al., 2022). Finally, mesozooplankton is evaluated using an observations-based machine learning product specifically developed for comparison with PISCES mesozooplankton fields (Clerc et al., 2024). No suitable gridded observational products were found for microzooplankton and large detrital POC concentrations.

In the case of sPOC, we used two estimates based on the same observational dataset. The approach of Galí et al. (2022), based on the variable ratio between POC and particulate backscattering, was used to estimate POC between 0–1000 m. For completeness, the method of Koestner et al. (2024), based on the relationship between POC, particulate backscattering and chlorophyll a, was used to estimate POC between the surface and the euphotic layer depth (below which gridded chlorophyll a data is not available). Hereafter, we refer to these sPOC estimates as G22 and K24, respectively. Additional details on these datasets are provided in Supplement Sect. S1.3.

Simulated export fluxes were evaluated separately using a subsampling approach designed to account for the sparse and uneven spatiotemporal coverage of the observations. Two complementary data compilations were used (Table 2): (i) Shallow POC export fluxes measured with the ²³⁴Th disequilibrium method, mostly at 100–150 m (Le Moigne et al., 2013); (ii) full water-column POC fluxes (Mouw et al., 2016), here restricted to sediment–trap estimates to avoid overlap with the ²³⁴Th dataset. Trap deployments longer than 1 month were discarded.

Table 2 Descriptions of the gridded observational dataset used for the model evaluation.

3 Results

In this section, an overview of model performance and simulated patterns (Sect. 3.1), is followed by more detailed descriptions of simulated POC dynamics, including vertical patterns in stocks, export and decay rates (Sect. 3.2), mesopelagic budgets of detrital POC at annual (Sect. 3.3) and monthly scales (Sect. 3.4), and POC export and transfer efficiency metrics (Sect. 3.5). The data provided in the following text and figures are available in two summary spreadsheets included in Supplement Sect. S4.

3.1 Model evaluation and overview of spatiotemporal patterns

3.1.1 Mixed Layer Depth (MLD)

Winter mixing depths decrease from the subpolar area (mean MLD of ∼ 240 m in March) to the Trans_Area (∼ 120 m in February) and the subtropical area (∼ 70 m in February). In the three regions, the spatially averaged monthly MLD from the simulation shows a strong temporal correlation with the observations (r> 0.95), low relative mean bias (RMB < 6 %) and realistic peak timing and amplitude (Fig. S3.1 and S3.2, Table S2.3).

Figure 2 The left column includes maps for the annual climatology of the five biogeochemical modelled variables evaluated in this study: NPP (a), surface carbon from phytoplankton biomass; diatoms (c) and miscellaneous phytoplankton (pico and nano-phytoplankton, e), carbon from mesozooplankton integrated over 0–200 m, (g) the small POC (particulate organic carbon, i) integrated over 0–200 m and (k) the small POC integrated over 200–1000 m. The right column shows the mean spatial bias of the simulated variables compared to observational datasets.

3.1.2 Net Primary Production (NPP)

The model reproduces the meridional NPP gradient, with the largest annual production in the Trans_Area (Fig. 2). In the SPNA and the Trans_Area, the seasonal cycles are well correlated with satellite estimates and capture the timing of the spring–summer production peak (Figs. 3a–b and S3). In the SPNA, however, NPP is systematically lower than satellite estimates year-round, with a −22 % mean deviation. The Trans_Area shows the smallest annual mean deviation (−7 %), which results from the compensation between moderate seasonal biases. In the STNA, the model underestimates NPP year-round and produces almost no seasonal variability, resulting in a weak temporal correlation with the satellite product (Fig. 3c).

3.1.3 Phytoplankton carbon biomass

Surface phytoplankton carbon is evaluated separately for diatoms (phydiat; Fig. 2c–d) and miscellaneous phytoplankton (phymisc; Fig. 2e–f). Modelled phydiat overestimates satellite microphytoplankton north of 40° N, whereas modelled phymisc underestimates satellite pico+nanophytoplankton in the SPNA and shows smaller deviations in the other regions. The opposite bias of the two simulated phytoplankton classes is also evident in seasonal patterns (Fig. 3d–f). As a result of these compensating biases, total phytoplankton presents −4 % relative bias year-round in the SPNA, 8 % in the Trans_Area and 13 % in STNA. Overall, the model reproduces the latitudinal shift from a strongly seasonal, diatom-rich regime in the SPNA to a weakly seasonal, small-phytoplankton-dominated regime in the STNA.

3.1.4 Mesozooplankton carbon biomass

Mesozooplankton carbon integrated over the upper 200 m displays annual mean patterns that resemble those of NPP and sPOC stocks (Fig. 2g–h). The relative mean bias is almost 30 % in the SPNA, about 90 % in the Trans_Area, and around 20 % in the STNA, similar to the results obtained by Clerc et al. (2024). Compared to the observational product, the model tends to exhibit a wider seasonal amplitude and, in the SPNA, a delayed peak, lagging diatoms by at least 1 month. In the observational product, mesozooplankton seasonality is better aligned with NPP, phytoplankton biomass and detrital stocks (Fig. 3a, d, g, j).

3.1.5 Small Particulate Organic Carbon (sPOC)

Annual mean sPOC stocks integrated over 0–200 m show similar large-scale patterns in the model and the G22 dataset (Fig. 2i–j), broadly resembling NPP patterns. The spatially averaged relative bias in the upper 200 m over the North Atlantic is −16 % (Fig. 2j), with smaller bias in the Trans_Area and stronger underestimation in the STNA (Fig. 3i). The 0–200 m sPOC seasonal cycles are closely phaselocked to NPP in the SPNA and Trans_Area and show a high temporal correlation with G22 (r=0.95 in SPNA and r = 0.97 in Trans_Area, both p < $\mathrm{1.0}\times {\mathrm{10}}^{-\mathrm{6}}$), although the model overestimates the seasonal amplitude in the SPNA (Fig. 3g–h).

Between 200 and 1000 m, model performance degrades. The relative mean bias in mesopelagic sPOC reaches ∼ 50 % in the Trans_Area and STNA, with widespread underestimation across most of the North Atlantic. An exception occurs in the western SPNA, which moderates the regional SPNA bias to -18 %. The SPNA also exhibits the widest seasonal amplitude in mesopelagic sPOC, which is reasonably captured by the model (Fig. 3j).

Comparison with vertical profiles (Fig. 4) complements this analysis and suggests that the model (i) underestimates the vertical gradient of sPOC in the SPNA epipelagic, (ii) overestimates the sPOC attenuation gradient in the upper mesopelagic of the Trans_Area, and (iii) underestimates the subsurface biomass maximum in the STNA, and simulates it too deep in the water column.

For completeness, modelled sPOC was also compared with estimates obtained using the algorithm of Koestner et al. (2024) in the epipelagic layer (Supplement Sect. S1.3). K24 sPOC exhibits seasonal patterns similar to G22 but is systematically higher by around 200 mmol C m⁻² (and up to ∼ 400 mmol C m⁻² in the SPNA winter, Fig. 3), implying a stronger underestimation of epipelagic sPOC by the model. Comparison between G22 and K24 did not indicate systematic overperformance of either product (Fig. 13 in Koestner et al., 2024). Accordingly, and because G22 was specifically designed to provide consistent estimates across epi- and mesopelagic, here we retain it as the primary reference, while using the comparison with K24 to illustrate the range of algorithm-dependent observational uncertainty. A more thorough evaluation of uncertainties in the conversion from particulate backscattering to sPOC, which is beyond the scope of our study, could further refine these findings.

Figure 3 Seasonal Cycle of the biogeochemical variables evaluated. Modelled variables are represented in solid lines, and the observations in dashed lines. Different colours represent different variables. Each column represents the three selected biomes: Subpolar North Atlantic (SPNA), Transitional Area (Trans_Area) and Subtropical North Atlantic (STNA). Each row shows a different seasonal cycle: NPP (a–c), Surface Carbon from phytoplankton (d–f, dark blue for diatoms and light blue for miscellaneous phytoplankton), Mesozooplankton (in pink) and Small Particulate Organic Carbon stock (sPOC- in yellow) over the 0–200 m depths (g–i) and sPOC stock over the 200–1000 m depths (in brown, j–l). Beige squares contain skill metrics model vs observations, and white squares include the mean and the standard deviation of the modelled variables.

Figure 4 Annual vertical profiles of different POC stocks (sPOC modelled and observed) and their components (first column). Vertical fluxes (second column) include gravitational sinking export fluxes only for sdetoc (F_GSsdetoc, in solid light purple) and ldetoc (F_GSldetoc, in solid dark purple), and advective (F_ADVZ), in dashed line style) and diffusive (F_DIFFZ, in dotted line style) vertical fluxes for detritus (sdetoc in light purple, and ldetoc in dark purple) and phyto-zooplankton (pkt in light green). The insets at the lower left part of the subplots in the second column contain the diffusive and advective vertical fluxes of sdetoc and ldetoc, and the advective fluxes of pkt, since their magnitudes are not comparable to the detrital gravitational sinking export fluxes and the vertical diffusive flux of pkt. Remineralisation of detritus rate (remsdetoc and remldetoc) and Primary production (PP) (third column). Green rectangles in the right column represent the productive layer (defined in Sect. 2.4); numbers at their bottom indicate the vertically integrated annual mean primary production (NPP) in that layer. Rows correspond to the SPNA (a, b, c), Trans_Area (d, e, f) and STNA (g, h, i) study regions.

3.1.6 Export fluxes at 100 and 150 m depth

To evaluate modelled POC export, we compared the monthly climatology of total vertical fluxes with ²³⁴Th- and sediment trap-based estimates from the global compilations of Le Moigne (2013) and Mouw et al. (2016a). Based on the spatial distribution of the measurements, we focus on the entire subpolar North Atlantic (SPNA) and two specific study sites (Fig. 5a–b). The Porcupine Abyssal Plain (PAP) lies near the boundary between the SPNA and Trans_Area regions in the eastern Atlantic. The other study site is located near the Trans_Area–STNA boundary in the Sargasso Sea, and includes the merged datasets from the Bermuda Atlantic Time Series and the Oceanic Flux Program (BATS/OFP). In SPNA and PAP, measurements are available only from May through August in scattered years, whereas BATS/OFP measurements are available year-round between 1985 and 2013.

In the SPNA (Figs. 5c, 6), the model captures the observed order of magnitude of the observations, which are nevertheless extremely variable at both 100 and 150 m. This variability, combined with spatial and temporal sparseness, precludes drawing clear conclusions on model behaviour. At the PAP site, modelled fluxes exceed observations by about ∼ 70 % at 100 m but are in good agreement with measurements at 150 m (bias = −21 %, Fig. 5d). At BATS/OFP, where observations are available only at 150 m, the model reproduces the weak seasonality but shows a mean positive bias of ∼ 90 %, largely driven by winter-spring months (Fig. 5e).

Figure 5 Modelled annual mean (1998–2019) vertical sinking fluxes (gravitational, diffusive and advective) at 100 m (a) and 150 m (b) are compared with in situ measurement compilations. Circles indicate ²³⁴Th-based POC export flux, and triangles, sediment trap-based POC export flux. Modelled monthly climatologies (1998–2019) of the vertical sinking fluxes are extracted for SPNA (c), PAP (d) and BATS/OFP study sites (e) and compared with available observations. Light purple lines in (c)–(e) represent export fluxes at 100 m, and dark purple lines represent export at 150 m. Solid lines indicate modelled fluxes. Dashed lines and ribbons represent the monthly mean and standard deviation of the observations.

3.1.7 Water column export fluxes

To complete the export flux evaluation, we subsampled model output at the depth, year, month and location of data from Mouw et al. (2016a). We used data from the upper 4000 m in the North Atlantic, of which 61 % was obtained at BATS/OFP, and computed model skill metrics for both non-binned data (n=3388) and data binned by layer and region (n=12). These comparisons and their corresponding metrics are shown in Fig. 6, suggesting that the model captures water-column export fluxes with small overall bias across two orders of magnitude when the data are binned. Indeed, non-binned data show substantial scatter, especially for small flux values. Crucially, different conclusions are obtained depending on whether in situ measurements are matched to modelled sinking POC fluxes or to total downward fluxes (sinking, diffusion and advection). In the first case, the comparison suggests an overall tendency for model underestimation, whereas the second case suggests a slight tendency for overestimation in the shallowest depth bin. Data scarcity hampers assessment in the epipelagic SPNA and the upper 500 m of the Trans_Area.

Figure 6 Comparison between modelled export fluxes and a compilation of in situ sediment trap measurements from Mouw et al. (2016a). Collocated data are compared directly (non-binned dataset) and after computing regional and vertical bin averages. The corresponding skill metrics are shown, along with a linear regression on binned data (purple line and envelope showing SE of the fit). Diagonal lines with different model:observation ratios are shown as visual guides. The same in situ fluxes are compared to modelled POC gravitational fluxes (left) and total fluxes (right). For this comparison, we assigned data from BATS/OFP to the STNA and data from PAP to the Trans_Area (see text and Fig. 5). Vertical boundaries were chosen to optimise data abundance in each layer while deviating minimally from the standard definitions of Z_prod and the mesopelagic layer.

Table 3 The upper part qualitatively summarises observed model biases or misfits. Arrows are used to indicate positive (↑) or negative (↓) model biases of large (2) or moderate (1) magnitude; indicates small bias, ∼ indicates temporal phase shift, and “?” indicates uncertain evaluation outcome usually due to data sparseness or too-large observational uncertainties. The lower part provides interpretations based on the comparison of various bias patterns. More speculative interpretations are highlighted in italics and are based on insights from the modelled mesopelagic budgets.

3.2 Simulated vertical profiles of POC pools and detrital POC decay

Evaluation results confirm that PISCES captures the main biogeochemical features of the North Atlantic basin, despite some season- and region-dependent biases (Table 3). To investigate how these patterns translate into POC vertical structure, this section documents the distribution of POC components and their corresponding fluxes and remineralisation processes from the epipelagic to the mesopelagic layer.

On an annual basis, total gravitational export at the base of the productive layer (Z_prod) is similar in the transition and subpolar areas (3.9–4.3 mmol C m⁻² d⁻¹), and much lower in the subtropical area (1.3 mmol C m⁻² d⁻¹) (Fig. 4). In addition, a substantial amount of POC is exported through turbulent diffusion (vertical mixing) of plankton in the SPNA (1.4 mmol C m⁻² d⁻¹) and the Trans_Area (0.8 mmol C m⁻² d⁻¹). Smaller exports arise from the diffusion of detrital POC, especially sdetoc, while vertical advection of particles is negligible at the regional scale examined here.

Vertical flux profiles vary widely across regions and POC tracers (Fig. 4). Gravitational sdetoc fluxes attenuate rapidly in the upper mesopelagic and more slowly below ∼ 400 m. By contrast, ldetoc shows deeper flux maxima, greater penetration through the upper mesopelagic, and gradual attenuation throughout the mesopelagic. Consequently, the proportion of export fluxes carried by sdetoc increases with depth in the mesopelagic. Diffusive plankton fluxes attenuate at shallower depths than gravitational ones. Yet, their penetration deepens towards the SPNA, consistent with deeper mixing in winter (Figs. S1 and S2).

Small detritus is the most abundant POC fraction in the mesopelagic layer, where it accounts for 52 % (SPNA) to 70 % (STNA) of the POC pool. The ldetoc fraction usually is ≤ 10 % of total POC and smaller than total plankton (Fig. 4a, d, g).

Vertical profiles of detritus decay rates also vary by detritus type and region (Fig. 4c, f and i). The integrated decay rate of detritus (sdetoc + ldetoc) between the surface and 1000 m is highest in the Trans_Area (13.9 mmol C m⁻² d⁻¹) compared to the other areas (8.4–8.6 mmol C m⁻² d⁻¹), mirroring NPP patterns. Small detritus drive most of the detritus decay, with higher sdetoc contributions in the epipelagic (85 %–92 %) than in the mesopelagic (53 %–61 %). To characterise distinct regimes, we calculated the percentage of POC degradation occurring in the epipelagic layer relative to the 0–1000 m integral. In SPNA, 58 % (sdetoc) and 28 % (ldetoc) of decay occurs in the epipelagic, concurrent with low temperatures, strong ldetoc export and shallow Z_prod. In contrast, the STNA shows 93 % (sdetoc) and 56 % (ldetoc) epipelagic decay, aligned with warm temperatures, reduced ldetoc export, and deep Z_prod.

Detritus decay rates result from the product of tracer concentration and the degradation rate constant k, which varies with depth and temperature (Fig. 7; see also Sect. 2.2). For sdetoc, 0 °C-normalized k typically ranges 0.010–0.023 d⁻¹ in the epipelagic and decreases sharply to ∼ 0.002 d⁻¹ at 1000 m; for ldetoc, k typically ranges 0.027–0.035 d⁻¹ in the epipelagic and 0.015–0.020 d⁻¹ at 1000 m. Thus, in PISCES, sdetoc becomes refractory at shallower depths than ldetoc, explaining the limited change in sdetoc concentration and sinking flux below ∼ 400 m compared to ldetoc (Fig. 4).

While the change in reactivity is best described by 0 °C-normalized k, the actual in-situ rates are differently affected by temperature across regions and depths (Fig. 7), which here we quantify using the enhancement factor with respect to k at 0 °C. In the epipelagic, mean thermal enhancement factors are 1.7, 3.2 and 4.5 in the SPNA, Trans_Area and STNA, respectively. Smaller thermal enhancement factors (1.5–2.5) are found in the mesopelagic. Thus, temperature control on k exacerbates the vertical gradients in reactivity at low latitudes.

Figure 7 Vertical profiles of the specific degradation rates of small and large detritus (sdetoc and ldetoc, respectively) at 0 degrees (k₀) and at in situ temperature (k_T) for the three regions.

3.3 Mesopelagic detrital POC budgets: annual-mean fluxes

Regional contrasts in POC vertical structure reflect distinct combinations of particle composition, export pathways, and remineralisation regimes (Figs. 4 and 7). To quantify how these differences arise in PISCES, we now examine the annual climatological budgets of detrital POC in the three regions of interest. We place particular emphasis on how regional differences in productive-layer depth, particle size structure, and export pathways shape mesopelagic budgets (Fig. 8).

For sdetoc, the relative contribution of gravitational inputs increases toward lower latitudes (Fig. 8a), from 36 % (SPNA) to 74 % (STNA). In contrast, the supply of sdetoc by the mesopelagic food web shows the opposite pattern and exceeds 50 % of total supply in the SPNA, driven primarily by mesozooplankton fragmentation of ldetoc (16 %), plankton mortalities (20 %), and microzooplankton non-assimilation (9 %). In the Trans_Area, plankton mortalities dominate biological supply (14 %), followed by non-assimilation (8 %), bacterial breakdown of ldetoc (7 %), and mesozooplankton fragmentation (6 %). In the STNA, bacterial breakdown of ldetoc (7 %) and abiotic DOC aggregation (12 %) contribute relatively more, whereas contributions from mortalities and microzooplankton non-assimilation are reduced. Removal of sdetoc is dominated by remineralisation in all regions (> 70 %), reaching 87 % in the STNA. Mesozooplankton removal via phagotrophic ingestion and flux feeding is also significant, peaking at 17 % in the Trans_Area, whereas microzooplankton removal remains minor (≤ 5 %).

For ldetoc, gravitational inputs dominate in all regions, with smaller contributions from plankton mortalities and mesozooplankton non-assimilation (Fig. 8b). As for sdetoc, the importance of these biological sources diminishes southward in favour of sinking. Removal of ldetoc is dominated by remineralisation, but at lower fractions (50 %–75 %) than sdetoc. In the SPNA, nearly half of ldetoc removal occurs through mesozooplankton flux feeding (30 %) and fragmentation into sdetoc (15 %). Mesozooplankton-mediated transformation rates decrease toward lower latitudes, especially fragmentation, which in the Trans_Area contributes half as much to sdetoc budgets as in the SPNA. Bacterial disaggregation of ldetoc is generally minor, reaching 7 % in the STNA.

Overall, both sdetoc- and ldetoc-specific processes jointly regulate the magnitude and evolution of regional detrital POC stocks in PISCES. Supply pathways are more heterogeneous for sdetoc, whereas removal pathways are more diverse for ldetoc. This asymmetry largely reflects ldetoc transformation by fragmentation and disaggregation processes, which act as sinks for ldetoc while simultaneously supplying sdetoc. The indirect contribution of plankton diffusion (mixing) also affects sdetoc differentially across regions (Fig. 4 b, e, i).

Figure 8 Climatological annual budget of small (a) and large (b) detrital POC over the mesopelagic layer (Z_prod−1000 m). To compare the contribution of each process across regions, the budget terms have been divided by the sum of the supply processes. Numbers on the right of each bar indicate the sum of the inputs (positive) and outputs (negative) for each detrital POC tracer. Vertical lines represent the temporal standard deviation of total monthly inputs and outputs after normalisation by the annual mean value.

3.4 Mesopelagic detrital POC budgets: seasonal cycle

Building on the annual-mean budgets presented in Sect. 3.3, we now examine how the processes controlling mesopelagic detrital POC vary over the seasonal cycle in PISCES. Monthly climatological budgets show that, consistent with the seasonality of physical forcing, temporal variability increases with latitude: it is minimal in the STNA but pronounced in the SPNA (Figs. 9, S5, S6 and vertical black lines in Fig. 8).

Gravitational sinking of large detritus during the high-latitude bloom and post-bloom period is the main seasonal event. The peak in gravitational ldetoc supply (June) is rapidly followed by maxima in ldetoc losses through remineralisation, mesozooplankton flux feeding and fragmentation. Enhanced mesozooplankton activity also generates secondary ldetoc sources via mortality and sloppy feeding. In contrast, sdetoc budgets exhibit a smoother seasonal succession, with sources shifting from mortalities and zooplankton non-assimilation in spring to ldetoc fragmentation in summer. Alongside dominant removal by remineralisation, sdetoc sinks transition gradually from microzooplankton ingestion to mesozooplankton ingestion and flux feeding through spring and summer.

Seasonal decomposition is also necessary to understand the role of vertical diffusive (mixing) fluxes and their contribution to detrital food webs. On an annual basis, vertical diffusion represents a significant net source of mesopelagic POC in the SPNA and the Trans_Area, whereas in the STNA, diffusive fluxes primarily redistribute POC within the epipelagic layer. Diffusive fluxes of POC are, however, largely driven by plankton and are therefore not included directly in detrital budgets. Instead, vertically mixed plankton enters detrital budgets through mortalities and sloppy feeding (Eqs. 1 and 2). In the SPNA, plankton diffusion supplies POC to the mesopelagic throughout the year and peaks in April (5 mmol m⁻² d⁻¹), when it is three times larger than total gravitational inputs (Supplement Sect. S4.2 file S4.2). Between December and April, when the MLD exceeds Z_prod, the ratio between [mortalities + non-assimilation] and plankton diffusion is 0.70 ± 0.27. Between May and October, this ratio increases to 7.0 ± 5.1. Thus, detrainment of plankton from the epipelagic layer dominates detritus supply (largely as sdetoc) before the spring bloom, whereas biological processing of sdetoc derived from large detritus fragmentation dominates during the post-bloom season. Similar, though weaker, patterns of diffusive supply occur in the Trans_Area.

Direct diffusion of detritus is comparatively small, with the largest contribution occurring in the Trans_Area (7 %). In the SPNA, diffusive fluxes of detritus switch sign over the seasonal cycle (Fig. 9e), resulting in a weak net annual supply (Fig. 8a). In late autumn and early winter, mixed-layer deepening entrains mesopelagic detritus into the epipelagic layer, whereas in spring, re-stratification promotes detrainment of detrital POC from the epipelagic to the mesopelagic layer.

Figure 9 Seasonal cycles of sdetoc (left) and ldetoc (right) budgets in the Subpolar North Atlantic (SPNA). Total tracer trends (a, b), sources (c, d), transports (e, f) and sinks (g, h) are shown separately. Note the different y-axis scales.

3.5 Export and transfer efficiency metrics

To synthesise the behaviour of the particulate C pump, we compare export and transfer efficiency metrics across regions (Fig. 10). Additionally, attenuation coefficients “b” (Martin et al., 1987) between Z_prod and 1000 m are reported in Fig. 11. All metrics are based on annual mean fluxes.

Considering gravitational fluxes only, PISCES-simulated export efficiency (EE) decreases with latitude, from 21 % in the SPNA to 8 % in the STNA. In contrast, TE varies non-monotonically with latitude and is lowest in the Trans_Area (13 %), intermediate in the SPNA (17 %), and highest in the subtropical area (23 %). As a result, overall particulate C pump efficiency at 1000 m, expressed as EExTE, is highest in the subpolar area (3.6 %) and lower in the transition (1.5 %) and subtropical (1.8 %) areas.

Given the substantial contribution of non-gravitational fluxes in NEMO–PISCES, especially plankton diffusion (Fig. 4), we re-computed EE and TE for the total POC vertical flux, including gravitational, diffusive and advective components. In the subpolar and transition areas, diffusive plankton fluxes are comparable in magnitude to gravitational ones at Z_prod but become negligible at 1000 m (Fig. 4). Therefore, their inclusion increases EE but decreases TE. In the SPNA, where the greatest non-gravitational fluxes occur, accounting for all vertical transport pathways increases EE by 8 % but reduces TE by 3 %, resulting in a modest increase in overall efficiency (EExTE) from 3.6 % to 4.0 %.

Figure 10 Export Efficiency vs Transfer Efficiency metrics in the three regions of the North Atlantic. Empty and filled symbols indicate, respectively, if calculations account for only detrital POC or total POC (detritus and plankton). Note that the two symbols overlap in the STNA. Symbol size is proportional to the NPP of each region. The isolines represent overall particulate C pump efficiency at 1000 m calculated as EExTE. See Sect. 2.4 for details.

4 Discussion

This study provides a comprehensive analysis of modelled POC dynamics in the upper 1000 m of the North Atlantic, linking regional and seasonal changes in primary production, POC stocks and export fluxes to the underlying mechanisms of detrital POC production, transport, transformation and decay. By comparing the simulation to available observations, our results shed light on how interactions between ecosystem structure and the physical environment together regulate epipelagic export and mesopelagic transformation of POC.

At the basin scale, and despite significant regional biases in some variables (Figs. 2–6), PISCES-simulated productivity and gravitational export fluxes are broadly consistent with observations and previous studies. As in satellite-based assessments (Henson et al., 2012; Siegel et al., 2023), maximum export fluxes (Fig. 5a, b) are displaced northward relative to NPP (Fig. 2a), reflecting regional differences in ecosystem structure and vertical mixing. Export efficiency increases toward subpolar waters (Fig. 10), concurrent with higher diatom abundance (Mouw et al., 2016) and an increasing contribution of fast-sinking large detritus (Fig. 11). Conversely, at lower latitudes, dominance of small organisms and detritus, together with warmer temperatures and deeper productive layers, enhances degradation within the epipelagic and limits export efficiency (Figs. 4d, 7c, 10). These large-scale patterns are consistent with early conceptual frameworks linking new production, community structure, and gravitational export (Eppley and Peterson, 1979), and ultimately reflect the joint influence of vertical mixing regimes (Sverdrup, 1953; Margalef, 1978) and temperature-dependent remineralisation (Cael and Follows, 2016).

Building on this context, the following sections examine in detail the biophysical controls on epipelagic export fluxes (Sect. 4.1), the interplay between zooplanktonic and bacterial cycling of mesopelagic POC (Sect. 4.2), and how modelled processes affect POC export and transfer efficiency metrics and relate to real-world observational constraints (Sect. 4.3). Since most POC fluxes and derived metrics cannot be comprehensively evaluated against gridded observational products, we compare them with sparse available observations whenever possible. To close the Discussion, we assess the limitations of the current PISCES model formulation, leveraging observed biases to outline potential pathways for improving the representation of mesopelagic POC cycling in models (Sect. 4.4).

4.1 POC inputs to the mesopelagic layer: POC sinking and vertical mixing

Evaluation against ²³⁴Th- and sediment trap–based estimates indicates that PISCES reproduces the order of magnitude and broad seasonality of epipelagic export fluxes in the North Atlantic, while exhibiting regionally distinct biases. In the SPNA and at PAP, model–data agreement generally improves with depth and toward summer. At BATS/OFP the model captures the weak seasonal cycle but overestimates export at 150 m, particularly in winter–spring. According to PISCES, large detritus drives 67 % of the gravitational flux at 100 m in the SPNA, exceeding the range (37 %–64 %) estimated in the western subpolar gyre by Wang and Fennel (2022). The latter study further suggests that the contribution of large detritus increases with depth, opposite to model estimates (Fig. 4). In the eastern subtropical Atlantic at 260 m, the contribution of slow-sinking (likely small) particles was > 60 % during most of the year (Alonso-González et al., 2010), roughly consistent with ∼ 50 % in PISCES at 200 m (Fig. 11). Model-data discrepancies possibly reflect a combination of observational uncertainties (Buesseler et al., 2007; Bishop et al., 2012) and limitations in how the model represents POC export using rigid relationships between size classes and sinking velocities (Jackson and Burd, 2015; Iversen and Lampitt, 2020; Cael et al., 2021)

Vertical gradients in simulated export are usually strongest near Z_prod (Palevsky and Doney, 2018; Buesseler et al., 2020) implying intense flux attenuation in the lower epipelagic and upper mesopelagic layers. Thus, uncertainties in the mathematical representation of export and attenuation processes, together with small spatiotemporal mismatches between simulations and observations, can lead to substantial model-data deviations. Importantly, it remains unclear whether observational export estimates from sediment traps and the ²³⁴Th technique isolate gravitational sinking alone (Fig. 5) or can also account for downward transport by diffusion and advection (and, if so, how efficiently). We therefore interpret shallow export comparisons as quantitative consistency checks rather than direct validation (Aumont et al., 2017), and analyse mesopelagic POC inputs in terms of the total simulated downward POC flux, integrating gravitational sinking of detritus with diffusive and advective transport of both detrital and living particles.

Compared to gravitational fluxes, the transport of POC by vertical mixing has received little attention until recent years. Classical formulations of the biological pump have assumed vertical particle diffusion to be negligible at large scales (Martin et al., 1987), but oscillations in mixed layer depth can effectively entrain and detrain particles (Gardner et al., 1995). In models, these transports are diagnosed as diffusive fluxes, which our results show to be quantitatively relevant at large scales, in agreement with recent research (Bellacicco et al., 2025). In the NEMO-PISCES simulation (Fig. 4), POC diffusion amounts to 37 % (25 %) of annual gravitational fluxes in the SPNA (Trans_Area), in line with independent estimates. For instance, Dall'Olmo et al. (2016) inferred that the POC mixed-layer pump could represent ∼ 23 % of annual mean gravitational fluxes at high latitudes. Lacour et al. (2019) used biogeochemical Argo float profiles to quantify springtime net POC detrainment at 4.6 mmol C m⁻² d⁻¹ (55 mg C m⁻² d⁻¹) in the SPNA. Our corresponding estimate for March–May (4.3±1.1 mmol C m⁻² d⁻¹) is strikingly similar, despite relying on a fundamentally different approach.

Crucially, PISCES further suggests that over 85 % of this diffusive POC flux is associated with living plankton rather than detritus, reflecting the sharper vertical gradient of plankton biomass and their enhanced transport during convection events (Galí et al., 2026). Altogether, gravitational and diffusive export fluxes exhibit different seasonality and deliver to the mesopelagic layer distinct particle mixtures in terms of size, sinking speed, and composition. Differences in supplied POC entail different transformation pathways in the mesopelagic zone (Figs. 7–9 and 11), which are examined in the next section.

Figure 11 Modelled annual budgets of detrital POC in the top km of the North Atlantic, focusing on mesopelagic transformations of small and large detritus. Epipelagic fluxes are shown more schematically, distinguishing the fraction of primary production (NPP) that flows through detrital particles and, within it, the fraction that results in gravitational (filled violet arrows) and diffusive (dotted purple lines) inputs to the mesopelagic layer. Additional particles enter the mesopelagic detrital pool via plankton vertical diffusion (dotted green lines) and subsequent conversion to detritus via mortality and zooplankton processing. Within the mesopelagic layer, blue and red arrows are used to distinguish sources and sinks, respectively, and violet arrows represent gravitational fluxes, whose narrowing represents flux attenuation. Arrow size is indicative of flux magnitude, but not to scale. All fluxes are in mmol C m⁻² d⁻¹, and, for simplicity, fluxes smaller than 0.05 have been omitted. The vertical attenuation coefficient “b” (Martin et al., 1987), computed between Z_prod and 1000 m (boxed italic numbers, unitless), is shown separately for gravitational fluxes of small and large detritus and the gravitational + diffusive fluxes together. The mesopelagic vertically-integrated stocks, in mmol C m⁻², are shown in bold italics. The modelled percentage of diatoms is indicated in green at the top of each panel: 42 % for SPNA, 18 % for Trans_Area and 7 % in the STNA.

4.2 Mesopelagic POC budgets: joint regulation by zooplankton and bacteria

Zooplankton plays a central role in detritus processing in PISCES: 50 %–60 % of the gravitational input flux transits through zooplankton at mid and high latitudes, compared to ∼ 20 % in the subtropical area (Figs. 8 and 11). Importantly, zooplankton processing does not equate to POC flux attenuation, as a substantial and variable fraction returns to the sinking detritus pool via fragmentation, faecal pellet production and mortality (Figs. 8 and 9). As a result, POC recirculates within a mesopelagic detrital loop before final removal, primarily through bacterial decomposition (Fig. 11). In this section, we examine how zooplankton transformations and microbial degradation jointly shape mesopelagic detritus budgets, and how these simulated rates compare with the limited available observations.

A stringent evaluation is provided by comparison with observations at the PAP site (see Fig. 2), where Giering et al. (2014) (G14) were able to constrain, within reasonable uncertainty, POC flux attenuation and mesopelagic metabolism. Using July–August PISCES means for the SPNA and the Trans_Area and acknowledging imperfect correspondence between model- and observations-derived rates, we find that PISCES reproduces within ±15 % the magnitude of gravitational inputs, total zooplankton detritivory and bacterial degradation (Table 4). However, the model overestimates large-detritus fragmentation and the contribution of small detritus to POC decay. This discrepancy may partly reflect differences in pool definitions, as suspended POC in G14 is not fully equivalent to PISCES sdetoc, which accounts for suspended and slow-sinking fractions (Galí et al., 2022). Overall, this comparison is reassuring regarding simulated flux magnitudes, and further work is warranted to explore POC dynamics at the PAP site.

In PISCES, flux feeding and fragmentation of large detritus by mesozooplankton critically influence the vertical transfer of gravitational fluxes (Fig. 8), consistent with previous assessments (Mayor et al., 2020). Seasonal export pulses enhance mesozooplankton biomass and the proportion of flux feeders (Sect. 2.2; Supplement file 1), progressively intensifying removal of large sinking aggregates by mesozooplankton as the season advances (Fig. 9). Fragmentation of ldetoc into sdetoc additionally slows sinking fluxes, favouring microzooplankton consumption and microbial degradation over vertical transfer. Together, these trophic feedbacks limit mesopelagic TE during peak export periods (Fig. 9).

Although the functional representation of mesozooplankton activities in PISCES reasonably captures the magnitude of some observed process rates (Table 4), it necessarily simplifies the diversity of mesopelagic zooplankton communities and feeding strategies. In particular, the model cannot explicitly resolve the wide range of behaviours spanning flux-feeding Rhizaria to particle-attached copepods and other taxa that interact with sinking particles in distinct ways (Stukel et al., 2019; Mayor et al., 2020; Lampitt et al., 2023; Laget et al., 2024). As a result, the diagnosed trophic feedbacks should be interpreted as an emergent, functional representation of zooplankton–particle interactions rather than a taxonomically explicit description of mesopelagic food-web dynamics.

Large-detritus fragmentation deserves particular attention, as PISCES was the only CMIP6 model that explicitly represented this process (Henson et al., 2022). Because the model prescribes enhanced fragmentation of silica-rich aggregates, fragmentation is maximal in the SPNA (Figs. 8 and 9; see also Fig. 3d). From May to July, zooplanktonic fragmentation accounts for 19 % of ldetoc flux attenuation in the SPNA, with bacterial disaggregation contributing an additional 5 %. By comparison, Briggs et al. (2020) attributed to fragmentation 49 ± 22 % of fast-sinking particle flux attenuation during intense export events in the SPNA and the Southern Ocean. This discrepancy likely reflects differences in process definitions and temporal averaging. In particular, the attenuation attributed to fragmentation by Briggs et al. (2020) better matches the 55 % of ldetoc attenuation explained by the sum of flux feeding, fragmentation and bacterial disaggregation in PISCES. Further observational (Briggs et al., 2020; Wang and Fennel, 2022) and modelling efforts are needed to constrain fragmentation rates and mechanisms.

Despite the prominent role of zooplankton, most flux attenuation in PISCES proceeds through microbial degradation (Figs. 8, 9 and 11). Owing to the variable reactivity formulation, slow-sinking particles are degraded at shallower depths than fast-sinking ones. Simulated remineralisation rate constants (Fig. 7) fall within the wide range of in situ estimates, although natural and methodological variability preclude a conclusive assessment (Benner and Amon, 2015 and references therein; Belcher et al., 2016 and references therein; García-Martín et al., 2021). A more robust assessment is obtained from vertically-integrated POC degradation rates. In the mesopelagic SPNA, decay rates in May–June (7 mmol C m⁻² d⁻¹) resemble the GEOVIDE cruise mean (6.5 mmol C m⁻² d⁻¹; Lemaitre et al. (2018), and summer estimates at PAP also agree within 10 % (Table 4). Unfortunately, we are not aware of suitable datasets in the STNA, where available estimates reflect the respiration of POC+DOC (e.g., Ono et al., 2001) and therefore exceed POC decay simulated by PISCES.

The contribution of different particle classes to community metabolism remains poorly constrained in both epipelagic (García-Martín et al., 2021) and mesopelagic waters (Table 4). In PISCES, fast-sinking aggregates retain higher intrinsic reactivity (Fig. 7), whereas the larger stock of small detritus compensates for its lower reactivity. As a result, both particle types support comparable fractions of mesopelagic metabolism across regions (Fig. 11). Perhaps counterintuitively, the strongest sdetoc removal, both in absolute (2.6 mmol C m⁻² d⁻¹) and relative terms (61 %), occurs in the subpolar area, fuelled by particle diffusion and ldetoc fragmentation inputs. Overall, these results reinforce the emerging view that suspended and slow-sinking particles play a central role in the mesopelagic carbon cycling (Alonso-González et al., 2010; Baker et al., 2017; García-Martín et al., 2021; Wang and Fennel, 2022; Baumas and Bizic, 2024).

Table 4 Mesopelagic POC budgets: PISCES vs. data from Giering et al. (2014) at Porcupine Abyssal Plain (PAP, Jul-Aug). Data extracted from their Fig. 2 and converted to mmol C m⁻² d⁻¹. These budgets are not closed: some terms are omitted due to limited comparability, and a seasonal imbalance exists (at least in the model simulation).

4.3 Carbon pump efficiency metrics: from models to real-world measurements

In PISCES, annual EE increases from 8 % to 29 % between subtropical and subpolar waters (Fig. 10). As found by Mouw et al. (2016), this gradient mirrors the latitudinal increase in diatom abundance (Fig. 2), which in the model rises from 7 % to 42 % (Fig. 10). Although the modelled latitudinal EE gradient falls within the range of current estimates (Villa-Alfageme et al., 2016; Siegel et al., 2023), EE is likely overestimated at mid and high latitudes because diatoms are overrepresented (Figs. 2 and 3; Table 3). Improved representation of phytoplankton functional groups and size structure in models is key to better constrain biological carbon export.

Our results further underscore the importance of including non-gravitational POC fluxes in carbon pump assessments. Neglecting these fluxes, especially those driven by vertical mixing, can substantially underestimate total POC export at high latitudes (Figs. 4 and 11). These unaccounted fluxes bias EE and TE metrics (Fig. 10) and affect interpretations of mesopelagic food-web metabolism and its seasonality (Figs. 9, S5 and S6) (Dall'Olmo et al., 2016; Lacour et al., 2019). Current model intercomparison frameworks (e.g. CMIP) rely on simulated gravitational fluxes at standard depths (Palevsky and Doney, 2018; Henson et al., 2022; Wilson et al., 2022). Thus, our findings suggest that non-gravitational POC transports, which are internally simulated by models, should be included in standard model intercomparison metrics.

Compared to the mechanisms controlling EE (see Sect. 4.1), the factors regulating mesopelagic TE remain more controversial. Mechanistic interpretations have alternately emphasised lower lability of exported POC in subtropical areas, leading to a latitudinal decrease in TE (Henson et al., 2012; Guidi et al., 2015), or reduced organic matter decay in cold mesopelagic waters, leading to a latitudinal increase in TE (DeVries and Weber, 2017). In PISCES, mesopelagic TE does not vary monotonically from subtropical to subpolar waters, nor with the depth of the productive layer, but instead exhibits a minimum at midlatitudes (Fig. 10). Lowest TE coincides with highest productivity (Figs. 10 and 11), broadly consistent with the climatological patterns reported by Mouw et al. (2016). This relationship also resembles the opposite interannual trends in NPP and export relative to mesopelagic transfer, reported by Lomas et al. (2010) in the Sargasso Sea. Nevertheless, the occurrence of maximum TE in the STNA and intermediate TE in the SPNA indicates that mesopelagic transfer is not controlled by productivity or productive-layer depth alone. As we shall see, these patterns cannot be attributed to a single dominant factor.

Examination of the mechanisms driving mesopelagic TE in PISCES is informative because the model represents the interplay between variable POC reactivity and temperature-dependent decay rates (Fig. 7), as well as zooplankton detritivory. First, exported POC – especially sdetoc – is least reactive in the STNA, reflecting deeper and warmer Z_prod (Fig. 4i), consistent with Henson et al. (2012) and Guidi et al. (2015). Second, mesopelagic POC decay rate constants also increase toward low latitudes, primarily due to temperature effects, consistent with Marsay et al. (2015) and DeVries and Weber (2017). Third, simulated detritivory rates decrease toward low latitudes (Fig. 11), but their relative importance peaks in the Trans_Area (Fig. 8), concurrent with the highest proportion of flux feeders (Supplement file S1). As a result of these interacting mechanisms, detritus turnover times with respect to both microbial and zooplanktonic removal (calculated as the quotient between stocks and removal rates; Fig. 11) are shortest in the Trans_Area, and much longer in the STNA. While these outcomes remain model-dependent and subject to known biases (Table 3), they help reconciling the lability- and temperature-based hypotheses, yielding a more nuanced understanding of mesopelagic TE variability (Marsay et al., 2015), and further highlight that zooplankton processes must be explicitly accounted for.

Previous studies showed that considering the seasonal cycle is key when linking EE and TE observations and models (Ceballos-Romero et al., 2016; de Melo Viríssimo et al., 2024). In this study, however, we refrain from reporting these metrics at intra-annual scale for the following reasons: (i) for EE, there is a regionally varying time lag between NPP, detritus accumulation, and export rates (Figs. 2–4); (ii) in the case of the TE metric, distortion arises from the gradual vertical propagation of sinking POC fluxes (Giering et al., 2017), and from the lagged response of food-web processes that drive flux attenuation. In the highly seasonal SPNA, for example, sources exceed sinks by 25 % in May, and sinks exceed sources by 23 % in August (Fig. 9). Much larger imbalances may occur when examining smaller regions and non-climatological data. In summary, seasonal departure from steady state is an inherent feature of model-simulated POC budgets (Oliver et al., 2025), and helps explain difficulties in interpreting short-term in situ measurements (Giering et al., 2017).

4.4 Study limitations and pathways for model improvement

In the previous sections, we have dissected the drivers of epipelagic export and mesopelagic POC budgets in PISCES. Here, we synthesise the resulting process-level understanding with the biases identified through model evaluation (Figs. 2–6) to provide a unified diagnosis of model performance across regions (Table 3) and outline pathways for model improvement.

In the subpolar area, simulated primary production and sPOC show small-to-moderate negative biases, while a clear bias in POC export cannot be identified (compare Figs. 5 and 6). Total phytoplankton biomass is reasonably simulated, but substantial positive biases in diatoms and mesozooplankton likely lead to an overestimation of large detritus export (Wang and Fennel, 2022). Correcting NPP underestimation would tend to increase sPOC stocks and exports. In contrast, correcting diatom overestimation would favour smaller classes of phytoplankton and detritus, lowering export efficiency. Addressing these opposing biases would likely improve model realism. In the mesopelagic SPNA, simulated sPOC stocks and export fluxes remain broadly consistent with what would be expected given epipelagic inputs (Table 3). The limited available observations of POC decay are reasonably reproduced (Lemaitre et al., 2018), but large detritus fragmentation may be underestimated (Briggs et al., 2020). Uncertainties in the model representation of fragmentation will be addressed through sensitivity analyses in forthcoming work (Orihuela-García et al., 2026).

In the midlatitude Trans_Area, epipelagic and mesopelagic biases differ in nature. Although diatoms are overrepresented at the expense of miscellaneous phytoplankton, NPP is accurately simulated, and no evident biases emerge in epipelagic sPOC stock and export fluxes, suggesting bias compensation. However, sPOC underestimation increases between epipelagic and upper mesopelagic waters, indicating excessive detritus removal. This is consistent with the overestimation of microbial POC degradation inferred from the comparison with PAP (Table 4) and would contribute to an underestimation of TE. Such a bias is not readily apparent from export flux comparisons alone, highlighting the value of evaluating models against diverse flux and stock datasets.

In the subtropical area, PISCES severely underestimates NPP, a bias that propagates to epipelagic sPOC stocks (Galí et al., 2022). Export production is instead overestimated, suggesting a strong positive bias in EE – at least according to BATS/OFP data at 150 m (Lomas et al., 2010). Below 500 m, POC fluxes are no longer overestimated (Fig. 6) despite the positive bias in shallow export (Fig. 5e). This suggests that mesopelagic TE is underestimated, consistent with the increasing underestimation of sPOC with depth, also found in the Trans_Area. Together, these results indicate offsetting biases between excessive EE in the epipelagic and insufficient TE in the mesopelagic. Increasing NPP would likely exacerbate the positive bias in shallow export. This bias would further propagate vertically if mesopelagic transfer were corrected upward. Thus, PISCES' representation of POC cycling in epi- and upper mesopelagic waters in oligotrophic areas may require revision.

Observed biases arise not only from model parameterisations but also from the omission of biological functional groups and/or processes. Several mechanisms known to influence mesopelagic carbon cycling are not explicitly represented in the model used here, including: zooplankton vertical migrations (Carr et al., 2008; Jónasdóttir et al., 2015; Gorgues et al., 2019); the diversity of zooplankton feeding strategies (Stemmann et al., 2004b, 2004a; Kiørboe, 2011; Mayor et al., 2020; Lampitt et al., 2023; Laget et al., 2024); and the distinct roles of free-living and particle-attached bacteria, carbon fixation by chemotrophic prokaryotes, and protist bacterivory (Arístegui et al., 2009; Herndl et al., 2023). For example, including zooplankton migration would likely enhance deep carbon export and mesopelagic sPOC stocks, while reducing NPP and epipelagic sPOC stocks (Aumont et al., 2018; Gorgues et al., 2019). These changes would entail region- and depth-dependent effects on model biases (Table 3). While incorporating missing processes might help bridge the gaps between the model and observations, it would increase model complexity and require extensive re-tuning to maintain consistency across coupled biogeochemical pathways.

5 Conclusions

In this study, we evaluated the representation of POC dynamics in the upper 1000 m of the North Atlantic in the NEMO4-PISCESv2_RC model using diverse observational datasets. Combined with a detailed examination of detrital POC budgets extracted from the model, this approach provided insights into the processes regulating POC distribution, vertical export and biological transformations in biogeochemically-contrasting regions. Below we summarise our main conclusions:

• PISCES-simulated gravitational fluxes increase with latitude (from ∼ 1.3 to 4.3 mmol C m⁻² d⁻¹ annually), whereas maximal primary production occurs in midlatitudes, in reasonable agreement with observations. Consequently, epipelagic export efficiency (EE) increases with latitude, mirroring the diatom fraction and the contribution of large, fast-sinking detrital aggregates to gravitational export. However, the model overestimates diatom biomass at mid and high latitudes, likely implying excessively large detritus export flux and removal pathways.

• In areas experiencing deep winter mixing, diffusive POC export fluxes (dominated by plankton biomass) peak in early spring, preceding the summertime peak in gravitational export. Model-derived diffusive fluxes are consistent with independent observational estimates, and supply ∼ 1.6 mmol C m⁻² d⁻¹ to the mesopelagic subpolar Atlantic annually. Neglecting non-gravitational fluxes substantially underestimates EE in this region (21 % vs. 29 %) and obscures the coupling between primary production and mesopelagic metabolism.

• In the subtropical area, bias compensation between underestimated primary production and overestimated POC export leads to overestimation of modelled EE (8 %). In the lower mesopelagic, however, the model reproduces observed gravitational fluxes, suggesting that model-diagnosed transfer efficiency (23 %) may be underestimated. Negative biases in mesopelagic small POC further indicate an increasing imbalance in modelled detritus supply and removal toward low latitudes. These biases warrant further examination, given the large contribution of oligotrophic regions to global POC export.

• Mesopelagic POC budgets are jointly regulated by zooplankton transformations and microbial degradation. A large fraction of the gravitational POC supply – up to 60 % in the subpolar region – transits through zooplankton. Much of this material is recycled via fragmentation, faecal pellet production, and mortality, rather than being immediately attenuated. This mesopelagic detrital loop modulates particle size, sinking speed, and residence time, thereby shaping microbial remineralisation pathways and rates that ultimately drive POC decay. Modelled mesopelagic decay rates agree within ±15 % with scarce observations at mid-to-high latitudes.

• In PISCES, mesopelagic transfer efficiency is lowest in midlatitudes and coincides with maximal productivity, consistent with some observational assessments. This pattern arises from interacting mechanisms: (i) increasing lability of exported POC toward higher latitudes; (ii) increasing temperature-driven mesopelagic decay rates toward low latitudes; (iii) seasonal trophic feedback between aggregate export and zooplankton, by which enhanced export of large aggregates triggers zooplankton detritivory and fragmentation, limiting transfer efficiency in productive regions. These modelled mechanisms reconcile apparently contradictory hypotheses and highlight the essential role of zooplankton in regulating mesopelagic POC transfer.

• Small slow-sinking detritus supports 50 %–61 % of mesopelagic POC decay across regions and drives 33 %–50 % of the export flux at 1000 m. Key modelled responses enabling these contributions are (i) diverse sources for small detritus – such as zooplankton fragmentation and bacterial disaggregation of large aggregates, plankton mortality, and sloppy feeding – and (ii) the refractory nature of small particles reaching the lower mesopelagic in PISCES. Further research is needed to clarify functional links between small and large POC fractions and their model representation.

• Overall, this study demonstrates the value of a mechanistic, budget-based approach to evaluating models of the biological carbon pump. Detailed POC budget analysis reveals model-specific drivers of export and transfer efficiency and provides a reproducible method to diagnose inter-model biases. This approach can guide model tuning and the development of parameterisations for missing processes. Ultimately, biases in modelled POC budgets can inform assessments of biogenic DIC sequestration fluxes.