基于系统动力学模型的淡水渔业捕捞时间和喂养策略的优化分析

刘梅　原居林　倪蒙　张雷鸣　顾志敏

摘要：从系统动力学角度出发，着眼于动态系统的整体统筹，综合考虑环境、生态和经济3方面要素，利用系统模拟软件Stella 9.1.3构建我国内陆淡水渔业池塘养殖系统动力学模型，再结合Berkeley Madonna优化软件，以水产养殖户总利润最大化为目标函数对模型进行优化。获取最佳捕捞时间和最优喂养方案，并反推此时对应的最佳养殖容量，同时利用扰动法对所建模型进行参数敏感性分析，以进一步有效提高养殖容量，指导养殖生产，实现淡水渔业经济生态互利共赢局面，以期实现我国内陆淡水养殖业的可持续发展。

关键词：最佳捕捞时间;最优喂养方案;池塘养殖系统;系统动力学;Stella模型软件

中图分类号：S964

文献标志码： A

文章编号：1002-1302（2020）02-0164-06

收稿日期：2018-08-17

作者简介：刘 梅（1988—），女，安徽临泉人，博士，助理研究员，主要从事养殖水环境研究。E-mail：liumei@zju.edu.cn。

通信作者：顾志敏，研究员，主要从事水产养殖研究。Tel：（0572）2043911;E-mail：guzhimin2006@163.com。

我国是世界淡水渔业养殖大国，养殖产量已连续多年稳居世界第一[1]，淡水渔业总产量约占国内渔业产量的40%[2-3]。除了内陆湖泊、水库等天然养殖区域外，池塘养殖在我国淡水渔业中占据主导地位。但随着池塘养殖业的快速发展，养殖区的环境污染问题也越来越受到人们的关注[4-6]。由于养殖过程中片面追求高产，结果不仅导致养殖水体对周边水体富营养化影响日趋加剧，而且造成养殖对象免疫力下降、病害频发，造成巨大經济损失[7-8]。这些问题都严重影响我国水产养殖业健康可持续发展。因此，为了缓解养殖水体水质污染，保证淡水水产养殖业的可持续发展，一方面需要优化养殖结构，另一方面须加强对养殖容量的研究。

关于养殖容量的概念，不同研究领域的学者有不同的理解，目前尚无完整统一的定义。1988年，Carver等将贝类养殖的养殖容量定义为产量最大化的同时对生长率不产生负面影响的放养密度，该定义只考虑产量，未考虑生态环境因素的影响，具有明显的缺陷[9]。董双林等把养殖容量定义为单位水体内在保护环境、节约资源和保证应有效益的各个方面都符合可持续发展要求的最大养殖量[10]。杨红生等把浅海贝类养殖业的经济效益与社会效益、生态效益结合起来，将养殖容量定义为对养殖海区的环境不会造成不利影响，又能保证养殖业可持续发展并有最大效益的最适产量[11]。刘剑昭等将养殖容量定义为特定的水域、单位水体养殖对象在不危害环境、保持生态系统相对稳定、保证经济效益最大，并且符合可持续发展要求条件下的最大产量[12]。可见，养殖容量的含义在不断地被充实、丰富，越来越趋于完善。

鱼塘生态系统是一种复杂的人工生态系统，受人类控制和影响非常大。若想实现单位水体养殖对象在不危害环境、保持生态系统相对稳定、保证经济效益最大的同时，符合可持续发展要求条件下的最大产量，只有控制鱼塘参数达到最佳时才能实现，如捕捞时间及饵料的输入等。可以看出，鱼塘生态系统涉及环境、生态和经济社会等多学科多领域，是一个复杂的系统，具有多层次、高阶非线性、动态性、自反馈性等特征[8，13]，一般数学方法难以对其进行量化描述和分析。鉴于系统动力学建模软件（structure thinking experimental learning laboratory with animation，Stella）具有处理动态性、非线性和高阶次复杂问题的功能，运用该软件建立综合考虑环境、生态和经济之间相互作用的养殖生态系统动力学模型，然后以水产养殖总利润最大化为目标函数，利用Berkeley Madonna（http：//www.berkeleymadonna.com）优化Stella所建池塘动态系统模型，以确定最佳捕捞时间及最优喂养策略，并反推此时对应的最佳养殖容量。然后推广到实际生产中，实现经济生态互利共赢局面，以期实现我国内陆淡水养殖业的可持续发展。

1 建模工具介绍

