徐 冰，劉 存，周東美，王玉軍

（1.中國科學(xué)院土壤環(huán)境與污染修復(fù)重點實驗室（南京土壤研究所），南京 210008；2.中國科學(xué)院大學(xué)，北京 100049）

2014年原環(huán)境環(huán)保部和原國土資源部發(fā)布的《全國土壤污染狀況調(diào)查公報》顯示，受快速的工業(yè)化和城市化進程影響，我國土壤污染點位超標(biāo)嚴(yán)重，其中主要污染物為重金屬和類金屬[1]。農(nóng)作物吸收污染土壤中的重金屬不僅會對農(nóng)作物本身造成減產(chǎn)等危害，還會通過食物鏈傳遞對人體健康造成不良影響。植物根表質(zhì)膜是保持農(nóng)作物內(nèi)部環(huán)境相對平衡穩(wěn)定的屏障，也是土壤溶液中自由態(tài)重金屬離子首先接觸、結(jié)合并進入細胞被農(nóng)作物吸收乃至產(chǎn)生毒性作用的主要通道[2]。因此，了解重金屬離子與植物根表質(zhì)膜之間的相互作用規(guī)律，是從根本上理解重金屬在土壤-植物體系中遷移機制的關(guān)鍵之一。目前，大多數(shù)研究是從宏觀上的生物配體模型[3-5]和分子生物學(xué)上單純的金屬通道蛋白[6-8]等角度探討植物根表質(zhì)膜表面重金屬離子吸附、吸收等過程以及其對植物毒性作用的影響。但重金屬離子與細胞質(zhì)膜相互作用的微觀過程，很難用實驗手段觀測到，因而缺乏分子水平上對其變化過程的描述。

細胞質(zhì)膜因其磷脂分子頭部的磷酸基團、膜蛋白上所帶的氨基酸殘基等基團[9-11]的存在，使得質(zhì)膜表面帶負(fù)電荷，引起植物根表質(zhì)膜-土壤溶液界面上電勢的變化，影響著重金屬離子在質(zhì)膜表面的吸附及土壤溶液中的分布。Kinraide等[12-13]將溶液化學(xué)中經(jīng)典的雙電層Gouy-Chapman-Stern模型（GCS模型）引入植物根表-土壤溶液體系，建立了表征植物根細胞膜表面離子活度的根表面絡(luò)合模型。Kinraide等[13-14]、Wang等[15-16]和Wang等[17]的實驗和模擬數(shù)據(jù)表明GCS模型可以用于量化 Cu2+、Cd2+、Ca2+、Mg2+等陽離子和砷酸根等陰離子對質(zhì)膜膜表面電勢（Ψ0）的影響和實際膜表面電荷分布情況。經(jīng)典GCS模型能夠成功地揭示膜表面電勢影響下不同離子競爭吸附的本質(zhì)規(guī)律，即細胞膜表面分布著負(fù)電和電中性等不同結(jié)合位點，不同離子在這些位點上有著不同的結(jié)合常數(shù)，離子在這些結(jié)合位點的競爭吸附改變了膜表面電勢，并控制離子在擴散層中的分布。但這一概念化數(shù)學(xué)模型為了簡化運算將細胞質(zhì)膜視為剛性平板電容，而實際的細胞質(zhì)膜是結(jié)構(gòu)組成復(fù)雜的兩性雙分子層，表面官能團較多，具有流動性等動態(tài)特性，因此很難通過實驗手段和這一經(jīng)典模型獲取分子水平上重金屬離子與質(zhì)膜相互作用的微觀動態(tài)和機制。分子動力學(xué)（Mo?lecular dynamics，MD）已被證明可以有效模擬細胞質(zhì)膜與不同分子、離子間的相互作用，可以有效地彌補常規(guī)實驗和結(jié)構(gòu)分析方法對分子水平上界面動態(tài)變化認(rèn)識的不足[18-22]。分子動力模擬可以在分子水平上描述重金屬離子在細胞質(zhì)膜表面的分布，以及不同官能團如羰基（C=O）、羥基（-OH）、磷酸基團等對結(jié)合構(gòu)型的影響；模擬得到的運動軌跡可以定量分析質(zhì)膜表面電勢的動態(tài)變化，以及其對溶液中離子濃度分布的影響，有助于深入了解分子尺度下重金屬離子在實際生物質(zhì)膜中吸附、穿透和吸收等過程。因此，本文采用分子動力模擬方法：（1）從分子層面解釋質(zhì)膜表面電勢對溶液中重金屬離子吸附和分布的影響；（2）與經(jīng)典雙電層GCS模型結(jié)合，探究不同陽離子或陰、陽離子共存情況下，質(zhì)膜結(jié)構(gòu)和表面電勢的變化規(guī)律，從而更好地描述質(zhì)膜表面電勢對溶液中重金屬離子吸附和分布影響的分子機制。

