中国生物工程杂志  2016, Vol. 36 Issue (1): 68-75

文章信息

张许, 丁健, 高鹏, 高敏杰, 贾禄强, 涂庭勇, 史仲平
ZHANG Xu, DING Jian, GAO Peng, GAO Min-jie, JIA Lu-qiang, TU Ting-yong, SHI Zhong-ping
基于差分进化算法的酿酒酵母分批补料培养在线自适应控制
Fed-batch Culture of Saccharomyces cerevisiae with Adaptive Control Based on Differential Evolution Algorithm
中国生物工程杂志, 2016, 36(1): 68-75
China Biotechnology, 2016, 36(1): 68-75
http://dx.doi.org/10.13523/j.cb.20160110

文章历史

收稿日期: 2015-09-09
修回日期: 2015-09-23
基于差分进化算法的酿酒酵母分批补料培养在线自适应控制
张许, 丁健, 高鹏, 高敏杰, 贾禄强, 涂庭勇, 史仲平     
江南大学生物工程学院 工业生物技术教育部重点实验室 无锡 214122
摘要: 酿酒酵母分批补料培养中,葡萄糖添加过量会导致乙醇大量积累,破坏细胞结构及功能,降低葡萄糖利用效率;葡萄糖添加不足会限制细胞生长。为解决这一矛盾,提出了一种基于差分进化算法的在线自适应控制策略,并利用计算机仿真方法对该策略、传统的间歇流加、分段恒速流加及PID控制策略的控制性能进行了研究和比较。结果表明,在该控制策略下,发酵液中的乙醇浓度能够被稳定地控制在1g/L的低水平,而细胞浓度却达到34.45g/L的高水平,比采用间歇流加、分段恒速流加及PID控制策略的批次分别提高了243%、18%和29%。由此可知,该自适应控制策略能够将葡萄糖流加速率控制在适宜水平,避免乙醇过量积累的同时保证细胞的快速增殖。
关键词: 酿酒酵母     乙醇积累     差分进化算法     系统辨识     自适应控制    
Fed-batch Culture of Saccharomyces cerevisiae with Adaptive Control Based on Differential Evolution Algorithm
ZHANG Xu, DING Jian, GAO Peng, GAO Min-jie, JIA Lu-qiang, TU Ting-yong, SHI Zhong-ping     
The Key Laboratory of Industrial Biotechnology, Ministry of Education, School of Biotechnology, Jiangnan University, Wuxi 214122, China
Abstract: In fed-batch culture of Saccharomyces cerevisiae, excessive glucose addition leads to much ethanol accumulation, destroying structure and function of cell and decreasing glucose utilization efficiency, while insufficient glucose addition limits cell growth. To solve this problem, a self-adaptive control strategy based on differential evolution algorithm was proposed. In addition, performances of the proposed strategy, traditional strategy, were tested and compared using computer simulation. As a result, under the proposed control strategy, ethanol concentration could be maintained at the low level of 1.0g/L, while the biomass concentration could reach to the high level of 34.45g/L, which was 243%, 18% and 29% higher than those under intermittent feed, stepped constant feed and PID control strategy, respectively. In conclusion, the proposed self-adaptive control strategy was capable of controlling glucose feed rate at proper level, and thus ensured the rapid growth of yeast when repressing ethanol accumulation.
Key words: Saccharomyces cerevisiae     Ethanol accumulation     Differential evolution algorithm     System identification     Adaptive control    

