Ecologists are increasingly using stable isotope analysis to inform questions about variation in resource and habitat use from the individual to community level. In this study we investigate data sets from 2 California sea otter (Enhydra lutris nereis) populations to illustrate the advantages and potential pitfalls of applying various statistical and quantitative approaches to isotopic data. We have subdivided these tools, or metrics, into 3 categories: IsoSpace metrics, stable isotope mixing models, and DietSpace metrics. IsoSpace metrics are used to quantify the spatial attributes of isotopic data that are typically presented in bivariate (e.g., δ13C versus δ15N) 2-dimensional space. We review IsoSpace metrics currently in use and present a technique by which uncertainty can be included to calculate the convex hull area of consumers or prey, or both. We then apply a Bayesian-based mixing model to quantify the proportion of potential dietary sources to the diet of each sea otter population and compare this to observational foraging data. Finally, we assess individual dietary specialization by comparing a previously published technique, variance components analysis, to 2 novel DietSpace metrics that are based on mixing model output. As the use of stable isotope analysis in ecology continues to grow, the field will need a set of quantitative tools for assessing isotopic variance at the individual to community level. Along with recent advances in Bayesian-based mixing models, we hope that the IsoSpace and DietSpace metrics described here will provide another set of interpretive tools for ecologists.
Stable isotope analysis has rapidly transitioned from a novel technique of limited interest to one of the most valuable tools in an ecologist's tool set, and the number of published papers that utilize isotopic approaches is growing at an exponential rate. Stable isotope analysis has been adopted by nearly every subdiscipline of ecology because it allows scientists to trace resources within and between animals, plants, and microbes, at scales ranging from the individual to the community level. Recent reviews of the subject summarize how isotopic data can be used to evaluate information on resource or habitat use, or both, that would help define an organism's niche in a traditional Grinnellian (Grinnell 1917), Eltonian (Elton 1927), or Hutchinsonian (Hutchinson 1957) sense. Some have even gone so far as to utilize the term isotopic niche (Flaherty and Ben-David 2010; Martínez del Rio et al. 2009a; Newsome et al. 2007) because isotopic axes provide information on the bionomic and scenopoetic aspects of the niche (Hutchinson 1978) and arguably can be as informative as many other environmental variables traditionally used to define niche hypervolumes. A clear distinction, however, must be made between the isotopic and realized niche space. Only with the conversion of isotopic data to numerical estimates of resource or habitat use, or both, via mixing models can traditional components of an organism's niche be evaluated with isotopic tools.
As one might expect of a tool that was adopted quickly and with little experimental groundwork, there are a number of unresolved issues concerning the use of stable isotopes as intrinsic ecological tools. In general, these issues center around ecophysiological and methodological considerations that in some (but not all) cases can have a major effect on how stable isotope data are interpreted. For example, a recent wave of controlled feeding experiments examined isotopic incorporation (or turnover) rates for a wide variety of tissues and taxonomic levels (Martínez del Rio et al. 2009b), which is essential for accurate interpretation of applied studies. Likewise, recent laboratory and field-based studies (Gaye-Siessegger et al. 2003; Newsome et al. 2010; Vander Zanden and Rasmussen 2001) have investigated potential mechanisms responsible for variation in trophic discrimination factors, or the difference in isotopic composition between a consumer's tissues and its diet (Δtissue−diet), that have important implications for examining diet or trophic level, or both, of individuals or populations, which is the most common use of stable isotope analysis in ecology. Methodological considerations include sample pretreatment, such as the proper preparation of tissues for hydrogen isotope (δD) analysis (Bowen et al. 2005) and whether consumer and prey tissues should be lipid-extracted prior to carbon isotope (δ13C) analysis (Newsome et al. 2010; Post 2002; Ricca et al. 2007). As the use of stable isotope analysis in ecology continues to grow, these issues and others currently unrecognized, present substantial challenges to the research community. In our opinion, these challenges represent an intriguing opportunity because they transcend the disciplines of ecology and physiology and require the need for both experimental and applied work.
Another methodological issue concerns the question of how isotopic data should be treated from a statistical and interpretive standpoint, which is quickly becoming a rich source of literature. Most of this work has focused on various forms of stable isotope mixing models that can be used to determine source proportions in consumer diets (Moore and Semmens 2008; Parnell et al. 2010; Phillips and Gregg 2003; Phillips and Koch 2002); on statistically based interpretation of δD data that is used to assess movement and migratory patterns (Farmer et al. 2008; Wunder and Norris 2008); on the use of single- versus multiple-compartment models for evaluating isotopic incorporation rates (Carleton et al. 2008; Cerling et al. 2007); and on the use of spatial metrics to characterize community-level variation in trophic structure across space and time (Layman et al. 2007; Turner et al. 2010). In this paper we address the 1st and 4th topics, and briefly summarize here the various forms of mixing models, the utility of spatial metrics, and the type of information these 2 approaches can provide; see the “Materials and Methods” for a more detailed description.
Initially, mixing models used a linear framework to determine a unique mathematical solution for the relative contributions of n + 1 prey sources using n isotope systems (e.g., δ13C or δ15N, or both). Later iterations of these models allowed users to compensate for differences in elemental concentrations among potential sources (Phillips and Koch 2002), which have important implications for the interpretation of isotopic data derived from omnivores that consume resources of varying quality (i.e., nitrogen content). For generalists that consume a wide variety of potential food sources, however, these models were sometimes difficult to use in practice. In response, Phillips and Gregg (2003) produced IsoSource, which generates a frequency distribution of potential proportions for greater than n + 1 sources when only using n isotope systems. In the past few years, Bayesian-based models such as MixSIR (Moore and Semmens 2008) and SIAR (Parnell et al. 2010) have been developed that build on the capabilities of IsoSource and allow users to include information regarding isotopic variation in potential sources and trophic discrimination factors, as well as input prior information on resource or habitat use obtained from other types of data (e.g., scat analysis or observation).
Mixing models are useful tools for converting isotopic data into a form that can be directly compared to traditional ecological information. In essence, these models provide estimates of trophic interaction strengths, or interaction distributions in the case of IsoSource or Bayesian-based approaches, between consumers and their prey, and therefore have greater ecological traction than comparisons of isotopic data presented in 2-dimensional bivariate space (e.g., δ13C versus δ15N). But even Bayesian-based mixing models can be cumbersome, especially when applied to generalist species with diverse diets, in situations with an even distribution of resources in the isotopic mixing space, or when making community-level comparisons of isotopic variation among consumers. In response, ecologists have started using spatial metrics to quantify various aspects of their data in the bivariate δ13C versus δ15N framework that is often used to present isotopic data (Layman et al. 2007, in press; Turner et al. 2010). Spatial metrics were originally used by paleontologists to study the evolution of morphospace in the fossil record (Foote 1990; Gould 1991), but were quickly adopted by ecologists to quantify differences in form and function among extant groups in an evolutionary context (Watters 1991; Winemiller 1991).
A recent concept paper by Layman et al. (2007) focused on how various spatial metrics, including convex hull area (CHA), nearest-neighbor distance (NND), and distance to centroid (DC) can be used to quantify community-wide measures of trophic structure. For example, the CHA of δ13C and δ15N values in bivariate space for a population (or species) can provide an estimate of dietary diversity; see the “Materials and Methods” for detailed descriptions of how these metrics are calculated and the kind of information they provide. A sharp critique of this approach (Hoeinghaus and Zeug 2008) highlighted the need to control for isotopic variation among potential sources available to consumers that occupy different habitats (Matthews and Mazumder 2004) and concluded that the conversion of isotopic data into source proportions via mixing models was the best way to address this problem (Newsome et al. 2007).
Here we show that when properly applied, spatial metrics are an intuitive set of tools for evaluating dietary variation, trophic structure, and habitat use at not only the level of the community, but also at the population and individual scale. We draw on existing spatial metrics, which we call IsoSpace metrics, and offer ways in which these approaches can control for isotopic variation among sources. We also offer a few novel dietary metrics, which we call DietSpace metrics, that are useful for assessing dietary specialization and complementary to other approaches used in the literature (e.g., variance components analysis). Our paper utilizes published isotopic data sets for 2 California sea otter (Enhydra lutris nereis) populations (Newsome et al. 2009, 2010) to highlight both the advantages and caveats of using these tools to interpret isotopic data. Although the sea otter data sets are particularly useful examples because they are paired with extensive observational data on diet composition, the tools described here could easily be applied to isotopic data from other taxonomic groups, including microbes and plants.
Materials and Methods
Insight into isotopic systems can be gained with the use of metrics that quantify aspects of isotopic bivariate space. IsoSpace is typically defined by δ13C and δ15N values of both the consumer and its potential resources (prey), and is often assumed to be equivalent to a consumer's resource isotopic niche space (Martínez del Rio et al. 2009a; Newsome et al. 2007). Because isotopic values are normalized to account for trophic discrimination, a consumer will fall within the isotopic space defined by those resources it has consumed (Fig. 1). The interior of this space is often referred to a prey or mixing space. Characteristics that are descriptive of both ecological and environmental processes that underlie these isotopic values can be quantified by Euclidean measurements of the various shapes and distances within this bivariate space. Although previous investigations have established frameworks for evaluating consumer-level differences when the isotopic values of resources do not change (Turner et al. 2010), such approaches have not incorporated isotopic variance in potential prey available to different populations that occupy different food webs. Here we briefly describe these measurements, as well as procedures by which variance can be incorporated. In contrast to a mixing model approach, note that consumer isotopic data do not have to be trophic-corrected prior to the calculation of the spatial metrics described below.
The CHA is the simplest area that can be drawn around the outermost coordinates that define the mixing space; interior sources are excluded and will not be referred to here. A mixing space can be subdivided into a series of triangular components (Fig. 1). The CHA can be calculated by the sum of the triangular areas within the convex hull; triangles are drawn between a single arbitrary source and each pair of nearest-neighbor sources that define the remaining hull. If vectors are labeled clockwise (P0–P5 in Fig. 1) the CHA can be calculated:
where ν;i and ν;i + 1 are vectors defined as in Fig. 1, and |ν;i × ν;i + 1| is the magnitude of the cross product of the 2 vectors (equation 1). Because potential sources in IsoSpace are typically distributions of values, it is important to take variance into account. Although more rigorous methods may be employed, here we use a simple algorithm to incorporate source variation into CHA measurements. To incorporate variance into CHA measurements, coordinates are randomly chosen from each unique source distribution; from this series the CHA is calculated. This process is iterated such that a representative sample of potential CHAs is quantified across all source distributions. For the mixing spaces that we evaluate here, 1 × 105 iterations is sufficient to accurately calculate the CHA, although the number of iterations required increases with greater isotopic variability of the convex hull, which increases the number of potential hull shapes. In principle, the CHA can be extrapolated to more than 2 dimensions (i.e., when using more than 2 isotopic tracers) such that the measurement would define a hypervolume, although this quickly becomes computationally expensive and is beyond the scope of this paper. We limit our discussion to a 2-dimensional IsoSpace for simplicity.
The NND is defined by the minimum Euclidean distance between an isotopic coordinate relative to all other coordinates in a set. In the 2-dimensional IsoSpace defined above, a coordinate is defined by its δ13C and δ15N values. We note that the Euclidean distance is not limited with respect to the number of isotopes (dimensions) that are used. NNDs can be calculated for a set of isotopic coordinates that comprise individuals within a population of a single species, or across multiple species. The NND provides information regarding the clustering of points within a set. If the NND of a set of coordinates has low variance, the isotopic values of a group are distributed across IsoSpace homogeneously; if NND values are highly variable, the set of isotopic values is heterogeneously distributed (e.g., clustered or overdispersed).
The DC metric is defined as the Euclidean distance between an isotopic coordinate and a predetermined central coordinate, or centroid. In IsoSpace, the centroid coordinate is defined by the average of each respective isotopic tracer (e.g., δ13C) across all sources. The distance of each source to the centroid provides information regarding the distribution of data points in bivariate space. Similar to the NND measurement, low variance of the DC implies a convex hull that is more circular in 2-dimensional space, and data points are distributed along the periphery of the convex hull. High variance in DC and NND implies that data points are more evenly distributed in bivariate space. The relative distribution of sources (i.e., prey) in bivariate space also has implications for the use of nonlinear mixing models because they tend to produce more defined results if potential sources lie on the periphery of the mixing space (see below). For the above metrics, we used a Welch's 2-way t-test to assess whether differences in mean values between sites were statistically significant.
Stable isotope mixing models
Mixing models are designed to determine the contributions of a given set of resources to a consumer's diet. Traditionally, such models have been limited to linear approaches, where the number of potential sources (i.e., prey) must be less than or equal to the number of isotopic tracers + 1; in such a situation, a unique analytical solution always exists. For example, a 3-source mixing space defined by 2 isotopic systems is given by the mass balance equations (equation 2):
where δ13C and δ15N are the respective isotopic values of carbon and nitrogen in δ notation (where δ = 1,000[(Rsample/Rstandard) − 1], and R = either 13C/12C or 15N/14N) for both the mix (m; i.e., the consumer) and sources (i; i.e., the resources or prey), whereas fi denotes the contribution of each source i to the mix. When n = 3, all proportional contributions can be solved analytically because there are 3 equations and 3 unknowns (Phillips et al. 2005).
Most ecological scenarios are more complex and involve many more sources than the number of isotopic tracers utilized in the study. In addition, there is both natural variability and error associated with isotopic measurements that cannot be included in the above framework. To address these issues, a numerical approximation called IsoSource was developed to determine proportional source contributions (Phillips and Gregg 2003). These numerical tools allowed the range of proportional contribution-to-diet values to be determined, even if the number of sources was larger than the number of isotopic tracers + 1, where no unique mathematical solution is possible. Other approaches, such as binning sources with similar isotopic values, can be utilized to further decrease the number of potential sources and increase the accuracy of the results (Phillips et al. 2005; Ward et al. 2011). Importantly, although ranges of contribution-to-diet values can be calculated with these numerical procedures, all values within the ranges are equally likely, thereby limiting interpretations when ranges are large.
Bayesian isotope mixing models were developed to cope with many potential sources, the uncertainties inherent in isotopic measurements, and the incorporation of prior knowledge. This approach results in an accurate quantification of the uncertainty that characterizes scenarios with many more sources than isotopic tracers being utilized, as well as both measurement and discrimination uncertainty. As such, the use of a Bayesian framework results in true posterior probability distributions of the potential contributions of each source to a mix.
Current Bayesian mixing models employ a sampling-importance-resampling or Markov chain Monte Carlo approach to determine the likelihood of potential source contributions to a mix. In general, for each source, a random proportional contribution vector is proposed (fq; where fi elements in the vector fq sum to unity). From this proposed vector, the mean and standard deviation are calculated, and the likelihood of the mixture, given these parameters, is determined:
where x represents the isotopic data describing the mix (a vector where each element is a set of isotopic measurements of the consumer), xkj is the value of the jth isotope of the kth element of the mix, μ ˆj is the mean of jth isotope of the proposal, and σ ˆj is the standard deviation of the jth isotope of the proposal (equation 3). The posterior probability is then calculated:
where L(x|fq) is as described above, P(fq) is the probability of the given proposal based on prior information, and the denominator is a normalizing constant (equation 4). According to the sampling-importance-resampling or Markov chain Monte Carlo algorithms, proposals are generated randomly, and a given proposal is accepted if the unnormalized posterior probability is higher than the previously proposed unnormalized posterior probability. As such, the most likely contribution is iteratively approached across the likelihood space. The final posterior probabilities are typically considered robust if there have been ≥1,000 accepted, unduplicated, contribution-to-diet proposals (see Moore and Semmens  for details). The output of Bayesian isotope mixing models appears similar to the simulations of numerical linear models; however, they are not equivalent. Values within a range given by Bayesian mixing models have associated probability densities, whereas values within the range given by linear mixing models all have the same probability (i.e., a uniform distribution). Recent advances in Bayesian mixing models have resulted in approaches that assess hierarchies of isotopic data (e.g., individuals within populations within species—Semmens et al. 2009), as well as sophisticated methods for binning sources (Ward et al. 2011).
For this study, we used previously published δ13C and δ15N data for 2 California sea otter (E. l. nereis) populations from San Nicolas Island (SNI) and Monterey Bay (MB) and their respective prey; see Newsome et al. (2009, 2010) for details on sampling strategy, tissue pretreatment, and isotopic analyses. Individual means and standard deviations were calculated from the subsampled whisker segments for each sea otter. Trophic discrimination factors of 2.0‰ and 3.5‰ were used for δ13C and δ15N, respectively (Newsome et al. 2010) in the mixing models for the SNI population. For MB, we used slightly different trophic discrimination factors of 2.5‰ and 3.5‰ for δ13C and δ15N, respectively (Newsome et al. 2009).
IsoSpace metrics were assessed for both SNI and MB sea otters, as well as their potential prey. We then calculated the proportional contribution of each prey to the diets of SNI and MB sea otters at the individual level with the Bayesian isotope mixing model MixSIR. Because the multiple measurements obtained from individual sea otter whiskers cannot be considered independent, a parameterized bootstrapping procedure was employed to incorporate within-individual variance into the model. By doing so, we were able to more accurately include variance associated with individual sea otter foraging behaviors while maintaining assumptions intrinsic to the model. Posterior distributions for the SNI and MB sea otter populations were obtained by bootstrapping and pooling individual sea otter MixSIR results such that each individual contributed equally to the final posterior probability distributions. Such an analysis can also be implemented with a hierarchical stable isotope mixing model (Semmens et al. 2009).
In its simplest incarnation, a consumer's diet is represented as a vector f = (f1, f2, …, fn) , such that each element of the vector (fi) represents the proportional contribution of a prey source (e.g., mixing model output), and n is the total number of prey. As stated previously, f must sum to unity. With the use of a Bayesian isotope mixing model, the diet of a consumer is defined by a posterior probability distribution, represented by a series of numerically calculated vectors, rather than a single vector. Accordingly, the distribution of each source quantifies the error associated with the isotopic measurements of both the consumer and its resources or the natural ecological variability of the consumer's diet, or both. Here we present a useful method by which to measure the degree of specialization for an individual consumer and compare groups of consumers (or populations). Note that in this context we define “specialization” in the classical sense of niche specialization, or the degree to which a consumer relies on a subset of prey, with respect to the total number of available prey, and distinguish this concept from individual diet specialization, which we discuss below.
The quantification of dietary specialization is particularly straightforward if the source contribution-to-diet probabilities are known. Here we establish a dietary Euclidean space with as many dimensions as prey, and define the centroid as an ultrageneralist consumer (Fig. 2A). For example, if a consumer's diet consists of n prey items, the centroid would be defined by the coordinate γ = (1/n,1/n,1/n,1/n,1/n), such that the consumer is an ultrageneralist, where every prey contributes equally. By contrast, we define a 2nd point in the Euclidean space that defines an ultraspecialist consumer by the coordinate φ; = (1,0,0,0,0), such that only 1 source is consumed. We note that for this metric it does not matter which element of the coordinate φ; has a value of 1 (e.g., (1,0,0,0,0) is equivalent to (0,1,0,0,0)). We can now define the distance from the ultrageneralist centroid of a consumer represented by a particular dietary coordinate f, relative to an ultraspecialist, and across all proposed contribution-to-diet vectors (equation 5), such that:
Thus the degree of dietary specialization at the population level (ε) varies between 0 and 1, where a value of 0 denotes the ultrageneralist consumer and a value of 1 denotes the ultraspecialist consumer (Fig. 2A). Because we have normalized this metric to the distance between the ultraspecialist and ultrageneralist, it is comparable across consumers with different numbers of potential prey, though it does not distinguish between consumers that specialize on different subsets of prey. When the diet of a consumer is quantified by probability distributions, ε also will be defined by a probability distribution, because each proposed dietary vector has an associated ε value. Here we use this metric to analyze output from a Bayesian isotope mixing model (MixSIR). Because mixing model results are expressed as a series of contribution to diet vectors, ε can be calculated for each vector independently such that a distribution of ε values is obtained. Furthermore, ε values for individuals can be pooled to obtain niche specialization values for a population.
It is often convenient to compare the dietary habits of 2 or more consumers, or an individual consumer with the mean dietary habits of its population. Again we borrow a well-known relationship from linear algebra such that pairwise dietary comparisons can be easily made. We note that this and similar metrics often have been used to compare species' niches and even dietary input, although, to our knowledge, it has not been used in the context of isotopic data or mixing model output (Bolnick et al. 2002; Kohn and Riggs 1982; Smith et al. 1990; Tinker et al. 2008). Because a consumer's diet can be thought of as a unique vector in diet-space, the angle that exists between 2 dietary vectors will define their relative similarity. As such, we define the dietary similarity index:
where f1 and f2 are vectors composed of proportional prey contribution-to-diet values for consumer 1 and 2, respectivey (equation 6; Fig. 2B). The similarity index is equal to the cosine of the angle between the vectors f1 and f2. As such, it can vary between 0 (exactly dissimilar) and 1 (exactly similar). As before, when the diets of consumers are quantified by distributions, the similarity index is itself a distribution, and naturally incorporates the uncertainty derived from mixing models. An alternative similarity metric that may be useful is the Bhattacharyya distance (Bhattacharyya 1943), which measures the distance between 2 discrete or continuous probability distributions. Here we employ the dietary similarity metric presented above to assess the similarity of individual otter diets, as quantified by MixSIR, to those of the whole population. A finding of low similarity would provide evidence for among-individual variation, or individual diet specialization (sensu Bolnick et al. 2003; Estes et al. 2003).
Observational dietary data
Data on foraging behavior and prey consumption by radiotagged sea otters were collected as described in Tinker et al. (2008). We restrict analysis to adult sea otters for which we recorded a minimum of 300 feeding dives over a 2-year period between 2003 and 2006 when data were collected: this resulted in a data set of 30,651 dives for 39 radiotagged study animals (11 at SNI and 28 at MB). We also assembled information on diameter–biomass relationships for each prey type (Oftedal et al. 2007). For each study animal, we estimated diet composition on the basis of consumed wet edible biomass using a Monte Carlo, resampling algorithm designed to account for uncertainty and biases inherent in the raw data (Dean et al. 2002; Tinker et al. 2008); see the Supplemental Material for more details (file can be found online at http://dx.doi.org/10.1644/11-MAMM-S-187.S1 (10.1644_11-MAMM-S-187.S1.pdf)).
Isotopic values of California sea otters and putative prey
The SNI sea otter population (n = 13) had mean (±SD) isotopic values of δ13C = −16.8 ‰ ± 0.6 ‰, δ15N = 14.9 ‰ ± 0.7 ‰. The MB sea otter population (n = 31) had mean (±SD) isotope values of δ13C = −14.5 ‰ ± 0.9 ‰, δ15N = 11.6 ‰ ± 0.8 ‰. Trophic discrimination factors were applied to measured sea otter isotope values in Fig. 3; see the “Materials and Methods” for actual discrimination factors used for each population. Isotope values were determined for the edible tissue of potential sea otter prey from SNI and MB (Newsome et al. 2009, 2010); see Table 1 for scientific names, sample sizes, mean isotope data, and [C]/[N] ratios of prey types. We were unable to obtain a permit to collect abalone at SNI and therefore used abalone data collected from the central California mainland coast from San Simeon to Monterey Bay (see Newsome et al. 2010).
To simplify the mixing space and increase the accuracy of our dietary estimates, we binned spiny lobsters with Cancer crabs, as well as sea urchins with northern kelp crabs in the SNI invertebrate community. The isotopic values of our resultant bins for SNI were: spiny lobsters + Cancer crabs: δ13C = −15.3‰ ± 1.0‰, δ15N = 14.9‰ ± 0.4‰; sea urchins + northern kelp crabs: δ13C = −14.4‰ ± 0.6‰, δ15N = 10.9‰ ± 0.8‰. Similarly, we grouped purple sea urchins with mussels in the MB invertebrate community. The isotopic values of this bin resulted in: purple sea urchins + mussels: δ13C = −17.2‰ ± 1.0‰, δ15N = 9.3‰ ± 0.5‰.
Nearest-neighbor distance (NND) values were determined for the mean δ13C and δ15N values of individuals within SNI and MB sea otter populations. NND values (±SD) were: NNDSNI otter = 0.43 ± 0.21, NNDMB otter = 1.27 ± 0.41. Similarly, the NND values of potential sea otter prey were: NNDSNI prey = 1.85 ± 1.19, NNDMB prey = 2.01 ± 0.84. DC measurements also were calculated for the mean δ13C and δ15N values of sea otters, their potential prey, and the respective centroids of each group, such that DCSNI otter = 0.92 ± 0.21, DCMB otter = 1.58 ± 0.43, DCSNI prey = 2.58 ± 0.77, and DCMB prey = 2.13 ± 0.90. Units for both NND and DC measurements are expressed as ‰. CHA measurements were calculated by 3 separate algorithms to assess the importance of covariance: 1) CHA using only the mean values of consumers and prey, respectively (cf. Layman et al. 2007); 2) CHA incorporating variance of δ13C and δ15N values of consumers and prey, respectively, which assumes independence of δ13C and δ15N values; and 3) CHA incorporating covariance of δ13C and δ15N values of consumer and prey, respectively. Units for each method are expressed as ‰2. CHA results from method 1 were: CHASNI otter = 1.75, CHAMB otter = 6.73, CHASNI prey = 15.80, and CHASNI otter = 10.01. CHA results from method 2 were: CHASNI otter = 1.86 ± 1.19, CHAMB otter = 6.73 ± 1.97, CHASNI prey = 15.81 ± 3.88, and CHAMB prey = 10.06 ± 4.30. Finally, CHA results from method 3 were: CHASNI otter = 1.75 ± 0.88, CHAMB otter = 6.59 ± 1.58, CHASNI prey = 15.93 ± 4.03, and CHAMB prey = 9.87 ± 4.22. If Megastraea snails are not included as potential prey for the SNI sea otter population, CHA estimates for the prey become: method 1 CHASNI prey = 5.24; method 2 CHASNI prey = 5.53 ± 3.4; and method 3 CHASNI prey = 5.67 ± 3.5.
Stable isotope mixing models
We used the Bayesian-based stable isotope mixing model MixSIR (version 1.0.4—Moore and Semmens 2008) to calculate posterior probability densities of the proportional contributions of prey to SNI and MB sea otter population diets. We used uninformative priors and only accepted results if there were ≥1,000 accepted draws. Dietary contributions for SNI and MB populations were calculated by bootstrapping and pooling MixSIR results for individual sea otters in both populations; the densities are displayed in Fig. 4. Median (1st quartile, 3rd quartile) contribution estimates for the SNI sea otter population are: spiny lobsters + Cancer crabs = 0.04 (0.01, 0.10); kelp crabs + sea urchins = 0.26 (0.13, 0.41); abalone = 0.20 (0.05, 0.40); Chlorostoma snails = 0.26 (0.01, 0.66); and Megastraea snails = 0.03 (0.01, 0.09). Median (1st quartile, 3rd quartile) contribution estimates for the MB sea otter population are: purple sea urchins + mussels = 0.02 (0.01, 0.05); kelp crabs = 0.38 (0.09, 0.61); abalone = 0.04 (0.01, 0.13); Cancer crabs = 0.17 (0.06, .028); Chlorostoma snails = 0.05 (0.02, 0.15); and clams = 0.04 (0.02, 0.11). Note that modeling results show a few cases where posterior probability distributions are multimodal or highly variable, such that the median is not an optimal descriptive statistic.
Niche specialization indexes were calculated from MixSIR results for each individual within populations SNI and MB; specialization values for the population were calculated by bootstrapping and pooling individual sea otter MixSIR results such that each individual contributed equally to the final population specialization values (Fig. 5). The SNI sea otter population had mean (±SD) specialization values of 0.53 ± 0.16, whereas MB sea otters had specialization values of 0.56 ± 0.15. Dietary similarity values (s) also were calculated for both SNI and MB sea otters, but because is a comparative measurement it can only be calculated for diet pairs. We calculated the dietary similarity of each otter individual relative to the dietary habits of its population for both SNI and MB, respectively (Fig. 5). The SNI sea otters had individual ∶ population dietary mean (±SD) similarity values of 0.65 ± 0.27 (dimensionless units), whereas MB sea otters had individual ∶ population dietary similarity values that were strongly bimodal. The mode corresponding to lower similarity values was 0.13, whereas the mode corresponding to higher similarity values was 0.96.
We used previously published data sets from 2 populations of California sea otters to identify questions regarding population and individual dietary variation that are informed by both a spatial metric and mixing model approach. We then examine evidence for individual dietary specialization in these populations using variance components analysis and 2 new indexes based on the output of Bayesian-based mixing models. We also discuss challenges and caveats with the use of these tools and offer new ways isotopic data may be utilized to identify differences in population- and individual-level dietary variation.
The isotopic patterns presented as δ13C versus δ15N biplots among the sea otter populations in Fig. 3 clearly suggest a difference in dietary variation at the population level. The amount of variation among individuals at SNI (Fig. 3A) is lower than at MB (Fig. 2B), and this distinction is accurately captured in the significantly larger CHA estimate (t156 = −19.4, P < 0.001) for the MB versus SNI population. The degree of isotopic variation among consumers in a population is fundamentally driven by the amount of variation among prey sources available to each population. Thus, spatial metrics as applied solely to consumer isotopic data with no consideration of variation among putative prey is problematic and can lead to flawed interpretations of dietary variation, individual specialization, and food-web structure (Hoeinghaus and Zeug 2008; Matthews and Mazumder 2004). To account for isotopic variation in prey available to each sea otter population, we calculated the CHA for the prey and present the percentage of this area occupied by each sea otter population (Fig. 3C). Despite removing a single prey species (Megastraea snails) at SNI that significantly reduced CHA estimates for prey at this locality (t190.9 = 20.9, P < 0.001), the sea otters occupied a significantly lower proportion of the prey space than at MB (t99.9 = −5.64, P < 0.001). This suggests that otters at SNI are using a smaller portion of the available niche space, and that interindividual dietary variation is larger at MB in comparison to SNI.
As with any data set, the existence of outliers can lead to an overestimation of the CHA of prey or consumer, or both, isotopic space. Calculation of CHA with and without outliers, as done here with Megastraea snails at SNI, is one way of assessing their impact. The DC metric is a 2nd useful method for assessing the degree of isotopic variation, but is less sensitive to the effects of outliers than CHA because it includes all individuals in a data set, not just those on the periphery that define the convex hull. For example, the mean (±SD) DC of MB sea otters (1.1 ± 0.5) was slightly, but not significantly (t96.6 = 1.87, P = 0.06), larger than that of SNI sea otters (0.8 ± 0.4). A 3rd measurement, the NND, is another informative method to assess the relative distribution of prey and consumer data in isotopic biplots. The mean (±SD) NND among sea otters is similar (P > 0.3) at SNI (0.3 ± 0.2) and MB (0.3 ± 0.1), providing support that the significant difference in CHA between the 2 populations is driven by population-level dietary variation rather than the presence of outliers in the MB data set.
A discussion of outliers and the relative distribution of data in bivariate space highlights differences in the nature of mixing spaces that can affect the utility of spatial metrics. Low variance in the mean NND and DC metrics is characteristic of an even distribution of data in bivariate space. In contrast, when most prey are distributed along the periphery of the mixing space, the existence of outliers can create a large degree of variance in NND and DC metrics. For comparison among populations or communities, a comparison of NND and DC variance (i.e., standard deviation) provides an easy way to determine whether differences in CHA result from outliers in bivariate space. The NND and DC metrics for the SNI and MB sea otters (Fig. 3) have similar standard deviations, suggesting that the distribution of data points in δ13C versus δ15N bivariate space is similar between the 2 populations. This supports the hypothesis that the significant difference in CHA between SNI and MB is indeed driven by a difference in dietary variation at the population level, and not because the CHA of 1 population (MB) is inflated because of outliers. This concept also has implications for the use of mixing models to interpret isotopic data (see below).
Lastly, the difference in sample sizes among individuals in a population or species in a community is a factor to consider when interpreting spatial metrics. For example, we only have data for 13 sea otters from SNI, but data for more than 30 individuals from MB. In this particular scenario, the small overall sea otter population at SNI (n = 30–40—Hatfield 2005) mediates the discrepancy in sample size, because the 13 individuals analyzed here represent a large portion (approximately 30–40%) of the total population in 2003 when vibrissae were collected. Although larger, the approximately 30 MB individuals represent a much smaller fraction (approximately 5%) of the sea otter population in MB. In situations where population sizes are unknown, yet sample sizes are uneven among populations or species, a bootstrap modeling approach can be used. Such an approach would randomly select x number of individuals from the larger of the 2 data sets, where x equals the total number of individuals analyzed in the smaller data set, and calculate a spatial metric (e.g., CHA or NND) several thousand times to provide a conservative estimate of the mean and variance for a subset of individuals in the population for which there are more data.
Stable isotope mixing model
When feasible, the use of stable isotope mixing models to convert isotopic data in resource proportion estimates provides the most useful ecological information that can be directly compared to traditional types of data. Mixing models are ideal for scenarios when trying to parse dietary information among 3 sources (i.e., prey) using 2 isotope systems, or between 2 sources using a single isotope system. As mentioned above, most ecological scenarios are much more complex than these ideal situations, and our example of California sea otters is no exception.
Observational data from SNI and MB show that sea otters consume more than 30 species of invertebrate prey, with substantial overlap in the prey species available to each population. For our isotopic study, we chose to analyze the 7 most important prey species for each population (Fig. 3), which based on observational data combine to represent >95% and approximately 90% of the prey consumed at SNI and MB, respectively (Table 2). In some situations, 2 prey types had similar mean δ13C and δ15N values and we chose to group them into a single source for the mixing models (Phillips et al. 2005). In some cases this strategy produced combinations of prey with similar ecological functions (e.g., Cancer crabs and spiny lobsters at SNI). Other combinations, however, did not include prey types with similar functions, such as purple sea urchins (macroalgae grazer–browser) and mussels (filter feeders) at MB. This overlap is likely related to the large variation in δ13C values of various types of macroalgae (e.g., brown versus red) previously reported from California kelp forest ecosystems (Hamilton et al. 2011; Page et al. 2008). In general, however, the isotopic patterns observed among kelp forest invertebrates in California (Fig. 3) conform to expectations based on the isotopic gradients associated with primary producers (i.e., macroalgae versus microalgae) and food-web structure (i.e., trophic level).
Mean carbon (δ13C) and nitrogen (δ15N) values, sample sizes (n), associated variance, [C]/[N] ratios, and population-level diet composition (±SD) based on observational data for sea otter prey groups from San Nicolas Island and Monterey Bay, California.
Diet composition for California sea otter populations at San Nicolas Island (SNI) and Monterey Bay (MB) based on observation and reported as the percentage of consumed organic biomass. Numbers in parentheses represent standard error of dietary contribution at the population level based on observation of n number of individuals.
We used a Bayesian-based mixing model (MixSIR version 1.0.4—Moore and Semmens 2008) to determine source proportions of the various prey types or groups in the diets of the 2 sea otter populations (Fig. 4). Again, the sea otter scenario presented here is unique because the overall performance of mixing models can be judged because we know a priori the relative contributions of prey types for these populations based on observational data (Table 1; Tinker et al. 2008). At SNI, mixing model results (Fig. 4) identify the most important prey items consumed by sea otters; however, the relative contributions of prey in the population's diet do not conform to observational data. For example, the median contribution for the sea urchin–kelp crab prey group was 26%, where observational data show that these 2 prey types combine to contribute 82% of diet (Table 2). The mixing model also suggests that abalone (20%) and Chlorostoma snails (26%) were more important dietary components than shown by observational data (1.8% and 3.3%, respectively). Other prey types consumed by this population have relatively minor median contributions, which is supported by the observational data (Table 2). This pattern, where only 2 or 3 prey types combine to represent the large majority of resources consumed by this population, also is supported by spatial metrics. Estimates of CHA, NND, and DC for the sea otters at SNI are small relative to those for the prey types available to the population (Fig. 3A). In addition, sea otters occupy a small proportion of the available prey mixing space (Fig. 3C).
Like the observational data, mixing model results for the MB sea otter population shows a more diverse and even contribution of prey that at SNI (Table 2). Similar to SNI, however, the model does a poor job of capturing the relative proportions of prey. For example, the mixing model shows that kelp crabs are the most important prey type (Fig. 4) and overestimates this prey's contribution to population diet with a median contribution of 39% versus 11% via observation. Furthermore, the model shows that the median contribution of sea urchins and mussels to the population's diet is only 2%, but observational data show that this prey is the 2nd most common prey consumed by MB sea otters, contributing approximately 14% to the population diet (Table 2). Lastly, Cancer crabs combine to contribute about 25% to diet based on observation, but the model shows that these prey have a median contribution of 17%.
The Bayesian framework is not a panacea. As with traditional linear mixing models, if the number of sources is too large relative to the number of isotope systems (e.g., δ13C) that are used, the answers provided by the model will be highly uncertain, because many source combinations may contribute to the observed mix. The opposite is also true—if all potential prey items are not included and residual variation is not explicitly estimated, results will be biased. Like any statistical model, if sources of uncertainty—individual or population variation, or variation in source estimates—are ignored or modeled inappropriately, estimates may be biased. Furthermore, if inaccurate trophic discrimination factors or variance in trophic discrimination factors are applied, erroneous posterior probabilities will result (Bond and Diamond 2011). Finally, the nature of mixing spaces can have a dramatic impact on the resultant distributions. For example, if a straight line can be drawn between the mixture (i.e., consumer) and 2 consecutive sources in the isotopic mixing space, results for those 2 potential sources will be highly variable, because many combinations of one, the other, or both sources could contribute to the consumer.
In the discussion of spatial metrics and mixing models above, we focus on aspects of population-level diet. The conventional calculation of a population's dietary breadth integrates prey selection across all individuals but ignores inter- and intraindividual variation in diet. A growing number of studies, however, show that individual dietary specialization (or individuality) is pervasive in many taxa and communities (Bolnick et al. 2003). The growing recognition of individuality in diet has important implications for population biology, food-web dynamics, and stability (Kondoh 2003), and may even contribute to interindividual variation in fitness (Annett and Pierotti 1999; Darimont et al. 2007) that may result in phenotypic diversification and speciation (Moodie et al. 2007; Svanback and Bolnick 2005).
Some recent studies (Lewis et al. 2006; Newsome et al. 2009, 2010) have shown that an isotopic approach is a reliable and cost-effective alternative to observational or gut content data (Estes et al. 2003; Tinker et al. 2008; Werner and Sherry 1987) that has been traditionally used to assess dietary specialization at the sex or individual level. Here we compare a previously published approach, variance component analysis (Newsome et al. 2009), with 2 novel indexes to assess individuality in the SNI and MB sea otter populations. The 2 latter indexes, specialization index (ε) and individual ∶ population similarity (s), are based on the output of the Bayesian-based mixing model (Fig. 4). In contrast, the variance components analysis only uses isotopic data from sea otters (not prey) to evaluate the within-individual components (WIC) and between-individual components (BIC) of diet.
The variance components results show there is nearly twice as much isotopic variance in the MB sea otter data set than from SNI (Fig. 5A), and that this pattern is driven by differences in the BIC (SNI = 0.77 versus MB = 1.40). Traditionally, the total niche width (TNW) of a population has been defined as the sum of the WIC and BIC of diet (Roughgarden 1972), and the ratio of the WIC to the TNW (WIC/TNW) of a population has traditionally been used to evaluate individual specialization (Bolnick et al. 2002). As WIC/TNW approaches 1, all individuals utilize the full spectrum of resources used by the population (i.e., all individuals are generalists), whereas a value close to 0 denotes that individuals are utilizing a small proportion of resources consumed at the population level. In other words, individual specialization increases as the WIC/TNW decreases. In a similar fashion, we can calculate a total isotopic niche width (TINW) for each population by summing the WIC and BIC in the isotopic data set (TINW = WIC + BIC). The mean (±SD) WIC/TINW for the MB population (0.33 ± 0.03) was lower than that of the SNI population (0.40 ± 0.08), suggesting a greater degree of individual specialization at MB versus SNI, which is also supported by observation (Table 2; Tinker et al. 2008).
As an alternative to the WIC/TINW index, we use our specialization metric (ε) to assess the degree of niche specialization for SNI and MB sea otter individuals. The benefit of this approach is that it does not require estimates of within-individual dietary variation, even though we do include within-individual variability in our analysis of sea otter diet. When isotope data are used to calculate ε, however, both predator and potential prey isotopic values are required to calculate the proportional contribution of each prey to a predator's diet. Because isotopic approaches are often of most use when quantifying species interactions that are difficult or impossible to observe directly, we suggest that this approach has merit. We observe that at both the population and individual level, the MB sea otters have higher specialization values than the SNI sea otters, which is qualitatively similar to the WIC/TINW metric (Fig. 5A).
An integration of the specialization metric with the measurement for individual ∶ population dietary similarity (s) may elucidate more information regarding sea otter diet (Fig. 5). Here we observe that most individuals at SNI have high similarity to the population mean (high s values), and range from being niche specialists to generalists. A few individuals exhibit slightly different diets, represented by much weaker peaks at lower values of s. By contrast, MB sea otters have similarity values that are strongly bimodal, a pattern associated with a high degree of individuality, with individuals falling into 1 of a number of potential specialist types.
The analysis of SNI and MB sea otters with our proposed DietSpace metrics presents a straightforward means to assess niche specialization and dietary similarity with mixing model output derived from isotopic measurements. Our results confirm the expectation that SNI sea otters are more similar to the population as a whole, whereas MB sea otters are more variable. Our results also confirm that the degree of niche specialization at the individual level is generally high for the MB sea otters, whereas the greater range in similarity values indicates that there is greater population dietary diversity at MB. A result that is less obvious when analyzing only isotopic data is that SNI sea otters span a range of niche specialization values and have low individuality (high s values), suggesting that individuals are less likely to include prey that is less preferred at the population level.
Future developments and words of caution
Stable isotope analysis is often attractive to ecologists because the results are numerical data and naturally lend themselves to descriptive statistics and quantitative manipulation. We foresee a strong future for quantitative tools that aid ecologists in interpreting isotopic data from the individual to community level. Future work, however, should look back at potential weaknesses in an isotopic approach; in particular ecologists should examine assumptions that underpin isotopic mixing models. Martínez del Rio et al. (2009b) demonstrated that inter- and intrataxon variation in isotope trophic discrimination factors can vary depending upon a number of factors. These discrimination factors are central in placing the consumer in the mixing space and affect some of the IsoSpace metrics and, to an unknown degree, the isotopic mixing models. Bond and Diamond (2011) found that mixing model results are highly contingent upon discrimination factors, but this likely depends upon the nature of the IsoSpace.
The similarity among functionally distinct prey groups at MB highlights an important caveat in using an isotopic approach to study variation in resource (or habitat) use. Isotopic data do not typically provide dietary composition data at the species level. Instead, isotopic variation within or among ecosystems is created by biochemical processes, such as the formation of a 3-carbon or 4-carbon sugar in the 1st step of photosynthesis or the transamination–deamination of amino acids in the tricarboxylic acid cycle that occurs during nutrient assimilation and amino acid synthesis. Biological processes that sort isotopes also can be driven by physicochemical variation in the environment (e.g., temperature), which governs the rate of biochemical reactions by which organisms assimilate nutrients, grow, and reproduce. Stable isotope analysis will only be useful when potential sources (prey) have distinct isotope values (but see Yeakel et al. 2011) and an ideal approach would be to combine stable isotope analysis with traditional diet proxies (e.g., scat analysis) that provide higher-resolution information on dietary diversity.
Finally, we look to the future of DietSpace metrics and other ways to assess individuality with isotope data. We anticipate that further quantitative descriptors relating dietary similarity or differences either among individuals or among populations will be useful in a number of ecological subfields. One potential example is a metric of DietSpace that uses hierarchical clustering to create a tree of dietary similarity. Such an approach could be used at the population level to examine animals that compete for resources. At the individual level, a dietary tree of individuals could be compared to genetic data to determine if foraging behavior is culturally transmitted from parents to offspring. As isotopic research continues and techniques improve in accuracy and rigor we expect that researcher creativity and ingenuity in exploring their data will expand as well. We welcome the further development of tools that seek to translate isotopic measurements into a language that can be understood by all ecologists.
Thanks to the Monterey Bay National Marine Sanctuary Foundation and the United States Marine Mammal Commission for funding the collection and processing of sea otter prey samples. We thank A. Green, M. Kenner, K. Miles, and J. Bodkin for assistance in prey collection and E. Snyder, C. Mancuso, W. Wurzel, E. Heil, and R. Harley for laboratory assistance. We thank A. C. Jakle, E. J. Ward, and L. Feliciti for constructive reviews. We thank A. Green for obtaining the necessary California Fish and Game permit for marine invertebrate collection and the United States Navy for permission to conduct research on San Nicolas Island. SDN was partially funded by the National Science Foundation (ATM-0502491, DIOS-0848028), Carnegie Institution of Washington, and the W. M. Keck Foundation (072000). S. D. Newsome and J. D. Yeakel contributed equally to this study.