系统动力学（system dynamics，简称SD）于1956年由美国麻省理工学院福瑞斯特教授开发创建。本研究采用Isee Systems公司的一款可视化系统动力学模拟软件Stella 9.0.1作为模型构建的平台，即图形导向的系统动力学模型发展软件，它与Vinsim以及Matlab的Simulink软件包相似，可以提供一个实现动态系统建模、仿真的集成环境[14-15]，动态展现和表达复杂系统和概念实际是如何进行运作的。Stella在研究复杂系统的行为和在处理高度非线性、高阶次、多变量、多重反馈问题方面具有优势，随着系统动力学理论和方法的不断发展与完善以及计算机技术的改进，目前已被广泛应用于规划学、生物学、生态学和环境科学等领域[16-18]。但是Stella软件并不能对所建动态系统模型的参数进行优化，以使系统达到所期望的最佳状态。

而Berkeley Madonna是一款优秀的动态系统建模分析优化软件，且与Stella兼容，可以运行优化Stella所建模拟系统。由著名的美国加利福尼亚大学伯克利分校研究人员开发，具体介绍详见网址http：//www.berkeleymadonna.com，其算法优秀，可以在几秒钟之内完成数百万个方程式的计算，可做参数校准（曲线拟合）和优化，有更多解决常微分方程问题的数值方法。因此，本研究利用该软件对所建鱼塘生态系统模型参数进行优化，以方便相关管理部门以及养殖户选取最佳养殖策略，获取经济生态共赢的局面。

2 鱼塘生态系统模型的构建

构建系统模型的指导思想是用尽可能简单的模型来描述系统结构，而构造高层结构图的指导思想则是简单、概括，即抓住系统的关键性流程，建立系统内的高层关系。由于鱼塘是一个汇聚了生物、社会、生态、经济等各个领域交互作用的动态系统，其包含的因素较多，为便于分析，假设这是一个完全人工饲养的池塘，鱼类生存的一个重要条件是在池塘中有足够的溶解氧。池塘中溶解氧的消耗主要分为2个部分，一部分是鱼的呼吸消耗，另一部分则是被用来分解死有机物质即碎屑[19]。碎屑是池塘生态系统中一个非常重要的因素，随着其浓度的增加，就会出现缺氧甚至厌氧条件[20]。当氧气浓度低于一定阈值时，鱼就会死光，因此本研究假设耗氧量随着碎屑浓度的增加而增加，简单地认为鱼的死亡是由高碎屑浓度导致的。基于以上假设分析，可以将鱼塘生态系统简单地划分为社会子系统、经济子系统、环境子系统3个子系统，每个子系统都有自己的结构特点和独特功能，其中1个子系统的输出是其他子系统的输入，子系统之间彼此联系，构成了模型的反馈结构，各子系统的构成及相互关系分析具体见下文。

2.1 系统结构图分析

建立系统动力学模型需要4个基本构造块：栈（stock）、流（flow）、转换器（converter）和连接器（connector）。栈用方框表示，代表事物的积累，即用来代表存储进入和流出栈的物质;流用带阀门的箭头表示，用来描述系统中的活动，连接到栈上的流会引起栈中存量的增加或减少，流上的箭头表示物质或非物质正向流动方向;转换器用圆圈表示，用以存储常量，并把它传输到其他模块中为其他变量所用;连接器用1条带有箭头的线段或虚线表示，作用是连接模型中的元素，与流不同的是，连接器传送的是信息。

由于研究区域生态问题的复杂性，因此，在建立系统流程图时，选择影响库区生态安全的核心因子，最终确定参数，建立系统动力学模型的系统流程，具体如图1所示。

2.2 系统分析及参数的确立

养殖塘生态系统主要有3个子系统，分别是池塘养殖子系统、社会子系统以及经济子系统，每个子系统的具体控制参数如下所列。

（1）池塘养殖子系统：池塘渔业资源生物量（Fish）与放渔量（Number）、生长率（Growth）死亡率（Mortality）、鱼质量（Weight）有关。其中，死亡率导致渔业资源生物量减少，而生长率及放渔量则构成了渔业资源量的补充量。模型中设置渔业资源生物量（Fish）和放渔量（Number）为状态变量，生长率（Growth）、死亡率（Mortality）为相应的流率变量。

（2）社会子系统：模型的生态部分将设置饵料生物资源（Feed）和碎屑量（Detritus）为状态变量，饵料喂养量（Feeding）、饵料损失量（Loss）和碎屑分解量（Decomp）為流率变量。

（3）经济子系统：模型该部分中，设置总利润Total_profit为状态变量。流率变量成本（Cost）、饵料价格（Feed_price）、收入（Revenue）、出售时间（Time_of_scale）为流量变量。其中收入（Revenue）表示流入，它会带来总利润的积累和增加，其余4个变量表示流出，会导致利润的减少。

2.3 模型方程分析

根据系统动力学流图中的参数关系，本研究构造了32个系统动力学方程，主要模型方程式如下：

d/dt（Fish）=Growth-Mortality;

