Southern Ocean silicon trap: Data-constrained estimates of regenerated silicic acid, trapping efficiencies, and global transport paths

Received 16 August 2013; revised 10 November 2013; accepted 7 December 2013; published 16 January 2014. [1] We analyze an optimized model of the global silicon cycle embedded in a dataassimilated steady ocean circulation. Biological uptake is modeled by conditionally restoring silicic acid in the euphotic zone to observed concentrations where the modeled concentrations exceed the observational climatology. An equivalent linear model is formulated to which Green-function-based transport diagnostics are applied. We find that the models’ opal export through 133 m depth is 166 6 24 Tmol Si/yr, with the Southern Ocean (SO) providing � 70% of this export, � 50% of which dissolves above 2000 m depth. The global-scale gradients of the opal dissolution rate are primarily meridional, while the global-scale gradients of phosphate remineralization are primarily vertical. The mean depth of the temperature-dependent silicic-acid regeneration reaches 2300 m in the SO, compared to 600 m for phosphate remineralization. Silicic acid is stripped out of the euphotic zone far more efficiently than phosphate, with only (34 6 5)% of the global silicic-acid inventory being preformed, compared to (61 6 7)% for phosphate. Subantarctic and tropical waters contribute most of the ocean’s regenerated silicic acid, while Antarctic waters provide most of the preformed silicic acid. About half of the global silicic-acid inventory is trapped in transport paths connecting successive SO utilizations, with silicic acid last utilized in the SO having only a (5 6 2)% chance of being next utilized outside the SO. This trapping depletes subantarctic mode waters of silicic acid relative to phosphate, which has a (44 6 2)% probability of escaping successive SO utilization.


Introduction
[2] The silicon cycle is a key component of the ocean's biological pump, with diatom phytoplankton accounting for as much as half of the ocean's primary productivity [e.g., Nelson et al., 1995] and a disproportionate fraction of organic matter export [e.g., Buesseler, 1998;Moore et al., 2004;Brzezinski et al., 2011]. The distribution and cycling of silicic acid (Si(OH) 4 ) are very different from those of the other macronutrients, phosphate (PO 32 4 ), and nitrate. The differences stem from the fact that while all phytoplankton require phosphate and nitrate to sustain their metabolism, silicic acid is only required by diatoms, radiolaria, silicoflagellates, and other siliceous plankton that form amorphous silica (opal) to build their cell walls or other supporting structures. Unlike phosphate and nitrate, opal is regenerated to silicic acid through temperaturedependent dissolution without passing through a pool of dissolved organic matter first. These differences in uptake and regeneration translate to profound differences in the global distribution and transport pathways of silicic acid as compared to the other macronutrients. Among other things, the fact that the dissolution of opal and the remineralization of organic matter have very different spatial distributions causes the ratio of silicic acid to nitrate or phosphate to vary from place to place. Therefore, a simple universally applicable ''Redfield'' stochiometric ratio between silicon and organic matter does not exist, in spite of the fact that silicon and nitrogen are taken up in a 1:1 ratio by healthy diatoms [e.g., Brzezinski, 1985;Sarmiento and Gruber, 2006]. In this work, we combine a very simple restoringtype model of the silicon cycle with a data-assimilated, annual-mean circulation model to estimate opal fluxes and Si(OH) 4 transport pathways, and to further elucidate and quantify the differences between the silicon and phosphate cycles.
[3] Considerable research has been devoted to the oceanic silicon cycle, underscoring its importance. Simple models of the silicon cycle have been constructed by Gnanadesikan [1999] to explore different opal dissolution schemes and to attempt to constrain the subgrid parameterizations of ocean circulation models. Sarmiento et al. [2004] use observations of the difference of silicic-acid and nitrate concentrations to trace out the global-scale nutrient transport by subantarctic mode waters (SAMW), and to argue that the observed globally low upper-ocean silicicacid concentrations are the result of silicic acid being efficiently stripped out of the source region of SAMW. In a follow-on study, Sarmiento et al. [2007] employ classic water mass analyses and an adjoint biogeochemical model [Schlitzer, 2007] to place constraints on opal export, opal flux attenuation, and deep opal dissolution rates, which serve as a reference point for our study here. A number of relatively complex ecosystem models have also been used to estimate opal export and the coupling between the silicon and carbon cycles [e.g., Aumont et al., 2003;Moore et al., 2004;Jin et al., 2006]. Dunne et al. [2007] estimate opal fluxes from satellite measurements of chlorophyll using empirical relationships and provide a useful summary of all the other published estimates at the time. Significant uncertainty in global opal export flux remains with estimates ranging from 70 to 180 Tmol/yr [Dunne et al., 2007, and references therein].
[4] Most studies to date have focused on model validation by comparison with observations, with much less emphasis on quantitative diagnostics of the transport pathways of silicic acid. Our paper is an effort toward understanding and quantifying how the global-scale pelagic silicon cycle (the silicon cycle in the open sea beyond the continental shelf) is shaped by the interaction between the biogenic particle flux and the advective-eddy-diffusive circulation of preformed and remineralized Si(OH) 4 . The key scientific questions we address here are: [5] 1. What are the opal fluxes inferred from an optimized restoring-type silicon-cycle model embedded in a data-assimilated circulation ?
[6] 2. What are the distributions of regenerated and preformed silicic acid in the ocean, how do these differ from those of phosphate, and what do they reveal about the pelagic silicon cycle ?
[7] 3. What are the dominant transport pathways of silicic acid, and what are the global teleconnections in biological productivity mediated by these paths?
[8] To answer these questions, we use a simple restoring-type model of the global silicon cycle in conjunction with a data-assimilated circulation to estimate opal export flux, opal dissolution rates, and the transport pathways of silicic acid. Efficient numerics afforded by a matrix formulation of the transport operator allow us to optimize the parameters of the silicon-cycling model. Applying novel transport diagnostics to the silicon cycle is possible by replacing the nonlinear restoring-type siliconcycling model with an equivalent linear model that has exactly the same solution as the nonlinear model. The lin-ear model allows-for the first time-the determination of the preformed and regenerated components of silicic acid as well as the number of transport paths per unit volume that silicic-acid molecules follow between successive biological utilizations. Throughout, we will contrast the cycling and transport of silicic acid with that of phosphate.

Circulation
[9] Our silicon-cycling model is embedded in the dataassimilated, steady circulation of Primeau et al. [2013], which is similar to that of DeVries and Primeau [2011], except that the horizontal resolution has been doubled to 2 3 2 , and CFC-11 and phosphate have been incorporated as tracer constraints in addition to temperature, salinity, and radiocarbon. The model of Primeau et al. [2013] jointly optimized both the circulation and a simple phosphorus cycling model, so that we have a circulation that is internally consistent with the pelagic phosphorus cycle. Throughout, we will use the data-assimilated phosphorus cycle as a point of reference for understanding the behavior of the silicon cycle. The data assimilation used the wind stress climatology of Trenberth et al. [1989] and specifies horizontal and vertical viscosities of A h 55310 4 m 2 s 21 and A v 510 24 m 2 s 21 . The circulation's advective-diffusive transport operator has fixed horizontal and vertical eddy diffusivities of 10 3 and 10 25 m 2 s 21 , respectively, typical of coarse-resolution models. The realism of the model's Atlantic meridional overturning circulation and of the meridional heat and freshwater fluxes was assessed at 4 3 4 resolution by DeVries and Primeau [2011]. The model's ability to capture transport in Antarctic Bottom Water, Antarctic Intermediate Water, SAMW, and North Atlantic Deep Water was quantified in terms of water mass fractions [DeVries and Primeau, 2011] and, at 2 3 2 resolution, in terms of the ability to simulate the observed distribution of CFC-12 in the Southern Ocean [Waugh et al., 2013].

