章明清，李 娟，張立成，姚寶全，張 華

（1. 福建省農業(yè)科學院土壤肥料研究所，福州 350013；2. 福建省農田建設與土壤肥料技術推廣站，福州 350003）

二次多項式函數(shù)是當前常用的計量施肥理論模型[1-6]，具有扎實的統(tǒng)計學理論基礎和計算簡便易行等優(yōu)點。但在施肥實踐中，二元、三元二次多項式肥效模型存在大量的非典型式[7-10]，至今仍然困擾著計量施肥的研究和應用。近年的研究顯示，二次多項式肥效模型及其他類似的多項式肥效模型存在設定偏誤、多重共線性和異方差等問題[2]，制約了普通最小二乘（OLS）的建模有效性。為此，針對二次多項式肥效模型設定偏誤，作者構建了非結構肥效模型[11-12]；針對三元二次多項式肥效模型的多重共線性和異方差危害，分別提出了主成分回歸（PCR）[13]和可行廣義最小二乘回歸（FGLS）[14]的建模技術，顯著提升了建模成功率。

迄今在三元肥效模型建模中，除了經(jīng)典的二次多項式肥效模型OLS建模法外，還有PCR和FGLS建模法，以及三元非結構肥效模型及其非線性最小二乘（NLS）建模法等。那么，針對具體作物的肥效試驗結果，該如何選擇這些模型和建模法以達到最佳建模效果呢？作者在探討PCR建模法時，曾發(fā)現(xiàn)“OLS建模法結合PCR建模法”的建模策略可使早稻的典型三元二次多項式肥效模型比例明顯高于單獨使用OLS回歸或者PCR回歸[13]的建模方法。本文在分析整合各三元肥效模型及其建模法的專業(yè)適用性基礎上，提出優(yōu)化建模策略，旨在盡可能提高典型式的比例和田間肥效試驗結果的建模成功率。

1 材料與方法

1.1 水稻和露地蔬菜氮磷鉀田間肥效試驗資料的收集整理

近十年來，福建省完成了眾多氮磷鉀“3414”設計的田間肥效試驗。在匯總這些試驗資料時，僅收集整理土肥技術力量較強的或項目負責人有能力嚴格把關的相關縣市完成的試驗結果。在此基礎上，根據(jù)“3414”試驗設計中，氮、磷、鉀各4個施肥水平及其試驗產量分別繪圖，若3個圖形均大致呈現(xiàn)拋物線型關系，則保留該試驗點資料，反之，則棄除該點試驗資料。截止2017年底，共收集整理了水稻和露地蔬菜田間試驗資料1 122個（表1）。供試土壤包括灰泥田、黃泥田、灰沙田、赤沙土等福建主要耕作土壤類型（土屬），供試土壤主要理化性狀（表2）采用常規(guī)方法測定[15]，處理（6）的N2、P2、K2施肥量及其試驗產量結果見表2。相關作物田間試驗設計、土樣采集和測定等，請參照文獻[16]。

表1 福建水稻和露地蔬菜“3414”設計試驗資料收集情況Table 1 Data collection of the field experiments in paddy fields and open vegetable gardens following the “3414” designing in Fujian Province

表2 水稻和蔬菜供試土壤主要理化性狀及其處理（6）施肥量和產量Table 2 Main physical and chemical properties of the tested soils，fertilizer application rates and yields of the 6th treatment

1.2 三元肥效模型的整合

1.2.1 非結構肥效模型 近年的水稻盆栽和田間肥效試驗均表明，施用單位養(yǎng)分的稻谷增產量與施肥量之間的關系是典型的指數(shù)函數(shù)關系[11]。在此基礎上，作者構建了一元非結構肥效模型[11]：

式中，Y為作物產量，kg·hm-2；X為施肥量，kg·hm-2；s0為土壤供肥當量，以N、P2O5、K2O養(yǎng)分形態(tài)計量，kg·hm-2；c為施肥對產量的效應系數(shù)；A表示基礎土壤對作物產量的生產能力。式（1）模型克服了一元二次多項式肥效模型的設定偏誤以及一次項和二次項回歸變量間高度線性相關的問題，具有更高的擬合精度。

