洞庭湖镉迁移转化的马尔可夫模型评价_张建新.pdf

第 26 卷第 4 期 2007 年 8 月 岩矿测试 ROCK AND MINERAL ANALYSIS Vol． 26，No． 4 August， 2007 文章编号 02545357 2007 04033904 洞庭湖镉迁移转化的马尔可夫模型评价 张建新，邢旭东 湖南省地质研究所，湖南 长沙410007 摘要 讨论了在现时泥沙、 水体、 鱼类、 藻类的含镉状态下洞庭湖的输入输出通量; 并以其为例， 在简要介绍马 尔可夫链数学方法的基础上， 建立了湖泊生态系统中镉迁移转化的马尔可夫链数学模型， 初步评价了当洞庭湖不 断有镉加入、 系统到达状态稳定时各种介质中的镉总量及镉的迁移转化趋势。 关键词 生态地球化学; 马尔可夫链; 湖泊生态系统; 镉通量 中图分类号 O211． 62 P595文献标识码 A 收稿日期 2007- 03- 12; 修订日期 2007- 04- 25 基金项目 国土资源地质大调查 湖南省洞庭湖区生态地球化学调查项目资助 1212010310304 作者简介 张建新 1958 ， 男， 湖南长沙市人， 教授级高工， 主要从事环境地球化学、 农业 生态 地学研究工作。 E- mail zjx0731163． com。 uation on Migration and Conversion of Cd in the Dongting Lake Based on Markov Model ZHANG Jian- xin，XING Xu- dong Geological Institute of Hunan Province，Changsha410007，China Abstract The and output fluxes of heavy metal cadmium are studied based on the data of cadmium distribution in media of lake sediment，water body， fishes， algae in the Dongting lake． On the basis of the results from this study， a Markov chain model for the prediction of migration and conversion of Cd in lake ecosystem is established． Using this model，a conclusion that Cd is continu- ously imported is obtained and Cd contents in different media in the lake are calculated when the system reaches the stable state． And migration and conversion trends of Cd are predicted． Key words ecology- geochemistry; Markov chain; lake ecosystem; cadmium flux 湖南的1 ∶ 25 万多目标地球化学调查发现洞庭湖区土 壤中的 Cd 含量高， 是我国土壤 Cd 丰度最高的区域。如何 正确评价其生态效应及发展趋势， 是目前生态地球化学评价 中亟待解决的方法问题。 马尔可夫链 Markov chain 作为一种有用的数学分析 工具， 目前已广泛应用于社会经济和自然科学各领域的预 测和评价中。其在我国地学研究中的应用始于 20 世纪 70 年代 ［1 －3 ］， 如运用马尔可夫模型 Markov model 查明火山 岩系的喷出顺序和侵入体形成的先后顺序， 划分矿床的成 矿期和成矿阶段、 揭示各个成矿阶段的空间分布等， 均与地 层剖面的马尔科夫过程模拟有关， 而地质体中元素转移趋 势 概率 的研究则鲜有涉及。 本文拟通过马尔可夫链的基本数学模型方法 ［4 ］的介 绍， 试对洞庭湖中 Cd 的生态地球化学循环建立马尔可夫 过程数学模型， 研究 Cd 在洞庭湖生态系统中的持续时间 和变化趋势， 为区域生态地球化学评价提供科学依据 ［5 ］。 1马尔可夫链的数学基础 1． 1马尔可夫链的转移矩阵 transition matrix 记事物在 t 时刻所处状态为 x t ， 称依赖于参数 t， 以 一定概率取值于某一状态的过程为随机过程 stochastic process 。当时间 t0时随机过程的状态为 i， 经 t t≥0 时 间， 状态达到 j， 过程的概率记作 pij t0， t0 t ， 若此概率仅 与时间 t0时的状态 i 有关， 而与 t0前所处的状态无关， 则 过程是无后效性的。无后效性的随机过程即为马尔可夫过 程 Markov process 。状态离散、 时间离散的马尔可夫过程 称为马尔可夫链。 在马尔可夫链理论中， 将从状态 i 经过 1 步转变成 j 状 态的概率 pij 1 称为转移概率 transition probability ， 记作 pij。让 i、 j 跑遍所有状态， pij ij 1， 2， ， n 排成的矩阵 933 ChaoXing P p11p12p1n p21p22p2n pn1pn2p           nn 1． 1 称为转移矩阵。它的第 i 行向量{ pil， pi2， ， pin} 分别代表 从 i 状态经一步达到所有状态的概率。显然， 它的和代表 必然事件， 概率等于1; 从状态 i 经过 k 步转变成 j 状态的转 移概率是矩阵乘积 pk的第 i 行 j 列元素值; 状态分配比率 为随机向量的事件， 经过 k 步转移后， 新状态分配比率向量 X k XPk。 1． 2吸收马尔可夫链 absorbing Markov chain 当马尔可夫链中状态 i 的转移概率满足两条件 ① pii1， ② pik 0 k≠i ， 则该状态为吸收状态 absorbing state 。通俗而言， 吸收状态就是不能离开的状态， 它从出现该 状态的所在行到自身的列的转移概率为1， 其余概率均为0。 吸收马尔可夫链是指至少存在一个吸收状态或从任何 状态经有限步终可到达一个吸收状态的马尔可夫链。 吸收马尔可夫链中的非吸收状态称为转移状态 transient state 。 对于吸收马尔可夫链， 从任何状态出发最终进入吸收 状态的概率为 1， 最终不能进入吸收状态的概率为 0。 可将吸收马尔可夫链的转移矩阵写成统一的规范形 式， 称为典范式 canonical P EO [] RQ 1． 2 其中分块矩阵 E 是 r r 阶单位矩阵; O 是 r k 阶零矩阵; R 是k r 阶矩阵， 表示从转移状态一步就达到吸收状态的概 率; Q 是 k k 阶矩阵， 表示从一个转移状态转入新的转移 状态 包括自身 的概率。 关于吸收马尔可夫链矩阵的运算有下述定理。 定理 1． 1 典范式转移矩阵的乘幕按矩阵分块运算有 以下结果。 Pt EO RtQ []t 1． 3 其中 Rt I Q Q2 Qt －1 R， I 是 k k 阶单位矩阵。 定理 1． 2 对于吸收马尔可夫链典范式的转移矩阵 ， ① 当 t≠0 时， Qt≠O; ② 矩阵 I －Q 可逆; ③ N I － Q －1 I Q Q2 ， 其中 N 是 k 阶方阵， 称为该吸收马尔可 夫链的基本矩阵 fundamental matrix 。 换言之， 基本矩阵 N 中的每个元素表示从一个非吸收 状态出发， 过程到达另一个非吸收状态的平均转移次数。进 一步推论， 就是在所有转移状态之间传递步数的数学期望值 是基本矩阵 N I －Q －1的第 i 行元素之和， 即在转移状态 之间传递的步数也就是进入每个转移状态步数的总和。 定理 1． 3 如果有多个吸收状态， bij表示从非吸收状态 i 出发被吸收状态 j 吸收的概率， 则不同吸收状态的概率分 配值 B NR。式中 R 是典范式中的吸收状态， N 是基本矩 阵。从而解决了对最终达到不同吸收状态的可能性或比率 的估计问题。 1． 3带输入的马尔可夫链 假设 Y ［ y1y2 yr］ 代表吸收状态， X ［ x1x2 xk］ 代表转移状态， Z ［ Y X］ 表示初始时每个状态的目标元素 如 Cd 的数量， F ［ f1f2 fr］表示每次补充的目标元素 的数量， Y t、 X t和 Z t分别是时间 t t 0， 1， 2， 时向量 Y、 X 和 Z 所处的状态值。 经 k 步后， Z k Z k －1P ［ O F］ 。 其中向量［ O F］的前 r 个分量属于吸收状态， 即零矩 阵; Z k［ Y kX k］ k 0， 1， 2， 。 称序列 Z 0Z 1Z 2为带输入的吸收马尔可夫链的 状态向量序列， 它由 r 维吸收部分和 k 维转移部分构成， F 为它的输入向量。 带输入的吸收马尔可夫链有如下定理。 定理1． 4 带有输入的吸收马尔可夫链， 输入向量为 F， 则其状态向量序列的转移向量部分 X k k 0， 1， 2， 存 在极限 FN， 其中 N 是吸收马尔可夫链的基本矩阵。该定 理说明随着时间的推移， 具有带输入的吸收马尔可夫序列 的最后状态 FN 与初始状态 X 无关， 仅与输入向量 F 有关; 初始状态只是暂时的， 它对系统的影响将逐渐湮灭; 而输入 向量却是长期持续的因素， 它的影响随时间增加而逐渐增 强， 直至完全控制系统。 2洞庭湖生态系统镉迁移转化的通量 2． 1洞庭湖泥沙输入输出及其镉含量 表 1 为洞庭湖历年输沙量变化情况 ［6 ］。其中“长江” 指经长江分洪河道 松滋、 太平、 藕池和调弦四口 ， “四水” 指湘、 资、 沅和澧四水。从表可见， 近 50 年来， 洞庭湖之长 江及四水源分段的泥沙输入、 输出量变化不大， 可视其为 基本恒定态。 表 1洞庭湖年输沙量统计 Table 1The statistic data of sand imput into Dongting lake per year 时间段/ 年份 m入湖/亿吨 长江湖南四水∑ m四口 m四水 m出湖/亿吨 城陵矶口 m入湖 m出湖 1951 ～19582． 20390． 43762． 46155． 040． 69793． 785 1959 ～19661． 90430． 28342． 18776． 720． 57793． 786 1967 ～19721． 44160． 40791． 47393． 530． 52403． 530 1973 ～19801． 10760． 36631． 47393． 020． 38323． 846 1981 ～19901． 09170． 23811． 32984． 590． 32094． 144 1991 ～19980． 73450． 23030． 96483． 190． 26273． 673 1951 ～19991． 39610． 33411． 73024． 180． 44743． 867 表 2 为各流域泥沙输入输出 ［7 ］及其 Cd 含量情况。 洞庭湖底泥 Cd 含量据鲍征宇等数据平均值为 2． 7 mg/kg。另洞庭湖面积以 2600 km2、 淤泥密度为 1． 3 t/m3、 淤泥平均厚 2 m 计， 可估算洞庭湖底泥中 Cd 总含量为 18． 252 t。 043 第 4 期 岩矿测试 http ∥ykcs． i3t． com． cn/ 2007 年 ChaoXing 年出湖泥沙 Cd 总含量 2． 412 t 与湖底泥中 Cd 总含 量 18． 252 t 的比值即为 Cd 的年输出率为 0． 13215。 2． 2洞庭湖镉年输入量计算 据鲍 征 宇 等报 道， 洞 庭 湖 水 体 Cd 平 均 含 量 为 0． 00024 mg/kg; 鱼类的 Cd 平均含量约为0． 3 mg/kg; 芦苇 Cd 平均含量为 2． 25 ～5． 13 mg/kg， 藻类 Cd 平均含量取芦 苇平均， 约为 3 mg/kg。 据调查 ［8 ］， 可假设长江干流与洞庭湖之间鱼类交流率 为 5; “四水” 内繁殖的家鱼苗对补充洞庭湖家鱼资源意 义不大， 交换量为零。据洞庭湖渔业资源总量和捕捞总 量 ［7 ］可估计出鱼类的捕捞率约为 10。 表 2洞庭湖各流域带入泥沙及 Cd 含量① Table 2The amount of sediment imput from each drainage area of Dongting lake and Cd content 水系 名称 湖南省 内流域 面积/ km2 水土 流失 面积/ km2 平均侵蚀 模数/［t/ km2a ］ 泥沙 总量/ 108t 年均泥沙 量/ 108 t/a 泥沙 Cd 平均含 量/ g/t 年输入 Cd 总量/ t/a 湘 水853831956433650． 15850． 003171． 90． 602 资 水26738715630340． 03860． 000770． 340． 026 沅 水510061195533290． 12220． 002440． 450． 110 澧 水15505481449090． 08640． 001730． 590． 102 四水总计 1786324348934350． 40570． 00811 0． 840 长江总计1． 20000． 024000． 431． 032 总量1． 73020． 032111． 872 出湖总量0． 44740． 00803． 0152． 412 ① 据表 1 和鲍征宇等资料整理， 其中出湖泥沙 Cd 平均含量鲍征宇数据 为 5． 5 g/t， 多目标地球化学调查在城陵矶口处的Cd 平均含量为 0． 53 g/t， 取两者平均 3． 015 g/t。 据文献， 每年湖周边制糖厂和造纸厂等向湖泊排废 2． 428 106t， 养殖投放化肥约 8 103～10 103t， 共 2． 44 106t。其中 Cd 含量以农田区化肥平均值 0． 116 2 g/t 计， 共约 0． 28 t。 因此， 由泥沙、 水、 藻和鱼四种介质每年向洞庭湖输入 的 Cd 量分别为 1． 872、 0． 28、 0 和 0 t。 3洞庭湖无镉输入和带镉输入的马尔可夫链模型 3． 1转移矩阵的建立 生态地球化学系统中元素含量的变化是无后效的变化 过程。农田生态系统中保留于土壤中、 溶解于水或吸收于 植物中 Cd 的循环只与其当前状态 总量及各相态的分配 有关， 而与进入土壤之前的历史过程无关， 且来源于母质母 岩风化和施肥等人为活动加入的 Cd， 其结果都一样。因 此， 可用马尔可夫过程数学模型来研究生态地球化学系统 中的元素循环。 在洞庭湖生态系统中， 鱼以摄取水中藻类而生存， 藻类 吸取水中的养分而生长。同时通过鱼的排泄和藻类的生长 又可把部分养分归还于水体。如果把泥沙、 水、 藻和鱼四者 中的 Cd 含量视为 4 个状态， 并且考虑鱼、 藻被捕捞和从城 陵矶口的流出过程， Cd 的流失和捕捞两个状态， 系统是具 有两个吸收状态的马尔可夫链。 在 2 节分析的基础上可建立洞庭湖与渔业有关的湖泊 生态系统中 Cd 的转移矩阵。 因洞庭湖为过水型湖泊， 换水快， 仅为 4． 7 d ［7 ］， 不严格 地可忽略Cd 在湖泊水体中对流、 扩散等的所需时间， 粗略地 视其为 “黑箱” ， 但考虑对流、 扩散等作用， 认为水体完全交 换所需时间应大于4．7 d。故泥沙中的 Cd 迁出系统的概率 可用泥沙中 Cd 的年输出率的10， 即以0．013 为估计值。 由于洞庭湖水 pH 值较大， 为弱碱性， 故取泥沙中 Cd 转移至水体的概率为 0． 1。 通过比较水体和藻类 Cd 含量可知， 藻类具有富集水体中 Cd 的能力， 可取其富集率为50; 由于藻类死亡和代谢作用， 可视20的 Cd 归还于泥沙和水中; 藻类被鱼摄取的概率为 40; 藻类中的 Cd 因溢出而迁移的概率可参考泥沙中的 Cd 的年输出率取10; 则仍存在于藻类中的 Cd 的概率为10。 据洞庭湖的捕捞率为 10; 溢出概率可参考河湖间鱼 类交流率、 泥沙和藻类的溢出概率， 取 5; 鱼的排泄将 Cd 归还于水中的概率为 20。 综合可得该系统 Cd 的吸收马尔可夫链转移矩阵于 表 3所示。 3． 2洞庭湖无镉输入及带镉输入的吸收马尔可夫链模型 据吸收马尔可夫链矩阵的定理 1． 2 和 1． 3 对洞庭湖生 态系统 Cd 的转移矩阵进行计算， 结果见表 3。 按照吸收马尔可夫理论， 上述系统随着时间的发展， 所 有状态中的 Cd 都将转移到吸收状态， 转移状态将趋向零， 意味着整个系统的 Cd 将净化。从基本矩阵 N 第一行之和 可知， 在无外界输入的前提下， 经过约 74． 274 年， 泥沙中的 Cd 将基本净化; Cd 在泥沙、 水、 藻和鱼四种状态的平均保存 时间分别是74．274、 39．837、 37．856 和25．621 年。 表 3洞庭湖生态系统 Cd 的转移矩阵 Table 3The transfer matrix of Cd in the Dongting lake ecosystem 初始状态 到达状态 捕捞溢出泥沙水藻鱼 捕捞100000 溢出010000 泥沙00． 0130． 9860． 00100 水00． 0002400． 499760． 500 藻00． 100． 200． 200． 100． 40 鱼0． 100． 0500． 2000． 65 在带有输入的情况下， 马尔可夫链具有极限状态。 据带输入的吸收马尔可夫链矩阵定理1． 4 计算上述转 143 第 4 期张建新等 洞庭湖镉迁移转化的马尔可夫模型评价第 26 卷 鲍征宇， 姚志刚， 张建新， 等． 洞庭湖湿地的生态地球化学 评价与研究报告．2006． 湖南省环境保护局． 洞庭湖区生态环境现状调查报告． 2002 24 －41． 杨忠芳， 余涛， 申志军， 等． 湖南洞庭湖地区农田生态系统 区域生态地球化学评价报告．2006 44 －45． ChaoXing 移矩阵， 结果见表 4。 F 为每年向洞庭湖输入的 Cd 向量; FN 表示状态稳定 时系统中泥沙、 水、 藻和鱼 4 种状态的 Cd 总含量。 由表 4 可见， 稳定状态时， 泥沙、 水、 藻和鱼 4 种状态的 Cd 总含量将分别稳定为 146． 631、 1． 627、 0． 904 和 1． 033 t; 系统中 Cd 的总含量将稳定为 150． 196 t。 表 4洞庭湖生态系统 Cd 的矩阵计算 Table 4Calculation of matrix of Cd in the Dongting lake ecosystem 无 输 入 向量类型泥沙水藻鱼Σ 转移向量矩阵 Q 0． 9860． 00100 00． 499760． 500 0． 200． 200． 100． 40 00． 2000． 65 基本矩阵 N 73． 6590． 2810． 1560． 17874． 274 31． 2203． 9342． 1852． 49839． 837 31． 2351． 9362． 1862． 49937． 856 17． 8402． 2481． 2494． 28425． 621 带 输 入 输入向量 F1． 8720． 28 00 转移向量极限 FN 146． 6311． 6270． 9041． 033150． 196 捕捞溢出捕捞溢出 吸收状态泥沙00． 013不同吸收状态0． 0180． 982 矩阵 R水00． 00024 间的分配概率0． 2500． 750 藻00． 10NR0． 2500． 750 鱼 0． 100． 05 0． 4280． 572 由于存在两个吸收状态， 因而事物的变化有两种结果。 这也正是生态地球化学系统评价中的两个关键问题 元素 的最终去向和生态效应。从 Cd 在不同吸收状态间的分配概 率 NR 的计算结果可知 因废水排放和施肥输入鱼体的 Cd， 最 终将有42．8将由于鱼类的捕捞而带出， 进入人类食物中， 另 外的57．2则通过城陵矶口流出洞庭湖， 进入长江。 4结论和讨论 1 马尔可夫链模型是生态地球化学评价的有效方 法， 对研究污染物或有益元素的持续时间、 生态效应和变化 趋势具有优势， 可为生态环境评价提供科学依据。 2 以洞庭湖生态系统中 Cd 在泥沙、 水体、 藻类和鱼 类的现有平均含量 2． 7、 0． 00024、 3 和0． 3 mg/kg 为基础， 可建立吸收马尔可夫链的无输入模型， 得知 Cd 在系统中 的最大保存时间为泥沙中的 74． 274 年。 3 从带输入的吸收马尔可夫链理论可知， 若改变系 统与外界的关系， 允许外界不断输入 Cd 元素， 随着时间推 移， 系统最终状态仅与输入量有关， 初始状态对系统影响将 逐渐消失， 而输入量的影响将随时间而逐渐增强， 直至完全 控制系统。以洞庭湖周边区现年向湖排放废水和施肥产生 的 Cd 的输入量计算， 当达到稳定状态时， 洞庭湖的泥沙、 水、 藻和鱼 4 种介质中的 Cd 含量分别稳定在146． 631 t、 1． 627 t、 0． 904 t 和 1． 033 t; 系统 Cd 总含量将稳定于 150． 196 t; 其中洞庭湖鱼体中的 Cd 将有 42． 8 会通过渔 业捕捞进入人类食物链中。 4 利用马尔可夫链模型解决生态地球化学问题， 关 键在于确定系统的转移概率。由元素含量及其比例来估算 转移概率的理论和方法， 尚待进一步研究。 致谢 本文是笔者就进行“湖南省洞庭湖区生态地球 化学调查” 工作时所遇实际问题而作的初步探讨， 文中的 分析数据和资料多来自于该项目及其子课题的成果， 特别 要感谢中国地质调查局奚小环、 湖南省国土资源厅龙服忠、 中国地质大学鲍征宇和杨忠芳等专家的指导与帮助。 5参考文献 ［ 1］於崇文． 揭示地质现象的本质与核心 地质作用与时空 结构［ J］． 地学前缘， 2000， 7 1 1 －12． ［ 2］於崇文． 数学地质的方法与应用 地质与化探工作中的多 元分析［M］ ． 北京 冶金工业出版社， 1980 459 －520． ［ 3］於崇文． 地质系统的复杂性 地质科学的基本问题 Ⅱ ［J］ ． 地球科学， 2003， 28 1 31 －40． ［ 4］徐克学． 生物数学［M］． 北京 科学出版社， 1999． 170 －214． ［ 5］张建新， 邢旭东， 鲁江， 等． 马尔可夫模型在洞庭湖区 Cd 的 生态地球化学评价中的应用［J］． 地学前缘， 2007， 14 3 204 －212． ［ 6］李景保， 钟赛香， 杨燕， 等． 泥沙沉积与围垦对洞庭湖生态系 统服务功能的影响［ J］ ． 中国生态农业学报， 2005， 13 2 179 －182． ［ 7］湖南省政协经济科技委员会． 三峡工程与洞庭湖关系研究 ［M］． 长沙 湖南科学技术出版， 2002 85 －102． ［ 8］三峡工程生态与环境监测系统信息管理中心． 渔业资源与 环境监测重点站站点介绍 http ∥www． tgenviron． org/fishre- source/introduction． html ［ EB/OL］  ． 仪器信息网“仪器专场” 介绍 “仪器专场” 是仪器信息网 www． instruemnt． conm． cn 的一个专门频道， 是把不同厂家的同类仪器以列表形式集 中展示， 按照仪器的产地、 仪器种类等方式进行细分， 以每 台仪器的单月访问量的综合指数 3I 指数 由高到低进行 排名。用户可以按照综合排名、 仪器产地、 仪器种类、 厂商 信用积分等方式查找， 可在专场中任选两台仪器进行横向 比较， 发现产地、 价格、 技术指标、 主要特点、 仪器简介方面 的差异， 选择最合适的仪器。 到目前为止， 本网已开设气相色谱、 液相色谱、 紫外分 光光度计等各类仪器专场近 400 个。每个仪器专场包括了 该类仪器的主要生产厂家、 耗材配件、 技术资料信息和论 坛。例如“气相色谱专场” http ∥gc． instrument． com． cn ， 参展厂商已达 42 家， 参展气相色谱仪器 122 台。业内 知名企业如 安捷伦、 岛津、 上分、 天美、 福立、 东西电子等都 已加入此专场。 仪器信息网 243 第 4 期 岩矿测试 http ∥ykcs． i3t． com． cn/ 2007 年 ChaoXing