1 材料與方法

1.1 模擬體系

模擬體系由磷脂雙分子層（圖1）和溶液相構(gòu)成。溶液分別由約 0.05 mol·L-1的 Na+、K+等一價陽離子，或0.025 mol·L-1的 Ca2+、Mg2+、Cd2+、Cu2+等二價陽離子和53 790個水分子組成：磷脂雙分子層由114個POPC（1-palmitoyl-2-oleoyl-sn-glycero-3-phospho?choline）和 14個 POPG（1-palmitoyl-2-oleoyl-snglycero-3-phosphoglycerol）組成[23-25]（圖2）。模擬體系設(shè)定的周期邊界條件PBC（Periodic boundary condi?tions）：盒子大小為64 ?×64 ?×170 ?，溶液相中另外補充適當(dāng)?shù)腃l-，以保持盒子整體的電中性。所有模擬計算由美國Sandia國家實驗室開發(fā)的開源性大規(guī)模原子分子并行模擬器LAMMPS（Large-scale atomic/molecular massively parallel simulator）軟件完成，計算體系中原子間相互作用力的力場選用CHARMM 36[26-27]，平衡后體系中磷脂雙分子層表面可變電荷密度（σ0）約為25.4 mC·m-2，與實驗觀測結(jié)果相符。

1.2 平衡方法

建立好的模型體系先用正則系綜NVT（N表示原子數(shù)，V表示盒子體積，T表示溫度，三者不變）平衡，初始溫度為300 K，平衡時間1.5 ns，再進行模擬退火0.8 ns，溫度由300 K逐步升至450 K再逐步降至300 K，最后采用等溫等壓系綜NPT（P表示壓力恒定），壓力為0.1 MPa，弛豫時間為1.5 ns。

1.3 模擬過程

平衡后的模型體系繼續(xù)選用非均向性的NPT系綜進行計算并收集數(shù)據(jù)分析，體系溫度設(shè)定為300 K，體系壓強為0.1 MPa，積分步長設(shè)置為2 fs，為保證體系總模擬時長可代表實際情況，模擬過程采用LAMMPS中的SHAKE算法約束53 790個水分子的振動（鍵長或鍵角）以達到延長時間步長的目的，體系總模擬時長為60 ns。

圖1 磷脂雙分子層結(jié)構(gòu)組分POPC和POPG分子式Figure 1 Phospholipid bilayer structure components POPC and POPG formula

1.4 GCS模型計算磷脂表面離子濃度分布

設(shè)定質(zhì)膜雙分子層中心為Z=0，GCS理論計算溶液中陽離子Mn+和Cl-在體系中的密度分布公式為[18]：

式中：nb表示溶液中金屬離子的濃度；κZ=（8πnbe2/εkBT）1/2為Debye-Hückel參數(shù)，其中ε表示水的介電常數(shù)，kB為玻爾茲曼常數(shù)，T為溫度；Φs表示Z=0處的膜表面靜電勢，由式（2）計算得出：

2 結(jié)果與討論