在式（1）模型中，當施肥量和土壤供肥當量均為零時，作物產量必等于零。因此，根據(jù)植物營養(yǎng)元素功能不可相互替代的原理，三元非結構肥效模型可由三個一元非結構肥效模型相乘導出[12]：

式中，N0、P0、K0分別表示供試土壤的氮、磷、鉀供肥當量，kg·hm-2；c1、c2、c3分別表示施用氮、磷、鉀養(yǎng)分的增產效應系數(shù)；N、P、K分別表示N、P2O5、K2O施肥量，kg·hm-2；其他代數(shù)符號的含義與式（1）相同。

數(shù)學理論分析表明，式（1）和式（2）模型在一定施肥量范圍內存在一個作物產量峰值，該峰值對應的施肥量即為最高產量施肥量。因此，根據(jù)微積分原理，令作物產量Y分別對氮、磷、鉀施肥量的導數(shù)等于零，得到最高產量施肥量的計算式；令作物產量Y分別對氮、磷、鉀施肥量的導數(shù)等于農產品和肥料價格倒數(shù)比，得到經(jīng)濟產量施肥量的計算式[11-12]。

式（1）和式（2）肥效模型均為非線性模型，而且不能直接進行線性化處理，模型參數(shù)估計需采用非線性最小二乘法[17]。假設非線性模型為Y=f（X，a），為求得參數(shù)a的估計值，可求解最小二乘問題：

其解?a作為參數(shù)a的估計值。非結構肥效模型的回歸顯著性檢驗與二次多項式肥效模型相似，但式（1）模型的回歸自由度為2，式（2）模型的回歸自由度為6。

1.2.2 多項式肥效模型 傳統(tǒng)上，假設施用單位養(yǎng)分的增產量與施肥量之間滿足線性關系，導出了一元二次多項式肥效模型[18]。根據(jù)數(shù)理統(tǒng)計學的線性加和性原理，進一步導出三元二次多項式肥效模型；有的學者根據(jù)研究需要，對二次多項式模型進行變形處理，形成了0.5次方和1.5次方多項式肥效模型[18]。顯然，這種多項式肥效模型的導出過程，雖然具有統(tǒng)計學理論支持，但缺乏植物營養(yǎng)學專業(yè)理論依據(jù)。

同理，將三元非結構肥效模型的組成項e-cN-c2P-c3K轉化為e-c1N×e-c2P×e-c3K。由于模型參數(shù)c1、c2、c3均在10-3級[12]，取泰勒級數(shù)展開式的前兩項，并忽略c1、c2、c3的兩兩乘積項以及c1c2c3乘積項和NPK三因子交互項等高次項，再令b0=AN0P0K0，b1=A（1-N0c1）P0K0，b2=A（1-P0c2）N0K0，b3=A（1-K0c3）N0P0，b4= -c1P0K0，b5=-c2N0K0，b6= -c3N0P0，b7=（1-N0c1）（1-P0c2）K0，b8=（1-N0c1）（ 1-K0c3）P0，b9=（1-P0c2）（1-K0c3）N0。則式（2）模型可轉化為：

這就是當前常用的三元二次多項式肥效模型，是三元非結構肥效模型的簡化式。但是，這種經(jīng)過簡化導出的二次多項式肥效模型，導致最高施肥量之前和最高施肥量之后的施肥效應是對稱關系，不符合生產實際。同時，導致強烈的多重共線性和異方差危害，制約了經(jīng)典最小二乘法回歸建模的有效性。

上述肥效模型整合過程可歸納為圖1所示。顯然，非結構肥效模型的導出過程具有較強的專業(yè)邏輯性，二次方、0.5和1.5次方多項式肥效模型導出過程的專業(yè)邏輯性較差。