酿酒酵母(Saccharomyces cerevisiae)是与人们生活息息相关的一种真核微生物,可用于生产多种产物,具有遗传背景清楚、不产内毒素、可对外源蛋白进行翻译后修饰、发酵周期短、操作方便等优点[1]。酿酒酵母培养过程具有明显的Crabtree效应,即当碳源浓度过高时,大量碳源流向乙醇合成途径,导致乙醇大量积累。积累乙醇的同时,伴随着大量CO2的排放,大大降低碳源的利用效率[2, 3],且高浓度的乙醇也会对细胞结构及功能造成不可逆转的损伤,抑制产物合成[4]。若减少碳源添加,可有效缓解乙醇积累,但过低的碳源浓度同样会影响酵母细胞的生长和产物的合成。利用基因工程手段敲除乙醇合成途径中的关键酶能够有效降低乙醇积累[5, 6],但此类工程菌对碳源浓度过于敏感,在发酵生产中的应用十分有限[7]。利用过程控制技术,合理调节碳源添加,在保证细胞正常生长的前提下抑制乙醇积累,是解决这一难题的有效途径之一[8, 9]。发酵过程具有高度的时变性、滞后性和非线性,利用控制参数恒定的常规反馈控制策略(如比例-积分-微分控制,即PID控制)对发酵过程实施控制,往往难以取得令人满意的效果[10, 11]。于是,能够自动实时调节控制参数的在线自适应控制策略被越来越多地用于发酵过程控制[12, 13]。差分进化算法(differential evolution,DE)是由Stom和Price[14]建立的一种全局优化算法,其原理简单、受控参数少、收敛速率快,已受到各领域研究人员的广泛关注[15, 16, 17]。在此基础上,本研究提出了基于DE算法的在线自适应控制策略,以期解决上述酿酒酵母培养过程中存在的难题,并进一步利用酿酒酵母培养过程模型[18]对不同葡萄糖流加策略的控制性能进行了计算仿真研究,比较了各葡萄糖流加策略下的酿酒酵母培养性能。

1 材料与方法 1.1 培养基及培养条件

平板培养基(g/L):葡萄糖20,酵母粉10,蛋白胨20,琼脂粉20。种子培养基(g/L):葡萄糖20,酵母粉8.5,(NH4)2SO4 1.3,MgSO4·7H2O 0.1,CaCl2·2H2O 0.06。发酵初始培养基(g/L):葡萄糖5,酵母粉12,(NH4)2SO4 2.5,MgSO4·7H2O 1.5,CaCl2·2H2O 0.1。流加培养基(g/L):葡萄糖500。

本研究酿酒酵母的培养在7L发酵罐内进行,初始体积为3L。pH和温度分别控制在30℃和6.0。通风比为1.3vvm,通过调节搅拌转速将溶解氧浓度维持在20%以上。

1.2 酿酒酵母培养过程模型

本文选用Takamatsu等[18]报道的酿酒酵母细胞生长与乙醇积累模型,并对其进行适当修改和参数的重新优化。修改和优化后的模型可由如下微分方程组表示。

S > 0.335 1时,a1 = 0,a2 = 1;当S ≤ 0.335 1时,a1 = 1,a2 = 0。X、S、Ce、V、F、Ki、S0、t、μ、νe、m、πe、Yx/e,Yx/s,Ye/s分别表示菌体浓度、葡萄糖浓度、乙醇浓度、发酵液体积、葡萄糖流加速率、菌体抑制常数、流加培养基的葡萄糖浓度、发酵时间、细胞比生长速率、乙醇比消耗速率、维持系数、乙醇比生成速率、细胞对乙醇的得率、细胞对葡萄糖的得率和乙醇对葡萄糖的得率。其中m、Yx/s、Yx/eYe/s都为定值;πeνe为非负数,分别由式(5)和式(6)计算得到;μ由式(7)计算得到,μmaxKs分别表示最大细胞比生长速率和底物饱和常数。模型中各参数取值如表 1所示。

表 1 模型参数及其取值 Table 1 Model parameters and their values
参数Yx/eYx/sYe/sS0 (g/L)Ki (g/L)Ks (g/L)m (/h)μmax (/h)
取值0.777 40.927 50.189 650080.034 00.098 30.175 8
1.3 间歇流加与分段恒速流加策略

间歇流加和分段恒速流加是两种十分常用的底物流加策略。本研究首先采用这两种流加策略控制葡萄糖流加,操作细节简述如下:①间歇流加,每次发酵液中乙醇浓度低于0.1g/L后,一次性添加0.2L葡萄糖流加培养基;②分段恒速流加,即在培养阶段的0 ~ 20h、20 ~ 35h和35 ~ 50h分别以0.02L/h、0.04L/h和0.06L/h的速率流加葡萄糖流加培养基。

1.4 PID控制

PID控制是一种发展最成熟、应用范围最广泛的反馈控制策略。PID控制包括比例、积分和微分3个环节,可由式(8)表示。

离散化后可以改写为式(9):

递推可得

两式相减得

将该策略用于本研究,根据乙醇浓度反馈调节葡萄糖流加速率F,则F可由式(12)计算得到,图 1为PID控制框图。

