www.biogeosciences.net/5/891/2008/ © Author(s) 2008. This work is distributed under the Creative Commons Attribution 3.0 License.

We present a new planktonic foraminifera model developed for the global ocean mixed-layer. The main purpose of the model is to explore the response of planktonic foraminifera to different boundary conditions in the geological past, and to quantify the seasonal bias in foraminifera-based paleoceanographic proxy records. This model is forced with hydrographic data and with biological information taken from an ecosystem model to predict monthly concentrations of the most common planktonic foraminifera species used in paleoceanography: N. pachyderma (sinistral and dextral varieties), G. bulloides, G. ruber (white variety) and G. sacculifer . The sensitivity of each species with respect to temperature (optimal temperature and range of tolerance) is derived from previous sediment-trap studies. Overall, the spatial distribution patterns of most of the species are in agreement with core-top data. N. pachyderma (sin.) is limited to polar regions, N. pachyderma (dex.) and G. bulloides are the most common species in high productivity zones, while G. ruber and G. sacculifer are more abundant in tropical and subtropical oligotrophic waters. For N. pachyderma (sin) and N. pachyderma (dex.), the season of maximum production coincides with that observed in sediment-trap records. Model and sediment-trap data for G. ruber and G. sacculifer show, in general, lower concentrations and less seasonal variability at all sites. A sensitivity experiment suggest that, within the temperature-tolerance range of a species, food availability may be the main parameter controlling its abundance.


Introduction
Planktonic foraminifera are widely used for paleoceanographic reconstructions.The spatial distribution of planktonic foraminifera species is controlled by physiological requirements, feeding preferences and temperature (e.g., Bé and Hamilton, 1967;Bé and Tolderlund, 1971).Shells of planktonic foraminifera extracted from marine sediments serve as an archive of chemical and physical signals that can be used to quantify past environmental conditions, such as temperature (e.g., Pflaumann et al., 1996;Malmgren et al., 2001), ocean stratification (e.g., Mulitza et al., 1997), atmospheric CO 2 concentration (Pearson and Palmer, 2000) and biological productivity (Kiefer, 1998).Past sea-surface temperatures can be estimated by either quantifying differences between modern and fossil species assemblages (e.g., CLIMAP, 1976;Pflaumann et al., 1996;Malmgren et al., 2001), or by analyzing the isotopic or trace-element composition of the calcite in the shell (e.g., Rohling and Cooke, 1999;Lea, 1999).In general, all estimation procedures are based on a correlation between modern environmental condition and assemblage composition or shell chemistry.
Seasonal changes in the flux of planktonic foraminifera are strongly influenced by environmental factors, such as seasurface temperature, the stratification of the water column, and food supply (e.g., Bijma et al., 1990;Ortiz et al., 1995;Watkins et al., 1996;Watkins and Mix, 1998;Eguchi et al., 1999;Schnack-Schiel et al., 2001;King and Howard, 2003a;Morey et al., 2005;Zarić et al., 2005).The seasonality of foraminiferal production is an important factor which has to be taken into account in paleoceanographic interpretations (e.g., Deuser and Ross, 1989;Wefer, 1989;Mulitza et al., 1998;Ganssen and Kroon, 2000;King and Howard, 2001; Published by Copernicus Publications on behalf of the European Geosciences Union.Pflauman et al., 2003;Waelbroeck et al., 2005).Any change in the timing of the seasonal maximum of foraminiferal flux may lead to a bias in estimated paleotemperature.Mulitza et al. (1998) have shown how temperature sensitivity can alter the proxy record in the sediment.Moreover, this differences in seasonality make reconstructed temperatures based on planktonic foraminifera assemblages difficult to compare with those derived by other sea-surface temperature proxies.For example, Niebler et al. (2003) suggested that discrepancies between temperature reconstructions based on foraminifera and alkenones might be due to different ecological and thus seasonal preferences of alkenone producing algae and planktonic foraminifera.Climate change could induce variations in the seasonal succession of the planktonic foraminifera and such variations need to be quantified to correctly interpret corresponding proxy-based reconstructions.
To study the seasonal variations of planktonic foraminifera species we have developed a numerical model for planktonic foraminifera at a species level.Previously, Zarić et al. (2006) developed an statistical model based on hydrographic and productivity data.In contrast to the model of Zarić et al. (2006), we present a dynamic model which, considering ecological processes, calculates the growth rate of the foraminifera population between successive time steps.This study shows model predictions for spatial and temporal distribution patterns of the five most important modern planktonic foraminifera used as SST proxies.

Model setup
The geographical distribution and population density of each planktonic foraminifera species depend on biotic (e.g., food, symbionts) and abiotic factors (e.g., light, temperature).To supply the foraminifera model with ecological information, we run the foraminifera module within an ecosystem model.

Ecosystem model
The employed marine ecosystem model (Moore et al., 2002a) is configured for the global mixed-layer of the ocean.It predicts the distribution of zooplankton, diatoms, diazotrophs and a generic group of phytoplankton (so-called "small phytoplankton").The model considers sinking and non-sinking detrital pools, and carries nitrate, ammonium, phosphate, iron and silicate as nutrients.
The ecosystem model is driven by hydrographic data that are derived from a general ocean circulation model and from climatologies.The forcing data include local processes of turbulent mixing, vertical velocity at the base of the mixed layer, and seasonal mixed-layer entrainment/detrainment. Horizontal advection is not included; thus, there is no lateral exchange between grid points.Since our main interest is the ecosystem in the mixed layer, processes below the surface layer are ignored.
Previously, this two-dimensional model has been validated against a diverse set of field observations from several JGOFS (Joint Global Ocean Flux Study) and historical time series locations (Moore et al., 2002a), satellite observations, and global nutrient climatologies (Moore et al., 2002b).The full list of model terms, parametrizations, resolution, equations and behavior in the global domain is described in detail in Moore et al. (2002a,b) and the code is available at http://usjgofs.whoi.edu/mzweb/jkmoore/areadme.html.