1.2.3 三元肥效模型的典型性判別 肥效模型典型性涉及邊際產量導數(shù)法推薦施肥的可靠性。由于農業(yè)生產條件的復雜性，根據(jù)田間肥效試驗結果建立的肥效模型，方程效應曲線或曲面的形狀多種多樣[7]。在通過統(tǒng)計顯著性檢驗的前提下，三元二次多項式肥效模型存在典型式和3種不同類型的非典型式[9]。若肥效模型同時滿足：（1）一次項系數(shù)的代數(shù)符號為正數(shù)，二次項系數(shù)的代數(shù)符號為負數(shù)；（2）肥效模型存在全局最高產量點；（3）邊際產量導數(shù)法得到的推薦施肥量均落在試驗設計施肥量范圍內，這種肥效模型符合植物營養(yǎng)學的一般肥效規(guī)律，稱為典型式，可用邊際產量導數(shù)法推薦施肥。反之，若3個條件中有任何一個條件不能滿足，則分別稱為肥效模型系數(shù)符號不合理的非典型式、肥效模型無最高產量點的非典型式和肥效模型推薦施肥量外推的非典型式，此時邊際產量導數(shù)法推薦施肥結果不可靠。如何判斷三元二次多項式肥效模型是否存在全局最高產量點？章明清等[9]總結提出了簡易方法。

三元非結構肥效模型在通過統(tǒng)計顯著性檢驗的前提下，也可能存在不同類型的模型。（1）若模型參數(shù)A、N0、P0、K0、c1、c2、c3均大于零，而且氮磷鉀推薦施肥量均落在試驗設計施肥量范圍內，此類模型滿足了植物營養(yǎng)學的一般肥效規(guī)律，稱為典型式，可用邊際產量導數(shù)法計算推薦施肥量；（2）若模型系數(shù)A、N0、P0、K0、c1、c2、c3有一個或一個以上的系數(shù)值為負數(shù)，該類模型稱為系數(shù)符號不合理的非典型式；（3）若模型參數(shù)均大于零，但邊際產量導數(shù)法的推薦施肥量有一個或一個以上落在試驗設計施肥量范圍外，該類模型稱為推薦施肥量外推的非典型式。因非結構肥效模型的數(shù)學結構特點，若模型參數(shù)均大于零，模型必有全局最高產量點，因而三元非結構肥效模型不存在無最高產量點的非典型式這種類型。

1.3 三元肥效模型建模與推薦施肥量計算的計算機實現(xiàn)

有許多數(shù)學軟件和統(tǒng)計分析軟件均能進行肥效模型的參數(shù)估計和顯著性檢驗，本文采用MATLAB R2015b軟件進行相關的編程和計算。其中，三元非結構肥效模型的參數(shù)估計調用nlinfit功能函數(shù)；三元二次多項式肥效模型OLS建模法則調用regress功能函數(shù)；三元二次多項式肥效模型的PCR建模法是調用pca功能函數(shù)進行主成分分析，然后根據(jù)提取的主成分得分矩陣與各處理的試驗產量，用regress功能函數(shù)進行回歸建模；三元二次多項式肥效模型FGLS建模法則調用fgls功能函數(shù)進行模型參數(shù)估計和統(tǒng)計顯著性檢驗。

上述4種建模法的具體數(shù)學原理、計算過程和MATLAB R2015b軟件功能函數(shù)的使用方法可參閱相關專著[19-21]。

2 結 果

2.1 三元肥效模型不同建模法的擬合效果

針對1 122個水稻和露地蔬菜氮磷鉀田間肥效試驗結果，分別采用三元二次多項式肥效模型OLS、PCR[13]和FGLS[14]建模法對各個試驗點進行回歸分析。表3的統(tǒng)計結果表明，OLS建模法的三元典型式比例平均僅占試驗點總數(shù)的19.8%，而PCR和FGLS建模法則分別提高至34.0%和27.1%。應用三元非結構肥效模型[12]及其NLS建模法，三元典型式平均比例達到試驗點總數(shù)的41.4%，較三元二次多項式肥效模型的OLS、PCR和FGLS的建模法分別提高21.6、7.4、14.3個百分點。

表3 水稻和蔬菜三元肥效模型不同建模法的擬合效果比較Table 3 Comparison of fitting between ternary fertilizer response models using different modeling methods for rice and vegetable

土壤肥力水平對典型三元肥效模型比例有重要影響。菜后稻是蔬菜地的輪作水稻，因土壤肥力水平普遍較高，OLS建模法的典型三元二次多項式肥效模型比例平均低至4.8%，即使采用PCR、FGLS建模法或非結構肥效模型，三元典型式比例也分別僅有28.6%、12.7%和27.0%。