式中,k为当前采样点;e为乙醇浓度偏差,即e(k)= CesetCe(k);kp、ki、kd为控制参数,由操作人员设定后,在控制过程中其取值不发生变化。

图 1 PID控制示意图 Fig. 1 Block diagram of PID control
1.5 基于DE算法的在线自适应控制

以PID控制为基础,本研究提出了基于DE优化算法的在线自适应控制策略,以下简写为DE-PID策略,操作步骤如下所述。

1.5.1 初始化

设定目标乙醇浓度Ceset和初始PID控制参数:kp、ki、kd。由于反映系统信息的数据不足,自适应控制策略无法实施,因此在采样周期小于4(k<4)时,按照上述标准PID控制策略确定该采样点的葡萄糖流加速率F(k)。

1.5.2 基于ARMA的系统辨识

k≥4时,即认为反映发酵系统特性的数据已经充足,可以利用自回归移动平均模型(autoregressive moving average model,ARMA)对酵母培养系统进行辨识。ARMA能够描述数据之间的有序性和相关性,通过数据内部之间的相互联系来识别系统的变化规律,常被用于时间序列分析[19, 20]。在酿酒酵母的培养过程中,乙醇浓度变化是一个连续过程,相邻采样时刻的数据之间具有极强的相关性,且乙醇浓度受葡萄糖流加速率影响。因此,本研究利用该模型描述系统变化规律并预测乙醇浓度。ARMA模型的标准形式如式(13)所示。其中乙醇浓度Ce为模型的输出;流加速率F为模型的输入;k为当前采样点;a1a2,…,anb1b2,…,bm为模型参数;e为模型的误差值。取n = m = 1,ARMA模型的标准形式可以简化为式(14)。

利用最小二乘法原理估计模型参数a1b1,该过程可表示为

采集乙醇浓度Ce(k),定义如下矩阵X和Y:

则模型参数θ = [a1,b1]T可以由以下公式计算得到

随着培养过程的进行,系统的动力学特性会不断发生变化,利用所有采样点的数据估计ARMA模型的参数,显然是不合理的。距离当前采样点越近的数据越能反映当前的系统特性,在参数估计过程中应当根据数据组的重要程度,对其乘以不同的权重系数:距离当前采样点越近的数据,其权重系数越大;距离当前采样点越远的数据,其权重系数越小。定义忘却因子λ(0<λ<1,本文取值为0.95),构建如式(18)所示的权重系数矩阵W,式(17)可以被修正为式(19)。

利用这一模型,可以在已知Ce(k)并且给定F(k)的情况下,预测k+1采样点的乙醇浓度Ce(k+1)。

1.5.3 目标函数的确定

k-1采样点的PID控制参数分别为kp(k-1),ki(k-1),kd(k-1),在k采样点按照以下公式对其更新:

采集Ce(k),并计算误差值,即e(k) = CesetCe(k)。又由于kp(k-1)、ki(k-1)、kd(k-1)、e(k-1)及e(k)已知,k采样点的葡萄糖流加速率可以看成是未知数p1p2p3的函数,表示如下:

X′=[-Ce(k) F(p1,p2,p3)],则k+1采样点乙醇浓度的预测值可由下式计算得到

再令X″=[-Ĉe(k+1) F(p1,p2,p3)],则k+2采样点乙醇浓度的预测值可由下式计算得到

因此,求解如式(24)所示的最优化问题,可以得到最优的p1p2p3,保证k+1采样点和k+2采样点的乙醇浓度始终维持在接近Ceset的水平。

1.5.4 利用DE算法在线优化PID控制参数

DE算法是一种智能型优化算法,可以用来求解上述最优化问题。求解过程如下所述。

(1)生成初始种群。设定种群大小NP=100,随机生成初始种群,种群中的个体可以由三维向量Pit=[pi,1t pi,2t pi,3t]表示。其中i∈[1,2,3,…,NP],表示该个体在种群中的索引号;t表示传代次数,在初始种群中t=0,其中的三个分量分别代表p1p2p3。计算每个个体所对应的适应度函数值G(Pit),并将适应度函数最小的个体定义为最优解Pbest

(2)变异、交叉及选择。对种群中的每个个体依次进行变异、交叉和选择操作。首先将P10作为操作的目标矢量,在种群中选取除目标矢量外的三个不同个体Pr10Pr20以及Pr30,按照式(25)生成变异矢量V = [v1,v2,v3],其中FDE为缩放因子,本研究中取值为0.8。