d/dt（Feed）=-Growth+Feeding-Loss;

d/dt（Detritus）=+Accum-Decomp;

d/dt（Total_profit t）=+Profit;

d/dt（Numbers）=-J3;

Growth=if Feed>0 then C_growth*Feed*Fish_W/（Feed+C_Hs） else 0;

Feeding=if（Time0） then C_feed else 0;

Loss=C_loss*Feed+Growth*0;

Accum=Loss+Mortality;

Decomp=C_decomp*Detritus;

Profit=Revenue-Cost;

Mortality=if（TIME>Time_of_sale+1）then Fish/DT else（C_mort+Detritus^4）/（C_mort_d^4+Detritus^4））*Fish;

C_feed=A*（TIME+B）^2+C;

Fish_Price=10+2*Weight;

Revenue=if Time>Time_of_sale AND Time>Time_of_sale+2 then Fish_Price*Fish else 0;

Cost=Feed_price*Feeding。

3 结果与分析

3.1 养殖塘系统模型的建立

首先选择一些喂养策略，使模型产生一些合理的结果。当A=0.001、B=-10、C=0.1或0.2时，得到图2、图3中的喂养方案和鱼及碎屑量的动态变化。由图2、图3可以看出，鱼类种群逐渐生长，直到某一点即第74天时，由饵料产生的碎屑量超过一定阈值时，导致鱼类的大量死亡。鱼类种群的崩溃，从而进一步增加了碎屑量。当碎屑分解率增加，碎屑量达到种群崩溃时的阈值也相应提高，则种群崩溃的时间也相应延长。说明系统模型模拟结果合理可靠，可以用来模拟分析养殖塘水生态系统。

3.2 养殖塘系统模型的优化

鱼塘生态系统动力学模型建好后，再以总利润最大为目标函数，用Berkeley Madonna软件优化模型中的喂养策略因子A、B、C以及捕捞时间 Time_of_sale，优化界面见图4。

用这个模型进行优化可以得到最大总收益为550万元，如图5所示。此时A=0.000 997 764，B=-13.517 6，C=0.171 658，Time_of_sale=75.401 6（图6）。捕捞的时间定在碎屑即将达到阈值并使鱼群处于灭绝的危险时。个别个体质量的增加与越来越多的鱼死亡相抵消，鱼的生物量迅速下降。因此选在第76天或者提前几天将养殖塘中的鱼进行捕捞后出售， 则能保证水产养殖户的最大经济效益。

3.3 模型的敏感性分析

参数灵敏度分析的目的是确定对模型预测结果影响较大的参数，找到这些参数后，预示着模拟预测时必须认真选取这些参数，如取值不当，会导致模拟失败[21-22]。本研究利用综合偏导数分别对生长率（C_growth）、死亡率（C_mortality）、饵料残留率（C_loss）和碎屑分解率（C_decomp）参数进行敏感性分析，结果显示其敏感度系数分别为0.46、-0.32、-0.35、0.58，其中绝对值越大，参数越灵敏。因此，碎屑分解率最为灵敏，其他依次为生长率、饵料残留率和死亡率。检验的方式则是保持其他参数不变，只改变被检验参数，以最灵敏的碎屑分解率为例，以幅度为±20%、±50%、±100%进行扰动分析，观察总利润Total_profit的变化。由图7可以看出，随着C_decomp的逐渐增加，即池塘新陈代谢加快，池塘碎屑残留量明显减少，即水环境得到改善，因而鱼的生物量增加，从而提高了养殖效益。

假如当其他参数未变，碎屑分解率从C_decomp=0.1增加到C_decomp=0.2时，即鱼塘碎屑分解加快，则优化结果就会发生很大变化。在碎屑分解率较高的情况下，碎屑积累得很慢，很难达到使鱼类死亡的临界值。因此优化只能试图使鱼类的体质量尽可能增加，而消耗最少的饲料，具体追踪优化过程如图8所示。而此时优化结果则与之前发生了很大变化，可以看出得到了更高的利润868（图9），而最佳捕捞时间则从放养的第75天变为现在的第83天。因此，喂食情景对模型中的碎屑分解率（C_decomp 参数）非常敏感。在实际生产中可以通过池塘底部微孔曝气方式增加溶解氧和添加微生态制剂等方式[23-24]来提高碎屑分解率，从而提高单位面积产量和经济效益。

4 结论