統(tǒng)計結果還顯示，二次多項式肥效模型OLS、PCR和FGLS建模法的無最高產量點非典型式的平均比例大致相當；然而，非結構肥效模型由于數(shù)學結構特點，無最高產量點的非典型式比例則為零。在未通過統(tǒng)計顯著性檢驗、模型系數(shù)不合理、推薦施肥量外推等模型種類方面，不同模型建模法之間的平均比例則差異較大。PCR建模法使未通過顯著性檢驗的模型比例從OLS建模法的28.7%下降至19.2%，而FGLS建模法則下降至零；PCR建模法使二次多項式肥效模型系數(shù)符號不合理的非典型式比例從OLS建模法的30.1%下降至16.4%，非結構肥效模型則進一步下降至4.6%，但FGLS建模法則提高至48.3%。PCR建模法使二次多項式肥效模型推薦施肥量外推的非典型式比例從OLS建模法的6.1%提高至14.4%，F(xiàn)GLS建模法則與OLS建模法相當，但非結構肥效模型推薦施肥量外推的非典型式比例平均達30.8%。

由此可見，二次多項式肥效模型經(jīng)典建模法的成功率明顯偏低，克服多重共線性或異方差危害的建模法均有利于提高典型式比例，而同時克服了模型設定偏誤和多重共線性危害的非結構肥效模型則較大程度地提高了建模成功率。因此，不同三元肥效模型及其建模法的適用性有明顯差別，合理選擇肥效模型及其參數(shù)估計方法對提高建模成功率具有重要價值。

2.2 供肥當量參數(shù)N0、P0和K0與對照區(qū)產量的關系

在建模成功率較高的三元非結構肥效模型中，N0、P0和K0分別表示土壤N、P2O5、K2O的供肥當量。為考察模型參數(shù)估計值的可靠性，本文對能得到典型三元非結構肥效模型的試驗點，根據(jù)水稻和露地蔬菜收獲時的氮、磷、鉀對照區(qū)產量，與對應試驗點的N0、P0和K0估計值進行線性回歸分析。表4的結果表明，除了土壤肥力水平較高的菜后稻的3個回歸方程和煙后稻的磷素回歸方程外，其他具有10個以上試驗點的供試作物，N0、P0和K0估計值與對應作物對照區(qū)產量之間均有統(tǒng)計顯著水平的線性正相關關系，顯示N0、P0和K0估計值較好地反映了稻田和菜地土壤氮磷鉀供肥潛力。因此，三元非結構肥效模型不僅具有較高的建模成功率，模型參數(shù)N0、P0和K0也具有較明確的物理意義。

表4 非結構肥效模型參數(shù)N0、P0和K0與對照區(qū)作物產量的相關性Table 4 Correlation analysis of model parameters of N0，P0 and K0 with crop yields in CK

2.3 三元肥效模型的四步建模法

在上述建模和統(tǒng)計過程中，作者發(fā)現(xiàn)三元二次多項式肥效模型能得到典型式的某些試驗點資料，并不能保證三元非結構肥效模型也一定能得到典型式，反之亦然。某些田間試驗結果可能僅適用于4種建模法中的一種或者幾種。從盡可能提高典型肥效模型比例的實用角度出發(fā)，考慮到非結構肥效模型具有較高的擬合精度、更寬的適用范圍[11-12]和較高的建模成功率（表3），同時較好地克服了二次多項式肥效模型推薦施肥量偏高[1，22]的問題，故將非結構肥效模型作為數(shù)據(jù)擬合的第一步；OLS建模法是具有最佳統(tǒng)計性能的經(jīng)典方法，因而將第一步建模不能得到典型式的所有試驗點，采用三元二次多項式肥效模型OLS建模法作為數(shù)據(jù)擬合的第二步；因PCR和FGLS建模法均為有偏估計[20]，PCR建模成功率高于FGLS建模法（表3），故將前兩步不能得到典型式的余下試驗點依次采用PCR和FGLS作為第三步和第四步的數(shù)據(jù)擬合技術依據(jù)。經(jīng)過計算機多次反復試驗比較，這種建模程序具有最佳的建模效果，簡稱為四步建模法，技術路線見圖2。

