Science in China Series B: Chemistry © 2009
www.scichina.com chem.scichina.com www.springerlink.com
SCIENCE IN CHINA PRESS
Springer
Modeling the depuration rates of polychlorinated biphenyls in two mussel species with theoretical molecular descriptors XU MingZhu, LIU XinHui†, WANG Liang, WU Dan, SUN Tao, YANG ZhiFeng & CUI BaoShan State Key Laboratory of Water Environment Simulation, School of Environment, Beijing Normal University, Beijing 100875, China
Using theoretical molecular descriptors as well as partial least squares (PLS) regression, two quantitative structure-activity relationship (QSAR) models were developed for depuration rate constants (kd) of polychlorinated biphenyls (PCBs) in two species of mussels, Perna viridis and Dreissena polymorpha. The 2 cross-validated Q cum (an indicator of fitting of goodness) values for the two models are 0.501 and 0.756, 2 and the standard deviation (SD) is 0.084 and 0.076, respectively. The achievement of satisfactory Q cum
and low SD values indicates good predictive ability and precision of the two models. The significant descriptors governing lgkd include polarizability (α), molecular volume (MV), molecular weight (Mw), molecular surface area (S), and total energy (TE). The key descriptors in the models reflect that van der Waals interactions play a dominant role in the depuration of PCBs. The depuration of PCBs in the two mussel species may be mainly attributed to the biota-water phase partitioning processes. quantitative structure-activity relationship (QSAR), depuration, polychlorinated biphenyls (PCBs), partial least squares (PLS)
1
Introduction
Since 1929, polychlorinated biphenyls (PCBs) have been largely produced and widely utilized in various commercial applications. They are ubiquitous, persistent and highly lipophilic environmental pollutants. Through various transportation routes, they may go into water bodies and eventually sink in the sediments, which may cause hazards to the aquatic animals. There is growing evidence that these compounds can accumulate in animal body tissues[1, 2]. By the biomagnification processes, they may accumulate in predators at the top of a food chain and cause a wide range of unfavorable results to both animals and humans[3, 4]. Monitoring programmes of PCBs pollution in water bodies have been deployed worldwide. The method using biological indicators to monitor trace toxic substances in the aquatic environment has been well established[5, 6]. Mussels are commonly preferred in biomoni-
toring programmes due to their wide geographical distribution, easy sampling, sessile lifestyle, long life-span, resistance to stress, high accumulation as well as low metabolizing rates to a large range of chemicals and the important role linking to human health. As surrogates for measuring the toxic substances in the aquatic environment, mussels must be calibrated to determine the time required to achieve chemical equilibrium with the environment. The time to achieve the steady state can be calculated from the depuration rate constants of chemicals in organisms. The depuration rate constant is one of the most important kinetic parameters when modeling contaminants cycling. It can be used to estimate uptake rate constants, bioconcentration factors, and the time to steady Received August 30, 2008; accepted November 11, 2008 doi: 10.1007/s11426-009-0042-y † Corresponding author (email:
[email protected]) Supported by the Major State Basic Research Development Program (973) of China (Grant No. 2006403303), the National Natural Science Foundation of China (Grant No. 20647003), and the Program for New Century Excellent Talents in University
Sci China Ser B-Chem | Aug. 2009 | vol. 52 | no. 8 | 1281-1286
state[7, 8]. Accordingly, it is of great significance to explore the depuration rate constants of hazardous PCBs. However, research focusing on the PCBs depuration in mussels is hampered by high cost, and the time-consuming processes. Fortunately, quantitative structure-activity relationships (QSAR) could provide an effective method to expediently interpret depuration in bivalves with low cost. The main idea of QSAR is to construct the correlation by expressing a special activity property of organic compounds in terms of appropriate molecular descriptors. Once the QSAR is constructed, it may clarify which aspect of the molecular structure affects molecular activity, and the activities of chemicals can be predicted in the absence of experimental data. QSAR has been widely used in the research of toxicity, affinities, bio-concentration factors, and physico-chemical properties of chemicals[9−12]. Moreover, due to its easy operation and powerful predictability, it is widely used in the field of ecological risk assessment[13]. Some theoretical molecular descriptors can not only clearly describe some defined molecular properties, but also are not restricted to closely related compounds. Using theoretical molecular descriptors, many robust QSAR models have been developed[14−16]. Amid various model developing techniques for QSAR, PLS is adopted for the present study because of its strong ability to analyze data with noisy, collinear and even incomplete variables in both the independent and dependent variables. Moreover, PLS has the desirable property that the precision of the model parameters can be improved with the increasing number of relevant variables and observations[17]. In this paper, the QSAR models for the PCBs depuration rates in two mussel species are developed using theoretical molecular descriptors as well as PLS regression. In addition, the depuration mechanisms of PCBs in the mussels are discussed and deduced.
2 Materials and methods 2.1
The data set
The depuration rate constants in the unit of d−1 were taken from refs. [7,18], and were converted into the form of lgkd for QSAR analysis. They are listed in Table 1. Based on the general first order eq. (1), depuration rate constants (kd) were determined from the slope of the line representing the changes in the log Ct vs. the time of depuration. k t log Ct = log Ct = 0 − d , (1) 2.303 1282
where Ct is the PCB concentration in the organism at time t, Ct=0 is the initial concentration at time zero for the depuration period. Both Ct and Ct=0 are lipid-normalized. 2.2
Descriptor generation
Theoretical molecular descriptors of PCBs were obtained using the MNDO method[19] contained in MOPAC 2000 of the CS Chem3D Ultra (Cambridge Soft Corp., 2000). The molecular structures were optimized using eigenvector following a geometry optimization procedure, in which the geometry optimization criterion of Gradient Norm was set at 0.1. A total of 12 molecular descriptors listed in Table 2 were adopted in the present study. Their corresponding values are available upon request. The molecular surface and molecular volume values were calculated by the COSMO method[20]. 2.3
PLS analysis
The QSAR models with theoretical molecular descriptors were developed using PLS regression, as implemented in the Simca-S package (Umetrics AB, Sweden). The conditions for the computation were based on the default options of the software. The model adequacy is mainly as2 sessed by the Qcum , the standard deviation, the correla-
tion coefficient between observed values and fitted values (r), the number of principle component (A) and the significance level (p). The obtained QSAR models with the 2 smaller SD values may be preferred, and Qcum > 0.5 is
commonly accepted as satisfactory result. SD =
n 2 1 ⎡ ⎤ ; − ∑ n − A − 1 i =1 ⎣⎢Y im,obs Y im,pred ⎦⎥
⎪⎧ ⎡ 2 = 1.0 − ⎨∏ ⎢ ∑∑ Yim,obs − Yim,pred Qcum ⎪⎩ ⎣ i m
(
)
(2)
⎪⎫ ⎥ SS ⎬ , (3) ⎪⎭a ⎦
2⎤
where Yim,obs and Yim,pred denote the observed and predicted lgkd values, respectively. i stands for different observations, m stands for different dependent variables (m = 1 for this study), SS is the residual sum of squares of the previous component, and a = 1, 2,…, A.
3 Results and discussion 3.1
Model fitting results
Based on the unscaled pseudo-regression coefficients of the independent variables and a constant transformed from PLS results, QSAR models 1 (eq. (4)) and 2 (eq. (5)) were obtained for Perna viridis and Dreissena
XU MingZhu et al. Sci China Ser B-Chem | Aug. 2009 | vol. 52 | no. 8 | 1281-1286
Table 1 Experimental and predictive lgkd values of PCBs PCBs
Model 1
Model 2 Pred. b)
Res. c)
Obs.a)
Pred. b)
Model 2
Obs.a)
Pred. b)
Res. c)
PCB18
−0.97
−0.85
−0.12
PCB74
PCB20
−0.90
−0.85
−0.05
PCB83
−1.08
−0.94
−0.14
PCB33
−0.93
−0.86
−0.07
PCB85
−0.99
−0.94
−0.05
PCB40
−0.80
−0.85
0.05
PCB87
−1.01
−1.05
0.04
PCB42
−0.91
−0.86
−0.05
−0.81
−0.90
0.09
PCB97
−0.97
−1.02
0.05
PCB44
−0.92
−0.86
−0.06
−1.00
−0.90
−0.10
PCB99
−0.95
−0.96
0.01
−1.04
−1.07
0.03
PCB47
−0.83
−0.88
0.05
PCB101
−0.93
−0.95
0.02
−1.00
−1.04
0.04
−1.00
−0.94
−0.06
PCB102
−0.90
−0.95
0.05
PCB105
−0.96
−0.96
0.00
−1.19
−1.04
−0.15
−0.97
−0.94
−0.03
−1.03
−0.99
−0.04
−1.17
−1.13
−0.04
−1.01
−1.10
0.09
−1.22
−1.24
0.02
−1.38
−1.23
−0.15
−1.15
−1.21
0.06
PCB49
Obs.a)
Model 1
PCBs
Res. c)
PCB51
−0.73
−0.86
0.13
PCB52
−0.90
−0.87
−0.03
PCB53
−0.74
−0.86
0.12
PCB118
−0.99
−1.01
0.02
PCB60
−0.93
−0.94
0.01
PCB128
−0.97
−1.01
0.04
PCB66
−0.92
−0.90
−0.02
PCB132
−1.10
−0.99
−0.11
−0.76
PCB64
PCB70
−0.94
−0.93
−0.01
PCB136
−0.92
−0.99
0.07
PCB138
−1.08
−1.01
−0.07
PCB141
−0.75
−1.02
0.27
PCB146
−0.97
−0.92
−1.00
0.16
0.03
PCB110
PCB129
PCB134
−0.95
−1.01
0.06
PCB177
−1.19
−1.07
−0.12
PCB178
−1.08
−1.07
−0.01
PCB180
−0.98
−1.08
0.10
−1.11
−1.13
0.02
−1.21
−1.13
−0.08
PCB183
Obs. a)
Pred. b)
Res. c)
−1.00
−1.05
0.05
PCB151
−0.97
−1.01
0.04
−1.08
−1.13
0.05
PCB185
PCB153
−1.10
−1.02
−0.08
−1.33
−1.18
−0.15
PCB187
PCB156
−0.98
−1.04
0.06
PCB194
−1.28
−1.32
0.04
PCB159
−0.96
−1.06
0.10
PCB195
−1.29
−1.29
0.00
PCB170
−0.99
−1.08
0.09
−1.27
−1.30
0.03
−1.24
−1.30
0.06
PCB203
−1.29
−1.33
0.04
PCB206
−1.34
−1.40
0.06
PCB200
PCB171
−1.24
−1.20
−0.04
PCB201
PCB172
−1.28
−1.22
−0.06
PCB202
−1.28
−1.20
−0.08
PCB173
−1.10
−1.08
−0.02
PCB174 PCB176
−1.09
−1.07
−1.17
−1.19
−1.23
−1.08
−0.09
−1.14
−0.05
−1.13
−0.10
−0.02
a) The observed value; b) the predicted value; c) the residual value between the observed and predicted value. Table 2 List of molecular structural descriptors of PCBs Descriptors
Definitions
Mw IP TE
molecular weight ionization potential
Ehomo Elumo QH+ QC− QCl+ μ α S MV
total energy energy of the highest occupied molecular orbit energy of the lowest unoccupied molecular orbit most positive net atomic charges on hydrogen atoms most negative net atomic charges on carbon atoms most positive net atomic charges on chlorine atoms dipole moment molecular polarizability cosmo molecular surface area cosmo molecular volume
Polymorpha, respectively. lg kd 1 = −1.13 × 10−1 − 2.94 × 10−3 α − 5.36 × 10−4 M w +5.42 × 10−5 TE + 9.19 × 10−2 Elumo ,
(4)
lg k d 2 = 3.94 × 10−1 − 3.79 × 10−3 α − 2.08 × 10−3 S − 1.97 × 10 −3 M V + 1.64 × 10 −1 Elumo −1
−
(5)
−3
− 7.79 × 10 QC + 9.53 × 10 μ .
The model fitting results are listed in Table 3. The two 2 models yield the Qcum values of 0.501 and 0.756, and
the conventional correlation coefficient r2 values of 0.526 and 0.778, respectively.
XU MingZhu et al. Sci China Ser B-Chem | Aug. 2009 | vol. 52 | no. 8 | 1281-1286
1283
Table 3 Models fitting results Model 1 2
N a) 39 31
RX2
A 1 1
(adj.)(cum)
b)
RY2 (adj.)(cum)
0.750 0.490
2 Qcum
c)
0.513 0.771
0.501 0.756
r2
SD
0.526 0.778
0.084 0.076
p −6
<10 <10−10
a) The number of compounds in the models; b) and c) stand for the cumulative variances of all the independent and dependent variables explained by all the extracted components.
One PLS component was selected for the two models. Model 1includes four descriptors, while there are six descriptors in model 2. In a PLS model, variable importance in the projection (VIP) is a parameter showing the influence on dependent variable of every independent variable. Terms with large VIP values are the most relevant for explaining dependent variable. The VIP values of descriptors included in the two models are listed in Table 4. Descriptors with VIP values larger than one were believed to be significant, which were only discussed in the present study.
estimation and/or model computation. The satisfactory 2 Qcum values as well as the low SD and p values for the
models validate precision of the two models. Accordingly, the models can be used to predict kd values of PCB congeners in corresponding species of mussels.
Table 4 The VIPs for the molecular structural descriptors in models Variables α Mw TE S
VIP values model 1 model 2 1.106 1.284 1.068 1.068 1.258
Variables MV QC− Elumo μ
VIP values model 1 model 2 1.254 0.391 0.706 0.945 0.386
For the models, plots of observed and predicted values are shown in Figures 1 and 2. According to the figures, most of the predicted lgkd values were within a 95% confidence interval except for PCB141 in model 1.
Figure 1 Observed versus predicted lgkd values for Model 1.
3.2 Model validation The model is validated by leave-one-out cross-validation. According to the model fitting results listed in Table 3, 2 values for the models obtained are larger than the Qcum
0.5, indicating good robustness and predictive ability of the models. r2 values for the two models are satisfactory. Comparing the models, it can be concluded that the quality of Model 1 is slightly lower than that of Mode 2 2 in terms of Qcum and r2 values. However, since kd values
for most PCBs in Perna viridis are still unknown. Model 1 can give rough estimation. Observing Table 1, it can be concluded that the biggest differences between the observed and predicted lgkd values for Perna viridis and Dreissena polymorpha are 0.27 and 0.16 respectively, indicating the possible errors of prediction in the scale of logarithms. The relative larger residual for PCB141 in model 1 may be due to the uncertainty of experimental 1284
Figure 2 Observed versus predicted lgkd values for Model 2.
3.3
Model interpretation
Examining Table 1 leads to the conclusion that Dreissena polymorpha demonstrated much slower PCBs de-
XU MingZhu et al. Sci China Ser B-Chem | Aug. 2009 | vol. 52 | no. 8 | 1281-1286
puration than that observed for Perna viridis. This may be caused by the dissimilarities in tissue compositions among different species of mussels, such as differences in tissue lipid contents believed to govern the animal partitioning capacity for hydrophobic chemicals[21,22]. Other factors, such as the gill surface to volume ratio, feeding states, and age, may also contribute to the difference in depuration rate[7,23,24]. It could be concluded that PCBs depuration rate constants are specific to different species of mussels. Thus it might be more reasonable to develop different models to predict kd values in different mussel species. It has been generally accepted that equilibrium partitioning is the major factor determining the uptake and release rates of some lipophilic pollutants in gill-breathing aquatic animals[19,21,22]. Various intermolecular interactions, such as van der Waals interactions (including dispersive interactions, dipole-dipole interactions, and dipole- induced dipole interactions), electrostatic interactions and hydrogen bonding, may govern the magnitude of partitioning property of chemicals. Examining the models leads to the observation that polarizability plays a significant role in PCBs depuration. It possesses the largest VIP values in the two models. The polarizability of a molecule is a measure of its ability to acquire an electric dipole moment and to respond to an electric field. It is a tensor quantity that relates changes in three-dimensional electron density distribution of the molecule to the strength of applied electric fields. α contributes negatively to lgkd in the models. This makes sense since PCB molecules with larger α values may have stronger intermolecular dipole-induced dipole interactions with lipid molecules considering that the polarizability of water molecules is lower compared with lipid molecules. The dipole- induced dipole interactions hinder the partitioning of PCBs into water phase. VIP values of the descriptor Mw as well as TE in 1
model 1 and VIP values of the descriptors S and MV in model 2 are larger than 1. This implies that the four descriptors are significant variables. All of the four descriptors generally reflect molecular size. Generally, Mw, S and MV are positively correlated to molecular size while TE is negatively related to molecular size. In the models, Mw, S and MV contribute negatively to lgkd while TE is positively correlated to lgkd. This is consistent with greater dispersive forces of larger molecules which facilitate their partitioning into the less polar biophase. Based on the descriptors, it can be concluded that intermolecular dispersive interactions contribute significantly to lgkd and the larger the PCB congeners are, the lower the corresponding lgkd values will be. Obviously, the depuration rate constants are closely related to dipole-induced dipole and dispersive interactions. Important descriptors involved in the models show that the main force dominating the depuration processes is van der Waals interactions. Both models provide additional support to the concept of equilibrium partitioning as the primary mechanism for depuration of the hazardous PCB congeners in mussels.
4 Conclusions Using theoretical molecular descriptors and PLS method, two optimal QSAR models were obtained for the depuration rate constants of PCB congeners in two mussel species. The descriptors involved in the QSAR models generally reflect that van der Waals interactions dominate the PCBs depuration kinetics in the two species of mussels. It is clear that the biota-water partitioning processes may play an important role in the depuration mechanisms. Using these models, the depuration rate constants of PCB congeners in Perna viridis and Dreissena polymorpha could be predicted conveniently and they can be applied to computing the time for steady state for PCBs in the two species of mussels.
Naito W, Jin J C, Kang Y S, Yamamuro M, Masunaga S, Nakanishi J.
mental and mechanistic considerations which support the develop-
Dynamics of PCDDs/DFs and coplanar-PCBs in an aquatic food chain
ment of toxic equivalency factors (TEFs). Cri Rev Toxicol, 1990,
of Tokyo Bay. Chemosphere, 2003, 53(4): 347-362 2
Wan Y, Hu J Y, Yang M, An L H, An W, Jin X H, Hattori T, Itoh M.
4
Characterization of trophic transfer for polychlorinated dibenzo-p-
Petrov E A. Toxicokinetics of PCDD, PCDF, and coplanar PCB
dioxins, dibenzofurans, non- and mono-ortho polychlorinated bi-
congeners in Baikal Seals, Pusa sibirica: Age-related accumulation,
phenyls in the marine food web of Bohai Bay, north China. Environ
maternal transfer, and hepatic sequestration. Environ Sci Technol,
Sci Technol, 2005, 39(8): 2417-2425 3
21(1): 51-58 Iwata H, Watanabe M, Okajima Y, Tanabe S, Amano M, Miyazaki N,
Safe S. Polychlorinated biphenyls (PCBs), dibenzo-p-dioxins (PCDDs), dibenzofurans (PCDFs) and related compounds: Environ-
5
2004, 38(18): 3505-3513 Galassi S, Bettinetti R, Neri M C, Jeannot R, Dagnac T, Bristeau S, Sakkas V, Albanis T, Boti V, Valsamaki T, Falandysz J,
XU MingZhu et al. Sci China Ser B-Chem | Aug. 2009 | vol. 52 | no. 8 | 1281-1286
1285
6
Schulte-Oehlmann U. A multispecies approach for monitoring per-
quantitative structure-activity relationship studies of alkyl (1- phen-
sistent toxic substances in the Gulf of Gdańsk (Baltic sea). Ecotoxicol
ylsulfonyl) cycloalkane-carboxylates. Chemosphere, 1997, 35(3):
Environ Saf, 2008, 69(1): 39-48 Riva C, Binelli A, Provini A. Evaluation of several priority pollutants
15
623-631 Amat L, Carbó-Dorca R. Simple linear QSAR models based on quan-
16
tum similarity measures. J Med Chem, 1999, 42(25): 5169-5180 Liu X H, Wang B, Huang Z, Han S K, Wang L S. Acute toxicity and
in zebra mussels (Dreissena polymorpha) in the largest Italian subalpine lakes. Environ Pollut, 2008, 151(3): 652-662 7
Morrison H, Yankovich T, Lazar R, Haffner G D. Elimination rate
quantitative structure-activity relationships of α-branched phenyl-
constants of 36 PCBs in zebra mussels (Dreissena polymorpha) and
sulfonyl acetates to Daphnia magna. Chemosphere, 2003, 50(3): 17
403-408 Svante W, Michael S, Lennart E. PLS-regression: a basic tool of
18
chemometrics. Chemometr Intell Lab, 2001, 58(1): 109-130 Tanabe S, Tatsukawa R. Mussels as bioindicators of PCB pollution: a
exposure dynamics in the Lake St. Clair-Lake Erie corridor. Can J Fish Aquat Sci, 1995, 52(12): 2574-2582 8
Uno S, Shiraishi H, Hatakeyama S, Otsuki A. Uptake and depuration kinetics and BCFs of several pesticides in three species of shellfish (Corbicula leana, Corbicula japonica, and Cipangopludina chinen-
case study on uptake and release of PCB isomers and congeners in
sis): Comparison between field and laboratory experiment. Aquat
green-lipped mussels (Perna viridis) in Hong Kong Waters. Environ
Toxicol, 1997, 39(1): 23-43 9
Dai J Y, Jin L J, Yao S C, Wang L S. Prediction of partition coefficient
19
and toxicity for benzaldehyde compounds by their capacity factors
method approximations and parameters. J Am Chem Soc, 1977, 20
99(15): 4899-4907 Connolly M L. Analytical molecular surface calculation. J Appl
21
Crystallogr, 1983, 16(11): 548-558. Hamelink J L, Waybrant R C, Ball R C. A proposal: exchange equi-
and various molecular descriptors. Chemosphere, 2001, 42(8): 899-907 10
Wei D B, Zhang A Q, Wu C D, Han S K, Wang L S. Progressive study and robustness test of QSAR model based on quantum chemical pa-
11
rameters for predicting BCF of selected polychlorinated organic
libria control the degree chlorinated hydrocarbons are biologically
compounds. Chemosphere, 2001, 44(6): 1421-1428
magnified in lentic environments. Trans Amer Fish Soc, 1971, 100(1): 22
207-214 Mackay D. Correlation of bioconcentration factors. Environ Sci
23
Technol, 1982, 16(5): 274-278 Baumard P, Budzinski H, Garrigues P. PAHs in Arcachon Bay, France:
Wang X D, Lin Z F, Yin D Q, Liu S S, Wang L S. 2D/3D-QSAR comparative study on mutagenicity of nitroaromatics. Sci China Ser B-Chem, 2005, 48(3): 246-252
12
Wang Z Y, Han X Y, Zhai Z C. QSPR to aqueous solubility (lgSw) of alkyl (1-phenylsulfonyl) cycloalkane-carboxylates using MLSER
Origin and biomonitoring with caged organisms. Mar Pollut Bull,
model and ab initio. Chemosphere, 2006, 62(3): 349-356 13
14
Pollut, 1987, 47(1): 41-62 Dewar M J S, Thiel W. Ground states of molecules 38 the MNDO
Chen J W, Li X H, Yu H Y, Wang Y N, Qiao X L. Progress and per-
24
1998, 36(8): 577-586 Bruner K A, Fisher S W, Landrum P F. The role of the zebra mussels,
spectives of quantitative structure-activity relationships used for
dreissena polymorpha, in contaminant cycling: Ⅱ zebra mussel
ecological risk assessment of toxic organic compounds. Sci China Ser
contaminant accumulation from algae and suspended particles, and
B-Chem, 2008, 51(7): 593-606
transfer to the benthic invertebrate, gammarus fasciatus. J Great
Chen J W, Wang L S. Using MTLSER model and am1 hamiltonian in
Lakes Res. 1994, 20(6): 735-750
1286
XU MingZhu et al. Sci China Ser B-Chem | Aug. 2009 | vol. 52 | no. 8 | 1281-1286