. Can organic matter flux profiles be diagnosed using remineralisation rates derived from observed tracers and modelled ocean transport rates

. The average depth in the ocean at which the majority of sinking organic matter particles remineralise is a fundamental parameter in the ocean’s role in regulating atmospheric CO 2 . Observed spatial patterns in sinking ﬂuxes and relationships between the ﬂuxes of different particles in the modern ocean have widely been used to invoke controlling mechanisms with important implications for CO 2 regulation. However, such analyses are limited by the sparse spatial sampling of the available sediment trap data. Here we explore whether model ocean circulation rates, in the form of a transport matrix, can be used to derive remineralisation rates and infer sinking particle ﬂux curves from the much more highly resolved observations of dissolved nutrient concentrations. Initially we show an example of the method using a transport matrix from the MITgcm model and demonstrate that there are a number of potential uncertainties associated with the method. We then use the Earth system model GENIE to generate a synthetic tracer data set to explore the method and its sensitivity to key sources of uncertainty arising from errors in the tracer observations and in the model circulation. We use a 54-member ensemble of different, but plausible, estimates of the modern circulation to explore


Introduction
Sediment trap studies show that the vertical flux of particulate organic carbon (POC) can be described empirically by a power-law curve (e.g. the "Martin curve": Martin et al., 1987; see Fig. 1a) where POC is rapidly remineralised in the upper water column (< 1000 m), leaving a small fraction (5-10 %) of POC that sinks to greater depths.The exponent of the Martin curve (b) reflects the distribution of POC with depth and whether POC is remineralised higher or lower in the water column.Recent studies using deep global sediment trap observations, observations of 234 Th and particle concentration data have suggested a highly heterogeneous pattern of flux profiles and the existence of a general latitudinal trend (Lutz et al., 2007;Honjo et al., 2008;Lam et al., 2011;Henson et al., 2012).POC fluxes in high-latitude regions decrease faster with depth -i.e.flux curves correspond to shal-Published by Copernicus Publications on behalf of the European Geosciences Union.lower remineralisation compared to those in low-latitude regions, where a greater proportion of POC is remineralised at deeper depths (Fig. 1a).However these patterns may be biased by a focus on deep ocean fluxes due to sampling biases in shallower traps (Buesseler et al., 2007).An analysis of neutrally buoyant sediment traps has since suggested that shallower POC fluxes have a latitudinal pattern inverse to that found in deep traps (Marsay et al., 2015).Analysis of the POC fluxes supports a number of mechanisms that explain these spatial patterns, including increased sinking velocity of particles via a "ballast" effect from minerals such as CaCO 3 (Armstrong et al., 2002;Klaas and Archer, 2002); aggregation of particles (Passow, 2004); surface ecosystem dynamics, such as the degree of nutrient recycling and blooms altering characteristics of the particles being exported (Lam et al., 2011;Le Moigne et al., 2012;Henson et al., 2012); and temperature-dependent remineralisation rates (Marsay et al., 2015).However, further progress has been hindered by the relatively low sampling density of long-term sediment trap deployments, particularly in regions such as the Southern Ocean (e.g.Le Moigne et al., 2012;Wilson et al., 2012).
Understanding the underlying reasons for the spatial patterns in remineralisation is a valuable step in understanding mechanisms driving the biological pump.This is key to understanding how the biological pump will respond to both past and current changes in climate (Honjo et al., 2014).Atmospheric carbon dioxide (CO 2 ) concentrations have been shown to be sensitive to changes in flux profiles when modelled, primarily via the effect of redistributing DIC in the modern ocean interior (Kwon et al., 2009).Atmospheric CO 2 is also potentially sensitive to changes in the POC flux to deep-sea sediments relative to fluxes of calcium carbonate (CaCO 3 ) over longer timescales (Archer and Maier-Reimer, 1994;Roth et al., 2014).The range of potential mechanisms results in a range of potential feedbacks for the biological pump with different magnitudes and directions that will respond to a range of different environmental factors.An additional source of uncertainty in modelling sensitivity studies is that only global changes in remineralisation depths, not spatially variable changes, have been considered.
A potential approach to increasing and enhancing the resolution of POC observations is to use climatological fields of dissolved nutrients to estimate POC remineralisation rates.The water column profile of remineralisation rates can be related to flux curves by the fact that the vertical profile of remineralisation rates is the first derivative of the vertical profile of fluxes (e.g.Suess, 1980;Martin et al., 1987) (Fig. 1a and  b).The global distribution of a biological nutrient, such as phosphate (PO 4 ), or dissolved inorganic carbon (DIC), results from the net action of the biological pump (uptake during photosynthesis at the ocean surface and subsequent remineralisation of organic matter in the interior) in combination with physical processes (e.g.air-sea gas exchange in the case of DIC) and other biological processes (such as denitrification) that are integrated through time via ocean circulation., where R z 0 = 1 b z 0 ; see Martin et al. (1987) for a full derivation.
In particular, apparent oxygen utilisation (AOU) has a long history of use as a measure of net organic matter remineralisation in the ocean interior.When combined with a tracer tracking the age of water masses, AOU can be converted to apparent oxygen utilisation rates (AOUR) and related to carbon remineralisation via stoichiometric ratios (e.g.Jenkins, 1982).Feely et al. (2004) compiled AOUR-derived profiles of organic carbon remineralisation rates for 10 regions in the Pacific, noting higher rates in the North Pacific relative to the South Pacific and similarities between regions with high rates and higher CaCO 3 fluxes.However, the in situ AOUR of a water parcel reflects the history of remineralisation over its whole trajectory (e.g.along an isopycnal), producing an average AOUR that is biased towards shallower regions, where remineralisation rates are higher.Therefore, when relating AOURs back to flux profiles by integrating them vertically, estimates are representative of large oceanic regions onlyi.e.there is a high degree of spatial smoothing (Sonnerup et al., 2013;Stanley et al., 2012).
An alternative to combining AOU and age tracers is to use the spatial gradients in tracers to separate out and quantify the change in a concentration of a tracer at any point from circulation only.Spatial gradients of a tracer along the trajectory of a water mass in the ocean interior reflect mixing with other water masses and processes such as the remineralisation of organic matter.Gradient-based approaches aim to solve for the effect of mixing by defining a water mass as the sum of mass fractions from different sources (e.g.Anderson and Sarmiento, 1994;Broecker et al., 1998;Gebbie and Huybers, 2010).A recent development of this method, the total matrix intercomparison (TMI) method described by Gebbie and Huybers (2010), solves for up to six mass fractions for each grid box in a 4 • × 4 • resolution with 33 vertical levels using temperature, salinity, δ 18 O, and nutrients with an additional source term reflecting organic matter remineralisation.The source term is related to the nutrients using stoichiometric ratios.The TMI method therefore produces a high-resolution field of remineralisation estimates.However, these terms do not reflect a rate -i.e.TMI predicts mol PO 4 remineralised in a grid box rather than mol PO 4 yr −1 .This is because the TMI method reflects the pathways of ocean transport but not the rates of transport (Gebbie and Huybers, 2010), although Gebbie and Huybers (2012) have shown that rates can be estimated when combining the TMI method with radiocarbon data.
Spatial gradients in tracers have also been used to diagnose export fluxes of calcium carbonate (Sarmiento et al., 2002) and opal (Sarmiento et al., 2004) from the surface.The method is based on taking the ratio between the vertical gradients of alkalinity and nitrate in the upper ocean, to reflect the ratio of CaCO 3 : organic carbon exported from the surface assuming a dominance of vertical transport over horizontal transport (Sarmiento et al., 2002).The ratio is converted to a flux using estimates of organic carbon export (e.g.Henson et al., 2011).However, a model-based assessment of the method suggests that biases could occur due to the remineralisation of dissolved organic matter (DOM) and strong meridional transport that violate the assumption of processes occurring primarily in the vertical dimension (Jin et al., 2006).This highlights an additional potential issue of using remineralisation rates to infer flux curves that has not been considered by previous studies.
Ocean circulation models offer the opportunity to estimate remineralisation rates of organic matter from tracer data by exploiting the calculated modelled transport rates to account for the effects of ocean circulation on tracers.The aim of this paper is to explore the feasibility of inferring flux profiles of particulate organic matter from remineralisation rates that have been derived from observed tracers using this method.A method and example of estimating remineralisation rates using transport matrices is first introduced.We identify potential sources of error for this method.We then derive a set of model experiments that are used as a synthetic data set with which to test the sensitivity of the approaches to various sources of error and explore potential constraints on estimates.Finally, we explore the uncertainties associated with the broader concept of inferring flux curves from remineralisation rates.