本研究演示了系统动力学方法在鱼塘生态养殖系统中的应用，可以看出，系统动力学方法直观性强，在宏观上将鱼塘养殖过程中的生态、经济、社会3个方面整合在一起，形成一个大的系统，融会贯通，且操作简便，使用者即便没有相当好的数学基础，也可以根据不同的使用目的改变模型，根据实际情况改变相应参数设置。另外，本研究是将系统动力学模型应用于淡水渔业水产养殖领域的一项探索性研究，模型中大多采用简化处理，部分参数采用估计值，而实际养殖系统中远远比本研究中建立的模型要复杂得多，尤其是参数的变动会在比较大的程度上影响结果，所以开展各项基础调查研究以及如何完善模型以提高其仿真的准确性将是今后工作的重点和难点。

参考文献：

[1]袁晓初，赵文武. 中国渔业统计年鉴（2015）[M]. 北京：中国农业出版社，2015.

[2]唐启升，丁晓明，刘世禄，等. 我国水产养殖业绿色、可持续发展保障措施与政策建议[J]. 中国渔业经济，2014，32（2）：5-11.

[3]宋 超，孟顺龙，范立民，等. 中国淡水池塘养殖面临的环境问题及对策[J]. 中国农学通报，2012，28（26）：89-92.

[4]鲍旭腾，徐 皓，张建华，等. 水产养殖面源污染控制的最佳管理实践[J]. 南方水产科学，2012，8（3）：79-86.

[5]贾敬德. 淡水渔业环境现状及保护对策[J]. 淡水渔业，2004，34（5）：59-61.

[6]贾敬德. 保护环境是渔业持续发展的当务之急[J]. 淡水渔业，1997，27（4）：26-29.

[7]Kenneth A，Bert B，Robert C，et al. Economic growth，carrying capacity，and the environment[J]. Science，1995，268（5210）：520-521.

[8]Jiang W M，Gibbs M T. Predicting the carrying capacity of bivalve shellfish culture using a steady，linear food web model[J]. Aquaculture，2005，244（1/2/3/4）：171-185.

[9]Carver C E A，Mallet A L. Assessing the carrying capacity of a coastal inlet in terms of mussel culture[J]. Aquaculture，1988，88：39-53.

[10]董雙林，李德尚，潘克厚. 论海水养殖的养殖容量[J]. 青岛海洋大学学报，1998，28（2）：253-258.

[11]杨红生，张福绥. 浅海筏式养殖系统贝类养殖容量研究进展[J]. 水产学报，1999，23（1）：84-90.

[12]刘剑昭，李德尚，董双林. 关于水产养殖容量的研究[J]. 海洋科学，2000，24（9）：33-35.

[13]王其藩. 系统动力学[M]. 北京：清华大学出版社，1988：25-27.

[14]杰拉尔德·温伯格. 系统化思维导论[M]. 张 佐，万起光，董 菁，译. 北京：清华大学出版社，2003：50-52.

[15]Ouyang Y，Zhang J E，Lin D，et al. A STELLA model for the estimation of atrazine runoff，leaching，adsorption，and degradation from an agricultural land[J]. Journal of Soils and Sediments，2010，10（2）：263-271.

[16]狄乾斌，徐东升，周乐萍. 基于STELLA软件的海洋经济可持续发展系统动力学模型研究[J]. 海洋开发与管理，2012，29（3）：90-94.

[17]于 涛，钟 非，贺 锋，等. 基于STELLA的循环水养殖系统池塘总氨氮动态模拟[J]. 农业机械学报，2013，44（7）：199-203.

[18]冯东溥，魏晓妹，降亚楠，等. 基于STELLA和气候变化情景的灌区农业供需水量模拟[J]. 农业工程学报，2015，31（6）：122-128.

[19]黄岁樑，臧常娟，杜胜蓝，等. pH、溶解氧、叶绿素a之间相关性研究Ⅰ：养殖水体[J]. 环境工程学报，2011，5（6）：1201-1208.

[20]孙 耀，陈聚法. 中国对虾养殖水体中溶解氧的动态收支平衡模式[J]. 水产学报，1999，23（4）：424-428.

[21]邱康俊，温华洋，王 根，等. 基于高斯随机扰动的潜在蒸散敏感性分析[J]. 中国农学通报，2015，31（26）：156-160.

[22]何 亮，侯英雨，赵 刚，等. 基于全局敏感性分析和贝叶斯方法的WOFOST作物模型参数优化[J]. 农业工程学报，2016，32（2）：169-179.

[23]顾海涛，刘兴国，何雅萍，等. 微孔曝气式增氧机的性能及应用效果[J]. 渔业现代化，2017，44（3）：25-28.

[24]李 敏，段登选，许国晶，等. 大薸-微生态制剂协同净化养殖池塘富营养化水体的效果[J]. 生态与农村环境学报，2015，31（1）：94-99.王宏玉，裴雪莹，王 涛，等. 杂交黄颡鱼黄优1号体质量和内脏指标的相关分析[J]. 江苏农业科学，2020，48（2）：170-176.