基于免疫遺傳算法估計Van-Genuchten方程參數

李 琳，黃 銳，彭化偉

(1.密山市富密灌區水利水保總站連珠山分站，黑龍江密山158300;2.大慶地區防洪工程管理處，黑龍江大慶163000)

基于免疫遺傳算法估計Van-Genuchten方程參數

李 琳1，黃 銳2，彭化偉2

(1.密山市富密灌區水利水保總站連珠山分站，黑龍江密山158300;2.大慶地區防洪工程管理處，黑龍江大慶163000)

土壤水分特征曲線是模擬土壤水分運動和溶質運移的一個重要參數，Van Genuchten方程是最常用的土壤水分特征曲線方程。建立了2種Van Genuchten方程參數的優化模型:基于非線性擬合函數lsqcurvefit與免疫遺傳算法共同改進的混合免疫遺傳算法(LIGA)和非線性擬合函數lsqcurvefit。結果表明:分別針對粉壤土、細砂土、礫石的吸濕和脫濕曲線得到了Van Genuchten方程的相關參數，兩種方法計算的θ值與實測的θ值擬合較好，相對誤差較小，LIGA的計算精度高于lsqcurvefit函數的計算精度，改進的LIGA算法收斂速度快，計算精度高，為該領域研究提供了新的思路。

遺傳算法;免疫算法;非線性擬合函數;土壤水分特征曲線;土壤水吸力;土壤含水量;VG模型

土壤水分特征曲線(soil water characteristic curve，SWCC)是土壤的一個基本水力參數，表示土壤水吸力和土壤含水量之間的關系，反映了土壤水的能量與數量的關系，可以間接地反映土壤中孔隙大小的分布［1-2］。土壤水分特征曲線受多種因素的影響，其數學模型比較復雜，未知參數多由經驗得到，且參數較多，大部分用于描述土水特征曲線的數學模型都是根據經驗、土體結構特征和曲線的形狀而建立起來的。準確可靠的土壤水分特征曲線模型是分析土壤水勢的高低與植物根系吸水的速率的變化規律，研究土壤水分和相應的吸力對作物生長的有效性，對于非充分灌溉具有重要的實用價值。土壤水分特征曲線的確定方法總的來說有直接法和參數法，直接法是指通過實驗方法直接測定土壤水分特征曲線的方法，用數學表達(經驗公式)來描述水分特征曲線為參數估計法。目前，在已經建立的眾多數學模型中，國內外普遍使用的描述土壤水分特征曲線的方程是Van-Genuchten模型(簡稱VG模型)［3-4］。根據不同土壤類型的的土壤含水率、土壤吸力的數據，運用參數估計法求解VG模型的較多，如廖林仙等［5］運用智能算法推求了VG方程參數;馬英杰等［6］運用阻尼最小二乘法，擬合了VG方程參數;彭建平等［7］對VG模型中的參數利用Matlab軟件中的非線性擬合函數來確定，模型計算取得了較好的效果;郭向紅等［8］利用混合遺傳算法對VG方程進行了求解;還有一些其他方法在推求VG方程的參數中取得了良好的效果。

混合免疫遺傳算法是非線性擬合函數(lsqcurvefit)和免疫遺傳優化算法(immune genetic algorithm，IGA)的耦合模型，具有免疫遺傳算法全局信息交換和局部深度搜索相結合的特點，利用lsqcurvefit函數的非線性擬合能力來捕獲變量間的非線性關系，進一步提高模型的局部搜索能力，加快向最優解的收斂速度。模型優化結果證明，該模型具有較強的非線性逼近能力和較強的魯棒性、較快的收斂速度、較高的搜索精度。本文將混合免疫遺傳算法引用于VG方程中，為方程參數的確定提供參考依據。

1 Van Genuchten方程

Van Genuchten公式描述壓力水頭與含水量的關系，即土壤水分能量與數量的關系，方程形式為:

式中:θ為體積含水率，cm3/cm3;θr為殘留含水率，cm3/cm3;θS為飽和含水率，cm3/cm3;h 為負壓，cmH2O;α，n，m 分別表示土壤水分特征曲線形狀的參數。式(1)中含有4個參數α，n(其中n中含有m)，θS，θr，根據實測的土壤含水率θ和土壤水吸力h資料，采用最小二乘法擬合結果的殘差平方和最小為目標，即式(2)，通過優化得到這4個參數:

式中:θi為實測含水率，cm3/cm3;hi為實測土壤含水率θi對應的土壤水吸力，cm;θ( hi，X)為根據式(1)計算的土壤含水率，cm3/cm3;X 為待求參數向量(α，n，θS，θr);N 為實測數據個數。

2 算法的構建

基于非線性擬合函數改進的混合免疫遺傳算法(Lsqcurvefit Immune genetic algorithm，LIGA)將遺傳算法［9］(genetic algorithm，GA)、免疫算法 (Immune Algorithm，IA)、非線性擬合函數lsqcurvefit結合構造的一種新的混合免疫遺傳算法。首先，在Matlab優化工具箱中，lsqcurvefit函數采用高斯—牛頓(Gauss-Newton)方法或者麥夸特(Levenberg-Marquardt)方法，用于求解最小二乘非線性數據擬合問題。混合免疫遺傳算法以lsqcurvefit函數的非線性求解作為種群的初始染色體的最優值，避免了算法開始階段搜索的隨機性與盲目性，使得混合免疫遺傳算法具有很強的非線性擬合能力與非線性逼近能力，群體具有的免疫機制對群體進行控制和調節，把目標函數和制約條件作為群體的抗原，保證生成的染色體直接與問題相關聯，收斂方向得以控制，對與抗原親和力高的個體進行記憶，并從中提取免疫疫苗對抗體群進行注射，保證抗體群在更新過程中的多樣性，促進求解的速度，提高算法的效率，并防止群體的退化，在很大程度上抑制混合免疫遺傳算法的未成熟收斂［10］。結合免疫算法、遺傳算法、lsqcurvefit函數的原理，下面給出具體實現步驟:

Step1:抗原輸入。輸入目標函數式(2)和各種約束作為混合免疫遺傳算法的抗原。

Step2:產生初始群體。生成初始抗體群N、促進記憶細胞庫N1、檢測記憶細胞庫N2。對初次應答，初始抗體N、N1、N2全部隨機產生，而對再次應答，則抗體群N借助免疫系統的記憶機制，部分初始抗體由免疫記憶單元獲取，即包含最優抗體促進記憶細胞庫，其余抗體隨機產生。

Step3:計算抗體適應度。群體的抗原和抗體V之間的親和力:axv=optv，其中optv表示抗原和抗體之間的匹配程度，本文用抗原和抗體之間的適應值函數fittness(v)來表示，即axv=fittness(v)。N個抗體的信息熵為:

式中:Hj(N)為N個抗體第j位的信息熵，Pij為N個抗體中的第j位為數值ki的概率?？贵wv和抗體w的親和力為:

抗體v的濃度cv為:

式中:Tac1為一個預先確定的親和力閾值;

Step4:免疫記憶細胞庫更新。根據一定比例，從抗體群N中選出親和力高的抗體，用它們替換促進記憶細胞庫N1中親和力低的抗體;將抗體群N中親和力低的抗體選入檢測記憶細胞庫N2，用它們替換檢測記憶細胞庫N2中親和力高的抗體。

Step5:抗體群的促進與抑制。當記憶群體N1中抗體濃度的最大值cmax低于抗體濃度閾值th時，記憶群體中的抗體處于多樣化，為未飽和狀態;否則抗體趣于一致化，為飽和狀態，濃度大的抗體將被淘汰，隨機產生的新個體代替被淘汰的個體。用檢測記憶群體N2去檢測N中是否含有已搜索過的抗體(即N2中記錄的抗體)，如果有就用隨機抗體取代它。

Step6:將促進記憶群體N1與檢測后的抗體群N相結合生成新的抗體種群，即父代群體y。

Step7:利用Matlab優化工具箱中的非線性擬合函數lsqcurvefit對目標函數進行優化，將優化的非線性結果作為初始個體最優值置于父代群體y中，提高各染色體對初始最優解的搜索能力，非線性擬合的數學模型［11］［12］:

根據輸入數據xdata和得到的輸出數據ydata找到與目標函數F(x，xdata)最佳的擬合參數向量x，以x作為群體中每條染色體的初始最優個體值。

Step8:對于父代群體y進行選擇操作、雜交操作和變異操作［9-13］。

Step9:終止條件是否滿足?如果滿足，結束迭代，否則轉向Step3。

3 實例分析