Estimating remineralisation rates using modelled ocean transport rates
Remineralisation rates can be calculated as the amount of tracer supply needed to maintain tracer observations at steady state once the effects of model ocean transport have been ac-counted for, i.e. transport divergence (Deutsch et al., 2007).
Remineralisation rates can then be used to estimate the vertical particulate flux curve.A practical method for this is to apply a transport matrix approach.A transport matrix, hereafter abbreviated to TM, is a representation of steady-state transport rates in the form of a sparse matrix that is derived empirically from an ocean circulation model (Khatiwala et al., 2005;Khatiwala, 2007).Using a TM is akin to the fixed transport rates implied in a box model (e.g."LOSCAR": Zeebe, 2012) but with the advantage that the rates are diagnosed from a dynamic ocean model with a much higher spatial resolution.
For every grid box in the model the TM defines a set of coefficients for neighbouring grid boxes that represent the change in any tracer due to ocean circulation during a single time step of the model (see Table 1 for an illustrative example).Using the TM in place of the model circulation, the evolution of a tracer in time is where A is the TM (unitless), c is a vector representation of the gridded three-dimensional tracer field (e.g.mol kg −1 ) with the superscript referring to the time-step index, and q is a vector representation of any non-circulation-related sources and/or sinks for a tracer during the time step (e.g.mol kg −1 ), for instance due to remineralisation.Assuming steady state (where c = c n+1 = c n ), the simplest approach to estimating remineralisation rates using Eq. ( 1) is to solve for q directly, given a TM and a steady-state tracer: where I is the identity matrix (see Table 1 for the differences between A and A − I).Applying Eq. ( 2) with an observed tracer field is the same as initialising an ocean circulation model with observed tracer concentrations and then stepping the model forward for one time step.q will therefore reflect the exact interior source/sink terms needed in one time step to maintain the steady-state tracer concentrations given the model transport rates.We refer to q as the interior source/sink (ISS) although it is not necessarily representing a mechanistic process.q reflects the amount of remineralisation predicted over the time step which we convert to an amount per year by multiplying by the number of time steps within a year.To infer particulate flux curves from estimates of the remineralisation rates, we fit a linear function to the logtransformed profiles of rates in the ocean interior vs. depth (Suess, 1980;Martin et al., 1987;Stanley et al., 2012).The gradient of the function can be related to the flux profile exponent by adding 1 (Fig. 1

Example using a general circulation model transport matrix
As an example of the approach, we use the annual average TM derived from a 2.8 (2) (Fig. 2).The resulting ISS estimates at shallower depths of the ocean interior (170 m) show some spatial patterns that are consistent with expectations of general export patterns, such as higher rates in the equatorial upwelling regions, and the subpolar regions of the Southern Ocean (Henson et al., 2011) (Fig. 2a).However, negative ISSs also exist, indicating a sink of PO 4 , such as in regions of the subpolar Southern Ocean and North Atlantic.At deeper depths (2030 m) the elevated values in the Southern Ocean match the shallower pattern, but overall there are fewer clear spatial features and a more random pattern of positive and negative ISSs (Fig. 2b).The existence of so many negative values gives rise to near-zero values when averaging over large spatial scales.
Linear functions were fitted to the log-transformed water column profiles of estimated remineralisation rates to infer particulate flux curve exponents (Fig. 2c).The reference depth (depth at which organic matter is exported from the surface) is taken as the edge of the second shallowest grid box (120 m).The significance of this choice is discussed later in Sect.4.2.Functions were only fitted where water columns containing > 3 grid boxes with positive rates.The global mean value for b is −0.77 ± 1.76(±1 SD), close to the value found by Martin et al. (1987), the mean (b = −0.7)found in a reanalysis of globally distributed sediment trap data (Primeau, 2006) and the global mean value (b = −0.639)found by Henson et al. (2012).The value is, however, smaller than values around −1.0 found when optimising ocean biogeochemical models (Kwon and Primeau, 2006).However, there is a wide range of values for b, including positive values, indicating an increasing flux with depth.There is also no clear pattern in the flux curves (Fig. 2).The negative ISSs and positive values for b suggest that the simple inversion of [PO 4 ] observations using this approach is susceptible to large errors that warrant further analysis to characterise and quantify potential uncertainties with this method when used with any modelled circulation.
Table 1.Example of using a transport matrix to calculate PO 4 remineralisation (µmol kg −1 yr −1 ) in one grid box from PO 4 concentrations (µmol kg −1 ) given in c.Grid boxes, taken from a row of the TM, are arbitrarily numbered, where 1 is the central grid box where the calculation is taking place.The example shows a simplified situation where there are six neighbouring grid boxes with their relation to the central grid box given by the directions in brackets.Each coefficient in the six boxes represents the flux of PO 4 into the central box from that grid box.The coefficient in the central grid box for A (see Eq. 1) represents the amount of tracer left in the central grid box after one time step, whilst the coefficient for A − I (see Eq. 2) is the flux out, equal to 0.9816-1.The sum of coefficients is shown underneath, with the estimated remineralisation (q in Eq. 2) in bold calculated as the sum of the element-wise multiplications of A − I and c.

Model description
To explore the errors when using modelled transport rates, we first derive a synthetic data set using the Earth system model "GENIE" (Ridgwell et al., 2007a).GENIE features a 3-D ocean circulation model coupled to a 2-D energymoisture balance model of the atmosphere and a dynamicthermodynamic sea-ice model (Edwards and Marsh, 2005).
In the configuration used here, the ocean model is nonseasonally forced and solved on a simplified 36 × 36 equalarea horizontal grid (10 • longitude by 3-15 • latitude) with eight vertical layers.The biogeochemical model is that described in Ridgwell et al. (2007a).
Our choice of GENIE over other possible models and available transport matrices reflects a number of considerations.The configuration of GENIE used here was derived using a set of ensemble parameters relevant to the physical circulation that were sampled from ranges to test the sensitivity of the model to assumptions about ocean circulation and find an optimal set of parameters (Annan et al., 2005).The ensemble members are therefore different realisations of the modern ocean circulation (see Supplement).The availability of this ensemble enables us to explore errors associated with uncertainty in model circulation states.The non-seasonal circulation, coarse resolution and integration speed of GENIE also facilitate the relatively easy and fast retrieval of TMs.This is the first time that transport matrices have been con-structed from an Earth system model.The relative simplicity of GENIE also keeps a focus on the methodological concept.It is worth noting that much of this could also be achieved by "coarse-graining" a TM derived from a higher resolution model, described in Khatiwala (2007), but without the availability of alternative circulation states.

Diagnosing transport matrices in GENIE
The method of Khatiwala et al. (2005) is adapted to diagnose a steady-state ocean circulation simulation in GE-NIE (see Supplement).The nth ocean grid box in GENIE is dyed with 1 mol kg −1 of an inert "colour" tracer.After the model is integrated for one time step, the resulting pattern of tracer is recovered, vectorised, and forms the nth column of a 6210 × 6210 sparse matrix.This is repeated for all 6210 ocean grid boxes in GENIE.Our method of diagnosing the TM in GENIE differs from that detailed by Khatiwala et al. (2005) in two ways.First, we do not use smoother basis vectors and instead use the simpler method of dyeing a single grid box.Secondly, each grid box is initialised only once and there is no averaging because the circulation in the eight-level version of GENIE is non-seasonally forced.The circulation can therefore be diagnosed during a single continuous simulation.The length of the time step at which the TM is diagnosed in GENIE is 0.01 year.The TM diagnosed here includes the effects of virtual salinity fluxes applied in the ocean model to account for changes in volume (Edwards and Marsh, 2005;Ridgwell et al., 2007a).To account for this, all tracer concentrations multiplied by the TM are first normalised by the corresponding salinity field and then converted back to concentration units afterwards (Ridgwell et al., 2007a).

Synthetic data sets
We use the biogeochemical model described in Ridgwell et al. (2007a) with the biogeochemical parameter values described in Ridgwell et al. (2007b) to produce a synthetic data set of tracers.In this, nutrients in the surface grid boxes are utilised by biological activity limited by [PO 4 ] according to a Michaelis-Menten-type limitation (a maximum rate of 1.96 µmol kg −1 year −1 and half-saturation constant of 0.22 µmol kg −1 ) and the ambient light levels (a linear limitation term).A fixed fraction of the uptake (66 %) is immediately transformed to DOM which can be transported by circulation and remineralised with a time constant of 1/0.5 year.The remaining fraction is exported (34 %) as particulate organic matter (POM), which remineralises instantaneously at depth according to the Martin curve with a global b value of −0.858.The use of a Martin curve keeps the discussion relevant to its use with sediment trap data and is the only difference between our model setup and that of Ridgwell et al. (2007a, b) pletely remineralised to maintain a closed system -i.e.there are no interactions with a (hypothetical) sediment.
[PO 4 ] is used as the tracer for inversion by the TM.An alternative tracer for this could be AOU (related to PO 4 via stoichiometric ratios), as this tracks only organic matter remineralisation (regenerated PO 4 ), whilst total phosphate also includes phosphate that has been advected from the surface (preformed PO 4 ).However, AOU is subject to assumptions about oxygen saturation (Ito et al., 2004;Duteil et al., 2013) and may also reflect non-biological processes (Dietze and Oschlies, 2005).To focus on errors deriving only from the method of inverting tracers, we choose to use [PO 4 ] over AOU.We compare our TM-derived estimates of PO 4 remineralisation with the actual model remineralisation of PO 4 in each grid box as diagnosed within the experiment run as mol PO 4 kg −1 year −1 .The [PO 4 ] field is the annual average taken from the last year of a 10 000-year spin-up.A corresponding TM is diagnosed at the end of the 10 000-year run.This synthetic data set is referred to as SYN.A second synthetic data set (SYN-NODOM) is also produced in the same way as the SYN data set except no DOM is produced to explore the effect of DOM when inferring flux curves from remineralisation rates.

Experiments
We design a number of experiments to explore the sensitivity of the approach to various sources of error (experiment names are indicated in brackets): 1. (TWIN) We first use the TM corresponding to the synthetic data set (SYN) to estimate remineralisation rates from the corresponding [PO 4 ] field as a proof of concept of the method (a twin test), i.e. to test that the modelled remineralisation rates can be recovered by applying the TM to the modelled [PO 4 ] field.

(ERR-OBS)
The effect of errors from the tracer observations themselves is simulated by calculating 100 random perturbations to the synthetic PO 4 concentrations (SYN) within one SE.For each grid box, the mean is taken as the synthetic PO 4 concentration (SYN), and one SE from the SE of World Ocean Atlas (WOA) [PO 4 ] observations regridded to the GENIE grid (Garcia et al., 2010).

(ERR-CIRC)
To explore the effect of circulation uncertainty, we diagnose 54 individual TMs from an existing perturbed physics ensemble (Annan et al., 2005).The ensemble is the result of tuning circulation parameters to fit modern temperature and salinity fields using an ensemble Kalman filter.Each ensemble member is spun up for 10 000 years, after which the TM is diagnosed and used to invert the synthetic [PO 4 ] field (SYN).The circulation parameters in Ridgwell et al. (2007a, b) are an average of the ensemble parameters, such that the 54 TMs offer a range of different ocean circulation states.Details of the ensemble parameters and comparisons against the standard configuration can be found in the Supplement.

(ERR-DOM)
We explore the effect of DOM when inferring particulate flux curves from remineralisation rates.As a comparison to the synthetic data set, we run an identical experiment but with no DOM created (SYN-NODOM) -i.e.all PO 4 is exported as particulate organic matter.
3 Uncertainty analysis results

Assessment of the transport matrix inversion method
We use the output from GENIE as a synthetic data set from which to assess the transport matrix inversion method and identify the sources and nature of the errors involved.Figure 3a and b show the [PO 4 ] field in GENIE at two depths, directly below the surface (290 m) and in the deep ocean (2106 m), with the corresponding annual average input of PO 4 from the continual remineralisation of sinking particles and dissolved organic matter at the same depths (Fig. 3c and  d).The higher remineralisation values calculated for single grid boxes occur where the remaining particulate flux is remineralised at the seafloor to maintain a closed system, i.e. to ensure there are no losses to sediments (Fig. 3d).We also show the inventory of phosphorus in DOM integrated over the ocean interior below 175 m (Fig. 3e).DOM has a visually similar pattern to shallow PO 4 remineralisation, because they are both linked by export production.The close resemblance and magnitude of the DOM at shallow depths may have implications for using remineralisation rates to infer remineralisation fluxes; these will be discussed in Sect.4.2.Finally, we note a correspondence between areas in the Southern Ocean with a large magnitude of remineralisation at shallow depths and locations where convection occurs in the model (Fig. 3f).This is consistent with high productivity driven by nutrients delivered to the surface via strong vertical mixing.
To demonstrate and test the method described in Sect.2, we first invert the model-generated synthetic [PO 4 ] field using the corresponding TM as per Eq.(2) (Fig. 4a).The interior source/sink term (ISS) calculated by inverting the synthetic [PO 4 ] field is consistent with the model-calculated remineralisation with minor deviations from the 1 : 1 ratio line (Fig. 4a).This demonstrates the success of the approach as the errors in the synthetic tracer field and circulation scheme are effectively reduced to near zero in this example.The errors between estimated and modelled remineralisation rates cluster around zero with a median proportion of error of 6.3 × 10 −5 (Fig. 4b).Minor errors likely represent small differences in the circulation at individual time steps that have been sampled by our method of sequentially dyeing grid boxes during a single run.Overall, however, the TM inversion of the synthetic [PO 4 ] estimates the PO 4 remineralisation rates very well (compare Fig. 4c and d with Fig. 3c  and d) with minor errors.This demonstrates that, in theory at least, modelled circulation rates in the form of a transport matrix can be used to successfully estimate remineralisation rates from a tracer field.

Sensitivity of inversions to sources of errors
Although remineralisation rates can be estimated by applying transport rates to a tracer field as shown above, there are several assumptions that will introduce error when this is applied to observations.In the following sections, we detail the results of experiments designed to explore these sources of error.

Error from observations
Error related to the 1 • World Ocean Atlas annual mean climatology (Garcia et al., 2010) will introduce some uncertainty in the TM inversion due to measurement errors and biases in the climatology itself.As a measure of how sensitive the TM inversion method is to these errors, we regrid the SE of annual [PO 4 ] observations onto the GENIE grid (Fig. 5).
We do not consider any additional uncertainty here that may arise through the re-gridding process (e.g.Kriest et al., 2002).The SEs in each grid box are used to produce an illustrative estimate of the uncertainty inherent in the observations.The SEs are highest in the coastal regions and high latitudes, and at shallower depths (compare Fig. 5a and b).We note that errors in the observations may be spatially correlated, which will warrant consideration in future work.The SEs are used to produce 100 versions of the synthetic [PO 4 ] field that have been perturbed within the observation uncertainty which are then inverted using the TM.The resulting variability in the PO 4 ISSs, as characterised by the SD, is relatively large compared to the ISS values themselves, around 1-3 orders of magnitude larger than the ISS values (see also Fig. 8).There are positive linear trends between the SD of the 100 perturbed observations in each grid box and the SD of the 100 corresponding ISS estimates in each grid box (Fig. 5c).Having grid boxes with higher uncertainty in the observation results in greater uncertainty in the ISS estimates.However, two clusters can be broadly defined in Fig. 5c, both with separate linear trends that correspond well with the size of the "flux-out" term of the TM (i.e. the flux of tracer associated only with ocean transport; see Table 1).The colour scale gives an indication of the relative size of the flux out for each grid box due to circulation.Where the central grid-box coefficient is larger, the uncertainty in ISS arising from the uncertainty in the observations is much more sensitive.This suggests that the ISS uncertainty is partly a function of the circulation diagnosed in the TM -i.e. if the observation uncertainty were fixed to a constant across all grid boxes, ISSs in some grid boxes will have greater uncertainty than others due to this effect.The flux-out term is largest in areas in our TM where convection occurs, i.e. the Southern Ocean and North Atlantic (Fig. 3f), because this is where the largest transport fluxes are in the model (Fig. 5d).

Error from circulation estimates
Another potential source of error when inverting nutrient observations arises from the use of a modelled circulation field that will inevitably have a somewhat poorly quantified relationship to the circulation of the real ocean.Using TMs representing plausible but different realisations of the modern ocean circulation from an ensemble to invert our synthetic data set, we can explore the effect of errors arising from uncertainties in circulation rates only.Figure 6a shows the mean and 1 SD of the ISS estimates generated when inverting the synthetic [PO 4 ] field with all 54 ensemble TMs.The values furthest from the 1 : 1 ratio line and those with the largest error bars are located in regions where convection occurs in the model (Fig. 6a and c).The strength of overturning varies within the ensemble (see Supplement), suggesting that this structural uncertainty in the model is a likely cause for the wide range of remineralisation estimates.An additional issue is that, unlike most grid boxes in our TM, where the spread of a tracer over one time step is limited to neighbouring grid boxes, convection increases the number of grid-box connec-tions in the vertical.This could have the effect of increasing the range of remineralisation estimates because there are more grid boxes.Figure 6b shows the same as Fig. 6a but with the convection-related values removed.Even in areas where there is no convection in the model, the range of ISS estimates from circulation uncertainty is larger relative to the range of remineralisation values.The range of errors arising from using different circulations are also larger at shallower depths compared to deeper depths in the water column (Fig. 6c and d).
To understand why different circulation estimates can have a large impact on ISSs, we explore the size of the PO 4 remineralisation flux in a grid box relative to the size of PO 4 flux from the modelled circulation.To illustrate this, we compare the steady-state circulation flux of PO 4 into a grid box with the remineralisation flux of PO 4 into each grid box from the synthetic run (the sum of these at steady state will equal the flux of PO 4 out of the box).Across the whole model ocean interior, the mean proportion of remineralisation flux to the total flux into each grid box is 0.005 ± 0.025 (±1 SD).The proportion is generally 2 orders of magnitude higher at shallower depths (290 m) than at depth (2106 m) (compare Fig. 7a and b), reflecting the decrease in remineralisation fluxes with depth, whilst circulation fluxes vary much less in magnitude in GENIE (Fig. 7c and d errors in the circulation fluxes are therefore likely to dominate over the remineralisation fluxes, leading to large errors when using modelled circulation rates.

Error comparison
To compare the magnitude of the various possible errors, we show the global mean synthetic PO 4 remineralisation profile with the global mean and median SDs of the ISSs calculated for each of the ERR-OBS and ERR-CIRC experiments -i.e. a SD is calculated for the experiments in each grid box, which are then averaged (Fig. 8).We use both the mean and median SD in this figure because the mean SDs of the ERR-CIRC experiments are skewed by the large variability in high-latitude regions (Fig. 6a; see also Supplement).Both sources of error are larger at shallower depth and mostly decrease in magnitude with depth although the mean circulation uncertainty increases below 2000 m.The magnitude of median uncertainty from the observations is much larger than from our circulation estimates although the two are similar magnitudes when the mean SD is used.
Despite similar magnitudes of uncertainty arising from both potential errors in the observations and from the model circulation field, the nature of the uncertainty is different.Uncertainty from the observations is higher in regions where observations are more uncertain, e.g.coastal areas in Fig. 5a and b, and in regions of the model where convection occurs (Fig. 5c).In contrast, the uncertainty arising from the model circulation field used will be systematic and dependent on where the model circulation is most different to the real ocean, e.g.Fig. 6.Although the patterns in the surface PO 4 ISSs from the MITgcm inversion are systematic (patches of positive and negative ISSs in the Southern Ocean: Fig. 2a), it is difficult to tell whether this reflects a systematic difference between observed and simulated circulation rates or a caveat associated with convection in the model and larger circulation fluxes in these areas.The relative smoothness of the ISS estimates may reflect the use of the climatology which has already been smoothed (Garcia et al., 2010).

Discussion
We have presented a straightforward method of using a steady-state model circulation, as represented by a transport matrix, to estimate organic matter remineralisation rates 0 .0 0 0 .0 1 0 .0 2 0 .0 3 0 .0 4 0 .0 5 0 x 1 0 -4 2 x 1 0 -4 5 x 1 0 -4  4 x 1 0 -4  3 x 1 0 -4  1 x 1 0 -  from a tracer climatology.Our main goal is to explore the feasibility of using this method to infer spatially explicit organic matter flux curves, aiding additional understanding of the biological pump in the modern ocean.Our results show that this method is associated with a number of significant sources of error that give rise to the spatial patterns and negative values seen in an example inversion using a circulation field from a coarse-resolution ocean model (Fig. 2).In the following sections we discuss potential directions for estimating remineralisation rates from tracer data using model circulation and additional considerations needed when using these to infer particle flux curves.
The sensitivity to errors in the observations is a result of the way that the transport matrix (TM) is constructed.A change in a tracer due to circulation in a model time step is relatively localised due to the finite speed of advection and diffusion in the model (Khatiwala et al., 2005).Therefore, the ISS estimates are sensitive to the resulting large coefficient in the central grid box (see Table 1).This will be a feature of all TMs constructed using the method of Khatiwala et al. (2005) regardless of the resolution or model used.Our results also illustrate the sensitivity of remineralisation estimates to differences between the model transport rates and actual transport rates in the ocean.Remineralisation fluxes of PO 4 in the synthetic tracer experiments are orders of magnitude smaller than fluxes of PO 4 from circulation.Model transport rates would therefore need to attain a high level of accuracy to minimise the effect of error on the solutions (Anderson and Sarmiento, 1994;Sarmiento et al., 2002).Even a data assimilated model, such as that produced by the ECCO ("Estimating the Circulation and Climate of the Ocean"; Stammer et al., 2004) consortium, designed to be a dynamically con-sistent estimate of ocean circulation over recent decades may still lead to large diagnosed flux errors.An additional consideration may also be the representation of a steady state in the TM, either because a model is not spun up fully or because there is seasonality.

Are there constraints on the circulation uncertainty?
The flipside of the magnitude and nature of the circulation control on the diagnosed remineralisation rates is that tracers with a steady-state constraint, where it is expected that there should be no significant sources or sinks at depth, could be used to estimate the magnitude of the circulation-based error.As an example, we show an ISS field generated when inverting the salinity field from our synthetic data set with the synthetic transport matrix (Fig. 9a and b).Salinity ISSs are randomly distributed around zero, consistent with the concept that there are no significant local sources or sinks in the ocean interior.In comparison, applying a different transport matrix, arbitrarily chosen from the ensemble of alternative modern circulations, results in distinct spatial patterns in the ISSs (Fig. 9c and d).Comparing the errors from the inversion of the synthetic [PO 4 ] field using the same TM shows that the two have some visible similarities (Fig. 9e and f).For example, there are correspondences in the deep South Atlantic (Fig. 9d and e) and the subtropical regions of the Pacific at 290 m (Fig. 9c and e).We do not expect an exact correspondence, but the visual similarities support the idea that they are related via errors in the model circulation.This suggests that conservative tracers, i.e. tracers that do not have sources or sinks in the ocean interior, could provide a constraint on the magnitude of error.For example, considering Eq. ( 2) but focussing on an individual grid box, the ISS (q i=1 ) is a function of the TM coefficients (M i ) and the tracer concentrations (C i ) as per the example in Table 1: Redefining the modelled circulation terms to reflect that the modelled circulation is a function of a "true" circulation term and an error term (M i = F i + i ) and substituting into Eq.( 3) and expanding, (4) For a conservative tracer at steady state, it is expected that N i=1 F i C i = 0.A significant departure from zero in q i=1 is likely to result from the error terms.This may provide a way forward to constrain the ISSs produced by the TM method described here.Such a method would be conceptually similar to the TMI method (Gebbie and Huybers, 2010).There are, however, a limited number of tracers available (e.g.temperature, salinity, δ 18 O, 14 C, CFCs) that could realistically be used to constrain model circulation errors.Of these tracers, CFCs are not yet present in the deep ocean and thus are not in steady state.In our TM, there are typically 15 gridbox connections (6 neighbouring grid boxes plus extra connections due to the Gent-McWilliams parameterisation) used for each calculation which would lead to an underdetermined solution, i.e.where the unknowns outnumber the constraints.It would be interesting to see whether a simplified TM with fewer coefficients, such as matching the method of Gebbie and Huybers (2010) by adopting a seven-point stencil (six neighbouring grid boxes), could use this approach.Even so, with a seven-point stencil, the solution may still be underdetermined.

Inferring flux curves from remineralisation rates
In the previous sections, we have shown that a simple approach to estimating remineralisation rates using modelled transport rates is sensitive to different sources of errors.Taking the next step, in the case that remineralisation rates could be estimated with some reliability, we explore the sensitivity of inferring flux curves by vertically integrating remineralisation rates in the presence of DOM.The remineralisation of DOM can occur away from where it was exported, affecting the assumption that the remineralisation rates reflect only vertical processes.To explore this, we infer flux curves using remineralisation rates from the synthetic data set (SYN) and a second run where no DOM is created (SYN-NODOM).To infer a power-law curve, a linear trend is fitted to the logtransformed remineralisation rates following previous studies (Berelson, 2001;Lam et al., 2011;Stanley et al., 2012).The gradient of the linear trend gives the value of the exponent for the remineralisation curve, which is converted to a flux curve by adding 1 (Martin et al., 1987;Stanley et al., 2012; see also Fig. 1).
The exponents from power-law curves, fitted to vertical PO 4 remineralisation profiles, when PO 4 is only exported as particulate organic matter, are all close to −1.9 (Fig. 10a).This corresponds to a flux curve exponent of −0.9, in good correspondence with the −0.858 used for the run.The presence of remineralisation from DOM (see Fig. 3e) lowers the value of the fitted exponent reflecting a shallower Martin curve (Fig. 10b).In our experiments, the remineralisation of DOM lowers the fitted exponent by as much as 0.6.This occurs because the remineralisation of DOM inflates the remineralisation rate in the shallower grid boxes relative to those in the water column below.In a few grid boxes in the North Atlantic the estimate of the flux curve exponent increases.This has a strong correspondence with the deepest convection (Fig. 3f) where DOM is transported deeper into the water column before being remineralised.This highlights that an additional source of remineralisation not restricted to vertical processes, such as the remineralisation of DOM, can alter the flux curve in an unpredictable way.The DOM bias in GENIE occurs predominantly in the high latitudes, where DOM is efficiently advected into the ocean interior.The range of exponents purely from this bias is of a similar magnitude and spatial distribution to the Martin curve exponents found in Henson et al. (2012).
A range of issues relevant to fitting flux curves are also relevant to the retrieval of reliable estimates of flux curves from modelled remineralisation rates.Previous studies have noted that fitted values of b can be sensitive to the estimate of the export depth (z 0 in Fig. 1) (Berelson, 2001;Primeau, 2006;Lam et al., 2011;Olli, 2015).Both Lam et al. (2011) and Olli (2015) consider that the export depth of particles may be deeper than the euphotic zone when particles are physically mixed deeper in the water column, as measured by the mixed layer depth (MLD).The model grids used here impose a fixed export depth on flux curve estimates which could be shallower than expected in the real ocean, such as in the high latitudes, where MLDs are deeper.For example, if the upper-most grid box is considered as the export depth when fitting Martin curves to the MITgcm example in Sect.2, rather than the two uppermost grid boxes (Fig. 2c), the global mean Martin curve exponent is −0.90, equivalent to shallower remineralisation.Developments in understanding the uncertainty in fitting Martin curves to data are therefore equally important to consider when fitting to estimated remineralisation rates (Primeau, 2006;Olli, 2015).

The way forward
Any method of inferring particulate organic matter flux curves from estimated remineralisation rates, whether using model transport rates or using observations (e.g.Stanley et al., 2012), will therefore need to take the remineralisation of DOM (as well as other processes not constrained to the vertical water column such as denitrification) into ac- count, particularly in regions where DOM is advected to deeper depths (Hansell et al., 2009).Ultimately, this additional source of uncertainty when inferring flux curves combined with the uncertainties associated with estimating remineralisation rates suggests this is a difficult task.An alternative approach that explicitly includes the effects of DOM remineralisation is the optimisation of ocean biogeochemical models.In this approach, a biogeochemical model is set up and optimal values of the parameters found, such as the depth of remineralisation, export production and the rate of DOM remineralisation, that result in the best fit between modelled and observed nutrient fields (e.g.Bacastow and Maier-Reimer, 1991;Schlitzer et al., 2002;Kwon and Primeau, 2006;Ridgwell et al., 2007a, b;Kriest et al., 2012;Yao and Schlitzer, 2013;Teng et al., 2014).For example, remineralisation depths could be varied in space and a solution found that best fits tracer observations.A first step towards this has been taken by Yao and Schlitzer (2013).Such an approach allows additional constraints such as export production magnitude to be taken into account as well as providing further evaluation on circulation models via biogeochemical tracers (Najjar et al., 2007).Transport matrices are a useful tool in this approach to help reduce lengthy multiple-thousand-year integration times for model fields to reach equilibrium as well as providing alternative methods to efficiently reach steady-state solutions (Kwon and Primeau, 2006).Ultimately, the approach of optimising biogeochemical models is complimentary to the approach described in this paper as they both rely on the combination between observations and models but in different ways, which could help facilitate a better understanding of biogeochemistry as well as of the models used to study it.

Conclusions
Profiles of remineralisation rates derived from ocean tracers offer a potential method to estimate high-resolution fields of flux curves that could supplement existing global sediment trap studies.The use of model transport rates offers one way of estimating remineralisation rates that could avoid the spatial averaging issues of combining AOU with age tracers.Using the transport rates in the form of a transport matrix is a first step towards this, but the simple application to a tracer such as [PO 4 ], to estimate steady-state remineralisation rates, is subject to potentially large uncertainties related to both uncertainties in the observations and the modelled transport rates.There are potential methods to constrain the uncertainty in the modelled transport rates, but these are currently limited by the number of tracers available.We also highlight that there significant uncertainties when integrating remineralisation rates vertically, a step necessary to estimate particle flux curves, associated with processes such as the remineralisation of dissolved organic matter.This last uncertainty has implications for deriving flux curves from databased estimates of remineralisation rates as well as modelbased ones.The spatial variability in POC fluxes observed in the modern ocean has important implications for our understanding of how the biological pump may have changed in the past and in the future.New approaches to estimating POC fluxes will help provide estimates in regions, such as the Southern Ocean, that are currently undersampled by sediment traps and key to testing existing mechanistic hypotheses.Exploring the spatial variability of POC fluxes in ocean biogeochemical models by finding a set of POC flux profiles that best fits observed tracers will help approach the uncertainties highlighted in this paper but also provide a quantitative analysis of the significance of spatially varying POC fluxes.
The Supplement related to this article is available online at doi:10.5194/bg-12-5547-2015-supplement.

Figure 1 .
Figure 1.The range of observed Martin curves and associated remineralisation rate profiles.(a) The mean (b = −0.639)and global range (b = −1.18 to −0.24) of Martin curves found by Henson et al. (2012) calculated for a unit flux and export depth (z 0 = 100 m) via F z = 1 z z 0

Figure 2 .
Figure 2. Example of using a GCM transport matrix to estimate PO 4 remineralisation rates.(a) The estimated PO 4 remineralisation rates generated using the MITgcm transport matrix at 170 and (b) 2030 m.(c) Power-law flux curve exponents (b) inferred from the estimated remineralisation rates by fitting a linear function to the log-transformed data as in Fig. 1b.

Figure 3 .
Figure 3.The synthetic tracer data set used for transport matrix inversions.(a) [PO 4 ] (µmol kg −1 ) at 290 and (b) 2106 m.(c) The total annual remineralisation flux of PO 4 (nmol kg −1 yr −1 ) at 290 and (d) 2106 m.(e) The water-column-integrated inventory of phosphorus in dissolved organic matter in the ocean interior (mol P).(f) A measure of the convection occurring in a water column in one model year (unitless).Higher values indicate the number of vertical levels involved in convection in the model, i.e. stronger convection.

Figure 4 .
Figure 4. Results from inverting the synthetic data set with its corresponding transport matrix.(a) The interior source/sink estimate for PO 4 when inverting the synthetic [PO 4 ] field with the corresponding transport matrix vs. the synthetic PO 4 remineralisation with a 1 : 1 ratio line.(b) The distribution of errors for the PO 4 interior source/sink estimates (50 bins sized 0.15 × 10 −5 ).(c) The interior source/sink estimate for PO 4 at 290 and (d) 2106 m.(e) Difference between inverse interior source/sink estimates and the synthetic remineralisation field at 290 and (f) 2106 m.

Figure 5 .
Figure 5. Assessment of the errors arising from the uncertainty in [PO 4 ] observations.(a) The SE of [PO 4 ] from the World Ocean Atlas 2009 (Garcia et al., 2010) 1 • climatology regridded to the GENIE grid at 290 and (b) 2106 m.(c) SD of all PO 4 interior source/sink estimates when the synthetic [PO 4 ] field is randomly perturbed within a normal distribution given by the SE of observations.The colour scale indicates the size of the "flux-out" term, where a larger value indicates relatively larger circulation fluxes in that grid box.(d) The SD of [PO 4 ] ISSs at 290 m.

Figure 6 .
Figure 6.Assessment of error arising from using circulation estimates.(a) Comparison of ISS estimates for each grid box from the 54 ensemble members against the synthetic data set remineralisation.Error bars are 1 SD around the mean.Red values indicate regions with convection.(b) As panel (a) but with the red values removed.Note the difference in scale.(c) The SD of PO 4 ISS errors (ISS -model) in µmol kg −1 yr −1 for 290 m.Values are shown on a log scale.(d) As panel (c) but at 2106 m.

Figure 7 .Figure 8 .
Figure 7.Comparison of inputs of PO 4 from remineralisation and circulation at steady state.(a) PO 4 remineralisation as a proportion of the total PO 4 flux into each grid box calculated using the synthetic tracer field at 290 m and (b) 2106 m.(c) The flux of PO 4 into each grid box from circulation only (µmol kg −1 yr −1 ) from the synthetic tracer fields at 290 and (d) 2106 m.

Figure 9 .
Figure 9. Inversion of salinity as a possible constraint on the uncertainty from using a modelled circulation.(a) Inversion of the salinity field from the synthetic data set using the corresponding transport matrix at 290 m (PSU yr −1 ) and at (b) 2106 m.(c) Inversion of the salinity field using a alternative transport matrix from the perturbed-physics ensemble (PSU yr −1 ) at 290 and (d) 2106 m.(e) The error of the synthetic PO 4 ISS (µmol kg −1 yr −1 ) (ISS -synthetic) using the same transport matrix in panels (c, d) at 290 and (f) 2106 m.

Figure 10 .
Figure10.Assessment of the uncertainty associated with dissolved organic matter when inferring flux profiles.Value of the exponent when fitting a power law to the water column remineralisation rates from (a) a model with remineralisation from only sinking particulate and (b) a model with particulate and dissolved organic matter.A value of −0.858 for the Martin curve was used in both models.All curves were fitted with an R 2 > 0.9.The exponent for the remineralisation curve is equivalent to −0.858−1.More negative values indicate a Martin curve that predicts shallower remineralisation in the water column.Hatched areas indicate where the water column contained too few boxes to fit a remineralisation curve (n < 3).