Opal Dissolution
[10] We model the vertical flux profile, F(z), of biogenic silica particles as being attenuated with depth z by a temperature-dependent dissolution length scale, Z(T), following Gnanadesikan [1999]. (We distinguish the flux-profile function F from the actual opal flux U, which is modeled in slightly more detail below.) Thus, the flux profile obeys @F @z 52 F ZðTÞ ; (1) where Z(T) is taken to have an Arrhenius-type dependence on absolute temperature, T, based on the measurements of Erez et al. [1982] and Hurd and Birdwhistell [1983] : ZðT Þ5 w f ae 2TE=T ; (2) with parameters a51:32310 16 day 21 , T E 511; 481 K (effectively an energy scale for dissolution), and terminal fall speed w f of 50 m/d.
[11] Gnanadesikan [1999] chose these values for a and T E so that the dissolution rate constant R5a exp ð2T E =T Þ HOLZER ET AL.: OCEAN SILICON CYCLE is 0.25 day 21 at 25 C and 0.01 day 21 at 2 C, broadly in agreement with observations at Bermuda  and with laboratory studies [Hurd and Birdwhistell, 1983]. However, a single rate constant that is representative of the distribution of silica particles with a wide range of sizes and shapes, and with a protective organic layer in various stages of decay, must be considered uncertain for any given temperature by at least a factor of two. We therefore treat w f , or equivalently w f /a, as an optimizable parameter below. We also tried jointly optimizing T E but the corresponding optimal solution turns out to lie at T E 50 (no temperature dependence) with a dissolution length scale of Z 5 3048 m. Although the corresponding opal export, particle-flux profiles, and silicic-acid field are not very different from what we get by keeping T E fixed, we reject the T E 50 solution on physical grounds: dissolution is bound to have some temperature dependence. We therefore optimize only w f in (2) and keep T E 5 11,481 K. The differences between an optimal temperature-independent dissolution length and the temperature-dependent case fall within the uncertainties estimated below.

Silicic-Acid Uptake
[12] We follow Najjar et al. [1992] and model the local uptake rate per unit volume, J ðvÞ, of silicic acid by conditionally restoring the silicic-acid concentration v to observations in the euphotic zone with a time scale of s530 days. The annual mean observed silicic-acid concentrations are taken from the 2009 World Ocean Atlas (WOA09) [Garcia et al., 2010]. Specifically, the uptake rate per unit volume is modeled as where z c 5 273.4 m is the vertical height coordinate of the base of the model euphotic zone. In the aphotic zone, there is no uptake so that J 5 0 for z < z c . In (3), the func-tionĤðx; x 0 Þ ½tanh ðx=x 0 Þ11=2 ensures that the uptake switches on if v exceeds the observed concentration v obs . The switch turns on over a finite concentration-scale v 0 , chosen here as v 0 510 25 mmol=m 3 ; for v 0 this small, the solution of the model defined below is completely insensitive to its precise value. The functionĤðx; x 0 Þ is simply a Heaviside step at x 5 0 that has been smoothed over a scale, x 0 , with lim x 0 !0Ĥ ðx; x 0 Þ5HðxÞ, where HðxÞ is the usual Heaviside function. The smoothing is necessary here to make J ðvÞ differentiable with respect to v so that Newton's method can be employed to find the model's steady state. The same form of a smoothed conditional uptake was used by Kwon and Primeau [2006] in order to apply Newton's method to a steady phosphate model.
[13] Note that the conditional restoring makes the uptake J ðvÞ a nonlinear function of v. The advantage of this restoring approach is that the uptake rate is constrained by the observational silicic-acid climatology. This implicitly accounts for light and micronutrient limitations, while retaining mathematical simplicity. In coastal regions less than two grid cells deep, and in the few grid boxes of the sea of Japan, we set J ðvÞ to zero to avoid unphysical nutrient trapping where the flow is under-resolved. We therefore only attempt to model the pelagic silicon cycle here.

Silicon Cycle
[14] The silicon-cycling model simply expresses the conservation of Si(OH) 4 mass in the presence of the sink J ðvÞ due to biological uptake and the source SðvÞ due to the dissolution of biogenic opal. We assume that the dissolution time scale is negligible compared to any transport time scales so that SðvÞ can be modeled as an instantaneous redistribution of the uptake rate J ðvÞ over the water column in accord with the flux profile (1). We therefore write SðvÞ5SJ ðvÞ, where the operator S accomplishes the desired redistribution. Physically, SJ ðvÞ represents the divergence of the biogenic particle flux, which is equated to the dissolution rate of opal, and hence production rate of silicic acid, per unit volume. In that sense S is a source operator. Thus, the silicon-cycling model is simply given by @v @t 1T v52J ðvÞ1SJ ðvÞ; where T is the data-assimilated advection-diffusion operator of the ocean, which is constructed to satisfy zero-flux boundary conditions at both the sea surface and at all solid boundaries.
[15] The operator S is formulated here in analogy with the remineralization of organic matter in other restoringtype nutrient models [Najjar et al., 1992]: a fraction r of the opal production is assumed to be dissolved in the euphotic zone, while a fraction 12r is distributed vertically throughout the aphotic zone with the divergence of the opal flux profile F(z). We thus have SJ ðvÞ where FðzÞ5Fðz c Þexp ½2 R z c z Z 21 ðz 0 Þdz 0 is the solution to the flux-profile equation (1), re-written for height coordinates, with Fðz c Þ5 R 0 z c J ðrÞdz, and z b is the vertical height coordinate of the sea floor. Note that S is a linear operator.
[16] The vertically uniform dissolution above z c is intended as a simple parameterization of the effect of mixed-layer turbulence on biogenic opal particles and rapid recycling within the euphotic zone. The fraction r dissolved in the euphotic zone is treated as another optimizable parameter. Note that, mathematically, introducing r is equivalent to multiplying the uptake J, and hence the inverse uptake time scale 1=s, by ð12rÞ and dissolving opal only below z c in the aphotic zone. Our model is therefore equivalent to replacing s in (3) by an optimizablê s5s=ð12rÞ, and setting, in (5) only, r to zero. However, because in the real ocean a finite fraction of opal dissolves in the euphotic zone [e.g., Tr eguer and De La Rocha, 2013], we prefer to regard s as fixed with r capturing rapid (instantaneous here) recycling in the euphotic zone.
[17] The Heaviside function in equation (5) ensures that opal sinking into the bottom layer is dissolved there; we do not allow for the permanent burial of opal in sediments (estimated at a few percent of the global export production Tr eguer and De La Rocha, 2013]) because we do not explicitly represent any external input of HOLZER ET AL.: OCEAN SILICON CYCLE silicate into the ocean (such as from rivers), and therefore wish to conserve the total amount of biogenic silicon cycled in our simple model. The loss of silicon to burial was similarly neglected by Gnanadesikan [1999] and its explicit inclusion has been shown to have a negligible effect on particle fluxes and regeneration rates .