2.1 質(zhì)膜表面的離子吸附構(gòu)型

圖2 植物根表質(zhì)膜與土壤溶液界面的GCS宏觀數(shù)學(xué)模型（a,示意圖修改自文獻[28]）和分子動力學(xué)模型（b）Figure 2 GCS model（a,Adapted from[28]）and molecular dynamics（b）for metal ion-lipid bilayer interactions

圖 3分別為不同溶液體系中Na+、K+、Ca2+、Mg2+、Cd2+、Cu2+等陽離子及陰離子Cl-在質(zhì)膜雙分子層表面的吸附模型。模擬結(jié)果表明，K+在質(zhì)膜表面吸附較弱（圖3b），而其余陽離子都傾向于與質(zhì)膜表面位點強烈結(jié)合[如羰基（C=O）與磷酸基（-OP）]。這是因為與Na+、Ca2+、Mg2+、Cu2+、Cd2+相比，K+的離子半徑更大，離子表面的電荷密度更小，與周圍水分子的結(jié)合能力更弱，水合半徑更小，進而KCl體系的離子強度也比其他體系低，而離子強度的大小又是影響細胞質(zhì)膜表面電勢絕對值的關(guān)鍵。因此K+吸附并進入細胞質(zhì)膜的自由能壘低于Na+以及二價金屬陽離子[29]。

除K+外，其他陽離子在結(jié)合到質(zhì)膜表面時，傾向于與質(zhì)膜頭部的極性基團羥基、磷酸氧和羰基結(jié)合（圖3內(nèi)插圖中內(nèi)容所示），主要是這些官能團的負(fù)電荷與陽離子發(fā)生靜電絡(luò)合作用。金屬陽離子與質(zhì)膜頭部基團的結(jié)合導(dǎo)致了一個與磷脂雙分子層自身方向相反偶極矩的生成。質(zhì)膜頭部帶有正電荷的膽堿基團，也因為與陽離子之間的排斥作用以及與陰離子的吸引作用發(fā)生一定的構(gòu)型改變，傾向于轉(zhuǎn)向水相或者與磷酸基團呈相反的方向。金屬陽離子吸附導(dǎo)致了分子構(gòu)型變化進而促使磷脂面積減小和質(zhì)膜厚度增加，促使質(zhì)膜更加緊密的包裹[19，30-31]。盡管自然狀態(tài)下根細胞質(zhì)膜呈負(fù)電性，通過靜電排斥作用導(dǎo)致陰離子Cl-通常分布在溶液中，隨著距離質(zhì)膜表面越遠，陽離子的濃度越低，陰離子Cl-的濃度越高。但圖3c的Ca2+體系吸附構(gòu)型還顯示出在Ca2+與質(zhì)膜表面的官能團結(jié)合時有Cl-伴隨，證明了溶液相的陰離子可以通過陽離子橋鍵結(jié)合的形式吸附到質(zhì)膜表面。

圖3 金屬離子在質(zhì)膜表面的吸附分布（內(nèi)插圖放大顯示了不同離子與質(zhì)膜官能團的典型結(jié)合方式）Figure 3 The adsorption configuration of Na+，K+，Ca2+，Mg2+，Cu2+and Cd2+on the surface of the plasma membrane，respectively，and the insets enlarged the coordination environment of cations

2.2 磷脂表面密度分布

