基于气象数据的洱海蓝藻水华驱动因子及预警研究

2016-11-30 05:23陈莉琼陈晓玲蔡明祥李海军余永明
关键词:水华蓝藻洱海

陈莉琼, 张 娇, 陈晓玲*, 蔡明祥, 李海军, 余永明

(1.武汉大学 测绘遥感信息工程国家重点实验室, 武汉 430079;2.江西师范大学 鄱阳湖湿地与流域研究教育部重点实验室, 南昌 330022;3.北京城建勘测设计研究院有限责任公司, 北京100101)



基于气象数据的洱海蓝藻水华驱动因子及预警研究

陈莉琼1,2, 张 娇1, 陈晓玲1,2*, 蔡明祥1, 李海军1, 余永明3

(1.武汉大学 测绘遥感信息工程国家重点实验室, 武汉 430079;2.江西师范大学 鄱阳湖湿地与流域研究教育部重点实验室, 南昌 330022;3.北京城建勘测设计研究院有限责任公司, 北京100101)

以2013年洱海夏秋季蓝藻水华遥感监测结果为基础,探讨影响蓝藻水华形成的主要气象要素,统计分析发现,洱海蓝藻水华多发生在阴雨天转晴之后,长时间强日照及较大的气温日较差是诱发水华发生的主要因素,同时低风速、低气压有利于蓝藻上浮聚集,对水华形成起到促进作用.在此基础上,构建以水华发生与否的二值变量为因变量,各气象因子为预测变量的Logistic蓝藻水华气象风险概率预测模型.结果表明,模型各项检验指标均表现良好,预测准确率达到87.5%,该模型能够用于分析影响洱海蓝藻水华形成的关键气象因素,同时证明了利用气象数据辅助实现水华监测、预报预警的可行性.

洱海; 蓝藻水华; 气象因素; Logistic回归模型

洱海(25°36′~25°58′N、100°06′~100°18′E)是云南省第2大淡水湖,水域面积约250 km2,平均水深10.6 m,最大水深21.3 m,蓄水量达28.8×108m3(水位1 966 m).自1990s以来,由于人类活动的影响,洱海水质不断下降,当前正处于中营养向富营养化湖泊转变的过渡阶段[1].2013年夏秋季,洱海再次发生大规模蓝藻水华,影响了湖泊水生态环境,威胁着沿湖居民饮用水安全.因此有效实现水华监测,分析水华发生规律,探究影响水华发生的驱动因素,对于预防和治理洱海蓝藻水华具有重要意义.

国内外对内陆湖泊蓝藻水华的相关研究表明[2-3],在湖泊富营养化使水体内营养盐满足蓝藻生长的前提下,气象与水文条件的驱动作用,是诱发水华蓝藻上浮集聚的关键.Ndong等[4]利用温度、风速、风向构建水华气象指数模型,分析了在湖泊氮磷比满足蓝藻繁殖的先决条件下,气象条件对水华的影响;Zhang等[5]对太湖蓝藻水华的研究表明,水华发生起始时间和持续时间同温度、日照呈显著负相关,而与风速和降水呈正相关;谢国清等[6]发现日照和风速是影响滇池水华发生的关键因素.不同湖泊的气候环境条件及湖泊特性差异导致影响水华发生的因素也不尽相同,洱海属高原湖泊,该区域光照充足,太阳辐射量高出位于平原地区的太湖20%左右[7],相比海拔相近的滇池,洱海湖区的平均气温也明显偏高[8];同太湖等浅水湖泊不同,洱海湖区水深最高达30余米,由此造成的水温分层效应也可能影响着藻类生长繁殖[9].从洱海特殊地域气候环境出发,分析洱海蓝藻水华成因,王芸[10]发现藻类密度同光照强度及气温呈正相关;蔡燕凤等[11]认为高温、低降水量和低风速有利于洱海藻类大量繁殖.但上述分析仅从气象单因素角度考虑,未能明确气象条件的协同作用对水华的影响,并且以藻类密度(叶绿素浓度)作为衡量水华发生的指标,不能准确反映以遥感监测蓝藻覆盖水面为主要特征的水华现象,实际在夏秋季水华易发阶段,湖区内叶绿素浓度(藻类密度) 往往维持在较高的水平,但水华却不一定发生[12-13],而本文直接以发生与否衡量水华,更加直观地表现气象因素同水华的密切关系,同时利用Logistic回归模型,构建以气象因子为预测变量的蓝藻水华气象风险概率模型,以概率的形式表现气象条件对于水华暴发的影响,实现了从简单定性分析到定量表示气象因素协同作用于蓝藻水华.