定义U = [u1,u2,u3]为试验矢量,利用式(26)计算U中各分量的值。其中j = 1,2,3;CR为交叉概率,本研究中取值为0.8;Rj为0 ~ 1之间均匀分布的随机数。RN为集合[1, 2, 3]中随机选取的一个索引,以保证试验矢量至少有一个分量由变异矢量贡献,否则试验矢量有可能与目标矢量相同而不能生成新个体。

选择操作则是将试验矢量U与目标矢量P10进行竞争,若U优于P10,则用U代替P10,进入下一代[式(27)]。

按照变异-交叉-选择的运算过程对初始种群中的每个个体进行传代,得到的新种群即为第一代种群,其中的个体可表示为Pi1,其中i∈[1,2,3,…,NP]。计算每个个体所对应的适应度函数G(Pi1),若min{G(Pi1)} < G(Pbest),则由min{G(Pi1)}对应的个体代替原有Pbest,否则Pbest维持不变。

(3)迭代。按照上述步骤进行多次迭代,若在连续8次迭代中,Pbest没有被更新,认为适应度函数已经收敛,此时的Pbest= [Pbest,1 Pbest,2 Pbest,3]则是式(24)所述问题的最优解。将p1=Pbest,1p2=Pbest,2p3=Pbest,3代入式(20)和式(21)更新PID控制参数并计算k采样点最优的葡萄糖流加速率F(k)。综上所述,基于DE算法的在线自适应控制策略可由图 2表示。

图 2 DE-PID控制策略示意图 Fig. 2 Block diagram of DE-PID control
2 结果与讨论 2.1 间歇流加与分段恒速流加策略下酿酒酵母培养性能

在不同葡萄糖流加策略下,利用四阶龙格-库塔法求解式(1) ~式(4)所述的微分方程组,对不同葡萄糖流加策略下的酿酒酵母培养性能进行计算机仿真研究。求解过程利用Matlab2010b编程实现。批次1采用发酵试验中最常用的间歇流加策略,当葡萄糖耗尽且副产物乙醇浓度低于0.1g/L时一次性添加0.2L流加培养基;批次2采用分段恒速流加策略,0~20h、20~35h和35~50h分别以0.02L/h、0.04L/h和0.06L/h添加流加培养基。批次1和批次2的葡萄糖流加速率如图 3(a)所示,在此流加速率条件下,以X=0.375g/L、S=5g/L、Ce=0.37g/L和V=3L为初始条件,求解酵母培养过程模型,结果如图 3(b)、(c)和图 4所示。

图 3 批次1及批次2葡萄糖流加操作曲线及培养性能 Fig. 3 Glucose feed rates and cultivation performances of run#1 and run#2
图 4 各批次细胞浓度变化 Fig. 4 Biomass concentration variation of each run

在批次1中,每次添加葡萄糖后,发酵液内的葡萄糖浓度偏高[图 3(b)],引发Crabtree效应,大量葡萄糖被用于合成乙醇[图 3(c)],同时排放大量CO2[3],只有少部分葡萄糖用于细胞生长。葡萄糖被迅速耗尽后,细胞长时间处于葡萄糖匮乏且乙醇过剩的环境中,以乙醇为替代碳源维持其生长。由于酵母细胞对乙醇的得率远远低于其对葡萄糖的得率(Yx/e<Yx/s表 1),导致细胞生长缓慢,培养结束时细胞浓度只达到10.42g dry cell weight/L (g DCW/L)(图 4)。批次2中,虽然葡萄糖浓度始终被控制在较低的水平,乙醇仍然持续积累,直至45.5h才开始降低,最终的细胞浓度达到29.22g DCW/L(图 4)。理想的酿酒酵母培养过程具有如下特征:乙醇积累浓度低,避免高浓度乙醇对细胞的损伤;细胞快速生长且在培养结束时达到较高浓度。按照此标准,批次1和批次2的培养性能较差。批次1的最高乙醇积累浓度达到5g/L以上,但最终的细胞浓度过低;批次2的细胞密度虽然较高,但是其最高乙醇积累浓度却达到16.91g/L。

2.2 PID控制策略下酿酒酵母培养性能