以文獻［14］中的粉壤土、細沙土、礫石3種介質的吸濕和脫濕曲線的實測資料為例，說明混合免疫遺傳算法在Van Genuchten方程參數擬合中的應用效果。表1為分別采用基于非線性擬合函數lsqcurvefit改進的LIGA和非線性擬合函數lsqcurvefit兩種方法得到的不同土壤類型的Van Genuchten方程參數結果。采用MATLAB7.0編程處理，選取選定父代初始群體規模為n=400，抗體濃度閾值th=0.85，交叉概率pc=0.8，變異概率pm=0.8，優秀個體數目選定20個，Van-Genuchten方程的估計參數取值范圍，α取值為0～1、n取值為1～10、θS取值為 0.2～0.7、θr取值為 0～0.2;通過MATLAB程序設計與其工具箱命令的條用實現非線性曲線擬合函數lsqcurvefit的計算。

表1 土壤水分特征曲線Van Genuchten方程參數計算結果

從表1中結果可以看出，LIGA與lsqcurvefit函數優化的參數估計結果比較接近，進而比較2種方法的求解擬合值與實測值之間的殘差平方和范數可知，對于目標函數值優化計算，LIGA優化結果計算結果θ值與實測的θ的殘差平方和范數小于lsqcurvefit函數的計算，可見改進的LIGA優化的參數結果優于lsqcurvefit函數優化的參數結果，并且這兩種方法的計算均達到了很高的精度。

4 結論

Van Genuchten方程是最常用的土壤水分特征曲線方程，方程中的參數擬合具有非線性特點，本文將免疫算法、遺傳算法和lsqcurvefit函數三者相結合構建新的混合免疫遺傳算法(LIGA)，并且嘗試性地將改進后的LIGA用于Van Genuchten方程參數擬合中，并進行了數值試驗，結果表明:改進的LIGA的計算精度比lsqcurvefit模型高，LIGA模型用于Van Genuchten方程參數求解問題中是可行的，既拓展了LIGA的應用領域，又為解決Van Genuchten方程參數擬合問題提供了新的思路和方法。

［1］趙愛輝，黃明斌，史竹葉.兩種土壤水分特征曲線間接推求方法對黃土的適應性評價［J］，農業工程學報，2008，24(9):11.

［2］盧小慧，靳孟貴，汪丙國.欒城農業生態系統試驗站土壤水分特征曲線分析［J］，中國農村水利水電，2006(12):30.

［3］李春友，任理，李保國.利用優化方法求解VanGenuchten方程參數［J］，水科學進展，2001，12(4):473-475.

［4］呂殿青，邵明安.非飽和土壤水力參數的模型及確定方法，應用生態學報［J］，應用生態學報，2004，15(1):63-165.

［5］廖林仙，邵孝侯，徐俊增.基于智能算法推求VanGenuchten方程的參數［J］，水利學報，2007(1):696-697.

［6］馬英杰，虎膽·吐馬爾拜，沈 冰.利用阻尼最小二乘法求解VanGenuchten方程參數［J］，農業工程學報，2005，21(8):179-180.

［7］王小華，賈克力，劉景輝，李立軍.VanGenuchten模型在土壤水分特征曲線擬合分析中的應用［J］，干旱地區農業研究，2009，27(2):179-180.

［8］郭向紅，孫西歡，馬娟娟.基于混合遺傳算法估計van Genuchten方程參數［J］，水科學進展，2009，20(5):677-680.

［9］付強，金菊良，梁川.基于實碼加速遺傳算法的投影尋蹤分類模型在水稻灌溉制度優化中的應用［J］，水利學報，2002(10):39-42.

［10］葛紅，毛宗源.免疫算法的改進［J］，計算機工程與應用，2002，14(47):47-48.

［11］彭建平，邵愛軍.用MatLab確定土壤水分特征曲線參數［J］，土壤，2007，39(3):433-436.

［12］唐家德.基于MATLAB的非線性曲線擬合［J］，計算機與現代化，2008(6):15-18.

［13］金菊良，楊曉華，丁晶.標準遺傳算法的改進方案—加速遺傳算法.系統工程理論與實踐［J］.2001(4):9-10.

［14］王金生，楊志峰，陳家軍，等.包氣帶土壤水分滯留特性研究［J］. 水利學報，2000(2):1-6.

S15

A

1007-7596(2012)02-0046-03

2011-12-19

李琳(1977-)，女，黑龍江密山人，工程師;董銳(1978-)女，黑龍江密山人，工程師;彭化偉(1976-)，男，黑龍江呼瑪人，工程師。。