Optimized Parameters
[18] We determine the steady state silicic-acid concentration v by discretizing the differential equation (4), organizing the transport operator T and source operator S into sparse matrices, and solving the resulting system of nonlinear equations using Newton's method. To ensure this can be done, the total amount of silicon in the ocean must be fixed because otherwise the system does not have a unique solution. (This issue does not arise when solving the time-dependent initial value problem for which the initial conditions determine the total mount of silicon.) We therefore add a small relaxational term 2c g ðv2vÞ to the right-hand side, where the ''geological'' restoring time scale c 21 g 510 6 years, and v is the observed average concentration of silicic acid from the WOA09. The solution is insensitive to the precise value of c g .
[19] Denoting the discretized matrix versions of the operators T and S by T and S, we therefore solve the steady state equation for the concentration field v, which is now a column vector whose entries map onto the discrete ocean grid boxes. (The matrix I is the identity matrix.) This equation is nonlinear in v because of the conditional uptake J ðvÞ, but is easily solved using Newton's method in typically less than 10 iterations.
[20] Optimal silicon-cycling parameters r and w f have been determined by minimizing a cost function equal to the volume-weighted squared difference between modeled and observed silicic-acid concentrations. The Matlab function fminsearch was used to locate the minimum of the cost function at the values of r50:80 and w f 540: m=d . The fidelity of the modeled silicic-acid concentration achieved with these optimal parameters is quantified in Figure 1 by the volume-weighted joint probability density function (pdf) of the modeled and observed fields. The joint pdf was estimated using the kernel density estimation method of Botev et al. [2010], modified to include grid-box-volume weights. (The joint pdf is simply the continuum limit of a volume-weighted histogram of the model-versusobservations scatter plot.) The figure shows that the modeled Si(OH) 4 field lies generally within 10% of the observed values, with no major systematic global-scale biases. The squared, volume-weighted Pearson correlation coefficient is r 2 50:95, while the volume-weighted rootmean-square (RMS) error is 12 mmol/m 3 , which is 13% of the global mean concentration. The model would likely fit the observations yet more closely if silicic acid had also been used to constrain the data-assimilated circulation, as has been done for phosphate, radiocarbon, and CFC-11  (the RMS error for PO 32 4 is only 3% of the global mean). Nevertheless, it is reassuring that the computationally much cheaper approach adopted here provides a reasonable match with observations.

Diagnostic Linear Silicon-Cycling Model
[21] To be able to apply Green-function-based transport diagnostics, it is necessary to construct a linear representation of the silicon cycle. Rather than linearize the siliconcycling model (4) and analyze small-amplitude perturbations, we build an equivalent linear diagnostic model using the same strategy that we demonstrated to be very useful in our analyses of the phosphate cycle Holzer and Primeau, 2013]. Specifically, we replace the nonlinear uptake J ðvðrÞÞ by the linear uptake cðrÞvðrÞ, where the rate coefficients cðrÞ are determined from the nonlinear solution as cðrÞ5J ðvðrÞÞ=vðrÞ and then prescribed for the equivalent linear problem. For steady state, the equivalent problem is thus where L5diagðcÞ.
[22] The key point is that the linear diagnostic model (7) and the nonlinear conditional restoring model (4), with (3), have exactly the same solution, that is, exactly the same silicic-acid concentration, export production, and regeneration fields. Because (7) is linear in v, it is amenable to Green-function-based transport diagnostics that allow us to quantify the preformed and regenerated components of silicic acid, as well as the detailed transport pathways and flow rates of the silicon cycle.

Opal Export and Flux
[23] For every column of water, the vertical flux of sinking opal through depth z must be balanced by the opal The shading corresponds to percentiles. The Nth percentile is defined such that N% of the distribution lies outside the N% contour. Contours with large percentiles thus correspond to high density. HOLZER ET AL.: OCEAN SILICON CYCLE dissolution rate integrated over the column below z because this dissolution is modeled as instantaneous with no horizontal opal fluxes. Hence, the export flux through depth z at every horizontal position (x, y) is given by where the S5SJ ðvÞ is the local rate of opal dissolution per unit volume defined in (5). Figure 2 shows a map of the opal export flux Uðx; y; z c Þ through the base of our model euphotic zone at z c 5273:4 m, together with its zonal integral. The global opal export is dominated by the Southern Ocean, with less intense export flux in the North Pacific and North Atlantic. The Weddell Sea has particularly intense opal export. Although the opal export in tropical regions is much weaker, the tropical export occurs over a much larger area and on a zonal integral both the tropical and Northern Hemisphere high latitudes have peak exports of 1 Tmol/yr/deg, while the Southern Ocean's peak export is about 5 times larger.
[24] The globally integrated opal export through z5 273:4 m is 171 6 31 Tmol/yr, while the opal export through z52133 m is 166 6 24 Tmol/yr, which compares with the estimate of Sarmiento et al. [2007] of 132 Tmol/yr. Uncertainties are simply estimated from the sensitivities with respect to the optimized parameters as described in Appendix A. Estimates of the global opal export from free-running ocean models range from 69 Tmol/yr [Moore et al., 2004] to 185 Tmol/yr [Heinze et al., 2003], with an estimate based on satellite-derived primary production being 101 6 35 Tmol/yr [Dunne et al., 2007]. It is important to keep in mind, however, that in the upper ocean the opal flux attenuates quickly with depth, and is therefore sensitive to the precise definition of the euphotic zone. It is also worth pointing out that the adjoint model used by Sarmiento et al. [2007] used a Martin-type powerlaw profile for the opal flux, with a spatially varying exponent. (When we use a Martin-type flux profile with a globally uniform, but optimized, exponent of b50:27 and a corresponding optimal value of r50:75, the total opal export through z5273:4 is 216 Tmol/yr. The corresponding Si(OH) 4 concentration has a 15% higher quadratic mismatch with observations, but the concentration and export fluxes are qualitatively similar to our solution with the temperature-dependent dissolution length scale.) [25] Figure 3 shows the area-integrated opal flux as a function of depth for each major ocean basin. The boundary of the Southern Ocean for these particular integrals was chosen as 36 S to facilitate comparison with the work of Sarmiento et al. [2007]. It is immediately evident that the Southern Ocean dominates the global export of opal, and that the Southern Ocean contributes to a degree far  Table 1). The other key feature of the opal flux profiles is their decay with depth, which is most rapid in the upper 500 m of the water column, and roughly linear below 1000 m.
[26] Following Sarmiento et al. [2007], we summarize the key information in Figure  Ocean, where opal export is negligible, was excluded from all these fractions for convenience and because the circulation of the Arctic Ocean is not well represented by our model.) These fractions are tabulated in Table 1, where we also provide the values obtained by Sarmiento et al. [2007] for comparison. Our analysis suggests that the Southern Ocean provides 70% of the global opal export, somewhat larger than the 57% in the study by Sarmiento et al. [2007]. Correspondingly, we estimate that the Pacific is responsible for only 18% of the global export, versus 30% in the study  by Sarmiento et al. [2007]. The Southern Ocean contributes nearly 3 times more ( prod i 52:67) to the global opal export than one would expect based on its area if the opal export were uniformly distributed over the global ocean.
[27] A simple metric of the attenuation of the opal flux due to dissolution is the ''preservation fraction'' f