PLAFOM
The planktonic foraminifera model determines the global distribution of the following 5 species: Neogloboquadrina pachyderma (sinistral and dextral coiling varieties), Globigerina bulloides, Globigerinoides ruber (white variety) and Globigerinoides sacculifer.These species have often been considered to be sensitive to sea-surface temperature, and therefore their assemblage can be used to estimate past seasurface temperatures by means of transfer functions.
Each species has a different food preference (Hemleben et al., 1989;Watkins et al., 1996;Schiebel et al., 1997;Watkins and Mix, 1998;Arnold and Parker, 1999).In general, spinose species prefer animal prey such as copepods (Spindler et al., 1984;Caron and Bé, 1984;Hemleben et al., 1989) while non-spinose species are largely herbivorous (Anderson et al., 1979;Spindler et al., 1984;Hemleben et al., 1985Hemleben et al., , 1989)), although in some specimens muscle tissue has been found in food vacuoles (Anderson et al., 1979;Hemleben et al., 1989).Many species also contain algal symbionts that may provide nutrition (Caron et al., 1981;Gastrich, 1987;Ortiz et al., 1995).On a seasonal scale, it is generally assumed that food is the predominant factor affecting the distribution of planktonic foraminifera under favorable temperatures (Ortiz et al., 1995).Planktonic foraminifera appear to respond to the redistribution of nutrients and phytoplankton very quickly, increasing in number of individuals within several days (Schiebel et al., 1995).Information about food availability is obtained from the ecosystem model.In the model, the food sources may be either zooplankton, small phytoplankton, diatoms or organic detritus.
For compatibility with the ecosystem model, the foraminifera model calculates foraminiferal abundance of each species via carbon biomass, the same as the ecosystem model.Since our study is directed to paleotemperature reconstructions, our main interest is in species relative abundances rather than in assessing the absolute biomass.
Accordingly, for each species the change in foraminifera concentration is calculated as follows: where F is the foraminifera carbon concentration, and GGE (gross growth efficiency) is the portion of grazed matter that is incorporated into foraminifera biomass, which we assume Biogeosciences, 5, 891-911, 2008 www.biogeosciences.net/5/891/2008/ to be constant regardless of the food source.TG and ML represent total grazing and mass loss, respectively.

