農(nóng)作物施肥分析

農(nóng)作物施肥分析

1問題重述某研究所為了研究氮肥(N).磷肥(P).鉀肥(K)三種肥料對于土豆和生菜的作用，分別對土豆和生菜兩種作物做了三組實驗，實驗中將每種肥料的施肥量分為10個水平，在考察其中一種肥料的施肥量對土豆或生菜的產(chǎn)量的影響時，其他兩種肥料的施肥量均在第七個水平上，實驗數(shù)據(jù)如下面表格所示，其中ha表示公頃,7表示噸，畑表示千克，試建立反映施肥量與產(chǎn)量關(guān)系的模型，并從應(yīng)用價值和如何改良等方面作出評價.施肥量與產(chǎn)量關(guān)系的實驗數(shù)據(jù):土豆：氮肥磷肥鉀肥施肥量伙g/加)產(chǎn)量(///??)施肥量伙g/加)產(chǎn)量(//ha)施肥量(kg/ha)產(chǎn)量{tlha}015.18033-46018.983421.362432-474727.356725.724936-069334.8610132.297337-9614038.5213534.039841-0418638.4420239.4514740-0927937.7325943.1519641-2637238.4333643.4624542-1746543.8740440.8329440-3655842-7747130.7534242-7365146.22生菜：氮肥磷肥鉀肥施肥量產(chǎn)量施肥量產(chǎn)量施肥量產(chǎn)量伙g/加)(///??)伙g/加)(//ha)(kg/Wha)06-39011.02伊15.75499-482812.704716.769812.465614.569316.8914714.338416.2714016.2419617.1011217.7518617.5629421.9416822.5927919.2039122.6422421.6337217.9748921.3428019.3446515.8458722.0733616.1255820.1168524.5339214.1165119.402模型假設(shè)八h忽略土壤現(xiàn)有肥力對農(nóng)作物生長的影響.假設(shè)在所有作物的生長環(huán)境的陽光，空氣、水.溫度變化都相同.假設(shè)施肥時肥料的質(zhì)量準(zhǔn)確.3問題分析:我們研究的是不同肥料的施肥量對農(nóng)作物產(chǎn)量的影響，就是要求出施肥量與作物產(chǎn)量之間的關(guān)系，冃前我們僅有單一施肥量與作物產(chǎn)量的數(shù)據(jù)，我們需要用線性回歸的方法擬合出每一種肥料對作物產(chǎn)量影響的方程，然后根據(jù)第七個水平的共同施肥量建立多種肥料對一種作物產(chǎn)量影響的方程，從而找到最正確施肥量使得作物的產(chǎn)量最大.4模型的建立和求解氮肥施肥量改變時對土豆產(chǎn)量的影響由現(xiàn)有數(shù)據(jù)建立氮肥施肥量對土豆產(chǎn)量的影響散點圖如下:4SO3-4003-358-3003-2505-20W-400500105 二8 , 3X3400500一N肥航肥吊圖1氮肥施肥■與土豆產(chǎn)■關(guān)系散點圖由圖1氮肥施肥■與土豆產(chǎn)■:關(guān)系散點圖我們發(fā)現(xiàn)氮肥的施肥量與土豆產(chǎn)量大致成二次關(guān)系，故我們建立如下模型：y二勺球+$*+%+￡N

現(xiàn)將上述模型線性化，算出爲(wèi)令值為；Vy那么上述模型可以簡化為：y="+bx+bxN}N32N在計算多元線性回歸方程時，假設(shè)回歸方程為：y=b+bx+bx+,??+bx+￡ ￡~N(0(7')0i] 22 {tf> ' /(4.1)其中?,b、，…，》,ci'都是X],?..,耳無關(guān)的未知參數(shù).設(shè)(X]&，?-時"(、，X??…，、,兒) (4.2)是一個樣本，和一元線性回歸的情況一樣，我們用最大似然估計法來估計參數(shù).即取&,&,???,￡使當(dāng)*二氐b嚴(yán)?！?bp=0時If(2=工5-瓜-6心If(2=工5-瓜-6心bp%$(4.3)到達(dá)最小.求Q分別關(guān)于%,b、，…,g的偏導(dǎo)數(shù)，并令它們等于零，得字二-2乞(H-%_b丹 bx)=0,pip“字二-逹。ODj/=|廠瓜-bjnbp%?*=o,j=\2“字二-逹。ODj/=|廠瓜-bjnbp%?*=o,j=\2…、p(4.4)化簡(4.4)式得加+勺工心+?&2+?r-1zzr-1n ,人遼"f=i+◎&“二？兒“ft丿M勺工兀+爲(wèi)工和2+???+勺￡護廠工勺兒

;=if=i n瓜工.5+b工巒、+b?H%兀2+.?.+"工痘=口』、r=l /=1i=li=li=\(4.5)我們將式(4.5)稱為正規(guī)方程組.為了求解的方便，將式(〉寫成矩陣形式.為此,引入矩陣：/-I11...1>1'XjY=■■x21...■■■■■■■■?■x2p…x嘰■兒./-]J-l2>汕于是〔4.5〕式可寫成(4.5)'XrXB=XY(4.5)'〔4.5）式兩邊左乘XX的逆矩陣〔X7%）-1〔設(shè)〔X「X〕存在）得到〔〕啲解.4bi=(Xbi=(XTxyiXTYtbp入入/K久+休￥]+..?切?二$作為Jt心1,兀二瓜+1価+..?+bpXpp元線性回歸方程方程:r14.7416經(jīng)計算得正規(guī)方程組的解為B=*b2那么求的原假設(shè)模型為川戲9715心+=(X'X)-lXrY=用爐來檢驗?zāi)Ｐ偷臄M合度，將所得模型與現(xiàn)有數(shù)據(jù)的散點繪制在同一張圖上得氮肥施肥量與土豆產(chǎn)量擬合圖形如下:50403020210対忙國豐102200300400500氮肥施肥址計算得疋擬合度良好，可以較好的反映氮肥施肥量與土豆產(chǎn)量的關(guān)系.磷肥施肥量改變時對土豆產(chǎn)量的影響由現(xiàn)有數(shù)據(jù)畫出磷肥施肥量與土豆產(chǎn)量關(guān)系的散點圖:40X)0-3sncr3250-.001COXJD2ODO333DCT1403CO圖2磷肥施肥量與土豆產(chǎn)?關(guān)系散點圖由圖2磷肥施肥量與土豆產(chǎn)量關(guān)系散點圖可以看出磷肥施肥母與土豆產(chǎn)量大致成二次關(guān)系，故可建立磷肥施肥量與土豆產(chǎn)量關(guān)系的模型如下：y=b+hx2+bxP 0 }P 2p為簡化計算，計算出尤令之為&,那么上述假設(shè)模型可化為：y=bix2+bx+b

Pp2p()? ■32.916■經(jīng)計算得正規(guī)方程組的解為B=八=(XtX)1XtY=入b那么假設(shè)模型可解得：?=-0.000138%;+0.07185%用疋來檢驗?zāi)Ｐ偷臄M合度，將所得模型與現(xiàn)有數(shù)據(jù)的散點繪制在同一張圖上得磷肥施肥量與土豆產(chǎn)量擬合圖形如下：50403020212+0.0719X+32.916R20100 200 300 400磷肥施肥雖經(jīng)計算得/?2=0.864，模型擬合程度良好，可以較好的反映磷肥施肥量與土豆產(chǎn)量之間的關(guān)系.鉀肥施肥量改變時對土豆產(chǎn)量的影響由現(xiàn)有鉀肥施肥量與土豆產(chǎn)量的數(shù)據(jù)畫出如下鉀肥施肥量與土豆產(chǎn)量關(guān)系散點圖50CONOCO-?0CCrco40000 6COCOK肥施肥圖3鉀肥施肥■與土豆產(chǎn)■關(guān)系散點圖由圖3鉀肥施肥量與土豆產(chǎn)量關(guān)系散點圖可以初步看出鉀肥的施肥量與土豆的產(chǎn)量大致成對數(shù)關(guān)系，我們因此假設(shè)土豆的產(chǎn)量與鉀肥的施肥量模型如T:為簡化計算令％=1咔先將全部的1咔計算出來用心來代換，貝I假設(shè)模型可化為：兒二a+bx+slnk求一元線性回歸的解法如下：模型假設(shè)y=a+bx+￡￡~7V(0,b()參數(shù)估計X=d+x+￡y—{a+bx〕=￡，勺?N〔0,b?〕i i i"ix=丈［牙-〔4+b?！常?=Q〔a,b〕i-i i-1翅仏.=0b〕daoQ〔a.h〕=0db

=0故方程有唯一解:故方程有唯一解:b= 弘4〕2

H一入一a=y-bx〔3〕得到回歸方程為：y=a+bx■A經(jīng)計算得那么假設(shè)模型為：丿兒=5.811%.+7.021?=因此可以得出原假設(shè)模型為：片=5-8111n兀5040:30:205040:30:2010o0 2 4y=03137x2+5.8579X+

R26 8 10 12鉀肥施肥量經(jīng)計算得模型/?2=0.865，因此模型擬合度良好，可以較好的反映鉀肥施肥量與上豆產(chǎn)量之間的關(guān)系.氮、磷、鉀肥施肥量共同改變時對土豆產(chǎn)量的影響氮肥、磷肥.鉀肥三種肥料對土豆的產(chǎn)量的影響模型組如下:yiV=-0.000344+0.19715XN+164<丹二一球859勺+yK由以上模型組可以假設(shè)土豆產(chǎn)量成〕與氮肥、磷肥、鉀肥三種肥料的施肥量的關(guān)系模型如下:>=%4+a>=%4+aNxN+bpXp+c|nX,+￡將30組數(shù)據(jù)匯總在一起，在沒有其他肥料的施肥量數(shù)據(jù)的地方填上對應(yīng)的肥料的第七水平的施肥量，分別計算出對、A-;、1BA;，然后做多元線性回歸分析，經(jīng)計算得:aN2a=0.1838Nbp〔<b二p因此解出假設(shè)模型為：134+0.1838".-147x;心+6.0667In無計算F來檢驗?zāi)Ｐ蛿M合度，得疋=0.953,接近1擬合度很好，可以很好的反映氮肥、磷肥、鉀肥的施肥量對土豆產(chǎn)量的影響.計算出最大產(chǎn)量時的氮磷肥三種肥料的施肥量如下表〔MATLAB代碼見附錄人表1氮確鉀肥施肥■對應(yīng)土豆產(chǎn)量表氮肥施肥量 磷肥施肥量 鉀肥施肥量 十石產(chǎn)量(kg/ha) (焙/ha) (kg/ha) (1/ha)293.6104 256.1250 651由表1氮磷鉀肥施肥量對應(yīng)土豆產(chǎn)量表得當(dāng)?shù)适┓柿繛?93.6104kg//"心磷肥施肥量256.1250kg/ha9鉀肥施肥量651kg/ha時土豆的產(chǎn)量最大值為46.2116t〃?n.氮肥施肥量改變時對生菜產(chǎn)量的影響由現(xiàn)有氮肥施肥量的改變對生菜產(chǎn)量的影響數(shù)據(jù)，建立氮肥施肥量與生菜產(chǎn)量關(guān)系的散點圖如下：2250-20X)0-秦10X)0-2250-20X)0-秦10X)0-00圖4氮肥施肥■與生菜產(chǎn)量關(guān)系散點圖由圖4氮肥施肥量與生菜產(chǎn)■關(guān)系散點圖可以看出生菜的產(chǎn)量與氮肥的施肥量大致成二次關(guān)系，故假設(shè)生菜的產(chǎn)量與氮肥施肥量之間的模型為：y生二+g為了簡化計算先計算出€用「代換，那么原模型可化為：y生二b&+E兀+%+&經(jīng)計算得正規(guī)方程組的解為：B=h\=(x!xy^XrY=0.10132.那么原假設(shè)模型解得：呪兀對模型進行疋檢驗，將所得模型與現(xiàn)有數(shù)據(jù)的散點繪制在同一張圖上得氮肥施肥量與生菜產(chǎn)量擬合圖形如下:2+0.1013X+2+0.1013X+1A2(R210100 200 300 400 500氮肥施肥址經(jīng)計算得/?2擬合度較好，可以較好的反映氮肥施肥量與生菜產(chǎn)量之間的關(guān)系.磷肥施肥■改變時對生菜產(chǎn)量的影響為了初步分析生菜產(chǎn)量在僅有磷肥改變時所受的影響，在現(xiàn)有數(shù)據(jù)的根底上建立磷肥施肥量與生菜產(chǎn)量關(guān)系的散點圖如下：圖S確肥施SDD圖S確肥施SDD-肥量與生菜.00200OT產(chǎn)■關(guān)系散點圖由圖5確肥施肥■與生菜產(chǎn)■關(guān)系散點圖初步可以看出磷月巴的施月巴量與生菜產(chǎn)量之間大致呈現(xiàn)二次關(guān)系，故建立如下模型：y生二%+隔+輕+￡計算出兀用「代換那么原假設(shè)模型可以化為如下多元線性回歸方程:

旳二b〔〕+b?+M+$經(jīng)計算得正規(guī)方程組的解為B=經(jīng)計算得正規(guī)方程組的解為B==(XTxyiXTY=0.060597.■■-■■-用F檢驗?zāi)Ｐ蛿M合度，將所得模型與現(xiàn)有數(shù)據(jù)的散點繪制在同一張圖上得磷肥施肥量與生菜產(chǎn)量擬合圖形如下:302520211050 y=y=-5E-O5X2R2?+0.0606X+經(jīng)計算的爐=0.959，模型擬合度較好，故可以較好的反映磷肥的施肥量與生菜產(chǎn)量之間的關(guān)系.200400600800200400600800磷肥施肥雖鉀肥施肥量改變時對生菜產(chǎn)量的影響由鉀肥施肥量與生菜產(chǎn)量關(guān)系數(shù)據(jù)，畫出鉀肥施肥量與生菜產(chǎn)量關(guān)系散點圖如下:21CO19￡O-1BCO-忡cor16CO-佔CO-.00圖6鉀肥施肥錄與生菜產(chǎn)■關(guān)系散點圖由圖6鉀肥施肥■與生菜產(chǎn)■關(guān)系散點圖可以看出第＞\組數(shù)據(jù)與其他數(shù)據(jù)存在較大差異，為了保證數(shù)據(jù)的準(zhǔn)確性，對本次實驗數(shù)據(jù)進行殘差檢驗，殘差分析圖如下：0ResidualPlot-4-00 200 300400500600700 圖7鉀肥施肥■對生菜產(chǎn)量影響的殘專分析圖從殘差圖可以看出數(shù)據(jù)的殘差離零點的遠(yuǎn)近，除了第八組數(shù)據(jù)的95熾差置信區(qū)間不包含零點，其他組數(shù)據(jù)的殘差95%置信區(qū)間均包含零點，那么第八組數(shù)據(jù)可視為異常點剔除.在剔除異常點后IB出鉀肥施肥量與生菜產(chǎn)量關(guān)系的修正散點圖如下：21CO-20CO-'19CO-產(chǎn)18CGIS17CC-1SCC-00W肥施JJ巴600000.11-200 400鉀肥施肥址6000.11-200 400鉀肥施肥址600800圖8鉀肥施肥■與股產(chǎn)*關(guān)系修止散點圖由圖8鉀肥施肥■與生菜產(chǎn)■關(guān)系修止散點圖可以看出鉀肥的施肥量與生菜的產(chǎn)量之間大致呈現(xiàn)二次關(guān)系，故建立模型如下：y生二bo+bA+E^+￡同理將岸計算出來，用S代換，進行多元線性回歸分析.%二(XTxylXr二(XTxylXrY=2那么原假設(shè)模型為：；1必用用檢驗?zāi)Ｐ偷臄M合度，將所得模型與現(xiàn)有數(shù)據(jù)的散點繪制在同一張圖上得鉀肥施肥量與生菜產(chǎn)量擬合圖形如下：y=-7E-06X2+0.0106X+15.836R2計算得F=0-830模型的擬合度良好，可以較好的反映鉀肥施肥量與生菜產(chǎn)量之間的關(guān)系.氮、磷'鉀肥施肥■共同改變時對生菜產(chǎn)量的影響由得氮肥、磷肥、鉀月歸種肥料對生菜的產(chǎn)量的影響模型組如下：畑；暫＜兒牛二—尤幵兒根據(jù)上述模型組，我們可以建立生菜產(chǎn)量〔y〕與氮肥〔心〕、磷肥〔丹〕、鉀肥〔?！持g的模型：y二%+碣+E叫+仇尤+*Xp+bX；+bX+￡5 6k-4.815425經(jīng)計算得正規(guī)方程組的解為：B=0.048748那么可以解出模型為：y=；M一兀+0.048748a；，兀用F檢驗?zāi)Ｐ蛿M合度，計算得/?2二接近于1,擬合度較高，可以較好的反映氮、磷、鉀三種肥料同時施肥時對生菜產(chǎn)量的影響.經(jīng)計算得出上述模型的最大產(chǎn)量時的氮磷鉀的施肥量如下表〔代碼見附錄〕：表2氮磷鉀肥施肥■對應(yīng)生菜產(chǎn)■表氮肥施肥量磷肥施肥量鉀肥施肥量生菜產(chǎn)量〔kgfha〕〔kg/ha〕Ckg/ha〕〔1/ha〕212.828364L4227457.685523.5586由表2氮磷鉀肥施肥量對應(yīng)生菜產(chǎn)量表可知當(dāng)?shù)适┓柿繛?13畑〃"心磷肥施肥量641kg/ha,鉀肥施肥量458kg/ha時生菜的產(chǎn)量最大值為23.561/ha(為了方便施肥，對小數(shù)點數(shù)據(jù)后進行四舍五入).5模型的評價與改良我們種植農(nóng)作物時，不能單一的只考慮髙產(chǎn)量，同時還需考慮作物的售價和肥料價格，以到達(dá)收益的最大化?目前各種作物價格以及肥料價格如下表表3作物及肥料價目衷種類 價格(元/kg)+石生菜氮肥磷肥鉀肥用叫分別表示土豆-生菜的利潤」、P,分別表示土豆.生菜的價格，Z八亞分別表示氮.磷、鉀三種肥料的單價.5.1土豆最正確施肥根據(jù)上述價格及方程那么可以建立土豆的利潤模型為：vv,=13現(xiàn)+-;心+6.0667Inx-174）門—V~xZ—xz5 V’ K / nppkk帶入單價得上述模型為（MATLAB代碼見附錄）:vv,=13琉+838心-147x；+0.0753?+6.0667Inx-74）-1.6x-0.7?-3-5忑Kn經(jīng)計算得：表4土豆最正確施肥計算結(jié)果表~~氮肥的施肥量 磷肥的施肥量 鉀肥的施肥量、 八八、 利潤(元/M)(kg/ha) (k/u/ha) (kg/ha)292.9003255.4611651.0000163430由表4土豆最正確施肥計算結(jié)果表可知道當(dāng)土豆的施肥量為氮肥293kg/ha,磷肥255kg/ha，鉀肥651kg/ha，可獲得最大利潤163430元/ha.生菜最正確施肥根據(jù)上述價格及方程那么可以建立生菜的利潤模型為：叫兀；+0.098706?—尤+0.048748?-0.000030#+5425)-16j-p—忑經(jīng)計算得(MATLAB代碼見附錄):附錄土豆計算最大產(chǎn)量的MATLAB代碼M文件：functionf=fun3(x)f=-(一0.000313*x(:L)*x(l)+0.1838*x(l)-0-000147*x(2)*x(2)+0.0753*x(2)+6?0667*log(x(3))-29.7174);主程序：x0=[l;l;l]；A=[100;010;001];b=[471;342;651];Aeq=[];beq=[]；VLB=[0;0;0];VUB=[];x0=[0;0;0];[