Opal Dissolution Rates
[28] Figure 4 compares the zonally averaged silicic-acid production rate per unit volume with the zonally averaged phosphate remineralization rate per unit volume in the data-assimilated phosphorus cycle of Primeau et al. [2013]. To facilitate a meaningful comparison between the silicicacid and phosphate regeneration patterns, we normalized the regeneration rates by their aphotic-zone (height coordinate z < z c ) volume-averaged values. Note that the contours of Figure 4 are spaced logarithmically. Apart from high regeneration rates in the euphotic zone for both Si(OH) 4 and PO 32 4 , the figure shows that in the aphotic zone the dominant gradients of PO 32 4 regeneration are vertical, while the dominant gradients of Si(OH) 4 regeneration are meridional. The high latitudes stand out as the primary regions of opal dissolution (equivalently silicic-acid production), with the Southern Ocean dominating. Furthermore, while PO 32 4 is regenerated primarily in the thermocline, high-latitude opal dissolution reaches much deeper into the water column, with first-order contributions from the deepest grid boxes (see also below).
[29] To summarize the vertical structure for each basin on a linear scale, Figure 5 shows the silicic-acid and phosphate regeneration rates, horizontally integrated over each basin and normalized by the aphotic-zone globally integrated regeneration rates. These profiles can be thought of as the regeneration rate per unit depth corresponding to the regeneration of one unit of nutrient per unit time in the aphotic zone. The profiles show that while the regeneration rate of silicic acid attenuates with depth through the PO 4 Figure 5. The horizontally area-integrated silicic-acid (blue) and phosphate (green) regeneration rates normalized by the rates globally volume integrated over the aphotic zone. The dashed horizontal lines indicate the aphotic volume-weighted mean regeneration depths for each basin. The dotted profiles show the horizontally integrated regeneration rates without the ''sediment effect,'' that is, without the contribution from the deepest box of each water column (but normalized by the volume integral of the total aphotic regeneration rate, including these contributions). Note the 3 times larger scale for the Southern Ocean.
HOLZER ET AL.: OCEAN SILICON CYCLE thermocline as expected from the temperature-dependent dissolution, below the thermocline the silicic-acid regeneration profiles fall off much more slowly with depth than the phosphate regeneration profiles, particularly in the Southern Ocean. The volume-weighted mean depth of the regeneration rates in the aphotic zone can be defined as z5 R zSd 3 r= R Sd 3 r, where the integrals with volume element d 3 r are performed over the aphotic volume of each basin. These mean depths are indicated as dashed lines on Figure 5. For phosphate, z 600 m in all basins, whereas for silicic acid z is significantly deeper in all basins: 2300 m in the Southern Ocean, 1700 m in the Pacific, 1400 m in the Atlantic, and 1100 m in the Indian Ocean.
[30] The nonmonotonic features of the opal dissolution rate profiles in Figure 5, particularly the increase with depth from 2.5 to 3.5 km, is due to dissolution in the deepest grid boxes. This is also evident in the zonal averages of Figure 4. Because the opal flux at 2000 m depth has attenuated only to 30-50% of its export value (Table 1), a significant portion of the opal flux arrives in the deepest grid boxes where it is dissolved according to equation (5). We find that 43% of the globally integrated aphotic opal dissolution occurs in the deepest grid boxes, compared to only 6% for phosphate. (The deepest grid boxes occupy 12% of the global aphotic volume.) The dissolution in the lowest grid boxes can be considered to occur at least partially in the sediments in the limit of zero residence time in the sediments. Because the remineralization-rate profile of phosphate attenuates more quickly with depth, this ''sediment'' effect is much less pronounced for phosphate, though still visible in Figures 4 and 5.
[31] The dotted curves of Figure 5 are the horizontally integrated regeneration rates with the contribution from the deepest grid boxes excluded (but still normalized by the total globally integrated aphotic regeneration rate) to remove the sediment effect. For silicic acid, the differences in the dotted curves among the basins show the effect of the temperature dependence of the opal dissolution. The Southern Ocean has the deepest reaching profile, consistent with dissolution lengths scales Z (equation (2)) that exceed 4000 m below 1 km depth because of cold polar waters. This deep-reaching profile also explains the large sediment effect for the Southern Ocean. The Pacific profile integrates shallower profiles associated with warm equatorial waters and deeper profiles associated with the colder waters of the North Pacific. The average Pacific Z still exceeds 4000 m below 2 km depth. The Atlantic has the warmest deep waters with the shortest average dissolution lengths at depth (exceeding 4000 m only below 3.7 km depth), and the Indian Ocean has the shortest average dissolution lengths above 1.3 km because of the predominantly warm waters of the upper Indian Ocean. As a consequence, the falloff with depth of the silicic-acid regeneration-rate profile is largest for the Indian Ocean, with the Atlantic profile being intermediate between those of the Indian and Pacific Oceans.
[32] To allow comparison with the dissolution rate estimates of Sarmiento et al. [2007], we computed both the volume averaged and volume-integrated deep dissolution rates for each basin below 2000 m depth. The rates are summarized in Table 2 and agree remarkably well with those of Sarmiento et al. [2007] for the Southern Ocean, Atlantic, and Indian Ocean. However, there is again a discrepancy in the Pacific, where our dissolution rate is 30% lower than that estimated by Sarmiento et al. [2007], as expected from our lower Pacific opal export production (Table 1).