批次3以乙醇浓度作为反馈指标,采用PID控制策略调节葡萄糖流加速率,经过遗传算法离线优化,控制参数kp、kikd分别被设定为0.001 7、0.059 9和0.022 4。为避免乙醇对细胞的毒害,同时保证葡萄糖的充足添加,将乙醇浓度设定于1g/L的较低水平,初始条件与批次1和批次2相同。培养过程中的葡萄糖流加速率、发酵液中葡萄糖和乙醇的浓度如图 5(a)、(c)、(e)所示。葡萄糖流加速率随培养过程的进行而不断提高,以适应细胞持续增殖的需求,但始终存在较大波动,导致葡萄糖浓度在0g/L与0.34g/L之间波动。受此影响,乙醇浓度波动随着培养的进行逐渐加剧,难以稳定地维持在设定值附近。虽然在流加葡萄糖期间乙醇浓度始终维持在2g/L以下,但由于培养状态在乙醇积累和葡萄糖匮乏之间频繁切换,而这两种状态都不利于细胞生长,最终的细胞浓度仅达到26.70 g DCW/L。批次3的最终细胞浓度虽然略低于批次2,但其乙醇浓度远远低于批次2。综合考虑多方面的因素,PID控制策略明显优于间歇流加和分段恒速流加策略。

图 5 批次3及批次4葡萄糖流加速率及培养性能 Fig. 5 Glucose feed rates and cultivation performances of run#3 and run#4
2.3 DE-PID控制策略下酿酒酵母培养性能

批次3中乙醇浓度控制性能不断恶化的原因是恒定的控制参数kp、kikd无法满足发酵系统不断变化的动力学特性。为解决这一问题,本研究提出DE-PID控制策略,利用DE算法在线优化控制参数。将该策略用于批次4,同样将乙醇浓度设定在1g/L,培养的初始条件也与批次3相同。控制参数kp、kikd的初始值分别设定为0.001 7、0.059 9和0.022 4。随着培养过程的进行,kp、kikd被不断调整(图 6),以适应培养系统不断变化的控制需求。葡萄糖流加速率、葡萄糖浓度和乙醇浓度分别如图 5(b)、(d)、(f)所示。葡萄糖流加速率随着培养过程的进行不断提高,同时在小范围内波动。这也使得葡萄糖浓度波动比批次3明显减小,乙醇浓度被稳定地控制在(1±0.2)g/L的范围内,葡萄糖匮乏(浓度为0g/L)的情况也极少出现。批次4最终的细胞浓度达到34.45g/L(图 4),高于其他各批次。由此可见,与其他各葡萄糖流加策略相比,DE-PID策略具备更好的性能,能够在酿酒酵母培养过程中同时实现细胞快速增殖和避免乙醇过量积累的培养目标。

图 6 批次4控制参数变化 Fig. 6 Control parametersvariation of run#4

综上所述,本研究提出的DE-PID控制策略具有如下优点:①通过ARMA模型在线辨识发酵系统的实时状态,使DE-PID控制策略能够很好地应对生物反应过程的滞后性、非线性和时变性特征;②DE-PID控制策略不依赖于被控对象的过程模型,可以被推广到许多其他的分批补料发酵过程中,还可以用于控制其他在线参数,如溶解氧浓度、呼吸商、CO2排放速率等;③DE算法运算速率比其他常规进化算法(如遗传算法)更快[14],使DE-PID控制策略更适合用于实际发酵过程的在线控制;④对于某些使用ARMA模型无法描述的复杂发酵过程,可以采用其他非线性智能模型代替ARMA模型进行系统辨识,这使得DE-PID控制策略能够应对更加复杂的培养过程。

3 结 论

本研究将DE算法与传统PID控制策略相结合,构建了用于酿酒酵母分批补料培养过程的在线自适应控制策略,并对其控制性能进行了计算机仿真研究。结果表明,该策略能够将乙醇浓度稳定地控制在设定值(1g/L)附近,同时确保细胞的快速增殖。在培养结束时,酵母细胞的最终浓度达到34.45g/L,比采用间歇流加、分段恒速流加和传统PID控制的批次分别提高了231%、18%和29%。