表5 的擬合結果表明，在1 122個試驗資料中，典型三元肥效模型占試驗點總數(shù)的平均比例達到57.5%，較單獨使用三元非結構肥效模型NLS建模法、三元二次多項式肥效模型OLS、PCR、FGLS建模法的典型式比例分別提高了16.1、37.7、23.5、30.4個百分點，大幅度提升了典型三元肥效模型比例和田間肥效試驗建模成功率。分析還表明，除了供試土壤肥力普遍較高的菜后稻典型式比例較低外，四步建模法在雙季稻、單季稻和露地蔬菜等作物間的典型式比例差異很小，顯示該建模程序具有較好的可靠性。

表5 四步建模法對水稻和露地蔬菜典型三元肥效模型比例的影響Table 5 Effect of the four-step modeling method on the proportion of typical ternary fertilizer response models for paddy fields and open vegetable gardens

3 討 論

3.1 肥效模型不同建模法的適用性

1 122個水稻和露地蔬菜氮磷鉀田間肥效試驗資料，采用三元非結構肥效模型NLS建模法以及三元二次多項式肥效模型OLS、PCR和FGLS建模法，三元典型式比例平均分別為41.4%、19.8%、34.0%和27.1%，顯示不同三元肥效模型及其建模法的建模成功率具有明顯的差異。

經(jīng)典回歸分析的OLS建模法具有理想的統(tǒng)計性能，但合理使用必需滿足10個基本假設條件[23]。在多項式肥效模型中，誤差項方差為常數(shù)、自變量之間不存在線性相關和模型被正確設定等三個假設條件通常得不到滿足，導致出現(xiàn)了異方差、多重共線性和模型設定偏誤等問題，制約了OLS有效性和統(tǒng)計檢驗可靠性[2]，成為當前二元、三元二次多項式肥效模型出現(xiàn)大量非典型式的主要技術原因。為提高建模精度和模型預測可靠性，PCR技術在土壤學中已得到普遍應用[24-26]。該法通過從試驗設計矩陣中提取互不相關的主成分，利用主成分得分值和試驗產量建立肥效模型，消除了多項式模型多重共線性危害，對提高肥效建模成功率具有重要作用[13]。FGLS建模法是消除異方差危害的有效方法，但在土壤學中的應用還較少[27]。在三元二次多項式肥效模型中，約有25%的肥效模型存在顯著水平的異方差[14]，導致建模結果發(fā)生異常。FGLS建模法是利用加權方法消除異方差，然后采用OLS法進行參數(shù)估計。因此，OLS建模法具有優(yōu)良統(tǒng)計性能，但受到諸多假設條件限制；PCR和FGLS建模法分別單獨克服多重共線性和異方差危害，且均屬于有偏估計[20]。

事實上，三元二次多項式肥效模型同時存在模型設定偏誤、多重共線性、異方差等問題，只有同時消除或緩解這些已知問題的危害，才能在較大程度上提高典型式比例。非結構肥效模型假設施用單位養(yǎng)分的增產量與施肥量之間滿足指數(shù)函數(shù)關系，克服了模型設定偏誤；該模型是非線性模型，且不能直接線性化處理，克服了多重共線性危害。非結構肥效模型同時較好地克服了多項式肥效模型的這些已知問題，是建模成功率得到較大提升（表3）的重要原因。

3.2 三元肥效模型的優(yōu)化建模策略

在表3的建模過程中，作者發(fā)現(xiàn)某些田間試驗結果要得到典型肥效模型，僅能適用上述4種建模法中的一種或者幾種，其余的方法達不到預期目的，原因是作物產量具有隨機性。三元二次多項式肥效模型和三元非結構肥效模型的數(shù)學形式不同，方程誤差項對回歸建模的影響程度不同。這種狀況就可能出現(xiàn)三元二次多項式肥效模型能得到典型式的某些試驗點，三元非結構肥效模型不能保證均得到典型式的結果，反之亦然。