Preformed Versus Regenerated Silicic Acid
[33] At an interior ocean point r, the preformed concentration of silicic acid is the concentration due to direct transport of silicic acid from the surface layer, taken here to be the euphotic zone, z > 273:4 m. The regenerated component of silicic acid is the concentration due to the dissolution of opal, analogous to the remineralized component of the other macro nutrients. We can thus write vðrÞ5v pre ðrÞ1v reg ðrÞ; and compute the concentrations v pre and v reg using our linear diagnostic silicon-cycling model as follows. Because all the silicic acid in the euphotic zone is preformed by definition, we can compute v pre as the response to enforcing the boundary condition v pre 5v in the euphotic zone, where v is the equilibrium solution of the silicon-cycling model. It is also interesting to compute the concentration of silicic acid v pre ðrjX i Þ that is due to the direct fluid transport from only region X i of the euphotic zone. The desired boundary condition for v pre ðrjX i Þ is vðrÞDðr; X i Þ, where Dðr; X i Þ51 if r lies in region X i and zero otherwise. This boundary condition essentially labels silicic acid in X i and removes this label as soon as any other euphotic region is contacted. A numerically convenient way to enforce this boundary condition is to force v pre ðrjX i Þ to vðrÞDðr; X i Þ in the euphotic zone by linear relaxation with a fast time scale, s a . We find that with s a 51 second, the problem is numerically stable with solutions that have converged to the s a 50 limit, within acceptable tolerance. The steady state preformed silicic-acid concentration can thus be determined by solving the linear problem Tv pre ðrjX i Þ52c a ½v pre ðrjX i Þ2vðrÞDðr; X i Þ; where c a diagðc a Þ, with c a ðrÞ51=s a in the euphotic zone, and c a ðrÞ50 in the aphotic zone. If the regions X i tile the global euphotic zone (without overlap), then the total concentration of preformed silicic acid is given v pre ðrÞ5 P i v pre ðrjX i Þ. The total preformed concentration can also directly be obtained as v pre ðrÞ5v pre ðrjXÞ, where X is the entire global euphotic zone.
[34] The regenerated component of silicic acid can be computed as a residual from (9) if the entire global euphotic zone X is considered. To compute the concentration of regenerated silicic acid v reg ðrjX i Þ that was last exported from X i , we mask the silicic-acid uptake field Lv with Dðr; X i Þ and enforce the boundary condition v reg 5 0 over the entire euphotic zone by fast relaxation. Thus, it follows from (7) that v reg ðrjX i Þ satisfies : Tv reg ðrjX i Þ5½2LvðrÞ1SLvðrÞDðr; X i Þ2c a v reg ðrjX i Þ; (11) where v on the right-hand side is again prescribed from the equilibrium solution (and we have neglected the small HOLZER ET AL.: OCEAN SILICON CYCLE geological restoring term of (7)). For regions X i that tile the global euphotic zone, the total concentration of regenerated silicic acid is given by v reg ðrÞ5 P i v reg ðrjX i Þ. In equation (11), the 2LvD term on the right-hand side could be omitted because we are interested in the limit of very large c a , which essentially forces v reg to zero in the euphotic zone regardless of the presence of any other sinks (or sources) in the euphotic zone. We include this term here so that it immediately follows that the sum of (10) and (11), summed over all regions of the euphotic zone, reduces to equation (7) for the total concentration.
[35] We now consider the global inventories of total, preformed, and regenerated silicic acid and phosphate regardless of surface origin (i.e., global origin, X i 5X), as well as partitioned according to the six X i regions defined in Figure 6. The same six regions were used in the analysis of the phosphorus cycle Holzer and Primeau, 2013]. The Southern Ocean is defined to be the ocean south of 38 S. (The sensitivity of key quantities computed for the Southern Ocean to the choice of its northern bounding latitude is quantified in Appendix B.) The SO euphotic zone is further divided into the SAAD (South of Antarctic Divergence) region and the NAAD (North of Antarctic Divergence) region by the latitude of the maximum Ekman divergence. In the SAAD region, surface waters and export production feed into deep and bottom waters, while in the NAAD region they feed into intermediate and mode waters.  and Holzer and Primeau [2013] used the nomenclature ANTZ and SANTZ for the SAAD and NAAD regions, respectively ; we use SAAD and NAAD here to avoid confusion with regions defined relative to ocean fronts.) The preformed and regenerated phosphate concentrations were computed using the analogues of (10) and (11) for the phosphate cycle [see also Primeau et al., 2013].
[36] The partitioning between preformed and regenerated nutrients in Figure 7 shows that silicon-utilizing organisms are much more efficient in stripping out silicic acid from the euphotic zone than they, and other plankton, are in stripping out phosphate: 66% of the silicic acid in the ocean is regenerated, compared to only 39% for phosphate. The bulk of the world ocean's silicic acid (preformed plus regenerated) was last in the euphotic zone in the Southern Ocean, split roughly equally between the SAAD and NAAD regions. However, only 33% of the silicic acid of SAAD origin is regenerated, while 85% of the silicic acid of NAAD origin is regenerated, which points to the fact that there is relatively little biological uptake of silicic acid in the SAAD region compared to the NAAD region. The SAAD region is the world ocean's dominant supplier of preformed silicic acid. The fact that phosphate is stripped out of the Southern Ocean euphotic zone far less efficiently than silicic acid is quantified by the fact that only 4% and 26% of the phosphate originating in the SAAD and NAAD regions, respectively, are regenerated. The sensitivities of the preformed and regenerated nutrient amounts from the NAAD region (and the region to its north) to the choice of the Southern Ocean boundary latitude is quantified in Appendix B.
[37] Figure 7 also shows that the latitudes between 38 S and 40 N are important sources of regenerated nutrients, with contributions of 23% and 69% to the global regenerated silicic-acid and phosphate pools, respectively. Most of the production between 38 S and 40 N occurs in regions of tropical upwelling, so that these latitudes make negligible contributions to the preformed pools of either nutrient. The North-Pacific and North-Atlantic euphotic zones are minor contributors of silicic acid, with a regenerated component of 85% and 54% for the North Pacific and North Atlantic, respectively. In our model, the Arctic Ocean makes negligible contributions to both the silicic-acid and phosphate inventories.
[38] Figure 8a shows zonal averages for each basin of the total silicic-acid concentrations and its preformed and regenerated components (regardless of origin region). Consistent with the SAAD region supplying most of the ocean's  preformed silicic acid, most of the preformed component lies in the deep oceans, with a minor excursion with Antarctic Intermediate Waters in the Atlantic sector. In the Atlantic sector, preformed silicic acid and regenerated silicic acid have a broadly similar structure, with the regenerated component penetrating further north at depth. This similarity points to opal exported in the NAAD region being dissolved in water masses ventilated in the SAAD region. Most of the regenerated silicic acid is found in the deep Pacific, with the highest concentrations in the deep North Pacific. Compared to regenerated phosphate, plotted in Figure 8b, the regenerated silicic acid is found much deeper in the ocean, with much lower relative concentrations near the surface. This is consistent with the remineralization of organic phosphate taking place higher in the water column than the dissolution of opal. That the deep North Pacific contains high concentrations of both regenerated nutrients, albeit with differing vertical profiles, is a reflection of the fact that the deep North Pacific is a holding pen of the world ocean's oldest waters Primeau, 2006, 2008;Primeau et al., 2013]. We note in passing that roughly 15% of the silicic-acid concentrations in the deep North-Pacific plume have recently been postulated to come from crustal fluid venting [Johnson et al., 2006]; such sources of silica external to the ocean are not included in our model. Overall, comparison between the PO 32 4 and Si(OH) 4 patterns of Figure 8 shows again that silicic acid is much more efficiently stripped out of the euphotic zone than phosphate. Correspondingly, the tongue of preformed nutrients in the Atlantic sector associated with mode and intermediate waters is strongly suppressed for silicic acid compared to phosphate.
[39] Given the different spatial distributions of the Si(OH) 4 and PO 32 4 regeneration rates per unit volume, it is not surprising that the silicic-acid versus phosphate scatterplot is far from linear, preventing the definition of a simple, globally applicable stoichiometric ratio. However, the spatial distribution of preformed phosphate originating in a sufficiently small region of the euphotic zone must be proportional to the corresponding preformed silicic-acid distribution, because the preformed component of these nutrients is shaped only by the ventilation and transport of water. A simple relationship between preformed silicic acid and preformed phosphate, regardless of surface origin, does however not exist because the different surface distributions of these nutrients weight their common transport paths differently. This point is illustrated by Figure 9, which shows the Si(OH) 4 versus PO 32 4 joint grid-volumeweighted pdfs both regardless of origin region (global origin) and for NAAD origin. These pdfs bin the tracer values over the entire ocean volume. The global panels show a more compact tracer-tracer relationship for the preformed component, but several mixing lines are superposed. When only nutrients that originated in the NAAD region are considered, the preformed component displays the expected tight linear relationship, while the regenerated component has a highly nonlinear joint pdf due to the differences in the interior regeneration rate distributions. The slope of the preformed-preformed correlation is simply set by the relative abundances of Si(OH) 4 and PO 32 4 in the euphotic zone and varies from region to region.