不同體系中陽離子、陰離子、水分子和質(zhì)膜雙分子頭部的羰基氧、羥基氧、磷酸氧和氮原子在質(zhì)膜表面密度分布情況如圖4表示。原子的位置對脂質(zhì)雙層結(jié)構(gòu)至關(guān)重要。分布峰之間的距離可以直接反映質(zhì)膜極性基團構(gòu)型的變化。以質(zhì)膜中心作為Z軸的原點（Z=0），圖中顯示的為質(zhì)膜雙分子層的一半，密度數(shù)值顯示的是質(zhì)膜兩側(cè)的分布均值。在膜-水界面的周圍，K+的密度分布沒有呈現(xiàn)出明顯的峰值（圖4b），這是因為K+半徑較大，電荷密度低，與負(fù)電基團親和性較弱[19，32]。Na+及二價金屬離子體系中（圖4a和圖4c至圖4f深藍線），陽離子在靠近質(zhì)膜表面的位置出現(xiàn)一個分布峰，其峰值出現(xiàn)的位置和位于脂質(zhì)的羰基（Mn+-Ocarb）和磷酸氧（Mn+-P）最大峰處的分布具有高度的相關(guān)性，羥基基團（Mn+-OH）密度分布的最大值與金屬離子概率密度的峰值重疊部分較少，說明羥基基團對吸附 Ca2+、Mg2+、Cd2+、Cu2+陽離子的貢獻比羰基和磷酸基團少，這一結(jié)論也驗證了不同金屬離子在質(zhì)膜表面吸附構(gòu)型的不同。同時我們觀察到氮原子（Mn+-N）與磷原子（Mn+-P）的密度分布相似，并且受陽離子結(jié)合的影響方式相似。這種構(gòu)型和電荷分布的相關(guān)變化是由于膽堿氮和陽離子之間的相互排斥作用以及陰離子和膽堿基團之間的吸引力引起的[19]。

圖4 不同體系中的金屬離子K+、Na+、Ca2+、Mg2+、Cu2+、Cd2+和Cl-（-Cl）、水分子（-H2O）、羰基氧（-Ocarb）、磷酸氧（-P）、質(zhì)膜頭部氮原子（-N）、羥基氧（-OH）等在質(zhì)膜表面的離子濃度分布Figure4 TheionconcentrationdistributiononthesurfaceoftheplasmamembraneofK+，Na+、Ca2+、Mg2+、Cu2+、Cd2+andCl-（-Cl），water molecule（-H2O），carbonyloxygen（-Ocarb），phosphoricacidoxygen（-P），plasmamembraneheadnitrogenatom（-N）andhydroxyoxygen（-OH）

除 K+以外，Na+和二價金屬離子 Ca2+、Mg2+、Cu2+、Cd2+都在水分子（Mn+-H2O）達到最大分布密度峰附近出現(xiàn)一個倒峰，水分子分布出現(xiàn)最大值的位置為質(zhì)膜吸附金屬陽離子緊密層（雙電層中的Stern層）的外側(cè)，二價金屬離子在從質(zhì)膜中心到體相溶液的分布過程中（距質(zhì)膜雙分子層中心距離由0到100 ?），先出現(xiàn)一個最高峰，隨后出現(xiàn)一個倒峰，然后趨于和伴隨陰離子相似的分布情況，但Na+沒有倒峰出現(xiàn)。Ca2+、Mg2+、Cu2+、Cd2+的倒峰值大小并不相同，這與離子半徑和表面電荷密度等特性有關(guān)。與陽離子分布不同，陰離子主要分布在質(zhì)膜的Stern層以外，由于靜電斥力的相互作用其分布密度與質(zhì)膜表面官能團沒有相關(guān)性，不同陽離子體系中的Cl-滲透到質(zhì)膜的能力也不同，表現(xiàn)為Mg2+＞Ca2+＞Na+＞K+[19]，圖4中顯示，在Ca2+和Mg2+等二價陽離子體系中，Cl-以與陽離子橋鍵結(jié)合的形式進入質(zhì)膜區(qū)域，而K+和Na+體系則沒有出現(xiàn)陰離子伴隨進入質(zhì)膜的情況。

2.3 磷脂表面離子濃度分布與GCS計算比較

圖5 不同體系中金屬陽離子和陰離子在質(zhì)膜表面的離子濃度分布和GCS計算比較Figure 5 The surface ion concentration distribution of metal cations and anion plasma membranes in different systems and compared with GCS model calculation