本文针对2013年洱海蓝藻水华暴发事件,借助遥感技术获取蓝藻水华信息,同时分析水华同期气温、日照、降水、气压等气象条件同水华的相关性,探究影响洱海蓝藻水华发生的主要气象因素,进而利用Logistic回归模型,构建洱海蓝藻水华气象风险概率模型,定量分析气象条件同蓝藻水华的关系,为提高蓝藻水华监测能力,科学预防治理水华提供理论依据.

1 数据和方法

1.1数据资料

本文选取2013年8~11月洱海地区Landsat和HJ-1B影像作为水华数据来源,经目视检验,共得到可用影像16景。基于ENVI软件,进行辐射定标、几何校正等预处理。同期气象数据获取自中国气象科学数据共享服务网,站点为大理站(25°42′N, 100°11′E),包括气温、气压、降水量、日照时数、风速风向等逐日气象数据。

同期气象数据下载自中国气象科学数据共享服务网(http://cdc.nmic.cn),站点为大理气象站(25°42′N, 100°11′E),包括气温、气压、降水量、日照时数、风速风向等逐日气象数据.

1.2研究方法

1.2.1水华信息提取及统计方法 大量蓝藻在水面聚集,形成水华,使水体表现出类似陆地植被的光谱特征,相比NDVI、EVI等植被指数,FAI方法加入短波红外波段,能消除复杂大气环境的干扰,结果更稳定:

FAI=(ρNIR-ρRED)+(ρRED-ρSWIR)×

(λNIR-λRED)/(λSWIR-λRED),

(1)

式中,ρRED,ρNIR,ρSWIR分别代表红光、近红外、短波红外波段的反射率;λRED,λNIR,λSWIR分别为对应波段的中心波长.

受到遥感数据时间分辨率和天气状况的限制,获取的水华时间信息不是连续的,参照太湖水华研究[15-16]及洱海历史水华资料[17-18],对洱海水华发生次数和持续时间进行界定,认为在气象条件比较稳定的时期内,相隔3~5 d内的水华日归为一次水华事件.例如遥感监测表明2013年11月2、4、6、8、12和14日发生水华,期间气象条件没有显著变化,故认定同属一次水华,持续时间为2~14 d.

1.2.2Logistic模型 在以蓝藻覆盖水面作为主要特征的水华分析中,叶绿素、藻类密度等并不能准确反映水华暴发,用二值变量“发生(1)”、“不发生(0)”表征水华是更为恰当的方法,此时分析气象因素对蓝藻水华的影响可以理解为分析连续气象自变量同二值因变量之间的关系.

Logistic回归模型作为一种非线性回归概率模型,通过转换问题为分析被预测变量的条件概率同预测变量之间的关系,建立二元因变量同连续型或离散型自变量之间的联系,在二分类因变量问题中有较多应用[19].

(2)

式中,ρ=P(y=1|x1,x2,…,xk)为事件发生概率.

利用SPSS软件,构建Logistic回归模型,选择Hosmer-Lemeshow检验模型拟合优度,并利用Wald统计量对回归系数进行检验,确定最优回归模型.

2 洱海蓝藻水华统计及气象条件分析

2.1洱海蓝藻水华信息统计

结合遥感监测及历史资料可知,2013年8~11月洱海共发生6次蓝藻水华,持续时间达26 d(图1),蓝藻水华最早出现在8月15日,最晚发生在11月25日;8月、10月水华持续时间较短,均为3 d,9月、11月分别发生两次水华,持续时间为5 d和14 d.北部和中部水域是蓝藻水华的多发地带,尤其是大型水华集中发生在中部水域,10、11月蓝藻水华分布区域最广,几乎覆盖整个洱海,最大覆盖面积分别达到28 km2、17 km2.

图1 2013年洱海蓝藻水华发生时间分布Fig.1 The time distribution of cyanobacterial blooms in Erhai Lake in 2013

2.2水华发生期气象条件分析

2.2.1气温 考虑到气温的季节性差异,水华发生时段内日平均气温和日最高气温变化幅度较大,在整个水华发生阶段,有整体下降的趋势,两者不能反映温度同水华的关系,但是在水华发生时,气温日较差普遍偏大,平均日较差为14.9℃,最大日较差为18.6℃(10月10日),气温日较差同水华的Pearson相关系数为0.58(在0.01水平上双侧显著相关),尤其是11月初连续多日气温日较差在15℃以上(图2(a)),这也可能是导致该次水华持续时间长达13 d的原因,因此较大的气温日较差对蓝藻水华爆发具有一定促进作用,故选取气温日较差作为气温指标分析对水华的影响.

图2 洱海蓝藻水华期间气象条件变化Fig.2 Variation of meteorological factors during the period of cyanobacterial bloom

2.2.2日照 从洱海蓝藻水华发生时日照变化可知(图2(b)),水华发生日(除有明显降雨的8月16日、17日和9月15日外),日照时数普遍偏长,平均日照时数为9.1 h,最高达10.4 h(9月14日).日照时数同水华的相关系数为0.53,剔除有降雨的时段,相关系数达到0.8(在0.01水平上双侧显著相关).同时发现在水华爆发前1~3日日照时数均较低,表明在水华发生前普遍会有阴雨天气出现,而天气转晴之后突增的强日照可能是诱发水华发生的原因.因此,日照对于水华发生有明显影响,强日照天气,尤其是天气由阴转晴,日照突增,有利于蓝藻上浮聚集,增大表层藻类浓度,从而形成水华.

2.2.3降雨 洱海水华发生时仅有8月16日、17日和9月15日、11月12日有少量降水(图2(b)),其余22日降雨量均为零,无降雨天数占水华总天数的85%.并且降雨并非发生在水华发生开始阶段,而是处于水华发生中期或末期,所以可以认为降雨不利于水华发生.

分析历次水华发生前降水情况,可知在水华发生前均有明显降雨或者是连续阴雨天气,如10月2~8日洱海出现连续降雨,雨量累积49.4 mm,持续降雨导致温度下降,日最高气温由25℃降低16℃,10月9日之后天气转晴,气温回升,最高气温超过25℃,次日洱海发生大面积水华.由此可知,降水导致气温降低,一定程度上抑制了蓝藻的生长,同时强降水对湖泊水体的搅动作用,影响藻类垂直分布,降低表层藻类浓度,抑制了水华的发生,但是当连续多日阴雨天气结束,天气转晴时,由此带来的气温升高,日照增加,则是蓝藻水华的诱发因素.

2.2.4气压 水华发生期间日均气压最高最低分别为807.7 mPa、798.5 mPa,变化幅度比较大,难以反映同水华间的关系,而比较水华发生前同水华发生当日气压变化,发现水华发生时,气压有明显的下降趋势,在全部26日水华中,有21日平均气压低于前3日气压值,占水华总天数的78%,气压同水华相关系数为-0.18,表明气压降低有利于水中藻类上浮,形成水华.

2.2.5风 风速过大,会对水体产生类似搅拌的作用,使水体中藻类垂直分布趋于均匀,从而不利于藻类上浮,所以水华形成初期,较低的风速是重要影响因素.

统计水华发生日平均风速,其中有22日小于2 m/s,占总天数的85%;最大风速在2.5~6.8 m/s之间,其中在2.5~4.5 m/s之间最为集中,共19 d,占总天数的73%.平均风速和最大风速均值分别为1.5 m/s和4.1 m/s,由此可知,风速较低时,对水体的扰动作用较小,有利于蓝藻上浮聚集.

3 洱海蓝藻水华气象预测模型

从众多气象变量中筛选影响水华发生的主要气象因子是建立气象预测模型的首要任务.洱海水华一般发生在夏秋季,在阴雨天转晴之后,强日照、气温日较差较大的天气下,当气压偏低,风速较低时,水华发生的风险会较高.基于上述分析,初步确定气温日较差、日照时数、平均风速、较前3日气压差、日降雨量5个参数指标作为自变量因子,以水华发生与否作为被预测变量,构建Logistic回归模型.

3.1多重共线性诊断

同多元线性回归类似,Logistic模型自变量也可能存在共线性问题.选择方差膨胀因子(VIF)对自变量因子做多重共线性诊断.结果表明,当自变量分别剔除气温日较差和日照时数,两者的VIF由6.1和5.9降至1.2以内,满足共线性检验条件,同时发现两者的Pearson相关系数为0.90,表明自变量因子间存在多重共线性是由于两者的强相关.不难理解,阴雨天,日照时数偏低,同时大气对太阳辐射的削弱作用导致白天增温较慢,而夜间大气保温导致降温慢,从而气温日较差较小;而晴天,日照时间较长,白天太阳辐射强烈,地面升温快,夜间地面辐射相对较强,降温较快,从而气温日较差较大.利用Logistic回归后向逐步选择剔除冗余变量,筛选显著因子,能够解决共线性问题,进而建立最佳回归方程.

3.2Logistic模型结果及分析

2013年8~11月期间122 d,共26 d发生蓝藻水华,采用交叉验证方法,将全部数据随机均分为5组,每次选取一个子集作为验证集,其余4组作为训练样本,依据最优预测精度得到Logistic回归模型.

按照逐步回归最终选取包括日照时数、平均风速、较前3日气压差、日降雨量在内的4个气象因子,在5组试验中,Logistic模型预测准确率最低为83.3%,最高为87.5%,充分说明了该模型在洱海水华气象预测中的适用性.根据最高预测准确率得到最佳回归方程:

ApΔ-1.141·W+0.171·P.

(3)

其中,ρ表示预测概率,S表示日照时数,ApΔ表示较前3日气压差,W表示平均风速,P表示日降雨量.

在模型拟合优度Hosmer-Lemeshow检验中,sig.=0.826>0.05,统计不显著,表明模型很好地拟合了数据;预测正确率为87.5%,表明模型具有较好的稳定性.采用Wald统计量检验自变量影响的显著程度,所有自变量均满足在a=0.05水平上的统计性显著,其中日照时数的Wald值为15.58明显高于其他气象要素,说明该要素对于蓝藻水华的影响最为显著,同单因素分析中水华同日照时数强相关相一致.

图3 Logistic模型概率预测结果验证Fig.3 Verification of Logistic forecast model

进一步分析模型预测结果(图3),可以看出,出现漏报的情况一般发生在水华发生后期,此时蓝藻水华发生不仅仅同气象条件有密切关系,水华前期发生规模和程度也对后期水华发展有很大影响,大规模的水华不可能在短时间内从水面消失,导致水华仍在继续,但是日照、气温等气象条件并不满足水华发生的一般规律.如8月中旬水华发生时,17日有11.9 mm降雨,日照时数为5.5 h,气温日较差也仅为10.9℃,并不具备发生水华的气象特征,但根据遥感监测结果,在洱海北部和中部均出现蓝藻水华,因此并不能单单以该日预测风险概率为16.4%,就判定为漏报.同样,对于结果中出现的共7日错报情况进行分析,除8月6日和9月23日外,其余错报情况均出现在水华发生前后5日以内,如11月1日预测概率为53.6%,次日监测到蓝藻水华,在缺少当日遥感监测数据的情况下,无法判断11月1日是否发生水华,而直接认定为错报显然是不合理的.

同时,分析模型预测结果同水华覆盖面积之间关系(图4),发现水华覆盖面积越大,对应预测概率也越大,如10月水华覆盖面积最大约28 km2,预测概率达到87.5%,而11月份水华规模有所减小,最大覆盖面积仅为17 km2,此时预测气象概率也降到67%.说明Logistic模型预测结果不仅能够表征气象条件是否满足水华发生要求,一定程度上还能够指示水华发生规模,对于大型水华预报预警具有重要意义.

图4 Logistic模型预测结果同水华覆盖面积的关系Fig.4 The relationship between the result of Logistic model and blooms area

4 讨论

本文以遥感监测结果“发生”、“不发生”的二值变量作为被预测值,利用Logistic回归模型,建立其同气象因子之间的联系,实现了基于气象数据分析水华发生风险概率估算,相比其他研究中以叶绿素浓度等作为表征水华的参量,该模型预测目标直观、明确,得到的气象风险概率,能够直接反映水华发生情况.

模型分析同气象单因子分析具有较好的一致性,验证了初步筛选关键气象因子结果的可信性.日照时数同气温日较差表现出较强相关,在水华形成初期起到关键作用,蓝藻细胞具有垂直迁移的特性,在光照充足的情况下能够上浮至水体表层,吸收光能进行光合作用,而此时较高的温度又能增加光合作用速率,是水华形成的关键影响因素.降雨、气压及风速均与水华呈负相关,降雨导致温度降低,日照减少,抑制了蓝藻生长;低气压,低风速有利于水体中已经大量存在的藻类细胞上浮,促进水华的发生.

在模型预测准确度分析中,发现错报、漏报现象较多出现在水华发生前后较短时间内,此时不能单单以当日预测结果作为评断标准,应该在考虑遥感监测确定水华起始时间的可靠性及分析历史水华数据基础上,加入合理先验知识,更加客观地评价模型预测结果.在水华监测预警中,人们往往关注对于水华起始阶段的确定,以便及时开展应急预防治理工作,可以考虑加入“在已知水华起始日前3日内预测结果有效”的判断条件,重新评估模型结果:如11月上旬水华事件中,1~3日预测概率分别为53.6%,25.2%,67.1%,虽然遥感监测结果表明水华始于2日,但依据上述设定,可认为1日预测结果在误差范围以内,不属于错报.对于水华期间的漏报情况,可以综合考虑前后几日内预测结果,若出现高风险概率情况,就不能直接做出无风险预报,而应当同样加强现场水华监测,根据其他环境条件做出综合判断.

将遥感技术应用于水华监测,相比其他现场调查等方法,更加快速、高效,并且能够获取湖泊历史水华信息,基于该数据得到的模型预测效果较好,表明该方法用来分析洱海蓝藻水华是可行的,但是还应当考虑到遥感数据受天气影响较大,导致获取的水华信息不可避免存在一定误差,所以应该考虑结合其他监测手段,以获取更为完整的水华监测数据,降低数据误差,提高模型精度.

在富营养化湖泊中,水体中氮磷等营养物质满足蓝藻生长需要,此时营养盐浓度不再是限制水华发生的主要因素,气象、水文等条件是诱发水华发生的关键,基于此建立的气象风险概率预测模型,重点分析气象因素对水华的影响,若将该模型应用到实际洱海水华预测中,还存在一定局限性,还需要综合分析水质、气象等环境条件的影响,进一步验证分析,完善预测模型.

5 结论

以2013年洱海蓝藻水华遥感监测结果为基础,分析了洱海湖区气象条件变化对水华暴发的影响,结果表明日照和气温同水华发生具有更显著的相关性,水华一般发生在强日照和气温日较差较大的情况下;虽然气压、风速对水华的影响较小,但在低气压、低风速的条件对水华的促进作用也不容忽视.利用Logistic回归建立气象风险概率预测模型,将气象因子同遥感监测水华结果建立直接联系,预测准确率达87.5%,该模型能够用于辅助蓝藻水华监测预报.针对模型存在的错误预测情况,下一步还应该考虑借助其他监测手段,获取更为完善准确的水华及其他监测数据,提高模型预测准确性,并应用于实际水华预测工作.

[1] 彭文启,王世岩,刘晓波. 洱海水质评价[J]. 中国水利水电科学研究院学报, 2005, 3(3): 192-198.

[2] 孔繁翔,高 光.大型浅水湖泊的蓝藻水华形成机理研究的思考[J].生态学报, 2005, 25(3): 589-595.

[3] 孔繁翔,马荣华,高俊峰, 等. 太湖蓝藻水华的预防, 预测和预警的理论与实践[J]. 湖泊科学, 2009, 21(3): 314-328.

[4] NDONG M, BIRD D, NGUYEN-QUANG T, et al. Estimating the risk of cyanobacterial occurrence using an index integrating meteorological factors: application to drinking water production[J]. Water Research, 2014, 56: 98-108.

[5] ZHANG M, DUAN H, SHI X, et al. Contributions of meteorology to the phenology of cyanobacterial blooms: implications for future climate change[J]. Water Research, 2012, 46(2): 442-452.

[6] 谢国清, 李 蒙, 鲁韦坤, 等. 滇池蓝藻水华光谱特征、遥感识别及暴发气象条件[J]. 湖泊科学, 2010, 22(3): 327-336.

[7] 陈 桥, 韩红娟, 翟水晶, 等. 太湖地区太阳辐射与水温的变化特征及其对叶绿素a的影响[J]. 环境科学学报, 2009, 29(1): 199-206.

[8] 李杰君. 洱海富营养化探析及防治建议[J]. 湖泊科学, 2001, 13(2): 187-192.

[9] 赵巧华, 孙绩华. 夏秋两季洱海、太湖表层混合层的深度变化特征及其机理分析[J]. 物理学报, 2013, 62(3): 1-9.

[10] 王 芸. 洱海夏秋季蓝藻种群动态变化及水华成因分析[J]. 大理学院学报(综合版), 2008, 7(12): 39-42.

[11] 蔡燕凤,张虎才,陈光杰,等. 洱海生态环境研究现状及存在问题[J]. 地球科学前沿, 2013(4):241-252.

[12] ONDERKA M. Correlations between several environmental factors affecting the bloom events of cyanobacteria in Liptovska Mara reservoir (Slovakia)—A simple regression model[J]. Ecological Modelling, 2007, 209(2): 412-416.

[13] 郁建桥, 吕学研. 太湖遥感可测性蓝藻水华发生条件分析[J]. 环境科学与技术, 2015(6): 015.

[14] HU C. A novel ocean color index to detect floating algae in the global oceans[J]. Remote Sensing of Environment, 2009, 113(10): 2118-2129.

[15] 马荣华,孔繁翔,段洪涛,等. 基于卫星遥感的太湖蓝藻水华时空分布规律认识[J]. 湖泊科学, 2008, 20(6): 687-694.

[16] DUAN H, MA R, XU X, et al. Two-decade reconstruction of algal blooms in China’s Lake Taihu[J]. Environmental Science & Technology, 2009, 43(10): 3522-3528.

[17] 吕兴菊, 朱 江, 孟 良. 洱海水华蓝藻多样性初步研究[J]. 环境科学导刊, 2010, 29(3): 32-35.

[18] 董云仙. 云南九大高原湖泊藻类研究进展[J]. 环境科学导刊, 2014(2): 001.

[19] LIU Y, WANG Z, GUO H, et al. Modelling the effect of weather conditions on cyanobacterial bloom outbreaks in Lake Dianchi: a rough decision-adjusted logistic regression model[J]. Environmental Modeling & Assessment, 2013, 18(2): 199-207.

Research on meteorological factors and Logistic prediction model of cyanobacterical blooms in Erhai Lake

CHEN Liqiong1,2, ZHANG Jiao1, CHEN Xiaoling1,2, CAI Mingxiang1, LI Haijun1, YU Yongming3

(1.State Key Laboratory of Information Engineering in Surveying, Mapping and Remote Sensing,Wuhan University,Wuhan 430079; 2.Key Laboratory of Poyang Lake Wetland and Watershed Research,Ministry of Education, Jiangxi Normal University, Nanchang 330022;3.Beijing Urban Construction Exploration & Surveying Design Research Institute CO.,Beijing 100101)

The main meteorological parameters of cyanobacterial blooms were explored based on the Remote Sensing monitoring results of Erhai Lake in 2013. Statistical analysis showed that the blooms mostly occurred after rainy days and the main cause were the strong sunshine and higher daily temperature variation. Lower wind speed and air pressure also contributed to the formation of blooms. Upon the above analysis, Logistic regression forecast model was established with dependent and independent variable of two-valued variable and meteorological factors, respectively. Results indicated that the Logistic prediction model performed well with forecast accuracy of 87.5%. The model is able to be used to analyse the key meteorological factors for bloom formation, and it also demonstrated the feasibility for monitoring and forecasting bloom by using the aided meteorological data.

Erhai Lake; cyanobacterial blooms; meteorological factors; Logistic model

2016-01-08.

国家水体污染控制与治理科技重大专项(2013ZX07105-005);国家自然科学基金项目(41301366);测绘地理信息公益性行业科研专项项目(201512026).

1000-1190(2016)04-0606-06

P334+.6;P237

A

*通讯联系人. E-mail: xiaoling.chen@whu.edu.cn.

猜你喜欢
水华蓝藻洱海
藻类水华控制技术及应用
河湖藻类水华应急治理决策研究
洱海月下
洱海,好美
南美白对虾养殖池塘蓝藻水华处理举措
南美白对虾养殖池塘蓝藻水华处理举措
针对八月高温蓝藻爆发的有效处理方案
洱海太湖石
爱上洱海,只需要这十个瞬间
可怕的蓝藻