Path Densities and Southern Ocean Trapping
[40] We now quantify the nutrient transport pathways from and to the Southern Ocean, and from and to the rest of the world oceans in terms of their path densities . The path density considered here is the number of nutrient paths per unit volume that connect last utilization in a specified region X i to subsequent utilization in a specified region X f (denoted as X i ! X f nutrient paths). Equivalently, the X i ! X f nutrient path density at point r is the concentration of nutrient at r that is in transit from last utilization in X i to next utilization in X f because each nutrient molecule undergoing such transport traces out one X i ! X f path. Here, we only consider the path densities regardless of the X i ! X f transit time. The nutrient path density will allow us to quantify how biogenic particle transport injects nutrients into different flow regimes of the advective-eddy-diffusive circulation, and how this differs between the phosphorus and silicon cycles.
[41] To be able to meaningfully compare the path densities of silicic acid and phosphate, we normalize the path density by the total amount of nutrient in the ocean. The normalized, residence-time integrated X i ! X f path density gðr; X i ! X f Þ is thus given by where vðrÞ is the local nutrient concentration, M v is the total number of nutrient moles in the ocean, f # ðrjX i Þ is the fraction of the nutrient at r last taken up biologically in region X i , and f " ðrjX f Þ is the fraction of nutrient at r that will next be taken up in region X f . It is important to recognize that the fractions f # ðrjX i Þ and f " ðrjX f Þ partition the nutrient at location r according to the region of last and next uptake, without regard to whether the nutrients at r are preformed or regenerated. Thus, for example, f # ðrjX i Þ includes nutrients last taken up in X i that subsequently passed through the euphotic zone, thereby becoming preformed en route to r. The fractions f # and f " are efficiently computed for steady flow using matrix inversion in the time-forward and time-reversed adjoint flows, respectively, as detailed by Holzer and Primeau [2013].
[42] The fractions f # and f " are normalized so that at every point The path density g therefore has the normalization : where the integral is over the global ocean volume, and it is assumed that the regions X i and X f each constitute a nonoverlapping tiling of the global euphotic zone, such as the six regions defined in Figure 6. Note that with this normalization the path density has dimensions of inverse volume (units of m 23 ). Further background and theory for the pathdensity transport diagnostic may be found in Primeau, 2006, 2008;Holzer, 2009].
[43] We compute the nutrient path densities for the sets of origin and destination regions fX i g5fX f g5 fSO; N38S g, where SO is the Southern Ocean south of 38 S and N38S is the rest of the world ocean. Figure 10 shows the zonal averages of the uptake-to-uptake path densities for each of the three ocean basins for all four X i ! X f combinations. Summing over all four combinations (columns of plots) and volume integrating over the global ocean yields unity in accord with (13). Figure 10a for the silicon cycle shows that the return paths from the SO back to the SO completely dominate the paths of silicic acid. In the Atlantic sector, the SO ! SO paths are confined largely to the Southern Ocean, but in the Pacific and Indian-Ocean basins, the SO ! SO paths reach far into the Northern Hemisphere with the abyssal circulation. The next most important Si(OH) 4 paths are those from the SO to N38S in the Pacific, and the N38S-to-N38S return paths in the North Pacific at thermocline depths. All other [44] For the same origin and destination regions, Figure  10b shows the uptake-to-uptake path densities for phosphate, which differ dramatically from those of silicate. The phosphate path densities are dominated by paths from outside the SO to the SO (N38S ! SO) in all ocean basins and by N38S ! N38S return paths in the North Pacific. The phosphate path densities also show trapping in SO ! SO paths, but these trapped paths account for a smaller portion of the total phosphate concentration than the N38S ! SO paths.
[45] The dominance of the SO ! SO uptake-to-uptake paths for silicic acid is due to the fact that silicic acid gets efficiently stripped out of the euphotic zone, shortcircuiting the transport of preformed silicic-acid out of the   Figure 10. The SO ! SO, SO ! N38S, N38S ! SO, and N38S ! N38S path densities (a) for silicicacid and (b) for phosphate, normalized by the total amount of each nutrient and zonally averaged for the Atlantic, Pacific, and Indian basins. These X i ! X f paths densities are the number of paths per unit volume connecting last utilization in region X i to next utilization in region X f . For each nutrient, the globally volume-integrated sum over all four X i ! X f cases gives unity.
SO in mode and intermediate waters. The SO-exported opal then dissolves on average at a depth of 2000 m in the deep overturning cell of the Pacific [e.g., Lumpkin and Speer, 2007] and in upwelling circumpolar deep water, which return the regenerated silicic acid to subsequent biological uptake in the SO euphotic zone, thereby trapping the silicate in SO ! SO uptake-to-uptake paths. Phosphate is much less efficiently stripped out of the SO euphotic zone, which produces only 20% of the regenerated global PO 23 4 pool compared to 64% of the regenerated global Si(OH) 4 pool (Figure 7). Consequently, phosphate has a larger preformed component that is available for lateral transport out of the Southern Ocean in mode and intermediate waters [Sarmiento et al., 2004;Palter et al., 2010;Holzer and Primeau, 2013]. The shallower mean phosphate-regeneration depth of 600 m retains more regenerated phosphate in the upper overturning cell of the SO [e.g., Speer et al., 2000], so that some regenerated phosphate is also available for lateral transport out of the Southern Ocean. Therefore, upwelling deep waters return less regenerated phosphate to the SO euphotic zone, where, compared to silicic acid, phosphate is much less likely to be utilized again, leading to relatively weak trapping of phosphate in SO ! SO uptake-to-uptake paths.
[46] We use two metrics to quantify the SO nutrient trapping. The metric most directly relating to Figure 10 is the fraction of the global nutrient inventory that is undergoing X i ! X f transport, which is given by We define the SO trapping metric R trap RðSO ! SOÞ: Table 3 collects R trap computed for the silicon and phosphorus cycles. Roughly 50% of the Si(OH) 4 inventory is trapped in SO ! SO paths between successive uptakes, compared to only 15% for phosphate.
[47] Another useful metric can be constructed by comparing SO ! SO with SO ! anywhere nutrient flow rates. The flow rate J ðX i ! X f Þ with which nutrients that were last utilized in X i are next taken up in X f can simply be computed as [48] We define the SO trapping efficiency as where X is the global euphotic zone, so that the denominator is the SO ! anywhere flow rate.
[49] We find that E trap 5ð9562Þ% for the silicon cycle, compared to E trap 5ð5662Þ% for the phosphorus cycle (Table 3). Because we model the silicon cycle as steady, this means that of every unit of silicon utilized in the SO, 0.95 units return to the SO euphotic zone for their next biological uptake. Thus, there is a probability p trap 5E trap that a silicon atom last utilized in the SO euphotic zone will next be utilized in the SO euphotic zone, and hence a probability p escp 512p trap that this atom will escape being again utilized in the SO to be next utilized in the ocean N38S, instead. The probability of n successive SO utilizations is therefore given by p n trap , from which we calculate the average number of successive SO uptakes as n5 P 1 n51 np n trap = P 1 n51 p n trap 51=p escp . For our optimized silicon cycle, p escp 512E trap 5ð562Þ%, so that we can infer that a silicic acid molecule last taken up in the SO experiences on average n52068 successive utilizations in the SO before being utilized N38S. By contrast, phosphate has a SO-utilization escape probability of (44 6 2)%, and hence undergoes, on average, only 2:360:1 successive utilizations in the SO.
[50] The ratio of R and J immediately provides a transport time scale. Specifically, the ratio C5M v R=J is the ''mean age on uptake'' introduced by Holzer and , that is, the transit time from last uptake in X i to next uptake in X f averaged over the population of nutrients being utilized per unit time in X f (here M v is the global number of nutrient moles as in equation (12)). For silicic acid, the mean age on uptake is 105671 and 8256189 years for SO ! SO and SO ! N38S paths, respectively. For phosphate, the corresponding mean ages on uptake are 394681 and 6006118 years. The shorter mean time between successive uptakes in the SO for silicic acid again points to more efficient utilization than is the case for phosphate: the longer phosphorus SO ! SO mean age on uptake may be thought of as the extra time required on average for phosphate to encounter the plankton that utilizes it. These time scales are the means of distributions with long exponential tails , and the centuries-long SO ! SO, and yet longer SO ! N38S, values of C reflect the fact that exported nutrients on average regenerate deep enough in the water column to access the slow abyssal circulation. The large relative uncertainty in the SO ! SO mean age on uptake reflects the fact that a large portion of the SO ! SO flow rate [ð5:63:Þ310 14 and ð1:160:2Þ310 12 mol=yr for silicon and phosphorus, respectively] is due to rapid recycling within the euphotic zone. This recycling is poorly constrained by the observed nutrient concentrations and the mean age on uptake is particularly sensitive to the fraction r of the production that is regenerated within the euphotic zone (equation (5)) when nutrient uptake is efficient.
[51] The fact that 50% of the ocean's silicic acid is trapped in SO ! SO paths, compared to only 15% for phosphate, implies that a smaller portion of the biological Si(OH) 4 uptake in the rest of the world oceans can be sustained by Si(OH) 4 last utilized in the Southern Ocean than is the case for phosphate. Following Holzer and Primeau [2013], we define the silicon cycle's productivity  (15) and (17).
HOLZER ET AL.: OCEAN SILICON CYCLE teleconnections with a given origin region as the fraction of the uptake that is sustained by nutrients last utilized in the origin region. The zonally averaged fractions of the production supported by nutrients last utilized in each of the six regions of Figure 6 are shown Figure 11a for silicic acid and in Figure 11b for phosphate.
[52] Figure 11 shows that, indeed, these teleconnections in the steady state silicon cycle are generally much weaker than those of the phosphorus cycle. For the silicon cycle, the dominant fraction of the uptake in each region comes from nutrients last utilized in the same region. Less than 10% of the silicic-acid uptake in the North Atlantic is due to Si(OH) 4 last utilized in the Southern Ocean, which compares to over 20% for the phosphorus cycle. Silicic acid last taken up at low latitudes between 38 S and 40 N only contributes a few percent to the Southern Ocean silicic-acid uptake, while nearly 40% of the Southern Ocean uptake of phosphate is fueled by phosphate last taken up between 38 S and 40 N. An interesting exception is the North Pacific. While most of the North-Pacific silicic-acid uptake is fueled by silicic acid last taken up in the North Pacific, and the North-Pacific portion fueled by silicic acid last utilized between 38 S and 40 N is much less than for phosphate, the fraction of the production sustained by silicic acid last taken up in the Southern Ocean is very similar to the phosphate case. This underscores that the North Pacific holds the majority of the ocean's silicic acid (Figure 8a), with about half of the silicic acid in the North Pacific thermocline that is fueling production north of 38 S having been last utilized in the Southern Ocean (see Figure 10, which shows comparable SO ! N38S and N38S ! N38S path densities in the upper North Pacific).