可惜的是，什么樣的試驗資料會造成模型設定偏誤或者多重共線性或者異方差對建模結果產生嚴重影響，目前在統(tǒng)計學上并無明確的答案[28]。考慮到農業(yè)生產條件的復雜性和多樣性，反映到作物施肥效應曲線或曲面也必然具有多樣性，期望通過一種模型或建模法來反映或概括這種多樣性是不現(xiàn)實的。因此，多種肥效模型及其建模法的綜合應用有其合理性。由于三元非結構肥效模型具有最高的建模成功率（表3），OLS建模法具有優(yōu)良的統(tǒng)計性能，PCR和FGLS能分別消除多重共線性和異方差危害。因此，從得到盡可能高的典型式比例的實用角度出發(fā)，多種模型或建模法綜合應用的四步建模法是一種現(xiàn)實可行的建模策略。根據(jù)該建模策略得到的推薦施肥量，在已完成的234個水稻和14個露地蔬菜大田示范中[29-30]，與習慣施肥相比，不同稻作間平均增產稻谷4.0%～12.5%，凈增收為875～2 616 yuan·hm-2；不同蔬菜種類間平均增產6.8%～10.6%，凈增收達692～3 834 yuan·hm-2。結果顯示，四步建模法得到的推薦施肥量具有較好的可靠性。

3.3 計量施肥模型構建的若干思考

肥效模型是實現(xiàn)計量施肥的主要技術手段。但在模型構建研究中，普遍存在一個認識誤區(qū)，即：施肥模型應包含盡可能多的變量，除施肥量外，還應包括氣候條件、土壤條件、作物品種、田間管理措施等，否則模型就不完善。實際上，計量模型考慮的變量過多，必然導致各變量之間提供的信息出現(xiàn)重疊，造成計量不準或出現(xiàn)異常，就如統(tǒng)計學上的多重共線性等問題一樣。影響施肥量的各種因素間大多數(shù)情況下是非線性關系，如何削弱或消除這種變量信息重疊對計量準確性的影響，目前尚無成熟和方便的方法可用。同時，作為研究對象的農作物具有生命特征，施肥量或其他變量與產量的關系是滿足統(tǒng)計學規(guī)律的隨機關系，既使將上述變量因子全部考慮進去也不可能變?yōu)榇_定性關系?？梢哉f，產生上述認識誤區(qū)的根源在于，研究者大多習慣于用確定性關系的經(jīng)典思維模式來考慮和評價具有隨機性質的施肥模型，而對現(xiàn)代統(tǒng)計思維模式的應用尚重視不夠。

作者認為，計量施肥模型至少可分為計量施肥基礎理論模型和區(qū)域推薦施肥應用模型兩個層次，前者是后者的基礎?；A理論模型的研究重點應在于探討如何建立普遍適用的作物施肥效應模型，甚至包括養(yǎng)分損失的施肥效應模型。這種模型應該簡潔美觀，既要滿足專業(yè)要求，又要符合統(tǒng)計學理論，就像理論生態(tài)學的邏輯斯蒂（Logistic）模型一樣[31]。因此，考慮的變量指標僅需要抓住那些具有共性特點的因子就可以了，例如，施肥量和土壤肥力等變量因子。

區(qū)域推薦施肥應用模型應在基礎理論模型的基礎上，結合服務區(qū)域的農業(yè)生產條件和區(qū)域特點，構建適合當?shù)厣a實際的區(qū)域推薦施肥系統(tǒng)。由于考慮變量太多會造成難以對模型進行精確的參數(shù)估計和施肥計量的問題，區(qū)域推薦施肥模型宜采用分區(qū)模式。例如，當區(qū)域小至縣域甚至村域，因為降水量、氣溫和作物管理水平等變量在區(qū)域內很可能無顯著差異，在構建區(qū)域推薦施肥模型時就可以不考慮了。因此，既使在應用層面上，評價區(qū)域推薦施肥模型的合理性和可靠性，關鍵也在于該計量模型是否抓住了服務區(qū)域內具有顯著水平差異的那些變量，并將這些變量納入考慮中，而不是全部變量。

基于上述思考，三元非結構肥效模型和三元二次多項式肥效模型均應屬于計量施肥的基礎理論模型范疇，可為區(qū)域推薦施肥提供基礎模型依據(jù)；三元非結構肥效模型的土壤供肥當量參數(shù)N0、P0、K0，與相應氮、磷、鉀對照區(qū)作物產量具有顯著水平的線性正相關（表4），為將“測土”與肥效模型緊密結合提供了一條可能途徑。

4 結 論

非結構肥效模型具有較強的專業(yè)邏輯性，其簡化式即為二次多項式肥效模型。四步建模法是提高典型三元肥效模型比例的優(yōu)化建模策略。