Environ. Sci. Technol. 2005, 39, 8395-8402
A Statistical Resampling Method To Calculate Biomagnification Factors Exemplified with Organochlorine Data from Herring (Clupea harengus) Muscle and Guillemot (Uria aalge) Egg from the Baltic Sea K A T R I N L U N D S T E D T - E N K E L , * ,†,‡ M A T S TYSKLIND,§ JOHAN TRYGG,| PETER SCHU ¨ LLER,+ LILLEMOR ASPLUND,@ ULLA E R I K S S O N , @ L I S B E T H H A¨ G G B E R G , @ TJELVAR ODSJO ¨ ,‡ MATS HJELMBERG,‡ M A T S O L S S O N , ‡,@ A N D J A N O ¨ RBERG† Environmental Toxicology, Department of Physiology and Developmental Biology, Evolutionary Biology Centre, Uppsala University, Norbyva¨gen 18A, SE-752 36, Sweden, and Contaminant Research Group, Swedish Museum of Natural History, P.O. Box 50007, SE-104 05 Stockholm, Sweden, Environmental Chemistry and Research Group for Chemometrics, Department of Chemistry, Umeå University, SE-901 87 Umeå, Sweden, InfiDyne AB, Murargatan 30A, SE-754 37 Uppsala, Sweden, and Institute of Applied Environmental Research, Stockholm University, SE-106 91 Stockholm, Sweden
A novel method for calculating biomagnification factors is presented and demonstrated using contaminant concentration data from the Swedish national monitoring program regarding organochlorine contaminants (OCs) in herring (Clupea harengus) muscle and guillemot (Uria aalge) egg, sampled from 1996 to 1999 from the Baltic Sea. With this randomly sampled ratios (RSR) method, biomagnification factors (BMFRSR) were generated and denoted with standard deviation (SD) as a measure of the variation. The BMFRSR were calculated by randomly selecting one guillemot egg out of a total of 29 and one herring out of a total of 74, and the ratio was determined between the concentration of a given OC in that egg and the concentration of the same OC in that herring. With the resampling technique, this was performed 50 000 times for any given OC, and from this new distribution of ratios, BMFRSR for each OC were calculated and given as geometric mean (GM) with GM standard deviation (GMSD) range, arithmetic mean (AM) with AMSD range, and minimum (BMFMIN) as well as maximum (BMFMAX) biomagnification factors. The 14 analyzed OCs were p,p′DDT and its metabolites p,p′DDE and p,p′DDD, polychlorinated biphenyls (PCB congeners: * Corresponding author phone: +46 18 4716498; fax: +46 18 518843; e-mail:
[email protected]. † Uppsala University. ‡ Swedish Museum of Natural History. § Environmental Chemistry, Department of Chemistry, Umeå University. | Research Group for Chemometrics, Department of Chemistry, Umeå University. + InfiDyne AB. @ Stockholm University. 10.1021/es048415y CCC: $30.25 Published on Web 09/21/2005
2005 American Chemical Society
CB28, CB52, CB101, CB118, CB138, CB153, and CB180), hexachlorocyclohexane isomers (R-, β-, and γHCH), and hexachlorobenzene (HCB). Multivariate data analysis (MVDA) methods, including principal components analysis (PCA), partial least squares regression (PLS), and PLS discriminant analyses (PLS-DA), were first used to extract information from the complex biological and chemical data generated from each individual animal. MVDA were used to model similarities/dissimilarities regarding species (PCA, PLS-DA), sample years (PLS), and sample location (PLS-DA) to give a deeper understanding of the data that the BMF modeling was based upon. Contaminants that biomagnify, that had BMFRSR significantly higher than one, were p,p′DDE, CB118, HCB, CB138, CB180, CB153, βHCH, and CB28. The contaminants that did not biomagnify were p,p′DDT, p,p′DDD, RHCH, CB101, and CB52. Eventual biomagnification for γHCH could not be determined. The BMFRSR for OCs present in herring muscle and guillemot egg showed a broad span with large variations for each contaminant. To be able to make reliable calculations of BMFs for different contaminants, we emphasize the importance of using data based upon large numbers of, as well as well-defined, individuals.
Introduction The aim of this paper is to introduce a novel method to calculate biomagnification factors (BMF), the randomly sampled ratios (RSR) method that denotes BMFs with an estimate of variation. With the method currently in use, BMFs are calculated as the ratio between the arithmetic mean (AM) or geometric mean (GM) concentration of a contaminant in the predator and the AM or GM concentration of the same contaminant in the prey. Therefore, these BMFs are always denoted as a single value with no estimate of the variation. There are many papers that deal with the issue of error/ uncertainty propagation regarding bioaccumulation of contaminants, for instance, to evaluate how model output changes with variations in input variables (1, 2) and/or how uncertainties in individual parameters can affect uncertainties in the results (3, 4). These papers use Bayesian methodologies (e.g., Monte Carlo or Marcov Chain Monte Carlo simulations) to generate new distributions of values for, for example, contaminant concentrations used as input/output values in the modeling. Often, such papers generate new distributions after a low number of analyzed samples from short-term experiments (i.e., the generated distributions are not based on the distribution of contaminants in the natural environment). The U.S. EPA recommends the use of techniques to clarify the impact/sources of uncertainty on the result of an ecological risk assessment but at the same time cautions that the poor execution of any method can obscure rather than facilitate understanding (5). In fact, Linkov and Burmistrov (6) found that modeler uncertainty and choices made by modelers contributed to as high as 7 orders of magnitude differences in model predictions. In a food-web of a specific ecosystem, nutrients are transported between abiotic and biotic compartments as well as between individual organisms within the biotic compartment. Contaminants are transported in a similar way forming a contaminant-web, where fluctuations in concentrations are caused by variations in biotic factors (e.g., species composition, age, and health status of individuals) as well as in abiotic factors (e.g., temperature, wind, and precipitation). VOL. 39, NO. 21, 2005 / ENVIRONMENTAL SCIENCE & TECHNOLOGY
9
8395
All together, this leads to variations in contaminant concentrations in any given individual at any given time (7-9). When modeling contaminant levels, transport, and biomagnification in an ecosystem, a thorough knowledge of variations in concentrations is important. Without this, modeling is difficult to perform at best, or glaringly faulty at worst (6, 8, 10-17). Biomagnification. As our present knowledge regarding quantitative data on uptake routes (e.g., gastrointestinal tract, respiratory organs, and integument) of contaminants in most animals in the wild is limited, this paper uses a wide definition of biomagnification (18). If the lipid normalized concentration of a substance is significantly higher in the predator than in the prey, regardless of how that higher concentration has been achieved, biomagnification of that substance has occurred. Properties of the contaminant essential for biomagnification, mainly bioavailability and persistence, are determined by a number of physical/chemical characteristics such as molecular size, structure (19-21), presence of reactive groups, lipophilicity, and water solubility. Other variables important for biomagnification are related to the animals: the predator and its prey (e.g., feeding habits and habitats, lipid characteristics and fat content, sex, age, and the animals’ capacity to biotransform/excrete contaminants). The fact that contaminants biomagnify to various degrees (22) ultimately leads to specific contaminant patterns in different animal species within the same ecosystem. The ultimate aim of our research is to study and model contaminants’ physical-chemical properties, including their structures, in relation to the varying contaminant patterns in different animal classes to gain an insight to biomagnification processes (i.e., factors favoring or counteracting biomagnification). To reach this goal, we first have to determine reliable BMFs with an estimate of the variation. In our present article, we have used large numbers of individually analyzed herring (Clupea harengus) muscle and guillemot (Uria aalge) egg from the Swedish monitoring program (23) carried out in the Baltic Sea to calculate biomagnification factors and to quantify the uncertainty in these calculations. The eggs are sampled from a bird colony whose feeding area covers the two sampling sites where the herring are caught; thus, the analyzed animals originate from the same contaminant-web; see Materials and Methods. Difficulties might arise with migratory birds where the birds’ exposure history is unknown for parts of the year (24). For guillemot from the Stora Karlso¨ colony, which is a stationary bird residing in the Baltic Sea all year-round, this problem does not arise. Only a few articles present calculations of BMFs using data from the biota of the Baltic Sea. On the basis of the levels of polychlorinated dibenzo-dioxins/furans (PCDD/F) in herring muscle and guillemot eggs, de Wit et al. reported a BMF of ∼16 for PCDD/F (25).
Materials and Methods Herring. Regarding transport of environmental contaminants in the Baltic Sea food-web, herring (C. harengus) is an important species; thus, a number of herring specimens are sampled each year from the Baltic Sea (23, 26-28). It is a planktivorous fish that feeds on selected pelagic zooplankton (29-31), and in turn, herring play an important role as food for piscivore animals such as other fish (e.g., cod (Gadus morhua) and salmon (Salmo salar)), birds (e.g., guillemot (U. aalge)), and mammals (e.g., grey seal (Halichoerus grypus)) and is an important human food source (32). Herring specimens used for this article were sampled from fish catches at two sites, Landsort (L) and Utla¨ngan (U) (23). L and U are situated to the north (L) and to the south (U) of the guillemot colony at Stora Karlso¨ (S), and both herring locations are within those birds’ large foraging area (33). The herring were 8396
9
ENVIRONMENTAL SCIENCE & TECHNOLOGY / VOL. 39, NO. 21, 2005
sampled autumn (October to November) of 1996, 1997, and 1998, and a total of 72 female fish (12 individuals from each site, each year) 3-4 years of age were later individually analyzed regarding their dorsal muscle concentration of a number of OCs. All collecting and sample preparation was carried out and recorded in a standardized manner (23, 34). (The sampling procedures as well as the resulting biological data regarding animal characteristics are presented in detail in the Supporting Information). Guillemot. Guillemot (U. aalge) is a bird species that resides year-round in the Baltic Sea; it is a piscivore that feeds mainly on the two clupeid fish species, herring and sprat (Sprattus sprattus) (23, 35-37). All guillemot eggs were collected from the guillemot colony at the island of Stora Karlso¨ (S) (23). At around 5 years of age, the guillemot female breeds for the first time, laying one large (∼110 g) egg per year. The high age at first egg laying, together with the fact that guillemot invest large quantities of lipids relative to their body weight of ∼1000 g to the one egg they lay each year, means that the contaminant composition of the egg mainly reflects the composition of maternal tissue rather than the female’s diet at time of yolk formation (38). The guillemot eggs were from spring (April to May) 1997, 1998, and 1999. A total of 30 guillemot eggs (10 eggs/year) were later individually analyzed regarding their (yolk and albumen) concentration of a number of OCs. One guillemot egg from 1999 was lost during sample preparation, so the total number of guillemot eggs analyzed for the 3 years were 29 (see Supporting Information). Analyzed OCs. Chemical analyses were carried out at an accredited laboratory (Institute for Applied Environmental Research (ITM) at Stockholm University), according to previously described procedures (39, 40). For both herring and guillemot samples, the following 14 OCs were analyzed: 2,2-bis(4-chlorophenyl)-1,1,1-trichloroethane compounds (p,p′DDT and its metabolites p,p′DDE and p,p′DDD); polychlorinated biphenyls (PCBs) with the congeners’ respective IUPAC number within parentheses: 2,4,4′-trichlorobipenyl (CB28), 2,2′,5,5′-tetrachlorobiphenyl (CB52), 2,2′,4,5,5′-pentachlorobiphenyl (CB101), 2,3′,4,4′,5-pentachlorobiphenyl (CB118), 2,2′,3,4,4′,5′-hexachlorobiphenyl (CB138), 2,2′,4,4′,5,5′hexachlorobiphenyl (CB153), and 2,2′,3,4,4′,5,5′-heptachlorobiphenyl (CB180); hexachlorocyclohexane isomers (R-, β-, and γHCH); and hexachlorobenzene (HCB). All compounds could be separated from one another, as well as from other PCBs of importance, except for CB138, where the interfering peak caused by CB163 means that CB138 can be overestimated by up to 20-30%. To facilitate reading, CB138 is henceforth written alone. (The chemical analyses as well as the obtained residue levels are presented in detail in the Supporting Information.) Missing Values. When calculating the various biomagnification factors, the level of quantification (LOQ) was used for a few samples with nonquantifiable levels of an OC. When possible for those samples, LOQ values determined from each individual chromatogram were used, denoted as LOQIND. For guillemot egg, this was not always possible, and then a general LOQ value was used. The LOQ was defined as 10 times the standard deviation of the measurements, and in practice, this was estimated by repetitive measurements of a sample containing analytes at a concentration close to the expected LOQ (41). These general LOQs were for p,p′DDT 50 ng/g lw, p,p′DDD 60 ng/g lw, R-HCH 30 ng/g lw, γ-HCB 40 ng/g lw, CB52 30 ng/g lw, and finally, CB101 35 ng/g lw. Statistics. For basic statistics regarding the biological variables and concentrations of OCs, the software GraphPad Prism 4.03 (42) was used. Multivariate data analyses (MVDA) (i.e., principal component analysis (PCA), partial least-squares regression (PLS), and PLS discriminant analyses (PLS-DA)) were performed using the software SIMCA-P 10.0.4 (43).
PCA was used to illustrate and to discern groupings in the data, a significance level of 0.05 was used, and centered and scaled data (to variance 1) were fitted to a principal components (PC) model (44). Values of the explained variation, R2, and predicted variation, Q2, were calculated, and R2 values >0.7 and Q2 values >0.4 denote an acceptable model when analyzing biological data (45). PLS was used to determine whether there were significant differences in pollution patterns between the 3 years, for herring and guillemot separately. PLS is an extension of multiple linear regressions combined with PCA to model the relationship between two matrixes, Y and X, that can both be multidimensional. This is performed by modeling X and Y separately and then relating their respective scores to each other (46, 47). PLS-DA was used to model contaminant patterns to discern differences between the two herring locations as well as between the two species. PLS-DA is used to determine whether an observation belongs to a preestablished class by modeling and quantifying the eventual discriminating variables that contribute to the class (e.g., X/Y matrixes) separation (48). Calculation of Biomagnification Factors (BMFs). To calculate BMFs, different methods were used. By the first method, the method currently in use, BMFs were calculated as the ratio between the geometric mean (GM) as well as the arithmetic mean (AM) concentration of a given OC in guillemot eggs for the 3 years 1997, 1998, and 1999 (n ) 29) and the GM as well as the AM concentration of the same OC in herring from both locations L and U for the 3 years 1996, 1997, and 1998 (n ) 72). This method gives single BMF values for each OC, which we denoted as BMFGM and BMFAM. With the new randomly sampled ratios (RSR) method, BMFs were calculated by randomly selecting one guillemot egg out of 29 and one herring out of 72, using the BMFRSR computer program especially designed for the task (49). This program uses a random sampling with a replacement technique where a single specimen is randomly sampled and can be sampled several times. The ratio between the concentrations of a given OC in that sampled guillemot egg and the same OC in that sampled herring muscle was determined. In this way, 50 000 iterations were performed, generating 50 000 observations of BMFRSR for each contaminant, and geometric mean (GM, BMFRSRGM) with geometric mean standard deviation (GMSD) as well as arithmetic mean (AM, BMFRSRAM) with arithmetic mean SD (AMSD) were calculated for this new distribution of ratios. BMFRSRGM as well as BMFRSRAM ranges ((1 SD) were then calculated. For BMFRSRGM, the lower limit is calculated as BMFRSRGM/GMSD and the upper limit as BMFRSRGM * GMSD. For BMFRSRAM, the lower limit is calculated as BMFRSRAM - AMSD and the upper limit as BMFRSRAM +AMSD. Standard deviation is a measure of the degree of dispersion of the data from the mean value. For approximately normally distributed data, as in this case after log transformation of the quotients, then about 68.27% of the values are within (1 SD around the mean, and about 95.45% of the values are within (2 SD around the mean. Also, minimum and maximum BMFs for all organochlorines were determined. These BMFs were calculated as the ratio between the guillemot egg with the lowest concentration of a given contaminant and the herring muscle with the highest concentration of the same contaminant (BMFMIN) and between the guillemot egg with the highest concentration and the herring muscle with the lowest concentration (BMFMAX ) (Figure 1). Number of Iterations. To discern variations in BMFRSR due to the numbers of iterations (number of randomly selected pairs), a series of different numbers of iterations (3, 5, 10, 50, 100, 500, 1000, 5000, 10 000, 50 000, and 100 000) was performed 10 times (e.g., 10 runs with 3 iterations, 10
FIGURE 1. Concentration on a logarithmic scale (ng/g of lipid weight) of CB153 in herring (C. harengus) muscle (n ) 72) from Landsort and Utla1 ngan and in guillemot (U. aalge) egg (n ) 29) from Stora Karlso1 , the Baltic Sea, 1996-1999. Biomagnification factors (BMFs) calculated as BMFMIN (3.0), BMFGM (24.7), and BMFMAX (158) (see Materials and Methods and Table 1). runs with 5 iterations...10 runs with 100 000 iterations). This was done using the CB153 concentration in herring muscle and guillemot egg. The coefficient of variation (CV%) between the 10 resulting BMFRSRGM from the 10 runs with the same number of iterations was calculated. This was also done in the same way for BMFRSRAM.
Results Biological Variables. The herring (n ) 72) had (GM with 95% CI) a fat content (F%) of 2.0% (1.8-2.3%), body weight (BW) of 33.7 g (32.0-35.5 g), and age of 3.7 years (3.6-3.8 years). The guillemot egg (n ) 29) had F% of 11.8% (11.312.3%) and total weight/egg of 110 g (107-113 g). (The measured biological variables are presented in detail in the Supporting Information.) Concentrations. The OC with the highest level in herring muscle was p,p′DDE with a concentration (GM with 95% CI) of 370 ng/g lw (320-430 ng/g lw), and this OC also showed the highest concentration in guillemot egg of 21 200 ng/g lw (19 400-23 200 ng/g lw). Among the OCs analyzed, the lowest concentration in herring muscle was found for CB28 that was 7.0 ng/g lw (6.2-8.0 ng/g lw), and here, 23 out of the total of 72 individuals had levels below the LOQ. In guillemot egg, p,p′DDT concentrations were below the LOQ (50 ng/g lw) in all samples. (The concentrations of organochlorines in herring muscle and guillemot eggs are presented in detail in the Supporting Information.) Multivariate Data Analysis. The PCA using the whole data set, with guillemot egg and herring muscle chemical variables (analyzed OCs) as X, gave a model (R2X ) 0.83, Q2 ) 0.74, two components), indicating that the species make up two separate groups (Figure 2a), meaning that the contaminant pattern differed between the two species (Figure 2b). PLS-DA with herring versus guillemot (R2X ) 0.78, R2Y ) 0.95, Q2 ) 0.94, two components) showed complete separation between the two species (score plot not shown, similar to Figure 2a), and the coefficient plot (Figure 3) showed that p,p′DDE, βHCH, HCB, CB28, CB118, CB138, CB153, and CB180 all had higher concentrations in guillemot egg, while p,p′DDD and CB101 had higher concentrations in herring muscle. Differences between Years? For herring muscle, the PLS model with the respective years as Y and biological as well as chemical variables as X revealed no significant differences in the variables between the 3 years (R2X ) 0.55, R2Y ) 0.23, Q2 ) 0.02, two nonsignificant components). For guillemot egg, the corresponding PLS model (R2X ) 0.36, R2Y ) 0.87, VOL. 39, NO. 21, 2005 / ENVIRONMENTAL SCIENCE & TECHNOLOGY
9
8397
FIGURE 2. Principal component analysis (PCA) based on concentrations of contaminants in guillemot (U. aalge) egg (n ) 29) from Stora Karlso1 (S) and herring (C. harengus) muscle (n ) 72) from Landsort and Utla1 ngan (L and U) from the Baltic Sea, 1996-1999. PCA model (R2X ) 0.83 and Q2 ) 0.74, two components). (a) Score plot and (b) loading plot. For abbreviations of contaminants, see Materials and Methods.
FIGURE 3. Coefficient plot with 95% CI for the respective variables from PLS-discriminant analysis (PLS-DA) based on concentrations of contaminants in guillemot (U. aalge) egg (n ) 29) from Stora Karlso1 (S) and herring (C. harengus) muscle (n ) 72) from Landsort and Utla1 ngan from the Baltic Sea, 1996-1999. PLS-DA model for guillemot vs herring (R2X ) 0.78, R2Y ) 0.95, Q2 ) 0.94, two components). For abbreviations, see Materials and Methods. Q2 ) 0.64, two components) showed that there were systematic and significant differences between the three years. The eggs from 1997 and 1998 had a significantly lower concentration of HCB than the eggs from 1999 and also a significantly higher concentration of other OCs. (PLS plots not shown, for concentrations, see Supporting Information.) Differences between Herring from Landsort and Utla1 ngan? The PLS-DA model (R2X ) 0.66, R2Y ) 0.65, Q2 ) 0.55, three components) based on herring data (biological as well as chemical variables) from the two locations L and U showed that the populations from the two locations overlap 8398
9
ENVIRONMENTAL SCIENCE & TECHNOLOGY / VOL. 39, NO. 21, 2005
FIGURE 4. PLS-discriminant analysis (PLS-DA) based on biological measurements and concentrations of contaminants in herring (C. harengus) muscle (n ) 72) from Landsort and Utla1 ngan (L and U) from the Baltic Sea, 1996-1998. PLS-DA model for L vs U (R2X ) 0.66, R2Y ) 0.65, Q2 ) 0.55, three components). (a) Score plot and (b) coefficient plot for herring from L with 95% CI for the respective variables, weight (BW), length (BL), condition factor (ConF), and lipid and dry matter content (F% and DM%, respectively). For abbreviations of contaminants, see Materials and Methods. to some small extent but that differences were discernible (Figure 4a). Herring from L were significantly older and had a lower body weight (BW) than fish from U (Figure 4b). The OCs with significantly higher concentrations in fish from L than in fish from U were R-HCH, β-HCH, CB101, and CB118 (Figure 4b). The p,p′DDT concentrations were significantly higher in fish from U than in fish from L (Figure 4b). Biomagnification Factors. The BMFs calculated using different methods are presented in Table 1. The OCs that biomagnify, with BMFRSRGM significantly higher than 1, were in the following decreasing order: p,p′DDE, CB118, HCB, CB138, CB180, CB153, βHCH, and CB28. The OCs that showed values below 1 for the calculated BMFRSRGM and did not biomagnify were p,p′DDD, p,p′DDT, RHCH, CB52, and CB101. Number of Iterations. When the series with increasing number of iterations (3, 5...100 000) was performed 10 times using the CB153 concentration in herring muscle and guillemot egg, the variation between the 10 resulting BMFRSRAM as well BMFRSRGM values from the 10 runs with same number of iterations showed similar patterns. The coefficient of variation (CV%) rapidly declined with an increase in number of iterations, from a CV% of 57% between the 10 runs with 3 iterations to 0.5% or below between the 10 runs with 10 000 iterations or more (Figure 5a). When using a higher number of iterations, the resulting BMFRSRGM converge to one consistent value (equal to BMFGM). Accordingly, when using a higher number of iterations, the resulting GMSD and AMSD converge to one consistent value (shown for GMSD in Figure 5b).
Discussion Modeling Biomagnification Processes. One important issue to consider when modeling biomagnification processes is
ΣDDT
24.6 20.9 24.6 1.97 12.5-48.4 30.5 20.0 10.4-50.5 3.77 156 yes
ΣPCB
p,pDDD 0.42 0.35 0.42 1.88 0.22-0.79 0.52 0.42 0.10-0.94 0.10 3.22 no
p,pDDE 57.1 47.4 57.1 1.98 28.8-113.1 71.3 49.2 22.2-120.5 9.30 403 yes
0.93 0.76 0.93 1.96 0.48-1.82 1.21 1.15 0.06-2.35 0.22 8.09 no
p,pDDT βHCH 18.6 18.7 18.6 1.47 12.7-27.3 20.1 8.78 11.3-28.9 7.17 84.4 yes
rHCH 1.37 1.35 1.38 1.27 1.08-1.75 1.42 0.34 1.08-1.75 0.52 2.63 no
2.10 2.08 2.10 1.16 1.82-2.43 2.12 0.32 1.81-2.44 1.11 3.24 m
γHCH 33.8 30.0 33.8 1.75 19.3-59.2 40.1 26.2 13.8-66.3 11.0 166 yes
HCB 12.5 10.3 12.5 1.74 7.18-21.8 14.2 6.65 7.54-20.8 0.79 47.3 yes
CB28 1.96 1.46 1.96 2.10 0.94-4.12 2.55 2.00 0.55-4.56 0.17 12.5 no
CB52 0.46 0.39 0.46 1.82 0.25-0.84 0.55 0.33 0.21-0.88 0.11 2.44 no
CB101 42.5 36.1 42.5 1.90 22.3-80.8 51.6 32.1 19.5-83.6 8.94 205 yes
CB118 29.9 25.3 30.0 1.98 15.1-59.4 37.1 23.9 13.2-61.0 3.65 166 yes
CB138
24.7 21.0 24.7 2.04 12.1-50.4 31.2 21.4 9.83-52.6 3.04 158 yes
CB153
27.8 21.8 27.8 2.25 12.3-62.6 37.4 29.5 7.89-67.0 2.24 287 yes
CB180
a The animal material was collected between 1996 and 1999. b BMF c GM; calculated with geometric mean (GM) concentrations for respective OC in Ua/Ch. BMFAM; calculated with arithmetic mean (AM) concentration for respective OC in Ua/Ch. d BMFRSRGM; calculated with the randomly sampled ratios (RSR) method using 50 000 iterations, resulting in GM ( GMSD. e GMSD; geometric mean standard deviation. f BMFRSRGM range; lower limit (BMFRSRGM /GMSD) and upper limit (BMFRSRGM * GMSD). g BMFRSRAM; calculated with RSR method using 50 000 iterations, resulting in AM ( AMSD. h AMSD; arithmetic mean standard deviation. i BMFRSRAM range; lower limit (BMFRSRAM - AMSD) and upper limit (BMFRSRAM +AMSD). j BMFMIN; minimum BMF, see Materials and Methods. k BMFMAX; maximum BMF, see Materials and Methods. l BM; yes ) biomagnification (BM) of contaminant with BMF significantly (p < 0.05) higher than 1 and no ) no BM of contaminant, BMF not significantly higher than 1. m Eventual BM cannot be determined, too high a level of quantification (LOQ) for γ-HCH in guillemot egg.
BMFRSRGMd GMSDe
36.0 31.3 36.0 1.86 BMFRSRGM rangef 19.4-66. 9 BMFRSRAMg 43.7 AMSDh 31.0 BMFRSRAM rangei 12.7-74. 7 j BMFMIN 6.94 BMFMAXk 328 BMl yes
BMFGMb BMFAMc
TABLE 1. Biomagnification Factors (BMF) (BMFGM, BMFAM, BMFRSRGM, BMFRSRAM, BMFMIN, and BMFMAX) as the Ratio of Concentration of Organochlorines between Herring (Clupea harengus) (Ch) Muscle, n ) 72, from Landsort and Utla1 ngan, and Guillemot (Uria aalge) (Ua) Egg, n ) 29, from Stora Karlso1 in the Baltic Sea, See Materials and Methodsa
FIGURE 5. Effect of increasing the number of iterations on the calculated biomagnification factor (BMFRSRGM) as well as on the size of the standard deviation (GMSD) based on CB153 concentration in herring (C. harengus) muscle and guillemot (U. aalge) egg from the Baltic Sea, 1996-1999. (a) Each dot represents the resulting value; BMFRSRGM after one run with 3, 5, 10...10 000 iterations/run. Horizontal line represent the BMFRSRGM for CB153 after 50 000 iterations. Increasing number of iterations leads to decreasing coefficients of variation, from ca. 60% (with 3 iterations) to 0.5% (with 10 000 iterations), see Materials and Methods. (b) Each triangle represents the resulting standard deviation for BMFRSRGM after one run with 3, 5, 10...10 000 iterations/run. Horizontal line represents the geometric mean standard deviation (GMSD) for CB153 after 50 000 iterations.
the wild animals’ diet as the exact composition is unknown. In the present study, the BMFs are based on data from herring, although the diet of guillemot is a mix, mainly of herring and sprat (23, 35-37). Herring and sprat are similar in morphology, feeding behavior, and ecology (e.g., the two species share the same fish shoals, have similar feed preferences as well as similar diurnal rhythms in vertical migration when feeding (30, 50-52)). A number of studies has shown that the contaminant levels in herring and sprat are remarkably similar (53-55). However, due to differences in growth rate, a sprat of a given size might be older and have higher concentrations of contaminants than a herring of the same size (55), and therefore, it is possible to overestimate BMFs for the step herring to guillemot if the guillemot diet was to a large extent comprised of old sprat. Ideally, for future
VOL. 39, NO. 21, 2005 / ENVIRONMENTAL SCIENCE & TECHNOLOGY
9
8399
research, we would therefore seek to calculate BMFs using contaminant data from the guillemots’ actual feed collected from the guillemots’ crop or stomach. The data from herring used in this article originate from muscle. We assume that the contaminant composition of herring muscle reflects the composition of the whole body burden. To our knowledge, a comparison between concentration in herring muscle and concentration in whole fish is not available in the literature. However, unpublished data from the Swedish National Food Administration (NFA) indicate that the concentrations of contaminants in herring and sprat muscle as compared to the concentration in the whole body are similar when expressed on a lipid weight basis. Further, the lipid-normalized concentration of polybrominated diphenyl ethers (PBDEs) was equal between herring muscle and liver (56). The data from guillemot used in this article originate from eggs. We assume that the contaminant composition of the guillemot egg mainly reflects the composition of maternal tissue rather than the female’s diet at the time of yolk formation (38). To our knowledge, any comparisons between the concentrations of contaminants in guillemot egg and that in guillemot muscle or the whole body are not available in the literature. The PLS-DA shows that there are differences between herring from the two locations regarding both biological as well as chemical variables (Figure 4). However, the guillemots at the Stora Karlso¨ (S) foraging area is large (33), so herring from L in the northern part of the foraging area, together with herring from U in the southern part, are likely to represent the fish feed, and thus the contaminants, to which S guillemot are exposed. By using animals from several locations within the predators foraging area, a better representation of the prey contaminant content is gained. As the concentration of any given contaminant in both prey and predator varies due to age, feed, location, as well as season of sampling, we deem that it is not possible to denote a true BMF for any contaminant based on data from wild living organisms. To make reliable comparisons with other calculations of BMFs possible, rigorously defined background information concerning predator and prey, as well as denoted variances for the BMFs, is needed. Use of LOQ. For some OCs predominantly from guillemot egg, we have used the LOQ or LOQIND for the BMF calculations (Table ST3 in Supporting Information). This procedure gives an overestimation of the calculated BMFs if the actual concentrations in guillemot eggs are lower than LOQ. This concerns p,p′DDD, p,p′DDT, RHCH, γHCH, CB52, and CB101. With one exception, γHCH, all of these OCs have BMFRSR values not significantly above 1. We can therefore conclude that because none of these OCs biomagnifies between herring muscle and guillemot egg, the overestimation of the concentration in guillemot egg has no significance. However, for γHCH, the high LOQ makes it impossible to determine whether γHCH is biomagnified or not, and a more sensitive analytical method is needed. Uncertainty and Variability. Recently, the U.S. EPA recommended probabilistic modeling approaches to characterize uncertainty and variability in ecological risk assessment (5). A method generating an estimate of the variation for the BMFs such as the RSR method as presented in this paper is in line with this recommendation. When dealing with error propagation of bioaccumulation factors, many existing models are based on one or a few analyzed contaminant/s from one single individual or pooled sample, or from a few individuals, and generation of new, larger distributions for the contaminant concentrations is done by using Monte Carlo techniques (1, 57) or Marcov Chain Monte Carlo techniques (4). The difference between BMFMIN and BMFMAX can be close to 100-fold (Figure 1 and Table 1) due 8400
9
ENVIRONMENTAL SCIENCE & TECHNOLOGY / VOL. 39, NO. 21, 2005
to the individual variation in concentrations between samples. Analyzed contaminant levels in a population of living animals tend to be skewed (Figure 1) with a few animals with more than 10 times higher contaminant levels as compared to the population mean (23). In this paper, we illustrate (Figure 5) that the accuracy of the BMF calculations is closely related to the number of observations. With low numbers of iterations, where 3 iterations simulate 3 individual samples, the calculated BMFRSRGM for CB153 could be lower than 10 or as high as ∼60. The high CV% between the runs with few iterations show that presence or the absence of outliers has a great impact on the resulting BMFs. This also demonstrates how uncertainties in the BMFs (e.g., BMFRSRGM with its GMSD) could be reduced with larger sampling programs. Further, by using a large number of iterations (g10 000), the impact of outliers decreases, and the CV% between the runs becomes as low as, or lower than, 0.5%. The result generated with a larger number of iterations is almost the same as the one obtained after a nonrandom all possibilities calculation that this study gave: 29 (guillemot) * 72 (herring) ) 2088 pairs. By computing the quotients for those pairs and calculating GM and GMSD, the determination of the most appropriate number of iterations is facilitated. Further, for our data, the geometric mean for the different BMFs, the old BMFGM and the new BMFRSRGM, does not differ from each other when we used 50 000 iterations or more. It is important to note that when performed correctly, these two values are identical. The main advantage with the RSR method is that it also gives an estimate of the variation. So, we conclude that when calculating BMFs using the RSR method, it is important to use the number of iterations that is necessary to get identical values for BMFGM and BMFRSRGM as well as reducing the CV% between runs to a minimum (