Summary and Conclusions
[53] We optimized a simple conditional restoring-type silicon-cycling model with a temperature-dependent opal dissolution parameterization embedded in a data-assimilated ocean circulation. An equivalent linear silicon-cycling model was constructed so that Green-function-based transport diagnostics could be used to quantify how the biogenic opal export interacts with the global advective-eddydiffusive circulation to shape the ocean's silicon cycle. Our main findings are as follows: [54] 1. The optimized conditional-restoring model produces an opal export through 133 m depth of 166 6 24 Tmol/yr. The Southern Ocean accounts for 70% of this export production, with 50% of the opal flux being preserved to 2000 m depth. The average specific dissolution rate below 2000 m in the Southern Ocean is estimated to be 0.36 6 0.07 lmol kg 21 yr 21 . On global scales, the dominant gradients of the Si(OH) 4 regeneration rates are meridional, in part because of opal dissolution that reaches as deep as the sediments, while the dominant gradients of the PO 32 4 regeneration rates are vertical. The mean depth of the phosphate regeneration rate is 600 m in all basins, while the corresponding mean depth of silicic-acid regeneration ranges from 2300 m in the Southern Ocean to 1100 m in the Indian Ocean.
[55] 2. Silicic acid is stripped out of the euphotic zone far more efficiently than phosphate. This has been quantified by the distribution of the preformed and regenerated nutrients, which is strikingly different for silicic acid than for phosphate: about (66 6 5)% of the global silicic-acid inventory is regenerated, compared to (39 6 7)% for the global phosphate inventory. The dominant regional contributors to the global regenerated silicic-acid inventory are the NAAD region of the Southern Ocean and the tropics, while the dominant contributor to the global preformed silicicacid inventory is the SAAD region of the Southern Ocean. By contrast, regenerated phosphate has its origin primarily in the low latitudes, with the NAAD region only making a secondary contribution.
[56] 3. The efficient utilization of silicic acid combined with its, on average, deep regeneration has profound consequences for the transport paths of silicic acid. Silicic acid is efficiently trapped in the Southern Ocean in the sense that nearly half of the global inventory undergoes a Southern-Ocean-to-Southern-Ocean round trip between successive utilizations. Nearly all the silicic acid last utilized in the Southern Ocean that reemerges in the Southern Ocean euphotic zone is again utilized there. We estimate that on average a silicic-acid molecule that was last utilized in the Southern Ocean has only a ð562Þ% chance of next being utilized north of the Southern Ocean (38 S), experiencing on average 20 6 8 successive utilizations within the Southern Ocean euphotic zone.
[57] 4. Southern Ocean trapping of phosphate in the current state of the ocean is much less pronounced. Phosphate last utilized in the Southern Ocean has a (44 6 2)% probability of escaping having its next utilization in the Southern Ocean, with only (15 6 1)% of the global phosphate inventory trapped in Southern-Ocean-to-Southern-Ocean paths between successive utilizations. As a result of the relatively weak Southern-Ocean trapping of phosphate and the strong trapping of silicic acid, subantarctic mode waters carry a significant amount of preformed phosphate from the Southern Ocean into the Northern Hemisphere [Sarmiento et al., 2004;Williams et al., 2006;Palter et al., 2010;Holzer and Primeau, 2013] but are largely depleted of silicic acid, consistent with the findings of Sarmiento et al. [2004Sarmiento et al. [ , 2007.
[58] 5. As a consequence of the efficient trapping of silicic acid in the Southern Ocean, productivity teleconnections with the Southern Ocean, defined as the fraction of local nutrient uptake fueled by nutrients last taken up in the Southern Ocean, are generally substantially weaker for silicic acid than for phosphate. The exception is the North Pacific, where silicic-acid teleconnections with the Southern Ocean are similar to those for phosphate. This is due to the fact that the bulk of the ocean's silicic acid lies in the North Pacific, with nearly half of the silicic acid in the North-Pacific thermocline having been last taken up in the Southern Ocean and being destined to fuel silicic-acid uptake north of 38 S.
[59] We emphasize that the teleconnections examined here quantify the unperturbed base state of the nutrient cycles. The teleconnections in the response to perturbations in silicic-acid uptake are potentially very different. (For the effects of idealized perturbations on Southern-Ocean phosphate uptake, see the work by Primeau et al. [2013] and Holzer and Primeau [2013].) For example, if Southern-Ocean silicic-acid uptake was to drop, there could be a substantial world-wide response because less silicon would be trapped in Southern-Ocean-to-Southern-Ocean paths, making more silicate available for biological production elsewhere. The response of the silicon cycle to decreased iron limitation in the Southern Ocean has been explored phenomenologically on the basis of isotope data , and in the context of a simple box model , to provide a possible explanation for carbon dioxide levels during glacial times. We plan to quantify the changes in export patterns, transport pathways, and the time scales of adjustment that result from perturbing silicic-acid uptake, in the context of our dataconstrained model in a future study.
[60] While our analysis has shed new light on how the ocean's silicon cycle is shaped by the interaction of the biological sources and sinks with the global circulation, it is important to keep in mind the limitations of our study. The lack of a seasonally varying circulation means that potentially important effects due to the covariance of seasonal diatom blooms and seasonal ventilation cannot be captured. The lack of seasonality also implies that there are unavoidable inconsistencies between the Si(OH) 4 distributions of the annually averaged model and the annually averaged observations. The effect of such inconsistencies is difficult to quantify without a seasonally varying model, but such inconsistencies may well contribute to our inability to reliably constrain the energy scale (temperature T E ) of the dissolution parameterization (2). It is also important to recall that our model applies only to the open ocean; biological production on the poorly resolved continental shelves and in coastal zones has not been modeled, although we expect that these regions, while sustaining high primary production, do not contribute strongly to the global export production. Another key simplification was to assume that the silicon-cycling parameters are spatially uniform. This assumption could in principle be relaxed at the expense of computational complexity.
[61] An interesting open question is what fraction of the carbon export production is due to siliceous plankton. One can empirically determine the ratio of the opal export flux to the total particulate organic carbon flux using observed ratios of the vertical silicic-acid and nitrate gradients across the euphotic zone, together with a bulk N:C ratio [e.g., Sarmiento et al., 2004;Dunne et al., 2007]. In the Southern Ocean these Si:N export ratios can exceed four. The determination of the portion of the organic carbon export carried HOLZER ET AL.: OCEAN SILICON CYCLE by silicate particles therefore not only requires a model of the Si:N ratio (or equivalently the Si:P ratio, assuming constant P:C) with which nutrients are utilized by diatoms stressed by iron and other limitations [e.g., Jin et al., 2006], but also a model of the average differential separation of organic matter from the cell's silicate structures in sinking detritus. Both of these aspects of the link between opal and carbon are beyond the scope of the simple restoring model explored here. In future work, we plan to extend our analysis to a seasonally cyclostationary system, and to explore systematically the coupling between the silicon, iron, and carbon cycles.

Appendix A: Sensitivity-Based Uncertainty Estimates
[62] The uncertainty dX in a quantity X was estimated from the linear sensitivities to the silicon-cycling model parameters r and w f evaluated at the optimal solution: ðdX Þ 2 ð@ r X drÞ 2 1ð@ wf X dw f Þ 2 : (A1) The sensitivities (partial derivatives) were estimated as finite differences of the nonlinear restoring model from the optimal solution. We took dr50:15 and dw f 520 m=day as reasonable uncertainties in the parameters themselves. The key quantities of the global opal export through 133 m depth, and the Southern-Ocean deep dissolution rate, have nearly the same relative sensitivities to a relative change in the parameters : r X @X @r 5277% and w f X @X @w f 527%.
[63] Uncertainties for quantities related to the phosphorus cycle were estimated as in the work of Holzer and Primeau [2013], using the sensitivities with respect to the phosphorus cycle's three parameters: the fraction r of the phosphate uptake allocated to dissolved organic phosphate production in the euphotic zone, the Martin power-law exponent b, and the inverse remineralization time scale j.

Appendix B: Sensitivity to Southern-Ocean Boundary Latitude
[64] Here, we quantify the sensitivity of key Southern-Ocean (SO) quantities to the choice of the northern bounding latitude, / SO , of the region we define to be the Southern Ocean (the union of the SAAD and NAAD regions in Figure 6). Figure B1 shows the amounts of preformed and regenerated silicic acid and phosphate that originate in the euphotic zone of the NAAD region (maximum Ekman divergence to / SO , blue), and in the ocean between / SO and 40 N (/ SO -to-40N , yellow), as a function of / SO . As / SO increases, the NAAD contribution increases and the / SO -to-40N contribution decreases by the same amount. In Figure B1, where the vertical width of the blue and yellow  Figure B1. The amount of preformed and regenerated nutrients that originated in the NAAD region of the Southern Ocean (vertical thickness of the blue band) and in the region between the Southern Ocean and 40 N (vertical thickness of the yellow band) as a function of the latitude separating these two regions. The left plots are for silicic acid, the right plots are for phosphate. The results presented in the main text are for the case where the northern boundary of the NAAD region is at 38 S (black dots).
bands represents the contributions from the NAAD and / SO -to-40N regions, this constraint is visible as the constant width of the sum of the blue and yellow bands. Because most of the silicic-acid uptake lies well south of 38 S, the NAAD// SO -to-40N partition of regenerated silicic acid is relatively insensitive to / SO at / SO 5 38 S. The partition is more sensitive for phosphate, because phosphate uptake  is more broadly distributed and the Southern-Ocean peak uptake is located further north than for silicic-acid. The relative sensitivities (X 21 dX =d/ SO , for quantity X) at 38 S for the silicon cycle are 0.5%/deg and 0.6%/deg, for the regenerated and preformed NAAD contributions, respectively, while the corresponding sensitivities for the phosphorus cycle are 4%/deg and 1%/deg.
[65] Figure B2 shows the trapping metrics E trap and R trap defined by equations (15) and (17) as a function of / SO . The sensitivity to / SO at 38 S is very small for the silicon cycle at 0.09%/deg for E trap and 0.4%/deg for R trap . The corresponding relative sensitivities for the phosphate cycle are again larger at 2%/deg and 5%/deg.  Figure B2. The trapping metrics E trap and R trap (defined in equations (15) and (17)) plotted as a function of the northern bounding latitude of our definition of the Southern Ocean. The left plot shows these metrics for silicic acid, and the right plot, for phosphate. The trapping metrics presented in the main text are for the case where the northern boundary of the Southern Ocean is defined to be 38 S. The values for 38 S are plotted with their error bars. HOLZER ET AL.: OCEAN SILICON CYCLE