MD模擬可以準(zhǔn)確預(yù)測不同體系表面電勢下的金屬離子濃度隨距質(zhì)膜雙分子層中心距離變化的分布情況，其結(jié)果和經(jīng)典的GCS理論模型計算結(jié)果具有一致性（圖5）。如圖5所示，通過GCS模型中的公式（1）計算得到的膜表面離子濃度分布規(guī)律為：與質(zhì)膜表面帶電性相反的金屬離子（帶正電），在質(zhì)膜表面的濃度整體呈隨距離減小而增大的變化趨勢（圖5黑線部分），這是由質(zhì)膜（負(fù)電）與金屬離子之間的靜電作用引起的；而體系中的陰離子Cl-在質(zhì)膜表面的濃度整體呈隨著距離的減小而減小的變化趨勢（圖5紅線部分）。這一變化趨勢與MD的模擬結(jié)果一致（圖5淺藍線與黑線重合，深藍線與紅線重合）。MD模擬結(jié)果表明，除KCl體系外（圖5a），Stern層內(nèi)金屬離子由于靜電吸引力的作用濃度明顯比Stern層外大，體系中Cl-在Stern層內(nèi)的濃度則相反。此外，CaCl2、Mg?Cl2、CuCl2和CdCl2等二價金屬離子體系中（圖5c至圖5e中淺藍線），MD計算結(jié)果還表明由Stern層外向Stern層內(nèi)靠近時，出現(xiàn)陽離子濃度先降低后顯著增加的過程，這是由于二價金屬離子與質(zhì)膜表面負(fù)電荷的官能團吸附強烈，吸附的正電荷超過了膜表面帶的負(fù)電荷量，使膜表面電勢由負(fù)變正，所以Stern層外優(yōu)先靜電吸引Cl-，排斥二價陽離子。通過比較MD模擬和GCS模型計算結(jié)果可以發(fā)現(xiàn)，MD模擬不但很好地吻合了GCS理論計算得到的Stern層外隨著Stern距離增加陰、陽離子濃度分布衰減趨勢明顯的變化特點，同時能夠預(yù)測離子在Stern層中的結(jié)合和分布情況。

但GCS理論計算離子在膜表面的電荷濃度公式是基于以下假設(shè)得到的：（1）體系中的離子為點電荷，沒有考慮實際的原子體積和半徑[33]；（2）GCS理論假設(shè)背景溶液為稀溶液，界面上水密度均為1 g·cm-3，不考慮界面相互作用的影響[34]。MD的模擬結(jié)果顯示實際平衡后的體系盒子中（圖2b）表面吸附離子濃度較高時，離子會占據(jù)相當(dāng)?shù)捏w積，且不同離子在體系中占據(jù)的體積大小并不相同，不能簡單視為點電荷；從圖5中可以看出，體系中的水分子由于與磷脂相互間的作用，在界面表面呈現(xiàn)出隨著距雙分子層中心距離的增加濃度呈由低到高的變化趨勢，而并不是均勻分布的，水分子在界面上的分布情況會影響到膜表面電勢的變化和離子的吸附濃度。MD模擬結(jié)果直觀而有力地彌補了GCS理論因簡化假設(shè)而造成的微觀描述缺陷，為拓展GCS模型的使用范圍和驗證MD模擬提供了良好的實驗數(shù)據(jù)。

3 結(jié)論

（1）MD模擬從分子水平上揭示了不同金屬陽離子與質(zhì)膜結(jié)合的構(gòu)型，發(fā)現(xiàn)質(zhì)膜表面電勢強烈影響著離子的分布，二價陽離子主要結(jié)合在質(zhì)膜頭部羰基基團和磷酸基團。

（2）MD模擬可以與GCS模型很好地吻合，并且能進一步準(zhǔn)確描述雙電層中的Stern層中陰、陽金屬離子的分布和結(jié)合位點，進而更好地為研究重金屬的農(nóng)作物吸收和毒性提供分子層面的機理解釋和直觀微觀圖像。