Growth (TG)
The growth rates are determined by available food using a modified form of Michaelis-Menton kinetics (Eq.2), where G max is the maximum grazing rate, g is the half saturation constant for grazing, α is the relative efficiency for grazing in relation to temperature (calibrated from relative abundances), C n represents the concentration of each type of food (diatoms (D), small phytoplankton (SP), zooplankton (Z) or detritus (DR)), and p is the preference for this food (assumed to be invariant in time).The values and units of all parameters are summarized in Table 1.Food requirements vary for the different foraminifera species.Many species of planktonic foraminifera consume a wide variety of zooplankton and phytoplankton prey, and they are capable of a reasonably flexible adaptation to varying trophic regimes.The food of N. pachyderma (sinistral and dextral varieties) consist almost exclusively phytoplankton, commonly diatoms (Hemleben et al., 1989).G. bulloides presents biological characteristics that place it on the border between spinose and non-spinose species; while most spinose species carry algal symbionts, G. bulloides does not (Gastrich, 1987;Hemleben et al., 1989;Schiebel et al., 1997).It is abundant in periods of high phytoplankton productivity (Prell and Curry, 1981;Reynolds and Thunell, 1985;Hemleben et al., 1989) and feeds on algal prey (Lee et al., 1966).G. bulloides is common in mid-latitude and subpolar waters, but it is also present in the subtropical waters of the Indian Ocean.It is generally more abundant in eutrophic waters with high phytoplankton productivity and for this reason it is commonly used as a productivity proxy (Hemleben et al., 1989;Sautter and Thunell, 1989;Ortiz et al., 1995;Guptha and Mohan, 1996;Watkins and Mix, 1998).G. ruber exhibits two varieties; a pink and a white form.The pink variety is limited to the Atlantic Ocean, and we have therefore only modeled the white variety.G. ruber (white) is a spinose species generally found in tropical to subtropical water masses.It hosts dinoflagellate endosymbionts, and feeds mostly on zooplankton, although it has lower zooplankton dependence than other spinose species (Hemleben et al., 1989).The characteristics of bearing spines, utilization of zooplankton prey and symbiotic association are typical of foraminifera adapted to oligotrophic waters.G. sacculifer is also a spinose species hosting dinoflagellate endosymbionts.Culture experiments with G. sacculifer confirm that it depends on zooplankton food (Bé et al., 1981).It is also adapted to low productivity areas, mainly the centers of the oceanic gyres.Watkins et al. (1996) suggested that the adaptation to oligotrophic waters is possible because these foraminifera obtain nutrition from their symbionts.However, the seasonal maximum abundance occurs when productivity in these regions is maximal.To account for adaptation to low productivity regions, we limited the growth of G. ruber and G. sacculifer to regions in which maximum nutrient and chlorophyll concentration does not exceed a threshold value.This is done multiplying "total grazing" (Eq.2) by a hyperbolic tangent function which, using maximum nitrate and chlorophyll concentration as input, identifies low productivity zones.
Maximum grazing rate for the foraminifera (G max ) varies with the food source.Zooplankton carbon concentration is generally much lower than phytoplankton carbon concentration.For this reason, when zooplankton is the food source G max , is set higher than if phytoplankton or detritus are the food source.Thus, under typical food availability conditions, carnivore species can grow as fast as herbivore species.
Based on the observation that most planktonic foraminifera distribution patterns are latitudinal and correlate with temperature, we assume that temperature is the most important physical parameter controlling the distribution of planktonic foraminifera.This is supported by the experimental work of Bijma et al. (1990), which showed evidence for direct temperature control over vital processes.These authors demonstrated that a correlation exists between in vitro temperature tolerance limits and the known natural limits of the species used in their experiments.The tolerance limits of most species are most likely progressive since a departure from optimal growth conditions causes a gradual reduction of vital processes (Arnold and Parker, 1999).Zarić et al. (2005Zarić et al. ( , 2006) ) compiled planktonic foraminiferal fluxes from sediment-trap observations across the World Ocean.They analyzed species sensitivity to temperature by relating fluxes and relative abundances of seven species to sea-surface temperature.Based on this work, we approximate the temperature relation with a normal distribution.Therefore each species exhibits an optimal SST and an SST tolerance range.The growth rate (Eq. 1) is limited by temperature through the parameter α (Eq.3).
The relationship with temperature assumes that the foraminifera concentration at any site is normally distributed, with an optimum temperature where the relative abundance is highest.Away from this optimum temperature the relative abundance decreases until a critical temperature beyond which the species does not occur.This pattern, with a central peak and symmetrical tails, can be approximated by Gaussian distribution (Eq.3).The value of α varies between 0 (out of limit of tolerance) and 1 (optimal temperature).
The parameter n is a arbitrary parameter that scales the values of α between 0 and 1. T opt and T s are the optimum and actual temperatures, respectively, and σ is the tolerance www.biogeosciences.net/5/891/2008/Biogeosciences, 5, 891-911, 2008 k=parameter for the range on temperature depending on the food availability.p(SP)=preference for grazing on small phytoplankton.p(D)=preference for grazing on diatoms.p(Z)=preference for grazing on zooplankton.p(DR)=preference for grazing on detritus.p′=preference for grazing when main food source is missing.G max (SP)=maximum grazing rate when grazing on small phytoplankton (per day).G max (D)=maximum grazing rate when grazing on diatoms (per day).G max (Z)=maximum grazing rate when grazing on zooplankton (per day).G max (DR)=maximum grazing rate when grazing on detritus (per day).g =half-saturation constant for grazing.GGE=portion of grazed matter that is incorporated into foraminifera biomass (Gross Growth Efficiency).pl=quadratic mortality rate coefficient.rl=respiration loss (per day).cl ij =effect of competition of the species i upon the species j .range of a species.Species with small σ are more sensitive to temperature.The values of all parameters for each species are summarized in Table 1.Of the five species, G. ruber (white) and G. sacculifer (both tropical species), together with N. pachyderma (sin.)exhibit the narrowest SST tolerance range.N. pachyderma (sin.) is absent above 23.7 • C ( Zarić et al., 2005).N. pachyderma (sin.) is a polar species and survives even within sea ice (Antarctic), where it feeds on diatoms (Dieckman et al., 1991;Spindler, 1996).N. pachyderma (dex.) and G. bulloides are present almost Biogeosciences, 5, 891-911, 2008 www.biogeosciences.net/5/891/2008/throughout the entire oceanic SST range; however, N. pachyderma (dex.)exhibits a clear preference for intermediate temperatures.For G. bulloides, temperature does not seem to be a controlling factor.It is generally more abundant in eutrophic waters with high phytoplankton productivity, and for this reason it is commonly used as a productivity proxy (Hemleben et al., 1989;Sautter and Thunell, 1989;Ortiz et al., 1995;Guptha and Mohan, 1996;Watkins and Mix, 1998).It has the second largest temperature tolerance, after N. pachyderma (dex.), and does not show a unimodal distribution when flux is plotted versus temperature ( Zarić et al., 2005).G. bulloides comprises at least six different genetic types and exhibits a polymodal distribution pattern (Darling et al., 1999(Darling et al., , 2000;;Stewart et al., 2001;Kucera and Darling, 2002;Darling et al., 2003).Zarić et al. (2005) showed that in the tropical Indian Ocean, G. bulloides is present at higher temperatures than in the Atlantic and Pacific Ocean.In this region, highest abundances of G. bulloides occur at SSTs at which Atlantic as well as Pacific samples show reduced fluxes.Since our study is applied at a global scale, the temperature calibration is based on the preferred temperatures of G. bulloides in the Pacific and Atlantic Ocean.In Eq. 5 we modified the normal distribution for N. pachyderma (dex.) and G. bulloides to accept wider limits under high food availability through the parameter k (see Table 1).

Mass loss (ML)
The mass loss (mortality) equation comprises of three terms representing losses due to natural death rate (respiration loss), predation by higher trophic levels and competition (Eqs.4-8). with Since our model does not include lateral advection, a minimum threshold is needed to preserve the foraminifera population over the winter at high latitudes or during periods with insufficient food supply in regions with high seasonal variability.We set the minimum foraminifera biomass at 0.01 mmolC/m 3 .When the populations reach this minimum level the mortality term is set to zero (Eq.6).Predators specialized on planktonic foraminifera are not known, and therefore, the mortality equation does not explicitly depend upon predator abundance.To represent predation, we choose a quadratic form which depends on foraminiferal biomass itself (Eq.5).This may be interpreted either as predation by a higher trophic level not being explicitly modeled (Steele and Henderson, 1992;Edwards and Yool, 2000).The parameter pl represents the quadratic mortality-rate coefficient, which is used to scale mass loss to grazing.From a bioenergetic perspective, predation is also temperature dependent.Food consumption rates typically increase with increasing temperature; therefore higher trophic levels will exert more predation pressure with increasing temperature (M.Peck, personal communication).The parameter b is used to scale the temperature function between 0 and 1.Note that T sk represents the absolute SST, and the maximum SST (T mk ) assumed in the model corresponds to 303.15 K (30 • C).Death rate refers to natural physiological biomass losses, including respiration (Eq.7).It is a linear term of 6% per day (rl), the same value used by Moore et al. (2002a) for zooplankton.
The presence and activity of one species influences negatively the resource availability for another species, leading to the assumption that competition occurs between different species of foraminifera inhabiting the same regions (Eq.8).In this equation, F i is the concentration of the foraminiferal species exerting competition, cl ij represents the maximum competition pressure of the species i upon the species j (varying from 0 to 1) and d is the e-folding constant, which controls the steepness of the Michaelis-Menton-type equation.

Standard model experiment: grid, forcing and boundary conditions
The foraminifera model is run within the ecosystem model for the global surface ocean, with a longitudinal resolution of 3.6 • , a varying latitudinal resolution (between 1-2 • , with higher resolution near the equator), and a temporal resolution of one month.This corresponds to the resolution of the underlying ecosystem model (Moore et al., 2002a,b).We used the same forcing as Moore et al. (2002a).Mixedlayer temperatures are taken from the World Ocean Atlas 1998 (Conkright et al., 1998), surface shortwave radiation from the ISCCP cloud-cover-corrected dataset (Bishop and Rossow, 1991;Rossow and Schiffer, 1991) and climatological mixed-layer depths from Monterey and Levitus (1997).The minimum mixed-layer depth is set at 25 m.The vertical velocity at the base of mixed layer is derived from the NCAR-3D ocean model (Gent et al., 1998).The turbulent exchange rate at the base of the mixed layer is set to a constant value of 0.15 m/day.Sea-ice coverage was obtained from the EOSDIS NSIDC satellite data (Cavalieri et al., 1990).Atmospheric iron flux was obtained from the dust deposition model study of Tegen andFung (1994, 1995).More details about the forcing can be found in Moore et al. (2002a).
Bottom boundary conditions are the same as for the zooplankton component of the ecosystem model.For all www.biogeosciences.net/5/891/2008/Biogeosciences, 5, 891-911, 2008 foraminifera species we assumed a uniform distribution inside the mixed layer, whereas below the mixed layer the concentration was calculated as a function of the surface concentration and the mixed-layer depth.When the mixed-layer is thin, the foraminiferal concentration below the mixed-layer is set to 75% of the surface concentration.With increasing mixed-layer depth, the concentration below decreases linearly, until reaching the value 0 at a mixed layer depth of 100 m.This is a realistic limit, as the maximum production of the species in question occurs within this depth range (Bé, 1977;Duplessy et al., 1981;Murray, 1991;Watkins and Mix, 1998).

Comparison to core-top data
Since our main interest is to understand the distribution of planktonic foraminifera at geological timescales, we used the Brown University Foraminiferal Database (Prell et al., 1999) to compare our model results with sedimentary faunal assemblages.This database contains core-top planktonic foraminifera counts from 1264 cores across the world ocean.We extended this database with the dataset by Pflaumann et al. (1996), which contains planktonic foraminifera counts for 738 surface sediment samples from the North and South Atlantic; and with another 57 core-top samples from the eastern Indian Ocean (Marchant et al., 1998).For comparison, the relative abundances were recalculated using only the five foraminifera species under consideration.Additionally, the number of individuals was transformed to biomass (mgC/m 3 ) to take into account the size differences of each species.For this transformation, we calculated the volume occupied by the cytoplasm approximating the shape of all the species to a sphere and assuming that all the volume is occupied by the cytoplasm.For the mean size of each species we used sediment-trap data from Peeters et al. (1999).We as-sumed the carbon content of the cytoplasm is 0.089 pgC/µm 3 (Michaels et al., 1995), the same in all species.
To assess the deviation between observed and modeled species distributions, we calculated the root mean squared error (RMSE).For this, the data from each each core-top sample was compared to the nearest model grid point.No averaging was applied to the core-top data.This is justified because the observational data base is identical for all species and our interest is only to test model performance for the five species.

Comparison to sediment-trap data
Several sediment-trap studies were used to compare measured and modeled foraminiferal fluxes (Table 2).Sediment traps show a high temporal resolution and record the flux continuously over several months or years.Because of the sinking speed of foraminiferal shells (150-1300 m/day depending on their weight and size; Takahashi and Bé, 1984), the sediment-trap samples are not significantly affected by dissolution, lateral advection or bioturbation, and therefore can be related directly to modern surface hydrography (e.g., Tedesco and Thunell, 2003;Marchant et al., 2004;Mohiuddin et al., 2004;Zarić et al., 2005).However, due to the short duration of the collecting periods those data may represent local processes of a particular year rather than a long term mean.We used the global database compiled by Zarić et al. (2005Zarić et al. ( , 2006)).This database contains planktonic foraminiferal fluxes calculated from various sediment-trap investigations across the world ocean.To compare the modeled and observed annual distribution of the different planktonic foraminifera species we used those datasets with a minimum collecting period of one year and at least monthly resolution.We extended the database of Zarić et al. (2005Zarić et al. ( , 2006) ) by adding trap data from the northwest Pacific (Oda and Yamasaki, 2005;Xu et al., 2005), Bering Sea (Asahi and Takahashi, 2007), South China Sea (Tian et al., 2005) and Arabian Sea (Schulz et al., 2002).Table 2 summarizes locations, details and references of the sediment-trap studies used in this study.Fig. 1 illustrates locations of the sediment traps.
Sediment-trap studies yield fluxes based on individual shells [ind.m −2 day −1 ] whereas the model provides concentrations [mmolC/m 3 ].To compare model output with observations, we assume that the flux through the water column is proportional to the surface concentration.The objective of our study is to detect relative changes in the seasonal distribution, rather than to assess foraminiferal biomass.

Sensitivity analysis of the parameters
To determine values for biological parameters is difficult as, unlike many chemical or physical parameters, they cannot strictly be regarded as constants.The free parameters have been tuned based on ecological knowledge about different species of foraminifera.In an attempt to assess the sensitivity Biogeosciences, 5, 891-911, 2008 www.biogeosciences.net/5/891/2008/ of the model to the chosen values, we developed a sensitivity analysis of the parameters.The procedure used was similar to other marine plankton models (e.g., Fasham et al., 1990).We kept the parameters that are common to the ecosystem model constant, and modified only the values chosen for the foraminifera module.We run the model with each parameter altered by half and twice the standard value respectively to determine which parameters have the most effect (Table 3).
In case of the parameter p (preference for the food type), we decreased the main food source was decreased by 50%, sharing this part between the other food sources (scaled rela-tive to the original p in experiment 1; and sharing equally between the other food sources during the experiment 2).
To test the effect of competiton we carried out three experiments: first, reducing experienced competition by 50%, a second experiment swithching it off completely and the last experiment increasing competition by 50%.For parameter choice, we compared modelled annual mean realtive abundances and core-top data.The sensitivity was quantified by calculating the change of RMSE between sensitivity experiment and standard run.Sensitivities of the species to each parameter are given in

Spatial distribution patterns
Modeled global distribution patterns of the five foraminifera species are shown together with the corresponding core-top data in Figs.2-6.The model results are expressed as relative abundances as derived from the biomass data.Relative abundances for core-top data consider only the five species included in the model.The global distribution of N. pachyderma (sin.)shows the lowest RMSE, around 9%, while for the remaining the species the error varies between 22% and 25%.
N. pachyderma (dex.) is typical of subpolar to transitional water masses.In the surface sediment samples, N. pachyderma (dex.)shows a very high relative abundance in the North Atlantic Ocean, the Benguela upwelling system, parts of the Southern Ocean and in the equatorial upwelling of the Pacific Ocean.It is also present, although at lower abundance, in the upwelling systems off northwest Africa.The model output shows very high concentrations in the Peru-Chile current and the eastern boundary upwelling systems, as well as south of Iceland, and moderate abundances at mid latitudes (Fig. 3).
Like N. pachyderma (dex.),G. bulloides typically occurs in subpolar and transitional water masses (Bradshaw, 1959;Tolderlund and Bé, 1971;Bé, 1977), and is also found in upwelling areas (Duplessy et al., 1981;Thunell and Reynolds, 1984;Hemleben et al., 1989).Temperature does not seem to be a controlling factor in the distribution of this species, although the exact relationship between environmental parameters and geographical distribution of G. bulloides may be masked by the fact that this species group comprises several distinct genotypes (Darling et al., 1999(Darling et al., , 2000;;Stewart et al., 2001;Kucera and Darling, 2002;Darling et al., 2003).Generally, the abundance of G. bulloides is related to high productivity areas (Prell and Curry, 1981;Bé et al., 1985;Hemleben et al., 1989;Giraudeau, 1993;Watkins and Mix, 1998;Zarić et al., 2005).G. bulloides shows a high relative abundance in the surface sediment samples of the North Atlantic Ocean, the upwelling systems off northwest and southwest Africa, the Southern Ocean, the northern Indian Ocean, and to a lesser extent, the upwelling region off Baja California.The model results show high concentrations of G. bulloides in the subpolar waters of both hemispheres, in the eastern Biogeosciences, 5, 891-911, 2008 www.biogeosciences.net/5/891/2008/Fig. 2. N. pachyderma (sin.)relative abundances (%) from core-top (left) foraminiferal assemblages (Pflaumann et al., 1996;Marchant et al., 1998;Prell et al., 1999) and model output (right).Relative abundances consider only the species included in the model.RMSE is 9%.
boundary currents of the southern hemisphere and in some locations of the Arabian Sea (Fig. 4).
The seafloor record shows high relative abundance of G. ruber (white) in the central North and South Atlantic as well as the South Pacific and less pronounced relative abundance in the South Indian Ocean up to 40 • S. The model output shows a similar pattern with high relative abundances in tropical and subtropical waters of the Atlantic, Pacific and Indian Oceans, and very low relative abundances in upwelling areas (Fig. 5).
G. sacculifer shows a clear preference for high temperatures (optimum of 28 • C) and is absent (or in concentrations ≤10%) below 23 • C ( Zarić et al., 2005).Core-top data show this species is limited to tropical waters, reflecting its narrow temperature tolerance (Fig. 6).The highest abundances occur in surface sediments from the equatorial Pacific and central Indian Ocean.The relative abundance in most of the core-top data from the upwelling region of the Arabian Sea is <10%.The annual mean distribution pattern of G. sacculifer in the model is limited to tropical waters, with highest concentrations in the equatorial Pacific.The model correctly simulates absence or low concentrations (<10%) in eastern upwelling systems as well as in the upwelling area of the Arabian Sea.

Temporal distribution patterns
We used several sediment-trap datasets to assess the modeled seasonal variations in foraminifera abundance.We limit the following comparison between the model output and the sediment-trap data to a few examples .
is forced with climatological data (i.e., long-term averages), and is therefore unable to reproduce interannual variability.
For that reason, we focused on the season with maximum production.In order to compare modeled and observed time series, we picked the season when maximum foraminiferal production occurs (Table 4).When sediment-trap were deployed for more than one year we considered the season in which most maxima occur.
Interannual variability of N. pachyderma (sin.) in all the locations is very high (Fig. 7), but the timing of the signal agrees between observed and predicted data.
The flux of N. pachyderma (dex.)increases during summer (July-October) in northern North Atlantic (Fig. 8a).The seasonal pattern of predicted concentrations corresponds well with the trap record.
Sediment-trap data located at Subantarctic Zone show an increase of the N. pachyderma (dex.)population during the summer (January-February).In accordance with the sediment-trap data, the model results also show the highest concentrations during the summer (Fig. 8b).
Examples for G. bulloides are shown in Fig. 9.At Walvis Ridge, the sediment-trap data reveals a strong seasonality, where maxima occurs in fall (September-November) and in spring (May-June).The model successfully captures this bimodal pattern, with the main bloom occurring in spring.The second example represents a station north of the Kuroshio current in the northwest Pacific (Fig. 9b).At this location the model predicts a small peak during winter (December-January) and the maximum during early summer (May-June).The time-series record also presents this bimodal pattern; nevertheless, model and sediment-trap show better correspondence during the second year.
Fig. 6.G. sacculifer relative abundances (%) from core-top (left) foraminiferal assemblages (Pflaumann et al., 1996;Marchant et al., 1998;Prell et al., 1999) and model output (right).Relative abundances consider only the species included in the model.RMSE is 23%.In general, G. ruber (white) show, less variability in the sediment-trap data (Fig. 10).Seasonal variations in the flux of G. ruber (white) off Cape Blanc, in the Canary Current, are shown for a period of four years (Fig. 10a).The first three years were characterized by a maximum in G. ruber flux during fall (September-October).During fall of 1991, however, the peak most likely occurred after the end of the sedimenttrap deployment.The model predicts a longer bloom (high concentrations from June to December), but the maximum in September coincides with sediment-trap data.The data recorded by the sediment-trap located at the western equatorial Atlantic does not show a clear pattern (Fig. 10b).The first sampling year is characterized by a bimodal pattern, with high fluxes in the austral summer and winter, whereas during the second year the winter bloom was missing.At this site, the model predicts a unimodal pattern with highest fluxes from September to October.
G. sacculifer shows low fluxes in all sediment-trap data used for model validation.In most of the sites where sediment traps were deployed, the model predicts very low concentrations (Fig. 11).
When the model reaches the threshold value set for the minimum population size (0.01 mmolC/m 3 ), the hydrographic component of the model starts to dominate over the population dynamic itself.Therefore, in most of the locations, the model output is too low for comparison.We present examples from the Bay of Bengal and the Arabian Sea (Fig. 11).At both sites, the model predicts little variability in the population of G. sacculifer.(Reynolds and Thunell, 1985;Sautter and Thunell, 1989;Wong et al., 1999); (b) Weddell Sea, 64.91 • S 2.55 • W, in the Southern Ocean (Donner and Wefer, 1994).

Spatio-temporal distribution pattern
We analyzed the model prediction for the temporal variation of G. bulloides in the North Atlantic.Fig. 12 shows the model output of G. bulloides concentrations throughout the year in the North Atlantic.The maximum concentrations occur around 40 • N during spring (March-April) and around 60 • N during summer (June-July), following the phytoplankton bloom in the model.

Sensitivity experiment: spatio-temporal distribution patterns with constant temperature
We ran the foraminifera module with a constant temperature of 12 • C everywhere to test the sensitivity of G. bulloides to other environmental parameters (mainly food availability).The chosen temperature corresponds to the optimal temperature of this species in the model.Fig. 13 shows the spatio-temporal distribution of G. bulloides in the North Atlantic for this experimental run.In general, absolute concentrations of G. bulloides are higher than in the standard run (Fig. 12).However, seasonal pattern does not change substantially from the standard run: During spring the model predicts the highest concentrations in the southern region (around 40 • ) while during summer, the bloom shifts to higher latitudes.Trull et al., 2001;King and Howard, 2003a,b).

Comparison with core-top data
In general, the global distribution patterns of foraminifera species predicted by the model are very close to those expected from core-top data.
The core-top data reflect the integrated flux through the water column, while our model reflects the situation in the mixed layer.As a consequence, some of the discrepancies between model and core-top distributions could be due to the different depth habitats of the species.However, the five species simulated in our model live for most of their life cycle in the upper part of the water-column, thus we expect only a small error at a global scale.In addition, fossil faunal assemblages may be altered by selective dissolution (Berger, 1968;Thunell and Honjo, 1981;Le and Thunell, 1996;Dittert and Henrich, 1999), and by displacement through subsurface currents or bioturbation processes (Bé, 1977;Bé and Hutson, 1977;Boltovskoy, 1994).Since we can not take into account any of these factors, these processes may explain some of the discrepancies between core-top data and the model results.
The global distribution pattern of N. pachyderma (sin.) is very similar to that in the core-top data (Fig. 2).Distinct N. pachyderma (sin.)genotypes were identified by Darling et al. (2004)  Those cryptic species seem to have different environmental preferences (Bauch et al., 2003).Accordingly, the temperature tolerance of the populations in the southern hemisphere is larger than that of North Atlantic population (Darling et al., 2006).For the calibration we did not take different genotypes into account.This causes a discrepancy at the edge of the distribution in the North Atlantic, where the temperature tolerance range is narrower.In fact, the model output agrees slightly better with the core-top data in the Southern Ocean, where the Root Mean Square Error is 8.7% compared to 9.3% in the northern hemisphere.This result reflects the fact that our parametrization is mainly based on the temperature tolerance of the southern population.However, the difference between hemispheres is not large, and treating the different genotypes of N. pachyderma (sin.) as a single ecological group in the model seems justified.
For N. pachyderma (dex.) the model was able to predict the high relative abundances (up to 90%) found in core-top data from the eastern equatorial Pacific and Benguela upwelling systems (Fig. 3).In contrast, in the upwelling system off NW Africa, the model predicts too high relative abundances of N. pachyderma (dex.), while those in North Atlantic (40 • -70 • N) were underestimated.A noticeable discrepancy arises in the equatorial Pacific, north of the eastern boundary upwelling region, where the model underestimates the relative abundance of N. pachyderma (dex.).In  (Fischer and Wefer, 1996;Fischer et al., 1996;Zarić et al., 2005); (b) western equatorial Atlantic, 7.51 • S 28.03 • W (Fischer and Wefer, 1996;Fischer et al., 1996;Zarić et al., 2005).(Guptha et al., 1997).this region, the surface temperature in the model is higher than that typical for N. pachyderma (dex.)(minimum temperatures above 22 • C).It is therefore possible that N. pachyderma (dex.) from these core-tops represent a population living below the mixed layer, as has been described in previous studies (Murray, 1991;Pujol and Grazzini, 1995; www.biogeosciences.net/5/891/2008/Biogeosciences, 5, 891-911, 2008 www.biogeosciences.net/5/891/2008/Biogeosciences, 5, 891-911, 2008 Kuroyanagi and Kawahata, 2004), or that they are expatriated specimens from the upwelling region.The model-generated global pattern of G. bulloides for the Atlantic, Pacific and Southern Oceans agrees well with core-top data (Fig. 4).The model, however, underestimates the abundance of G. bulloides in the northern Indian Ocean.This underestimation could be due to the different G. bulloides genotypes.The two warm water types are found mainly in tropical/subtropical regions, whereas cold water types are found in transitional to subpolar waters (Kucera and Darling, 2002).Zarić et al. (2005) studied the sensitivity of several planktonic foraminifera species to sea-surface temperature and concluded that the population of G. bulloides present in the tropical Indian Ocean comprises mainly the warm-water genotype.Since the parametrization of the model is done at a global scale without specifically considering the warm water type, the increased relative abundance of warm-water G. bulloides in the tropical Indian Ocean cannot be captured by the model.The high concentrations simulated in the Arabian Sea are due to an unrealistically high phytoplankton concentration in the ecosystem model.In addition, the modification of the normal distribution by the introduction of the fooddependent relation (through the parameter k in Eq. 3) allows G. bulloides to grow in tropical waters.The model underestimates relative abundances of G. bulloides in the upwelling regions off northwest Africa, Peru-Chile and Benguela as a result of the overestimation of N. pachyderma (dex.)concentration.This could be due to the higher turnover-rate of G. bulloides in comparison to that of N. pachyderma (dex.), which is not included in this model.
The simulated global distribution pattern of G. ruber (white) is in good agreement with the core-top data.Only in the northern Indian Ocean are abundances somewhat overestimated (Fig. 5).This may be due to the underestimation of G. bulloides in the model as we compare relative abundances.
Both core-top data and model output show that the distribution of G. sacculifer is limited to tropical areas (Fig. 6).The model favorably captures the distribution patterns in the Atlantic and eastern Pacific Oceans.The predicted relative abundances of G. sacculifer in the Indian and western Pacific Oceans are underestimated in the model.The observed distribution shows a wider spatial range than in the model.This could be due to competition exerted by G. ruber.The abundance of this species is overestimated in the mentioned regions and therefore competition is exerting a stronger influence on the other species.

Comparison with sediment-trap data
The comparison between model predictions and sedimenttrap records bears several difficulties.Most of time series only represent short sampling periods (single or a few years).Sediment traps that were deployed for more than one year show pronounced interannual variability that is not captured by the model due to the climatological forcing.Moreover, the majority of sediment-traps have been deployed close to the coast (Fig. 1), where environmental conditions and ecological successions are more complex than in the open ocean.The lack of sediment-trap data in open-ocean settings is an obstacle for a global scale comparison.Nevertheless, in general, the seasonal pattern of species concentrations in the model is similar to the sediment-trap records (Table . 4).Simulated variations of N. pachyderma (sin.)correlate well with sediment-trap data.In 80% of the cases, the preadicted season of maximum production coincides with observational data.In the model N. pachyderma (sin.)lives in polar/subpolar waters, and the maximum production of this species occurs during a short period in summer, together with the phytoplankton bloom.For N. pachyderma (dex.), the season of maximum production varies between locations.However, model prediction and sediment-trap data coincide in 75% of the cases (Table 4).For the stations in Peru-Chile current and northwest Pacific, model prediction and sedimenttrap data differ substantially.However, in the former location, the sampling period is limited to a single year and the data shows an irregular pattern without a distinct seasonal peak, whereas in the northwest Pacific the interannual variability in sediment-trap data is very high.It is likely that at these particular locations the sediment traps do not reflect the mean long-term flux pattern.
For most of the locations, the season of maximum production of G. bulloides simulated by the model does not coincide with the observations.At Ocean Station PAPA (northeast Pacific) the modeled seasonal peaks are delayed with respect to sediment-trap data, whereas in the Peru-Chile current and off Cape Blanc, the peak occurs too early.At the Walvis Ridge, close to the coastal upwelling zone, the model predicts successfully the bimodal pattern (Fig. 9a).However, when considering the absolute flux maximum in each year, the season when it occurs varies throughout the sediment-trap record.At Ryukyu Islands (northwest Pacific) two sediment traps were deployed during the same period, but maximum production seasons recorded in both traps are different.The differences observed in the sediment-trap data highlights difficulties in comparing model predictions and observational data.
The season of maximum production of G. ruber (white) in the model corresponds to the recorded data in 6 of 8 stations.However, in the western Atlantic, the model produces highest concentrations too early in the year (Table 4).On the other hand, the variations in the sediment-trap data are very small and seem to occur randomly.
Only very few data are available to compare seasonal variations of G. sacculifer.Model and sediment-trap data show, in general, lower concentrations than the other four species and little seasonal variability in all sites.This is not surprising considering that G. sacculifer is limited to tropical waters, with small seasonal temperature range.Biogeosciences, 5, 891-911, 2008 www.biogeosciences.net/5/891/2008/

Sensitivity analysis
In an attempt to assess the sensitivity of the model to the chosen parameter values, we performed a sensitivity analysis of the parameters.The procedure used was similar to other marine plankton models (e.g., Fasham et al., 1990).We kept the parameters that are common to the ecosystem model constant, and modified only the values chosen for the foraminifera module.We run the model with each parameter altered by half and twice the standard value respectively to determine which parameters have the most effect (Table 3).
The sensitivity was quantified by calculating the change of RMSE between the sensitivity experiment and the standard run.The results (Table 3) indicate that none of the parameters lead to uniform changes for all species.The model seems to be more robust for G. sacculifer than for other species.In several experiments, the error between model and core-top data decreases.This occurs because the standard parametrization is based on ecological data compiled from literature rather than "tuned" to obtain a better fit.Removing competition generates a general increase of RMSE.Not surprisingly, the temperature tolerance range (σ ) seems to be the most sensitive parameter.

Model experiment with constant mixed-layer temperature
When the foraminifera model is run with a constant temperature of 12 • C, G. bulloides in the North Atlantic still showed highest concentrations at low latitudes during spring and maximum concentrations at higher latitudes in June, linked to the seasonal migration of the phytoplankton bloom (Fig. 13).This indicates that temperature is not the only controlling factor, but that food supply plays an important role in the temporal distribution pattern of this species.The experiment confirms the results of Ganssen and Kroon (2000), who from isotopic studies on North Atlantic surface sediments, concluded that G. bulloides reflects temperatures of a northward migrating spring bloom.

Summary and conclusions
A global model has been developed that predicts monthly planktonic foraminifera concentrations for N. pachyderma (sin.),N. pachyderma (dex.),G. bulloides, G. ruber (white) and G. sacculifer.It is a nonlinear dynamic model simulating growth rate of foraminifera populations using information from an underlying ecosystem model (Moore et al., 2002a).
The model aims at predicting the distribution of planktonic foraminifera at geological timescales.Overall, the global distribution patterns of the predicted species are similar to core-top data.
Modeled seasonal variations overall agree with sedimenttrap records for most of the locations, although the compar-ison is hampered by interannual variability not captured by the model.
A sensitivity experiment using a constant temperature of 12 • C indicates that food availability (primary production in the case of G. bulloides) is an important factor controlling the distribution of some species.
Our model provides a tool that will contribute to better assessing how changing environmental conditions in the geological past affected the distribution of foraminifera in space and time.
Quantitative data and a better knowledge of ecological process from laboratory and field studies are essential for further improvement of the current model.Results may also be improved by including additional information, such as different classes of zooplankton, or by explicitly resolving depth.

Fig. 1 .
Fig. 1.Locations of the sediment-trap stations used to compare measured and modeled foraminiferal fluxes.See Table2for details.
Fig. 11.Comparison of measured fluxes of G. sacculifer in sediment-traps (orange bars) vs. modeled abundances (blue lines).Symbols and layout of the graphs are the same as in Fig. 7. (a) central Bay of Bengal, 13.15 • N 89.35 • E (Guptha et al., 1997).

Fig. 12 .Fig. 13 .
Fig. 12. Modeled monthly concentrations of G. bulloides in the North Atlantic in the standard run.

Table 3 .
Sensitivity analysis to the parameter values: Reduction of parameter values and resulting change of RMSE between the model and core-top relative abundances (RMSE sensitivity experiment minus RMSE standard run) Experiment Parameter change N. pachyderma (sin.)N. pachyderma (dex.)G.bulloides G.ruber (white) G. sacculifer

Table 4 .
Season with maximum foraminiferal production at each site in Fig.1.Empty cells denote if species is absent either in sediment-trap data or the model output.