参考文献
[1] 邱玲. 乙肝表面抗原在酿酒酵母中的异源表达及优化. 无锡:江南大学,生物工程学院, 2014. Qiu L. Heterologous Experssion and Optimization of HBsAg in Saccharomyces cerevisiae . Wuxi : Jiangnan University, School of Biotechnology,2014.
[2] 王萍, 王泽建, 张嗣良. 生理代谢参数RQ在指导发酵过程优化中的应用. 中国生物工程杂志, 2013, 33(2) : 88-95. Wang P, Wang Z J, Zhang S L. Fermentation optimization directed by physiological parameter respiratory quotient. China Biotechnol, 2013, 33(2) : 88-95.
[3] Zigova J. Effect of RQ and pre-seed conditions on biomass and galactosyl transferase production during fed-batch culture of S. cerevisiae BT150. J Biotechnol, 2000, 80(1) : 55-62.
[4] Gibson B R, Lawrence S J, Leclaire J P R, et al. Yeast responses to stresses associated with industrial brewery handling. FEMS Microbiol Rev, 2007, 31(5) : 535-569.
[5] Tokuhiro K, Ishida N, Nagamori E, et al. Double mutation of the PDC1 and ADH1 genes improves lactate production in the yeast Saccharomyces cerevisiae expressing the bovine lactate dehydrogenase gene. Appl Microbiol Biotechnol, 2009, 82(5) : 883-890.
[6] Ida Y, Hirasawa T, Furusawa C, et al. Utilization of Saccharomyces cerevisiae recombinant strain incapable of both ethanol and glycerol biosynthesis for anaerobic bioproduction. Appl Microbiol Biotechnol, 2013, 97(11) : 4811-4819.
[7] Oud B, Flores C L, Gancedo C, et al. An internal deletion in MTH1 enables growth on glucose of pyruvate-decarboxylase negative, non-fermentative Saccharomyces cerevisiae . Microb Cell Fact, 2012, 11:113.
[8] Arndt M, Hitzmann B. Kalman filter based glucose control at small set points during fed-batch cultivation of Saccharomyces cerevisiae . Biotechnol Prog, 2004, 20(1) : 377-383.
[9] Paalme T, Kevvai K, Vilbaste A, et al. Uptake and accumulation of B-group vitamers in Saccharomyces cerevisiae in ethanol-stat fed-batch culture. World J Microbiol Biotechnol, 2014, 30(9) : 2351-2359.
[10] Hitit Z Y, Boyacioglu H, Ozyurt B, et al. Self-tuning GMV control of glucose concentration in fed-batch baker's yeast production. Appl Biochem Biotechnol, 2014, 172(8) : 3761-3775.
[11] Rhee I L, Ritzka A, Scheper T. On-line monitoring and control of substrate concentrations in biological processes by flow injection analysis systems. Biotechnol Bioproc Eng, 2004, 9(3) : 156-165.
[12] Valentinotti S, Srinivasan B, Holmberg U, et al. Optimal operation of fed-batch fermentations via adaptive control of overflow metabolite. Control Eng Pract, 2003, 11(6) : 665-674.
[13] Vovsik J, Hrncirik P, Náhlík J, et al. Adaptive control of Saccharomyces cerevisiae yeasts fed-batch cultivations. Chem Biochem Eng Q, 2013, 27(3) : 297-306.
[14] Storn R, Price K. Differential evolution - a simple and efficient heuristic for global optimization over continuous spaces. J Global Optim, 1997, 11(4) : 341-359.
[15] Ponsich A, Coello Coello C A. Differential evolution performances for the solution of mixed-integer constrained process engineering problems. Appl Soft Comput, 2011, 11(1) : 399-409.
[16] Mallipeddi R, Suganthan P N, Pan Q K, et al. Differential evolution algorithm with ensemble of parameters and mutation strategies. Appl Soft Comput, 2011, 11(2) : 1679-1696.
[17] Piotrowski A P. Differential evolution algorithms applied to neural network training suffer from stagnation. Appl Soft Comput, 2014, 21 : 382-406.
[18] Takamatsu T, Shioya S, Okada Y, et al. Profile control scheme in a bakers' yeast fed-batch culture. Biotechnol Bioeng, 1985, 27(12) : 1675-1686.
[19] Can I, Tosunoglu F, Kahya E. Daily streamflow modelling using autoregressive moving average and artificial neural networks models: case study of Coruh basin, Turkey. Water Environ J, 2012, 26(4) : 567-576.
[20] Liu Y C, Wu J, Liu Y, et al. Analyzing effects of climate change on streamflow in a glacier mountain catchment using an ARMA model. Quatern Int, 2015, 358 : 137-145.