專利名稱:作物性能分析的方法和系統的制作方法
技術領域:
本發明涉及一個根據廣闊區域中試驗數據評價農作物品種性能的過程。本發明的一個過程根據一個統計混合效應模型,采用空間估計和空間預測來對比品種性能。
背景技術:
在一個農作物新品種的培育中,要對該品種和與之競爭的其它品種采集性能數據。這些性能數據包括與指定農作物有關的多種農藝學特性的測量值;例如,對于Zea mays,要測量其糧食產量、糧食濕度和作物倒伏。
評價一個新的農作物品種(下文中稱為“品種”)潛在的商業價值時,要將其農藝學性能與其它品種的農藝學性能進行對比。其它對照品種包括,培育該品種的公司中已商業化和未商業化的品種,以及競爭對手公司的商業化品種。注意,對現有的商業化品種也會進行這種同樣類型的評價,以確定它們是應當保留在市場上還是應當被培育中的、更新的品種替代。
新品種和對照品種的農藝學性能數據來自多個試驗點。這些試驗點通常廣泛分布于試驗所包括品種的適應區域。覆蓋這些試驗點的適應區域往往很大,其范圍達數百平方英里。例如,一個新的Zea mays培育品種的試驗區可能從衣阿華州西部到密執安州東部,從威斯康星州中部到伊利諾伊州南部。
由于試驗項目中的差異,對于一個指定的品種和它的競爭者,數據往往相當“不平衡”,指的是并非指定組中的所有品種都出現在所有的試驗點上。考慮只有一對品種的試驗數據,也就是新品種和單一的競爭者,某些試驗點兩個品種都有,而其它點只有這兩個品種之一。
對這些性能數據進行分析,以便確定與對照品種相比,新品種具有足夠大性能優勢的地理范圍,從而證明它應當引入這些地區的市場。在理想情況下,在整個試驗區域中,與所有對照品種相比,正在培育的品種都具有顯著的性能優勢。然而在某些情況下,一個品種可能僅在一個地區具有性能優勢,但是它仍然可以在該地區中滿足相當大的市場需求。所以,不僅要在整個試驗區域將指定品種與其它對照品種相比確定性能特征,而且也要在各個地區內這樣做,這一點很重要。
不同地理位置或地區兩個品種相對性能的差異來自所謂的“環境影響的基因型”(Sprague & Eberhart,1976)。環境影響的基因型是由不同品種對環境條件的反應差異造成的。這些環境條件可能包括如晝長、溫度、土壤濕度和病蟲害壓力。注意,術語“環境”可以指指定一年的不同位置,或者指定位置的不同年度,或者位置與年度的某種組合。
Bradley等人(1988)介紹了品種性能評價的傳統統計分析方法。這些傳統方法通常是基于“位置匹配的”數據,也就是,對于一個指定品種相對于一個對照品種,在分析中使用的數據僅僅來自兩個品種共同種植的試驗點。采用一次配對的t試驗,對兩個品種沒有性能差異的假設進行檢驗。此外,為了獲得指定地區相對性能的有關結論,分析中只使用該地區內試驗點的數據。
對于位置匹配的數據采用t試驗的傳統分析效率很低,至少有五個理由。第一,它沒有使用全部數據;它使用的數據僅僅來自兩個品種共同種植的試驗點。第二,為了進行地區性能對比,它不使用所關注地區以外的鄰近區域的數據。第三,它不利用與所關注之性能特征有關的共變量(如灌溉或土系)來幫助解釋和預測這些差異。第四,在分析模型中它只使用年內的數據;如果一個模型使用多年度的數據,可以獲得更穩健的結論。第五,它是基于一個不正確的假設一個試驗點的觀測結果獨立于其它點的結果。
有兩個概括的理由,說明以上列舉的低效為什么限制了廣闊區域中的農作物評價。第一是數據的使用效率。要作出結論時,實驗人員自然希望使用全部數據。在當前的情況下,這包括一個品種在沒有種植其它品種之處的數據,包括與所關注地區類似區域的數據,包括有關共變量的數據,還包括跨年度的數據。統計方法應當爭取充分利用所有可用的信息。傳統分析受到限制的第二個理由是它基于數據獨立的經典假設。品種性能數據幾乎毫無例外地不滿足這種獨立性假設,使統計的結論無效,當數據并不真正表明差異時,往往使人推斷出存在品種差異。
現有文獻綜述現有文獻的簡要綜述也強調現有技術中潛在的不足,本發明爭取解決這些問題。
對農場地塊的實驗應用空間統計的領域中,關系最密切的論文之一、可能是最值得關注的工作是由Zimmerman和Harville(1991)完成的。作者認為觀測結果是一個隨機田的表現,引入了所謂的隨機田線性模型(RFLM)。在這個模型中,趨勢是由一個平均結構所模擬的,而小尺度的相關性通過空間自相關結構模擬。參數估計是通過一種似然方法進行。通過實際數據分析,作者試圖表明他們的模型優于最近鄰域分析(NNA)方法。注意,他們的研究專門針對小區域估計的情形,空間相關性的范圍限于一個試驗點。
在空間預測中使用共變量的情況下,另一篇值得關注的論文是由Gotway和Hartford(1996)完成的。對糧食產量作為一個共變量的數據,通過應用協克里金法預測土壤硝酸鹽含量,作者提出在空間預測中使用輔助或次要變量。將他們的方法應用于一個試驗點的數據后,他們表明了其方法的效益超過更傳統的外表漂移方法。他們研究的范圍也是僅限于點上的預測。
涉及多地點產量試驗的近期論文之一是由Cullis等人(1998)完成的。在這篇論文中,作者提出了一個多環境初代品種試驗中空間分析的方法。這個方法使用最優線性無偏預測(BLUP)求取基因型效應和環境影響的基因型效應,使用REML求取空間參數和方差分量。不過,該方法對每個試驗單獨模擬協方差結構,也就是不考慮試驗之間的相關性。
在廣闊的陸地區域將地質統計學應用于土壤化學研究的領域中,Yost,Uehara和Fox(1982a,1982b)是農業科學研究的先驅。在兩篇相繼發表的論文中,他們論述了在夏威夷島上進行土壤化學空間預測的研究結果。不過,他們的研究限于只使用克里金方法,而不考慮任何協克里金的方式。此外,在他們的研究中,除了使用選定的土壤化學變量以外,沒有使用其它類型的變量。在土壤科學領域值得一提的另一項類似卻相對更近期的工作是由Ovalles和Collins(1988)完成的。在整個佛羅里達西北部大約覆蓋380公里×100公里的區域中,作者使用通用克里金方法研究選定土壤性質的空間變化,記載的自相關范圍大約為40km。在他們的研究中沒有試圖進行任何空間估計或協克里金。
涉及農作物產量空間預測的論文之一是由Bhatti,Mulla和Frazier(1991)完成的。作者使用一塊面積大約為655米×366米的實驗田,研究小麥產量與土壤有機質和土壤磷含量的關系。預測產量時使用了克里金和協克里金。然而,研究限于單一的實驗田,所以不包含任何模擬大尺度趨勢的方式,而該趨勢在廣闊區域的試驗中往往存在。農作物產量空間分析領域的另一篇論文是由Brownie,Bowman和Burton(1993)完成的。在這篇論文中,作者對比了三種可選的空間方法趨勢分析、最近鄰域Papadakis分析和相關誤差分析,以研究玉米(Zea mays)和大豆(Glycine max)產量數據的空間變化。與農作物產量的其它現有研究的情況相同,他們的研究也是僅在點上的,也就是不考慮各實驗點之間的空間相關性。事實上,作者在其論文的結論中也注意到,雖然分析中結合了來自多點的數據,文獻中還沒有農作物產量數據的跨點空間分析。
有關農作物產量空間分析的其它論文在文獻中也是存在的。Bhatti,Mulla,Kooehler和Gurmani(1991)使用了半變量圖來識別農作物產量的空間自相關性。他們通過研究應用NNA前后的半變量圖,顯示了NNA消除空間可變性的效果。他們研究的范圍也是僅限于點上分析,空間自相關的最大范圍大約是20米。此外,他們的研究中不包括空間估計或預測的方式。Wu等人(1997)就谷類培育試驗中消除產量數據的空間變化,對比了變量誤差一階差分(FD-EV)法(Besag和Kempton,1986)與更傳統的Papadakis最近鄰域法和經典的隨機化完整區塊(RCB)分析。不過,他們的方法不需要預先指定的趨勢和空間自相關結構模型。另外,他們的研究僅限于點上空間變化。
Stroup,Baenziger和Mulitze(1994)使用了若干育種區的數據,對于有效地消除空間可變性產生的干擾,就處理效果比較,對比了傳統的RCB分析、NNA和隨機場線性模型分析(Zimmerman andHarville,1991)。自然,他們的空間相關性研究限于每個育種區之內的數據。
Gotway和Stroup(1997)的一篇論文的獨特之處在于,作者把廣義線性模型理論拓展到包括離散和分類的空間變量的空間估計和預測。他們將其拓展應用于兩個數據集,一個有關作物蟲害,另一個有關雜草計量。不過,與所有其它研究的情況相同,他們研究的范圍及其應用限于每個實驗點之內的數據。
在以上幾段中,回顧了與本發明有關領域中現有的文獻。這些領域中研究工作的總體主要可以劃分為兩個范疇(1)地質統計學和空間統計學應用于土壤科學領域的大區域和小區域試驗,(2)地質統計學和空間統計學應用于農田內小區域試驗的農作物反應分析。現有文獻缺乏在大區域試驗情況下農作物反應分析(如Zea mays的糧食產量)領域的研究,其中空間相關性的范圍要擴展到超出單個的實驗點。
在農作物品種的性能評價中,根本問題是在多種環境下,也就是覆蓋多個試驗點的廣闊地理范圍內得出結論。在文獻綜述中也應當注意,現有文獻的任何研究中都沒有涉及空間估計和空間預測這二者的使用。相反,正如下面將要詳細討論的,本發明涉及解決品種評價問題的一種新方法,它不僅是基于多環境、大區域試驗,它還具有兩個部分來回答兩個明顯的問題(1)一個估計部分,回答的問題是一個品種的長期平均性能或兩個品種的長期平均性能差異,(2)一個預測部分,回答的問題是在若干指定的時間(年度)點上和/或若干指定的空間(地理)點上,一個品種的性能或兩個品種的性能差異。通過通用克里金方法和通用區塊克里金方法,本發明背后的方法考慮了大尺度的趨勢,而且通過協克里金方法容易使用共變量。現有文獻中的工作沒有把以上的全部功能結合在一個統一的方法中,以便研究農作物品種的性能。
以上列舉的問題并非面面俱到,對于目前已知的作物分析技術中有損于有效性的許多問題,只是涉及了一部分。其它值得關注的問題也可能存在;然而,上述問題應當足以表明現有技術中采用的方法總而言之還是不能令人滿意。
發明內容
本發明的若干實施例采用一種稱為線性混合模型的統計模型和地質統計方法,由廣闊區域的試驗數據評價某個農作物品種的性能。評價農作物性能可以通過測量商業上的重要特性,比如(但不限于)產量、糧食濕度或作物倒伏。混合模型能夠使用共變量,比如年度、土壤類型、灌溉等等,作為固定效應,而且能夠使用隨機效應,比如試驗位置,這樣做有助于解釋農作物的性能。固定效應和隨機效應不能解釋的剩余變化采用地質統計方法模擬。地質統計模型考慮了數據中的空間自相關性,能夠獲得有效的置信區間,以評價估計和預測中的不確定性。
本發明的實施例具有兩個明顯的部分空間估計和空間預測。估計部分可用于(1)估計一個品種在一個指定空間位置的平均性能(點估計),然后使用這些點估計結果構造(a)一個品種在廣闊區域(如多國或多州的集合)的平均性能曲面(跨越若干年度),以及(b)兩個品種的性能在廣闊區域(如多國或多州采集)的平均差異曲面(跨越若干年度),(2)估計一個品種在一個規定地理區域比如一個市場區(性能的區塊估計)跨越若干年度的平均性能,以及(3)估計兩個品種之間在一個規定地理區域比如一個市場區(相對性能的區塊估計)跨越若干年度的平均性能差異。
預測部分可用于(1)預測一個品種在一個指定年度一個指定空間位置的平均性能(點預測),然后使用這些點預測結果構造(a)一個品種在廣闊區域(如多國或多州的集合)每一年度的性能預測曲面,(b)兩個品種在廣闊區域(如多國或多州的集合)每一年度的性能差異預測曲面,(2)預測一個品種在一個規定地理區域比如一個市場區(區塊性能預測)每一年度的平均性能,以及(3)預測兩個品種之間在一個規定地理區域比如一個市場區(區塊相對性能預測)每一年度的平均性能差異。
決定選用新培育的、與其它品種相比具有最優性能的品種時,需要如上所述的品種性能估計和預測。對于一個候選品種是否進入商業化狀態需要作出決定,如果進入,將它置于適當的地區,使之與競爭對手相比表現良好。
從營銷的角度看,需要進行一個品種相對性能的定量評價有兩個明顯的理由,首先,為了向市場引入一個新品種而將該品種推進到商業化狀態,其次,為了決定由一個相對性能表現更好的新品種取代一個現有的商業化品種。
以上的決定所依據的最重要的準則之一,是一個品種的期望長期性能優勢。評價品種在市場中的商業價值時需要這一點,其中度量一個品種在市場中的價值要跨越它在市場中的整個壽命。通過評價一個品種在一個指定位置上或一個指定地理區域中的長期相對性能,本發明的估計部分提供了這個問題的答案。
在一個指定地點上或一個指定地理區域中,本發明的預測部分可以評價在一個指定年度,一個品種的性能以及品種之間的性能差異。確定一個品種跨越空間和/或時間情況下性能或相對性能的一致性時,需要這種特定空間和/或特定時間的性能評價。對于一個商業化品種來說,性能的一致性——在作物培育術語中稱為“品種穩定性”——是一個極為必要的屬性。
在廣闊區域的試驗中,本發明優于傳統方法的長處包括以下幾點。第一,它并不要求數據是位置匹配的(數據集內的每個位置都包含兩個品種的觀測結果),而是使用全部試驗點的數據,只要該點至少具有所研究的品種之一。第二,為了進行地區性能評價,不限制它只使用本地區內的數據,而是也使用周圍地區的試驗點上的數據。第三,它易于加入共變量的信息,這些信息與所關注的主要特性有關(例如當產量是所關注的主要特性時,有關土系的數據)。第四,它使用的模型加入年度作為一個模型因素,因此可以容納多年度的數據。第五,它并不假設不同試驗點的觀測結果具有獨立性,而是利用不同試驗點之間的空間相關性以提供改善的統計結論。
本發明的一個實施例中的步驟可歸納如下1.構建農作物在廣闊區域中性能的一個數據庫,包括試驗點的空間坐標、每個試驗點上品種性能的測量值和每個試驗點上有關共變量的數據。
2.選擇兩個品種進行對比,例如,一個正準備推進到商業化狀態的品種“頭領”和一個已經在市場中的品種“競爭者”。
3.通過目測以及采用基于帽子矩陣的統計試驗和Cook距離方法,從兩個品種的數據中排除異常值。
4.為線性混合模型選擇一個空間協方差模型。
5.估計空間協方差參數以及線性混合模型中其它固定效應和隨機效應的參數。
6.使用估計的參數(a)估計跨年度平均值或按年預測,(b)計算在地理區域中連續曲面或區塊平均值,(c)計算每個品種性能特征或性能特征的差異,以及(d)求取估計結果或預測結果的標準差。
7.使用估計結果和預測結果及其標準差評價品種性能,例如評價頭領品種的相對性能,以決定將其推進到商業化狀態。
一方面,本發明是一種在廣闊區域中評價農作物品種性能的方法,它使用的線性混合模型包含地質統計部分,還包括固定效應、隨機效應和共變量的參數。“廣闊區域”僅僅表明一個試驗品種假定的區域。例如,“廣闊區域”可能是指一個試驗品種的適應區域,它可能包括許多試驗點。“農作物品種”是指培育一個指定的植物種類以及業內周知的任何其它用途。構建一個廣闊區域的數據庫,其中包括一個或多個農作物品種的試驗點空間坐標以及一個或多個農作物品種的性能特征值。“數據庫”是指數據的任何集合。例如,一個“數據庫”可能指的是可搜索數據的一個電子形式的集合。“性能特征值”是指與特定品種相關聯的任何受關注的農藝學特征。比如“性能特征值”可能指的是業內周知特征的任何數字。例如,它可能指的是比如Zea mays的糧食產量、糧食濕度和作物倒伏。根據廣闊區域數據庫中的數據擬合線性混合模型,從而估計固定效應、隨機效應和共變量的參數。在一個或多個指定空間位置中的每一個,使用參數估計結果,估計農作物品種的長期期望性能。“估計長期期望性能”意為在多個試驗期間(如年度)估計一個品種農藝學特征的一個數值,以及業內周知的該詞組的正常用途。例如,它可能指的是對于該品種從產生直到它退出市場的整個時間周期,估計一個特征的一個數值,不過它也可以理解為該時間周期可隨估計過程的需要而改變,而且它可能比從產生到退出的時間短得多。
在其它方面,數據庫進一步可能也包括共變量數據。估計長期期望性能可能包括估計該農作物品種與另一個農作物品種之間的長期期望性能差異。估計這些參數的過程可能包括約束最大似然方法。估計長期期望性能的過程可能包括廣義最小平方方法。本方法可能也包括在參數估計之前,使用一種倍率或Cook距離的方法排除數據。本方法可能也包括計算與長期期望性能相關聯的一個標準差。本方法可能也包括形成長期期望性能的一個輸出。該輸出可能包括文本輸出。該輸出可能包括圖形輸出。該圖形輸出可能包括一個等值線圖,表示長期期望性能的連續曲面。
另一方面,本發明是一種在廣闊區域中評價農作物品種性能的方法,它使用的線性混合模型包含地質統計部分,還包括固定效應、隨機效應和共變量的參數。構建一個廣闊區域的數據庫,其中包括一個或多個農作物品種的試驗點空間坐標、試驗點所處的地理區域以及一個或多個農作物品種的性能特征值。根據廣闊區域數據庫中的數據擬合線性混合模型,從而估計固定效應、隨機效應和共變量的參數。在一個或多個指定地理區域中的每一個,使用參數估計結果,估計農作物品種的長期期望性能。
另一方面,本發明是一種在廣闊區域中評價農作物品種性能的方法,它使用的線性混合模型包含地質統計部分,還包括固定效應、隨機效應和共變量的參數。構建一個廣闊區域的數據庫,其中包括一個或多個農作物品種的試驗點空間坐標以及一個或多個農作物品種的性能特征值。根據廣闊區域數據庫中的數據擬合線性混合模型,從而估計固定效應、隨機效應和共變量的參數。在一個或多個指定空間位置中的每一個,在一個或多個指定期間中的每一個,使用參數估計結果,預測農作物品種的平均性能。“預測平均性能”僅僅表明預測一個品種的一個農藝學特征的一個數值以及業內周知的該詞組的正常用途。
在其它方面,估計這些參數的過程可能包括約束最大似然方法。數據庫可能也包括共變量數據。共變量數據可能包括一個反應變量。預測平均性能的過程可能包括通用協克里金方法。共變量數據可能包括固定效應。預測平均性能的過程可能包括通用克里金方法。預測平均性能的過程可能包括預測該農作物品種與另一個農作物品種之間的平均性能差異。本方法可能也包括在參數估計之前,使用一種倍率或Cook距離的方法排除數據。本方法可能也包括計算與預測平均性能相關聯的一個標準差。本方法可能也包括形成預測平均性能的一個輸出。
另一方面,本發明是一種在廣闊區域中評價農作物品種性能的方法,它使用的線性混合模型包含地質統計部分,還包括固定效應、隨機效應和共變量的參數。構建一個廣闊區域的數據庫,其中包括一個或多個農作物品種的試驗點空間坐標、試驗點所處的地理區域以及一個或多個農作物品種的性能特征值。根據廣闊區域數據庫中的數據擬合線性混合模型,從而估計固定效應、隨機效應和共變量的參數。在一個或多個指定地理區域中的每一個,在一個或多個指定期間中的每一個,使用參數估計結果,預測農作物品種的平均性能。
另一方面,本發明是一種雜交品種培育方法,培育一個雜交品種,獲得該雜交品種和對照雜交品種的性能數據。采用廣義最小平方法使一個三次多項式曲面擬合到每個雜交品種的性能數據上,使用球面變量圖模擬殘差方差。新的雜交品種與對照雜交品種的性能進行對比。
另一方面,本發明是一種系統,包括一臺計算機和一種程序。該程序在該計算機上執行并包括以下程序代碼采用廣義最小平方法使一個三次多項式曲面擬合到每個雜交品種的性能數據上;使用球面變量圖模擬殘差方差;以及新的雜交品種與對照雜交品種的性能進行對比。
本公開文件的優點應當理解為,預測、區塊預測、估計和點估計能夠以任何數目的各種排列進行組合,以獲得有價值的性能評價。所有這些組合都不脫離本發明的范圍。
附圖簡要說明本發明采用舉例說明,但不受附圖的限制。本專利的文件包含至少一張彩色圖。專利和商標事務所在收到請求和必需的費用之后,將會提供本專利的拷貝,包括彩色附圖。
圖1是一個等值線圖,顯示1996年兩個商業化雜交品種——DK512和DK527——產量差異的預測曲面。圖中顯示了三條產量差異等值線(-5、0和+5蒲式耳/英畝)和四條顯著級別等值線(產量差異±5%和±20%),也顯示了兩個雜交品種的試驗點以及市場區的邊界和名稱。
圖2是一個等值區域圖,顯示1996年在若干地理區域中的每一個區域內,兩個商業化雜交品種——DK512和DK527——產量差異的區塊預測。等值區域圖上的圖例既表明了產量差異的級別(-3、0和+3蒲式耳/英畝),又表明了差異是否達到5%的顯著級別。
圖3是一個等值區域圖,顯示1996年在若干地理區域中的每一個區域內,兩個商業化雜交品種——DK512和DK527——產量差異的簡單均值。等值區域圖上的圖例既表明了產量差異的級別(-3、0和+3蒲式耳/英畝),又表明了差異是否達到5%的顯著級別。
圖4、圖5、圖6表示的信息分別與圖1、圖2、圖3相同,只是時間為1997年。
圖7、圖8、圖9表示的信息分別與圖1、圖2、圖3相同,只是時間為1998年。
圖10是一個等值線圖,顯示1996年、1997年和1998年兩個商業化雜交品種——DK512和DK527——產量差異的跨年度估計曲面。圖中顯示了三條產量差異等值線(-5、0和+5蒲式耳/英畝)和四條顯著級別等值線(產量差異±5%和±20%),也顯示了兩個雜交品種的試驗點以及市場區的邊界和名稱。
圖11是一個等值區域圖,顯示1996年、1997年和1998年在若干地理區域中的每一個區域內,兩個商業化雜交品種——DK512和DK527——產量差異的跨年度區塊估計。等值區域圖上的圖例既表明了產量差異的級別(-3、0和+3蒲式耳/英畝),又表明了差異是否達到5%的顯著級別。
圖12是一個等值區域圖,顯示1996年、1997年和1998年在若干地理區域中的每一個區域內,兩個商業化雜交品種——DK512和DK527——產量差異的跨年度簡單均值。等值區域圖上的圖例既表明了產量差異的級別(-3、0和+3蒲式耳/英畝),又表明了差異是否達到5%的顯著級別。
示范實施例的說明現在將要詳細介紹在此公開之發明的優選實施例,它的實例在附圖中展示。
本發明涉及以上在背景技術中論述的傳統統計分析方法中的所有五個限制。第一,它并不要求位置匹配的數據,而是使用全部試驗點的數據,只要該點至少具有試驗中的品種之一。第二,為了進行地區性能對比,不限制它只使用本地區內的數據,而是也使用周圍地區的試驗點上的數據。第三,它易于使用共變量。第四,它使用的模型加入年度作為一個模型因素,因此可以容納多年度的數據。第五,它并不假設不同試驗點的觀測結果具有獨立性,而是利用不同試驗點之間的空間相關性以提供改善的統計結論。
期望品種性能數據具有空間相關性是相當自然的,因為品種具有對環境條件的敏感反應。然而,為了對品種性能數據尤其是產量數據中空間相關的程度作出定量評價,對來自Monsanto(DEKALB)玉米研究機構的Zea mays產量數據進行了一項研究。分析了三年中(1994-1996)包括許多Zea mays品種的歷史產量數據。分析結果表明,大約80%品種的數據具有很強的空間自相關性,其證明是擬合的協方差模型中很大的正范圍(Bhatti,Mulla,Kooehler and Gurmani,1991)。這一結果受到一個1994年和1995年Zea mays產量數據混合模型分析的支持,那是一項似然比試驗,分別對每個品種試驗空間自相關的程度。平均說來,75%的品種在其數據中表現出顯著的自相關性。
本發明使用一種統計模型,它既能進行空間估計,又能進行空間預測,這使它在品種性能評價中具有超過傳統方法的某些明顯的優點。空間估計能夠評價品種的長期平均性能,也能夠評價不同品種的長期平均性能差異。這種信息的價值在于確定一個指定的品種相對于其它競爭品種是否具有某種長期的性能優勢。可以在特定的位置確定,也可以在更大的區域中確定。這一點對于評價該品種在市場上的商業價值是必需的,因為度量一個品種的價值是要跨越它在市場上的整個壽命。空間預測能夠在指定的時間點(年度)和/或指定的空間點(地理位置)或聚集在更大區域的空間點上,評價品種的性能和性能差異。在確定一個品種跨越空間和/或時間情況下性能的一致性時,這種特定空間和/或特定時間的性能評價是有用的。對于一個商業化品種來說,性能的一致性——在作物培育術語中稱為品種穩定性——是一個極為必要的屬性。采用品種性能的混合模型分析,可以解決這種穩定性問題(Littell等人,1996,250頁)。
本發明在對比品種性能的方法中引入了地質統計模型。地質統計模型使用從一種地貌的不同位置獲取的樣本點,并由這些樣本點產生(插值出)一個連續曲面。這些樣本點為某種現象的測量結果,比如空氣中的污染物質含量或地質巖芯的礦物級別。據信地質統計學的起源有數個來源。使用協方差的插值——也稱為最優線性無偏預測(BLUP)——來自Wold(1938)、Kolmogorov(1941)、Wiener(1949)、Gandin(1959)和Goldberger(1962)。使用變量圖的插值歸功于Gandin(1959,1963)和Matheron(1962,1969)。Cressie(1990)對地質統計學的起源給出了更多細節。地質統計學最初應用于礦業領域,隨后應用的領域諸如氣象學和環境科學。近年來它已經應用于農業有關的領域,比如土壤性質和農作物產量的模擬。
混合效應模型“混合效應模型”形成了理解本發明的背景,所以應當進行更為詳細的討論,以向讀者提供與在此討論的實施例有關的背景信息。考慮以下普通模型
y=Xθ+Dr+ε, [1]式中y為反應向量,θ為未知固定效應參數的向量,r為未知隨機效應參數的向量,X和D為已知結構的矩陣,ε為未知隨機誤差分量的向量。因為這種模型具有固定和隨機兩種效應參數,所以稱為混合效應模型,或者更簡單地稱為“混合模型”。假設r和ε為正態分布,其參數為Erϵ=00,---[2]]]>以及varrϵ=G00Σ.----[3]]]>在以下段落中,將要說明地質統計模型與這種混合模型的關系。地質統計學中的混合模型地質統計學假設的一個空間線性模型由下式給定Z(s)=μ(s)+ε(s),[4]式中向量s包含空間坐標(如經緯度),表示隨機變量Z(s)的空間位置。[4]式由兩部分組成一個確定部分μ(s)和一個隨機部分ε(s)。對于兩個不同的點s和u,隨機變量ε(s)和ε(u)是空間相關的。因為這種相關性存在于同一隨機過程ε(·)內的各項之間,在術語中稱為“自相關”。為了在統計模擬中進行復制,假設這種空間自相關僅僅依賴于分離s和u的距離和方向,而與它們的確切位置無關。這種性質在術語中稱為“平穩性”。
在地質統計學中最普通的模型是[4]式的一種特殊情況,確定曲面μ(s)不隨位置s而變化。這種模型由下式給定Z(s)=μ+ε(s). [5][4]式和[5]式的模型都可以用于預測,術語“預測”對應于隨機量的推斷(Cressie,1993,15頁)。這里,推斷一詞是指隨機量的預測以及與該預測相關聯的不確定性。在[4]式的模型中,確定曲面μ(s)隨位置變化,對于某個位置s0,Z(s0)的預測在術語中稱為“通用克里金”。在[5]式的模型中,模型的確定部分是一個常數,不隨位置變化,Z(s0)的預測在術語中稱為“普通克里金”。
式和[5]式的模型也可以用于估計,術語“估計”是指固定或確定量的推斷(Cressie,1993,13頁)。這里,推斷是指固定量的估計以及與該估計相關聯的不確定性。注意,在[4]式的模型中估計指的是確定μ(s),而在[5]式的模型中估計指的是確定μ。
以上關于預測和估計的介紹涉及在若干點上的推斷;也就是在特定點上某些數量的預測或估計。預測和估計也可以應用于空間聚集或稱為區域的級別;例如在整個州。在這種情況下,在某個區域A中,我們或者預測隨機變量Z(s)的平均值,或者估計確定趨勢曲面μ(s)的平均值。這兩種運算在術語中分別稱為“區塊預測”Z(A)≡∫AZ(s)ds/|A|和“區塊估計”μ(A)≡∫Aμ(s)ds/|A|,式中|A|≡∫A1ds為A的面積。
式的模型可以寫為一個普通線性模型,z=Xθ+ε [6]表示一個品種所有隨機變量的聯合期望。這里,隨機誤差的方差-協方差矩陣給定為var(ε)=∑。注意,[6]式的模型可視為[1]式模型的一種特殊情況。它可以推廣為容納多變量的情況(例如兩個品種),寫為z1z2=X10102X2θ1θ2+ϵ1ϵ2,----[7]]]>式中var(ϵ)=varϵ1ϵ2≡Σ=Σ11Σ12Σ12Σ22,----[8]]]>它仍然是[1]式模型的一種特殊情況。例如,若是第i個品種有ni個觀測結果,i=1,2,那么[7]式中的向量zi和εi為ni維的向量,i=1,2。矩陣0i,i=1,2是元素為零的矩陣。0i和Xi的行數是ni,i=1,2,而這兩個矩陣的列數和向量θi,i=1,2,的長度取決于模型中使用的參數數目。
在[7]式的模型中,對第i個品種預測Zi(s0),i=1,2時,這個過程稱為簡單情況下的協克里金Xθ=1n10n10n21n2μ1μ2,----[9]]]>式中1ni和0ni分別為ni個1和0的向量,i=1,2。對于比[7]式的模型更一般的情況,這個過程稱為通用協克里金。當聯合預測多個品種時——由向量z(s0)表示,這個過程稱為多變量空間預測(Ver Hoef和Cressie,1993)。
協方差的地質統計模型到目前為止,我們還沒有指定隨機誤差協方差矩陣var(ε)=∑的結構。確實,使本方法成為“地質統計學”的,是使用空間坐標指定矩陣∑的結構。令一個位置的空間坐標包含在2×1的向量s中。令第i個品種(如第一個品種)的第j個觀測結果位于s=(xs,ys),并令第k個品種(如第二個品種)的第m個觀測結果位于u=(xu,yu)。那么第i個品種的第j個觀測結果與第k個品種的第m個觀測結果之間的協方差可以由一個地質統計協方差模型來模擬。這種協方差模型的一個實例是球面模型(Cressie,1993,61頁),由下式給出 式中h=u-s,‖h‖是向量h的長度。本質上‖h‖是從s到u的歐幾里得距離,等于(xs-xu)2+(ys-yu)2]]>。在地質統計學中,參數ik稱為門檻,而ρik稱為射程。如果i=k,[10]式就定義了自相關,否則它就定義了互相關。從物理上說,門檻參數表示方差,而射程表示一個距離,超過后相關性下降至零。
必須注意確保方差-協方差矩陣var(ε)=∑是正定的。在單一品種的情況下,使用[10]式保證了∑是正定的。例如,注意到對于自相關,每個i的φii和ρii都大于零,足以保證正定性。然而對于互相關,要保證正定性,對φii、φkk和φik就有所限制(Goovaerts,1977,108頁)。傳統方法中使用了一類稱為共區域化法的模型,以保證矩陣∑是正定的。Ver Hoef和Barry(1998)已經給出了模擬互相關的其它方式。注意到本實施例中∑12=∑21=0,當使用[10]式時,∑11和∑22中每個i的φii和ρii都大于零,顯而易見[8]式中的方差-協方差矩陣是正定的。
現在我們要介紹一類共區域化法的模型,作為本發明和本實例的細節。對于一個指定的整數m,m=1,2,……,M,讓我們開始為第i個變量構建一個過程如下,i=1,……,K,Yim(s)=1-rim2Wim(s)+rimW0m(s),----[11]]]>式中-1≤rim≤1,對于m=1,2,……,M,Wim是另一個空間過程,對于所有i≠k(i,k=0,1,……,K),具有以下協方差性質cov[Wim(s),Wkm(s+h)]=0,以及對于i=0,1,……,K,
fm(h;γm)≡cov[Wim(s),Wim(s+h)],在這里,函數fm(h;γm)是某種無尺度(fm(0;γm)=1)的模型,作為協方差模型是有效的,也就是所有可能的協方差矩陣都是正定的。例如,按照[10]式定義、但是沒有尺度參數φik的球面模型,作為協方差模型是有效的,如下式給出 式中γm=ρm。對于fm(h;γm),另一個顯而易見的例子是無尺度的“熔核效應”,對于所有h,都有fm(h;γm)=1。熔核效應通常是由于誤差和/或可變性的采樣尺度小于觀測數據集中最小的試驗點距離造成的。大多數地質統計模型對于某種其它的模型都包括一種熔核效應,比如以上[12]式中的球面模型,當h=0時,這種效應造成原點處的不連續。
使用[11]式和fm(h;γm)的定義,我們有 下面,對于i=1,2,……,K,我們通過定義以下求和,構建共區域化法的線性模型,ϵi(s)=Σm=1MφimYim(s),----[14]]]>注意,[14]式的構建提供了一個有效的空間隨機過程。換句話說,對于i=1,2,……,K;m=1,2,……,M,只要分量過程Yim(s)存在,對于i=1,2,……,K,過程εi(s)也就存在。還要注意,對于m=1,2,……,M,過程Yim(s)可能具有不同的基本協方差模型(例如球面的、指數的等等)。
式的構建為構建εi(s)提供了以下協方差模型(共區域化法的模型) 注意,當K=2時,參數r1m和r2m并不能各自確定,但是它們的積是可以確定的,我們簡單記為rm≡r1mr2m。現在只剩下估計以上[15]式協方差模型中的參數,以及固定效應θ。我們在下一節討論這些問題。
約束最大似然對于[1]式中論述的混合模型,[6]式中的地質統計模型是一種特殊情況,必須使用數據來估計模型的參數。一般說來,參數可以分為兩部分[1]式中給出的固定效應θ,以及[8]式給出的∑中包含的協方差參數。
讓我們把∑對參數向量的依賴記為∑()。例如,可能是[12]式中球面協方差模型的參數,或者[15]式中共區域化法模型中的參數。一個普通的估計方法是使用約束(也稱為殘差)最大似然(REML)方法估計θ和。
我們現在簡要介紹REML估計方法,它是由Patterson和Thompson(1971,1974)研究的。假設數據服從一個多變量正態分布N(Xθ,∑()),傳統的最大似然估計方法能夠用于估計θ和,只要在θ和的參數空間中,使雙倍負對數似然,L(θ,)=Nlog(2π)+1og|∑()|+(y-Xθ)′∑()-1(y-Xθ), [16]同時對于θ和達到最小。在的條件下,θ的一個解由下式給出 它是一個眾所周知的、θ的廣義最小平方估計量。把[17]式代入[16]式,協方差參數的似然函數變為 眾所周知,使[18]式中的L()最小化而得出的協方差參數的最大似然估計是有偏的(Mardia和Marshall,1984)。REML估計是最大似然估計的一種改進,它校正了的最大似然估計結果的偏置。對于最小化 +log| X′∑()X|就獲得了REML估計量,式中d為X的秩。[19]式需要一個數值解法。關于空間協方差結構的更多REML細節可以參見Cressie(1993,91頁)。[19]式最小化的計算細節可以參見Harville(1977),Zimmerman(1989)和Wolfinger等人(1994)。雖然研究REML時假設數據服從一個多變量正態分布,Heyde(1994)和Cressie和Lahiri(1996)說明了[19]式也是一個估計方程(Godambe,1991);所以即使數據并不服從一個多變量正態分布,的REML估計也是無偏的。
參考了以上討論,現在就可以詳細討論本發明的實施例了。
模型對于單一品種或者品種對比,為了在集合的點級別或某種區塊級別進行估計和預測,研究了一種統計模型。在品種培育的過程中,品種的對比是最重要的目的,幾乎不可能發生實驗者一次只考慮一個品種或者擬合單品種模型的情況。注意,如果對雜交品種各自的數據擬合單品種模型,就不可能進行品種對比的統計試驗。相反,使用一個雙品種或多品種的模型是最普通的情況,而且在品種對比的同時,也能夠得出對單一品種的推斷(估計結果和預測結果)。所以,詳細考慮這種一般模型。此處考慮的模型是一種統計模型,既包含固定效應,又包含隨機效應,稱為混合效應模型。該模型為y=Xθ+Dr+z+ε. [20][20]式模型中的向量y包含了所關注的主要變量(如糧食產量),可能還有其它伴隨變量。作為這些伴隨變量的一個例子,考慮生長季節中的平均降水量。提供這些降水量數據的氣象臺,處于測量所關注變量的試驗點所在大區域。向y引入這些其它變量的觀測結果的理由,是使用這些變量與預測和估計中所關注的變量之間的協方差。在預測的情況下,讀者可以參考前面有關協克里金的討論。
在一個實施例中,y可能包含恰好兩個品種中一個變量的反應,后面將要如此假設。注意,這個實施例仍然表示向量y的一般形式,因為它能夠使用兩個品種的該變量反應之間的協方差。在這個實施例中,向量y寫為y=(y′1,y′2)′,式中y′i=(Yi1,Yi2,..., ),i=1,2,包括第一個品種的n1個觀測結果和第二個品種的n2個觀測結果。換句話說,一共有來自兩個品種的N≡n1+n2個觀測結果。
式模型中的矩陣X是固定效應的設計矩陣,X=X100X2,----[21]]]>式中 i=1,2。在[22]式中,列向量vij,i=1,2,j=1,2,……,ni包含用于模擬趨勢曲面的固定效應。這些向量包括第j個觀測結果的空間坐標函數,j=1,2,……,ni。例如,第j個觀測結果的x坐標可能是是緯度,而y坐標可能是該觀測結果的經度。也可以使用其它坐標系統,比如通用橫墨卡托(UTM)。如果x坐標和y坐標記為xj和yj,i=1,2,那么對于第i個品種的第j個觀測結果,Vij的實施例包含一個三次多項式內的項,vij=(1,xj,yj,xj2,yj2,xjyj,xj3,yj3,xjyj2,xj2yj)′.----[23]]]>在[22]式中,對于第i個品種的第j個觀測結果,列向量xij,i=1,2,j=1,2,……,ni包含其它固定的共變量效應。這些共變量可能包括某些因素,比如有沒有灌溉、年度的影響等等。所以固定效應參數向量θ可以寫為θ=α1β1α2β2,----[24]]]>式中子向量αi是一個ai×1的固定效應向量,用于多項式v′ijαi形成的、反映第i個品種的空間趨勢曲面,βi是一個bi×1的固定效應向量,反映第i個品種的共變量,i=1,2。令θ的維數為d×1,也就是d=a1+b1+a2+b2。在[23]式給出的實施例中,對于i=1,2,αi是一個10×1的向量。
在[20]式的模型中,r、z和ε是隨機效應的向量。向量r包含位置的隨機效應。如果有nc個位置兩個品種都有觀測結果,那么就只有N-nc個不同的位置,r就變成(N-nc)×1的向量。假設E(r)=0以及var(r)=σ02Ir,----[25]]]>式中Ir為(N-nc)×(N-nc)的單位矩陣。在[20]式的模型中,D為一個N×(N-nc)的關聯矩陣,反映位置的隨機效應。關聯矩陣D的行中所有元素都置為0,每行只有一個1表明一個特定位置與該觀測結果的關聯(例如,參見Searle et al.,1992,138-139頁)。此外,一列中出現兩個1表明在一個位置兩個品種都有。
在[20]式的模型中,z表示一個空間隨機過程。向量z可能包含兩個品種的空間相關隨機變量;z=(z′1,z′2)′式中z′i=(Zi1,Zi2,..., ),i=1,2。假設E(z)=0以及cov(z)=Σ11Σ12Σ21Σ22,----[26]]]>具有某種空間自相關和互相關的變量圖/協方差模式。令一個位置的空間坐標包含在2X1的向量s中。令第i個品種的第j個觀測結果位于s=(xs,ys),并令第k個品種的第m個觀測結果位于u=(xu,yu)。那么第i個品種的第j個觀測結果與第k個品種的第m個觀測結果之間協方差的一個實施例是一個球面協方差模型(Cressie,1993,61頁),由下式給出 式中h=u-s,‖h‖是向量h的長度。本質上‖h‖是從s到u的歐幾里得距離,(xs-xu)2+(ys-yu)2.]]>在地質統計學中,參數φik稱為門檻,而ρik稱為射程。如果i=k,[10]式就定義了自相關,否則它就定義了互相關。
必須注意確保方差-協方差矩陣[26]是正定的。傳統方法中使用了一類稱為共區域化法的模型,以保證這個矩陣是正定的,但是這些方法相當受限。注意到對于自相關,每個i的φii和ρii都大于零,足以保證正定性。然而對于互相關,要保證正定性,對φii、φkk和φik就有所限制(Goovaerts,1977,108頁)。Ver Hoef and Barry(1998)已經給出了模擬互相關的其它方式。注意到一個實施例中∑12=∑21=0,當使用[10]式時,∑11和∑22中每個i的φii和ρii都大于零,顯而易見[26]式中的方差-協方差矩陣是正定的。
在[20]式的模型中,ε是一個包含獨立隨機誤差的向量。向量ε可能包含兩個品種的隨機變量;ε=(ε′1,ε′2)′式中ε′i=[ε′i1,εi2,..., ],i=1,2。假設E(ε)=0以及cov(ϵ)=σ12I100σ22I2,----[28]]]>式中Ii為ni×ni的單位矩陣,i=1,2。還假設r、z和ε相互獨立。
根據[20]、[25]、[26]和[28]式,以及r、z和ε相互獨立的事實,可知反應向量y的方差-協方差矩陣由下式給出Σ=σ02DD′+Σ11Σ12Σ21Σ22+σ12I100σ22I2.----[29]]]>注意,矩陣DD’僅包含0和1,1在對角線上。在DD’中偏離對角線的位置出現1表明在某些位置兩個品種一起進行了試驗。
干擾診斷在擬合[20]式的模型以便得出推斷之前,很重要的是定位有干擾的數據點以及評價它們對模型的沖擊。在進行模型擬合之前,最好排除這些數據點。注意,對于數據的空間分布進行一次目測,即可排除一個點或一簇點,因為在空間上脫離其它觀測結果。在我們的過程中,除了對數據的目測,我們還使用干擾的以下兩種目標測度,以甄別異常或有干擾的觀測結果。第一種干擾測度稱為倍率。對于混合效應模型最一般的情況,Chrestensen et al.(1992)定義了所謂的“帽子”矩陣,H=X(X′∑-1X)-1X′[30]作為倍率的一種測度。如果hjj>2t/N,那么第j個數據點就被視為一個高倍率點,式中hjj為H的第j個對角線元素,t為H的跡,N為數據點的總數目。在一個實施例中,當∑是一個單位矩陣時,帽子矩陣退化為以下形式(Hoaglin and Welsh,1978)
H=X(X′X)-1X′. [31]在這種情況下,如果hjj>2d/N,那么第j個數據點就被視為一個高倍率點(Montgomery and Peck,1992,182頁),正如前面的定義,式中d為θ中參數的數目。在我們的過程中使用的第二種干擾測度為對于第j個觀測結果定義的Cook距離Dj(Cook,1977)。對于混合效應模型最一般的情況,這個距離由下式給出,Dj=(θ^(j)-θ^)′X′Σ-1X(θ^(j)-θ^)d,----[32]]]>式中 為排除第j個觀測結果后θ的一個估計結果,∑為估計的方差-協方差矩陣。Dj>1的觀測結果通常視為有干擾的(Montgomery and Peck,1992,182頁)。在一個實施例中,當∑是一個單位矩陣時,為了甄別有干擾的觀測結果,采用以下表達式計算Cook距離Dj=(θ^(j)-θ^)′X′X(θ^(j)-θ^)dσ^2,---[33]]]>如前面一樣,Dj>1的觀測結果視為有干擾的。在我們的過程中為了進行干擾診斷,我們使用[31]式和[33]式,以及一個沒有任何共變量效應的模型。
為了篩選出有干擾的觀測結果,既使用H又使用Dj的理由如下。H的計算完全根據數據點的空間位置(而不是數值)。所以H對空間的異常敏感,可以用于篩選空間異常。相反,Cook距離Dj則是部分根據空間位置,部分根據觀測值。所以,Dj可以甄別的數據點具有異常的高低值但是不必是空間異常。
混合模型參數的估計介紹了[20]式給出的數據模型之后,下一步應當涉及參數估計的問題。一般說來,這些參數可以劃分為兩個范疇[24]式給出的固定效應θ和[29]式給出的∑中包含的協方差參數。讓我們把∑對于參數≡( φ11,φ12,φ21,φ22,ρ11,ρ12,ρ21,ρ22)′的依賴記為∑()。在一個優選實施例中,約束最大似然(REML)方法用于估計參數θ和。采用REML估計出∑之后,利用廣義最小平方估計θ,由Searle(1971,87頁)給出如下θ^=(X′Σ-1X)-1X′Σ-1y,----[34]]]>正如前面的定義,式中y為全部選定數據點的向量。
參數估計結果用于兩個不同的過程估計和預測。在兩個過程中,所關注的都是或者在單一的位置,或者是一個地理區域的平均。前者稱為點估計(預測),后者稱為區塊估計(預測)。注意,點估計結果也用于構建估計趨勢曲面,也就是對于覆蓋著數據的網格,在每一個網格點上進行點估計。同樣,對于覆蓋著數據空間范圍的網格,在所有網格點上進行點預測,就可以構建預測曲面。
在以下兩個小節中,我們將要考慮通過擬合[20]式的模型以進行估計和預測。注意,在各個區域進行區塊估計或區塊預測之前,使用所有選定的數據點(原始數據集減去有干擾的觀測結果)擬合[20]式的模型只進行了一次。由這個模型估計出的參數然后用于各個點和區域均值的估計和預測。
區塊和點估計假定目標是估計固定空間效應趨勢曲面的均值。對于一個位置s=(x,y)’和第i個品種,i=1,2,讓我們定義一個空間坐標的向量值函數wi(s)=(wi1(s),wi2(s),······,wiai(s))′.---[35]]]>注意,在[23]式給出的實施例中,我們有a1=a2=10和wi(s)=(1,x,y,x2,y2,xy,x3,y3,xy2,x2y)′。所以,在某個位置s的趨勢曲面值由w′i(s)α1給出。不過,往往受到關注的是某個地理區域A的趨勢曲面平均值。如上所述,估計A上的這個平均值稱為區塊估計,對于i=1,2,區塊估計結果由下式給出μi(A)=xi0′βi+1|A|∫Awi′(s)αids,----[36]]]>式中|A|為區域A的面積,對于i=1,2,向量xi0包含一組對應于區域A的共變量數值。假設共變量xi0與趨勢曲面沒有相互作用;指定xi0的數值僅僅在相應的區域內使趨勢曲面升降而不影響其形狀。
在[20]式中模型的一個實施例中,因為某些理由,比如沒有共變量數據或缺少一個適當的共變量,[22]式中的共變量效應xij可以完全從模型中消除。在這種情況下,對應于共變量的x′i0βi,項就會從[36]式中消失。注意,從[36]式的積分中提出x′i0βi項有助于區分共變量與趨勢曲面分量。這些共變量可能是空間的,也可能是非空間的,可能是觀測的,也可能是用戶指定的。它們的數值可能是連續的,也可能是離散的。這些共變量類型的某些實例是有規律的。作為離散非空間共變量的一個實例,考慮年度效應,它不是空間可積的。當一個共變量是空間的,并且數值是連續的,向量xi0中的相應分量就表示一個積分平均結果,與[36]式中的積分相似。對位于s的第i個品種,一個空間共變量的數值記為gi(s),并令此共變量對應于[22]式中向量xij的第k個分量。那么xi0的第k個分量會是 ,式中gi(u)假設是黎曼可積的。作為空間共變量觀測值的一個實例,考慮土系。在所有點上此共變量都必須已知。幸運的是所有主要的農業區域都有土系地圖,并對任何區域(例如縣)都按多種土系進行了分類。這意味著在所有可能的空間位置,都能獲得作為一個共變量的土系數據。作為用戶指定的空間共變量的一個實例,考慮灌溉,它可以作為二進制值函數(1≡灌溉的,2≡非灌溉的)對待。例如,用戶可能關注在100%灌溉情況下的一個區域中,對兩個品種之間的差異進行估計(預測)。如果作為一個共變量,灌溉對應于向量xij的第l個分量,那么這就意味著xi0的第l個分量指定的數值為1表示灌溉條件。
在一個實施例中有恰好兩個品種的一個反應變量,對于區域A中兩個品種的區塊平均,讓我們將其差異定義為μ0(A)。那么μ0(A)由下式給出,μ0(A)≡μ1(A)-μ2(A)=x10′β1-x20′β2+1|A|∫A[w1′(s)α1-w2′(s)α2]ds.---[37]]]>現在,讓我們定義一個集合SA={psp∈A},其中的點{sp}或者是位于覆蓋在A上的規則網格上,或者是在A之內隨機選取的,SA中有nA個元素。讓我們再定義λ1(sp)≡w1(sp)x1000,λ2(sp)≡00w2(sp)x20,andλ0(sp)≡w1(sp)x10-w2(sp)-x20.----[38]]]>那么,近似[36]式和[37]式中的積分,并使用[24]式中θ的定義,我們有μi(A)≈xi0′βi+1nAΣp∈SAwi′(sp)αi=1nAΣp∈SA[λi(sp)]′θ,---[39]]]>i=1,2,以及μ0(A)≈x10′β1-x20′β2+1nAΣp∈SA[w1′(sp)α1-w2′(sp)α2]=1nAΣp∈SA[λ0(sp)]′θ.---[40]]]>在[39]式和[40]式中代入[34]式中θ的廣義最小平方估計量,我們就得到了μi(A)的估計量μ^1(A)≈1nAΣp∈SA[λi(sp)]′(X′Σ-1X)-1X′Σ-1y,---[41]]]>i=0,1,2。 的方差為var(μ^1(A))≈1nA2Σp∈SAΣq∈SA[λi(sp)]′(X′Σ-1X)-1[λi(sp)],---(42)]]>其標準差為se(μ^1(A))=var[μ^i(A)]]]>i=0,1,2。
現在讓我們考慮點估計的情況。目的是估計以下函數,它是[36]式的空間情況,其中區域A現在退化為一個位置(點)。在s0的點估計為μi(s0)=x′i0βi+w′i(s0)αi.[43]正如區塊估計時的情況,共變量xi0的數值僅僅通過使點估計結果或增或減,在[43]式中引入了一個附加的效應。然而,由于這種附加效應的幅度和方向在不同的點可能有變化,模型中有了這些共變量之后,采用點估計結果構建的整個趨勢曲面的形狀可能發生實質性的變化。
作為[41]式中給出的區塊估計量的一種特殊情況,μi(s0)的估計量可以由下式獲得μ^i(s0)=[λi(s0)]′(X′Σ-1X)-1X′Σ-1y,----[44]]]>i=0,1,2。類似于[42]式中的區塊估計結果,可以獲得 標的估計方差,由下式給出,var[μ^1(s0)]=[λi(s0)]′(X′Σ-1X)-1[λi(s0)].----[45]]]>點估計結果的標準差為se[μ^i(s0)]=var[μ^1(s0)]]]>i=0,1,2。
通用區塊和點協克里金讓我們首先考慮區塊預測,它也稱為區塊協克里金,正如在本發明的背景技術中的討論。這里的目標是在某個指定的地理區域,估計反應曲面的平均值,它是以下三個分量的和固定的空間效應的趨勢曲面、固定的協方差效應和自相關的隨機曲面。定義隨機函數Zi(s),其中集合{Zi(s);s∈A};i=1,2形成了一個隨機曲面,可以進行Cressie(1993,107頁)定義的積分。令Z0(s)≡Z1(s)-Z2(s)為一個隨機函數,其中{Z0(s);s∈A}形成了一個隨機曲面,如果Z1(s)和Z2(s)可積,它就可積。現在,對于指定的、xi0,i=1,2給出的共變量數值集合,讓我們定義函數Si(A)≡xi0′βi+1|A|∫A[wi′(s)αi+Zi(s)]ds,----[46]]]>i=1,2,為第i個品種在區域A上平均實現值。正如區塊估計時的情況(參見[36]式后的討論),[46]式中的x′i0β1項可以容納所有類型的共變量,并表示對應于共變量xi0的一個平均的結果。
這里的目的是根據數據的某種線性組合,比如說ω′iy預測區塊均值S1(A)、S2(A)和S0(A)≡S1(A)-S2(A)。向量ωi,i=0,1,2的估計是由一種稱為最優線性無偏預測的方法完成的,預測量是數據的一個線性函數,是無偏的,并且在所有線性無偏預測量中具有最小的均方預測誤差。
讓我們首先假設[20]式給出的線性模型。為了預測區塊均值Si(A),i=0,1,2,通過最小化E[Si(A)-ω′iy]2, [47]可以獲得最優線性無偏預測量(BLUP),此時需要無偏條件E[Si(A)]=E(ω′iy), [48]i=0,1,2。現在,使用[20]式的模型和[46]式,i=1,2時,[48]式可以寫為1|A|∫Awi′(s)αids+x01′β1=ωi′Xθ,----[49]]]>而i=0時,[48]式可以寫為1|A|∫Aw1′(s)α1ds+x01′β1-1|A|∫Aw2′(s)α2ds-x02′β2=ω0′Xθ.---[50]]]>一般說來,使用[35]式中定義的向量值函數wi(s),i=0,1,2,[49]式和[50]式可以結合成單一的方程u′iθ=ω′iXθ, [51]i=0,1,2,式中 為了滿足無偏的條件,對于每一個可能的θ,[51]式都必須成立,因此有u′i=ω′iX,[53]i=0,1,2。[52]式中的向量可以近似為 式中SA已經在前面定義過。現在,使用[20]式、[46]式和[52]式,均方預測誤差[47]式可以寫為E[Si(A)-ωi′y]2=E[ui′θ+1|A|∫AZi(s)ds-ωi′(Xθ+Dr+z+ϵ)]2----[55]]]>i=0,1,2。
讓我們定義N向量δ=(δ11,δ12,......, ,δ21,δ22,......, )′δ=Dr+z+ε,并令Z1(A)≡1|A|∫AZi(s)ds,]]>i=0,1,2。然后,使用[51]式,[55]式可以寫為E[Si(A)-ω′iy]2=E[Zi(A)-ω′iδ]2, [56]i=0,1,2。求取[56]式中的數學期望需要求取cov[Zi(A),δ],它可以近似為 i=0,1,2。現在,由于cov(z,ε)=0和cov(z,r)=0,[57]式中的向量可以寫為 求取[56]式中的數學期望還需要求取Ci(A,A)≡Var[Zi(A)],它可以近似為C1(A,A)≈1nA2Σp∈SAΣq∈SAcov[Zi(sp),Zi(sq)],----[59]]]>i=0,1,2。現在,由于Z1(A)≡1|A|∫AZi(s)ds,]]>它可以近似為1nAΣp∈SAZi(sp),]]>我們注意到E[Zi2(A)]=E[1nAΣp∈SAZi(sp)]2]]>=1nA2Σp∈SAΣq∈SAE[Zi(sp)Z(sp)]]]>=1nA2Σp∈SAΣq∈SAcov[Zi(sp),Z1(sq)]]]>=Ci(A,A)因為對于任何p∈SA以及i=0,1,2,E[Zi(Sp)]=0。所以,[56]式中的均方預測誤差可以寫為E[Zi(A)-ωi′δ]2=E[Zi2(A)-2Zi(A)ωi′δ+ωi′δδ′ωi]]]>=E[Zi2(A)]-2ωi′E[Zi(A)δ]+ωi′E(δδ′)ωi----[60]]]>=Ci(A,A)-2ω′ici+ω′i∑ωi,i=0,1,2。現在,為了滿足[53]式中給出的無偏條件,增加一個拉格朗日乘子mi,[60]式可以寫為L(ωi,mi)≡Ci(A,A)-2ω′ici+ω′i∑ωi-2(u′i-ω′iX)mi, [61]i=0,1,2。[61]式對ωi和mi求偏導,有∂L(ω1,m1)∂ωi′=-2c1+2Σωi+2Xmi,----[62]]]>以及∂L(ωi,mi)∂mi′=-2(ui-X′ωi),----[63]]]>i=0,1,2。
式和[63]式中的這兩個偏導數置為0,然后聯立求解,使[60]式中的均方預測誤差在[53]式中給出的無偏條件下最小化,并產生區域A中品種性能(或性能差異)的最優線性無偏預測量(BLUP)如下。區塊預測(協克里金)方程為ΣXX′0ω1m1=ciui,----[64]]]>i=0,1,2。對方程[64]求解ωi由下式給出,∑ωi+Xmi=ciωi=∑-1(ci-Xmi), [65]i=0,1,2。所以,區域A中品種性能(或性能差異)的預測量由下式給出,S^1(A)=ωi′y,----[66]]]>i=0,1,2。把從方程[65]解出的ωi代入[60]式的表達式中,就提供了均方預測誤差或稱為預測方差的表達式如下,var[S^i(A)]=C1(A,A)-2ωi′ci+ωi′Σ[Σ-1(ci-Xm1)]]]>=Ci(A,A)-2ω′ici+ω′ici-ω′iXmi=Ci(A,A)-ω′ici-ω′iXmi[67]=Ci(A,A)-ω′i(ci+Xmi),i=0,1,2。把從方程[65]解出的ωi代入從矩陣方程[64]式下半部獲得的方程X’ωi=ui中,就提供了拉格朗日乘子mi的表達式如下,X′ωi=ui X′∑-1(ci-Xmi)=ui X′∑-1ci-X′∑-1Xmi=ui[68] mi=(X′∑-1X)-1(X′∑-1ci-ui),i=0,1,2。
我們現在考慮點預測或稱為點協克里金的情況。注意,方程[46]的一種特殊情況是區域A是單一的點,比如說s0,也就是,Si(s0)≡x′i0βi+w′i(s0)αi+Zi(s0),i=1,2。現在的目標是根據數據的某種線性組合,比如說 ,在點s0進行點預測,也就是,預測S1(s0)、S2(s0)和S0(s0)≡S1(s0)-S2(s0)。換句話說,點預測量由下式給出,S^i(s0)=ωi*′y,----[69]]]>i=0,1,2,式中協克里金權向量, 需要確定。與方程[65]類似,權向量 由下式給出,ω1*=Σ-1(ci*-Xmi*),---[70]]]>式中向量 類似于為了[58]式中的區塊預測而定義的向量ci,拉格朗日乘子向量 類似于為了[68]式中的區塊預測而定義的向量mi。所以向量 由下式給出, 以及mi*=(X′Σ-1X)-1(X′Σ-1ci*-ui*),----[72]]]>i=0,1,2。在方程[72]中,向量 ,i=0,1,2類似于為了區塊預測的情況而在方程[52]中給出的向量ui,i=0,1,2,現在定義為 最后,正如在方程[67]中一樣,獲得點預測量 的預測方差為var[S^i(s0)]=Ci(s0,s0)-ωi*′(ci*+Xmi*).---[74]]]>以下包括的實例是為了說明本發明的特定實施例。本領域的技術人員應當理解,在以下的實例中公開的技術是發明人發現的,在本發明的實踐中表現良好,因此可視為構成了其實踐的特定模式。不過,利用本公開文件的利益,本領域的技術人員將會理解,在這些公開的特定實施例中可以進行許多改變,依然可以獲得相象或相似的結果,而不脫離本發明的實質和范圍。閱讀以下實例時應當參考圖1-圖12,在以上的附圖簡要說明中已有小結。
實例一兩個Zea mays品種在廣闊區域試驗中的產量對比我們首先說明本發明如何用于兩個Zea mays品種在廣闊區域的產量對比。這兩個品種稱為DK512和DK527,使用從1996年到1998年三年的廣闊區域試驗數據進行對比。
本過程的第一步包括使用帽子矩陣[31]式和Cook距離[33]式整理數據。對位置也進行了目測,以排除任何異常的空間簇。整理之后,對于DK512和DK527這兩個雜交品種,跨年度的樣本數目分別為n1=445,n2=389。兩個雜交品種的數據合并之后,樣本數目為N=834。按年度的樣本數目由下表給出。
表1
本過程的下一步包括擬合[20]式中給出的混合模型。對于兩個雜交品種,[20]式中的固定效應Xθ都是由[23]式中給出,但是沒有截距項1。因此,[24]式中向量θ的β1和β2分量,每個的維數都是9×1。在向量θ中,α1和α2分量是3×1的子向量,表示1996年到1998年三年中每一年的截距參數。所以,[22]式中的Xii=1,2為ni’12,[21]式中的X為N’24矩陣。對于[26]式中的方差協方差矩陣,我們使用[15]式中給出的共區域化法模型,而且M=2,其中f1(h;γm)和f2(h;γm)都是[12]式中給出的球面模型。[28]式中給出了對每個變量帶有一種單獨的熔核效應的模型。使用的共區域化法模型帶有熔核效應,以及如方程[25]之前的段落所介紹的一種位置匹配的隨機效應,所有數據的方差-協方差模型由[29]式給出。
下一步是估計固定效應參數和協方差參數。使[19]式中給出的REML方程最小化,可以估計方差-協方差參數。估計出的參數由下表給出。
表2
獲得了參數估計結果之后,使用該模型進行兩個雜交品種之間產量差異的空間預測和空間估計。對于單獨的地理位置(點預測和估計)和地理區域(區塊預測和估計)都要進行產量差異的空間預測和空間估計。除了這些空間計算之外,還使用傳統的簡單均值法估計地理區域的平均產量差異。(在這個實例中選用的地理區域表示推銷商業化品種所使用的區域。)回顧一下,估計的目標是確定兩個品種產量差異的長期平均。預測的目標是在特定的年度和/或在特定的位置確定兩個品種的平均產量差異。因為這個理由,空間估計是以跨年度的方式進行的,而空間預測是以按年度的方式進行的。為了提供與傳統均值方法對比的基礎,對按年度和跨年度的分析,都進行了簡單均值分析。(跨年度簡單均值方法是基于所有三年的匯總數據。)使用表2中的參數估計結果和[38]式中的λ0(sp),以及x10=x0=(1/3,1/3,1/3)¢,從[41]式獲得市場區的跨年度平均產量差異的區塊估計結果。(簡單地設置x10=x0=(1/3,1/3,1/3)¢,在估計中為所有三個年度提供了相同的加權系數。)使用[41]式獲得區塊估計結果時,為每個市場區{sp}∈A產生一個致密、規則的點網格,其中A是市場區。另外,從[41]式中估計出的數值,以其標準差——[42]式的平方根——劃分時,假設具有標準正態分布;稱這為z統計。如果z統計小于-1.96或者大于1.96,對于兩個品種之間無差異的虛假設,就假定其顯著程度在α=0.05的級別。
使用表2中的參數估計結果和[54]式中的u0,從[66]式獲得多個市場區按年度平均產量差異的區塊預測。在使用[54]式時,對于1996、1997和1998年,x10和x20都分別設置為(1,0,0)¢、(0,1,0)¢和(0,0,1)¢。對于每個市場區,在使用致密、規則的點網格方面,以及為了確定統計顯著程度而計算和使用z統計方面,區塊預測過程都類似于區塊估計過程。在區塊預測中計算z統計時,從[66]式中預測出的數值,以其標準差——[67]式的平方根——劃分時,假設具有標準正態分布。
下一步我們將考慮兩個品種之間產量差異的點估計和點預測。點估計和點預測都是在25×25的矩形網格上進行的,該網格覆蓋著采集每個雜交品種數據的區域。在25×25網格上的點估計結果和點預測結果共同形成了產量差異的趨勢曲面,它可以用于研究雜交品種試驗的地理區域中產量差異的變化。
在以上討論的25×25網格上,跨年度的點估計使用[44]式計算產量差異μ0(s0)。為了確定統計顯著程度,已經發現,拒絕一個真實的、兩個品種之間無差異的虛假設的概率為2{1-Φ(|μ^0(s0)/var[μ^0(s0)]|)},]]>式中 由[45]式給出,Φ(·)為一個標準正態隨機變量的累積分布函數。
在以上討論的25×25網格上,按年度的點預測使用[69]式計算產量差異S0(s0)。為了確定統計顯著程度,已經發現,拒絕一個真實的、兩個品種之間無差異的虛假設的概率為2{1-Φ(|S^0(s0)/var[S^0(s0)]|)},]]>式中 由[74]式給出,如同前面的情況,Φ(·)為一個標準正態隨機變量的累積分布函數。
以上分析的全部結果以圖形方式展示在圖1至圖12中。這些圖分為兩種類型,一種是等值線圖,另一種是等值區域圖。等值區域圖是一種二維彩色地圖,表示一個變量數值在一個區域的變化。
在等值線圖中,三個級別(-5、0和+5蒲式耳/英畝)的產量差異和四個級別(產量差異±5%和±20%)的統計顯著程度都用等值線表示。在等值線圖中還顯示了兩個雜交品種的試驗點和市場區的邊界和名稱。
在等值線圖中,實等值線表示產量差異,虛等值線表示顯著程度級別。如同任何等值線圖,各條線的關系表示同一變量(如產量差異)的不同級別,用于確定曲面的大體形態。例如,對于產量差異,一個跨越紅色、藍色和綠色實線的橫切面(以此順序),就是在產量差異曲面上向上移動,分別是從-5至0至+5蒲式耳/英畝。
對于產量差異和顯著程度級別而言,等值線之間的關系用于在產量差異和統計顯著程度的多種級別,確定產量差異發生顯著變化的地理區域。例如,一條綠虛線圍繞的區域落在一條綠實線圍繞的區域中,后者又落在一條藍實線圍繞的區域中,表明該區域具有大于+5蒲式耳/英畝的、統計意義上顯著的產量差異(在顯著程度級別0.05)。
另一種類型的圖形,也就是等值區域圖,對于若干地理區域(市場區)中的每一個,顯示這兩個雜交品種在其中的平均產量差異。等值區域圖上的圖例既表明產量差異的級別(小于-5、-5和+5之間以及大于+5蒲式耳/英畝),又表明差異是否達到5%級別的統計顯著程度。
對于空間預測方法(僅僅按年度,1996年到1998年)和空間估計方法(僅僅跨年度)產生等值區域圖和等值線圖。等值線圖不適合簡單均值方法,所以只產生了等值區域圖(既有按年度也有跨年度)。為了保持一致,如果某個市場區沒有一個試驗點種植了兩個雜交品種,該市場區的產量差異就不會展示在等值區域圖上(回顧一下,簡單均值方法需要位置匹配的數據。)為了對比傳統的簡單均值方法和空間方法,簡單均值方法的按年度地圖可以與按年度的區塊預測地圖對比,簡單均值方法的跨年度地圖可以與區塊估計方法的跨年度地圖對比。
實例二有效性研究為了表明本發明在農作物品種性能評價的應用條件下的功效,進行了幾項有效性研究。對于本發明的兩個組件估計和預測,進行了分組的有效性研究。回顧以上討論的實例,這兩個組件中的每一個又包括兩個部分(1)點估計(或預測)部分,用于構建帶有等值線的估計(或預測)曲面地圖;(2)區塊估計(或區塊預測)部分,在實例中通過映射每個區域的平均值,用于構建市場區性能地圖。
我們以對預測的有效性研究開始。首先,我們討論點預測或點克里金的有效性研究。這項研究的基礎是互證有效(Cressie,1993,101頁),使用Monsanto玉米研究系統對于兩個DEKALB玉米雜交品種,DK512和DK527,從1996年到1998年的實際數據。對于兩個雜交品種估計模型的預測質量,互證有效用作一種診斷檢驗。在互證有效的過程中,基本的思路是逐個地每次刪除一個數據,然后使用其它數據預測已排除的數據。對于估計和預測的地質統計學方法,這是一種常見的評價方式。在進行互證有效時,834個數據中的每一個,都是一次排除一個,然后使用我們的過程預測已排除的數據。注意,這里我們將要預測隨機變量Yi(s0);i=1,2,它是從數據集中排除的觀測值之一。我們不想預測光滑值Si(s0);i=1,2,它具有與排除的熔核效應相關聯的誤差項。我們使用了協克里金預測量,使用Cressie(1993,138頁)介紹的協方差。
在互證有效的研究中,計算了以下五種統計量。1.均方預測誤差(MSPE),以量 的平均值計算,其中 是已排除數據的預測值,而Yi(s0)是觀測值。
2.平均估計預測方差(Cressie,1993,140頁,3.2.52式)。
3.80%預測區間包含真實值次數的比例。
4.標準化偏置[Yi(s0)-Y^i(s0)]/se(Y^i(s0)),]]>,(Cressie,1993,102頁,2.6.15式)。
5.標準化MSPE,以量[Y1(s0)-Y^i(s0)]2/var(Y^1(s0)).]]>的平均值的平方根計算。在構建預測區間時,預測區間的下限計算為Y^i(s0)-1.28*se(Y^1(s0)),]]>Y^i(s0)+1.28*se(Y^i(s0)),]]>上限計算為 。這里原文53頁19行2是預測標準差,它計算為估計預測方差(Cressie,1993,140頁,3.2.52式)的平方根。
跨越所有三個年度歸納了以上的統計量。我們現在按照支持本發明的需要和實際服從的規律,展示互證有效研究的結果。
a.MSPE(統計量1)和平均估計預測方差(統計量2)應當相互接近。這兩個統計量的互證有效研究結果分別為269.1和297.9。
b.預測區間(統計量3)應當接近名義(80%)覆蓋。研究中的平均覆蓋為82.7%。
c.標準化偏置(統計量4)應當接近0。研究中的標準化偏置計算結果為-0.021。
d.標準化MSPE(統計量5)應當接近1。研究中的標準化MSPE為1.049。
我們可以看到,互證有效的結果非常接近理論上需要的值。換句話說,該有效性結果支持以下說法作為點預測工具,在[20]式中提出的模型符合我們的期望。
對于本發明中區塊克里金組件的有效性研究,也是基于玉米雜交品種產量的實際數據分析。這種有效性研究的基礎是把這個組件與傳統的方法對比。回顧一下,在傳統的方法中,來自一個區域的所有位置匹配的、并排種植的品種差異數據的簡單均值,用于預測該區域的均值。我們將稱這種方法為均值方法。有效性研究包括的步驟如下1.選擇在國內有充足重疊數據的一對玉米雜交品種。
2.對于這對雜交品種,搜集其全部位置匹配的、并排種植的差異數據。這將被稱為“原始種群數據”。
3.在原始種群數據的空間范圍之內構建一個矩形。我們將稱這個矩形為“預測矩形”,我們將在其中構建平均預測值。
4.使用預測矩形內部的所有原始種群數據的均值表示該矩形的“真實均值”。
5.隨機選擇原始種群數據的50%。這表示可用于分析的“樣本數據”。
6.使用落在預測矩形內部的樣本數據的均值表示用均值方法預測的、該矩形的均值。
7.使用(預測矩形內外的)全部樣本數據,應用本發明的通用區塊克里金組件預測該預測矩形的均值。這表示區塊預測的數值。
8.在步驟4的真實均值與步驟6和步驟7的兩個預測均值之中的每一個之間,計算平方差異。
9.對于幾對雜交品種重復以上所有步驟。
10.對于所有的雜交品種對,使用步驟8的平方差異,對均值方法和區塊克里金方法計算均方預測誤差(MSPE)。
在這項有效性研究中,使用了Monsanto(DEKALB)研究系統的12對玉米雜交品種的數據。研究結果表明,區塊克里金方法的MSPE比均值方法的MSPE低36%。注意,MSPE是預測誤差的一種直接度量,所以,較低的MSPE對應于較高的預測精確度。
下一步,我們討論對本發明的估計組件進行的有效性研究。這項分析是基于一項廣泛的仿真研究,它以一種獨特的方式,既使用真實數據,又使用仿真數據,正如我們下面的介紹。首先回顧一下,本發明的估計組件是基于廣義最小平方(GLS)過程,它是對更傳統的普通最小平方(OLS)過程的一種改進。與GLS過程不同,OLS忽略數據中的空間自相關,僅僅模擬大尺度的趨勢。這項有效性研究的計劃,是要對比三種截然不同的方法(a)GLS方法、(b)OLS方法和(c)傳統的均值方法,在品種對比的情況下,后者僅僅包括有關的位置匹配的、并排種植的兩個品種的差異均值的計算。這項有效性研究的主要目的之一,是要了解先從均值方法到OLS方法,再從OLS方法到GLS方法,是否有增加的利益。
下面提供了這項有效性研究過程各步驟的介紹。
1.(從1994年、1995年或1996年中)選擇一年,再選擇相對成熟程度(RM)接近的一對玉米雜交品種。不考慮位置匹配,搜集這兩個雜交品種的全部研究試驗數據。
2.整理數據,排除空間異常值,并用一個矩形框圍繞整理后的數據。我們將稱之為數據矩形。這將表示可用于分析的樣本數據。
3.在數據矩形中構建另一個內部矩形,在東西方向和南北方向都覆蓋1/3到2/3的數據范圍。我們將命名這個矩形為估計矩形。這就是我們將要進行估計的矩形。
4.使用數據矩形內的全部觀測結果,對每個雜交品種分別擬合一個三次多項式趨勢曲面。擬合這個模型是通過GLS完成的,先擬合一個球面協方差模型(參見背景技術一節)到每個雜交品種的殘差協方差,然后使用兩個品種共同的估計熔核效應。注意,在這兩個雜交品種之間,并未假設相互的協方差。
5.在估計矩形中對兩個品種及其差異(也就是兩個擬合趨勢曲面的差異),計算趨勢曲面均值的估計結果。我們將稱它們為原始估計結果。
6.對于歐幾里得距離0.01和0.05(分別對應于大約40和200英里),計算每個雜交品種的自相關和兩個雜交品種之間的互相關。
7.使用對應趨勢曲面的估計參數和協方差參數以及步驟4的熔核效應,在數據矩形的0.01網格上,對每個雜交品種模擬數據。
8.在x方向和y方向都采用隨機的振幅和頻率產生兩種正弦波。對于第一種正弦波,隨機產生之頻率的均值和標準差分別為100和5,而對于第二種正弦波它們分別是17和3。同樣,對于第一種正弦波,隨機產生之振幅的均值和標準差分別為4和2,而對于第二種正弦波它們分別是3和1。(注意,在模擬的數據加入這些正弦波,意味著引入某種不易被三次模型擬合的系統性變化。現在的數據結果等于三次曲面加正弦波加自相關后的隨機誤差。)9.在估計矩形中,使用(通過步驟4中的GLS擬合出的)三次曲面參數和(步驟8中的)正弦波參數,計算兩個雜交品種均值的“真實”值及其差異。這些真實的均值將用于在以下步驟中計算估計的偏置和置信區間的覆蓋范圍。
10.獨立地對每個雜交品種隨機地去除40%的數據點,在網格點上產生兩個雜交品種的人工不匹配。
11.在估計矩形內部,僅在兩種模擬值都有的網格點上,使用模擬值計算雜交品種的平均值及其差異。這將稱為均值估計結果。
12.對于估計矩形內外的全部樣本值(也就是使用整個數據矩形),通過擬合三次多項式(兩個雜交品種為不同的多項式),以及假設兩個雜交品種的殘差相互獨立但是殘差的方差不同,使用OLS估計兩個雜交品種(在估計矩形上)的均值及其差異。這就提供了OLS估計結果。
13.使用估計矩形內外的全部樣本值(也就是使用整個數據矩形),通過擬合三次多項式(兩個雜交品種的不同)和一個帶有熔核效應和殘差互相關參數的球面協方差結構(兩個雜交品種的不同),使用GLS估計兩個雜交品種(在估計矩形上)的均值及其差異。這就提供了GLS估計結果。
14.對于使用(i)均值估計結果、(ii)OLS估計結果和(iii)GLS估計結果的真實均值及其相應的標準差構建80%置信區間。標注置信區間是否覆蓋步驟9中獲得的真實均值。對以上三個估計結果中的每一個,還要計算對于真實均值的偏置。
15.重復步驟7到步驟14一百次。這是模擬的數目。
16.使用所有100次模擬的估計結果,對于(a)均值的均值估計結果、(b)均值的OLS估計結果和(c)均值的GLS估計結果計算標準差。通過采用步驟14中獲得的偏置(包括所有100次模擬)的平均值,還要計算三個估計結果中每一個的平均偏置。
17.對于均值、OLS和GLS方法中的每一種,計算(在100次模擬中)置信區間覆蓋真實均值次數的百分比。
18.對于覆蓋寬廣的RM范圍的27對雜交品種和對于1994年、1995年和1996年三年的所有數據,重復以上的所有步驟。
注如上所述,注意到這種有效性研究以一個獨特的方式結合了實際數據和模擬結果,首先使用每對玉米雜交品種的實際數據獲得趨勢和相關參數,然后使用這些參數產生模擬研究所用的數據。注意,在這項有效性研究中不能僅僅使用實際數據的理由是,這個過程需要“真實”均值來計算覆蓋范圍的百分比和偏置,如步驟14中的介紹。
我們現在歸納這項有效性研究的結果。
1.以標準差度量提高了精度,從均值方法到GLS方法為37%。這種提高在從均值到OLS和從OLS再到GLS之間等分。
2.GLS方法的80%置信區間覆蓋真實均值次數的百分比(77%)非常接近名義值,而對于均值和OLS方法則比80%低很多(均值為46%,OLS方法為56%)。
3.對于均值和OLS方法,這項研究揭示了真實標準差的總體欠估計,這意味著對無差異的虛假設有太多的假拒絕,也就是品種之間確實沒有差異時,結論卻是品種差異存在。
結論我們在這里介紹過的所有這些有效性研究,共同表明了本發明在農作物品種性能的空間估計和空間預測方面的功效。
參考文獻這里專門引用下列參考文獻和說明書中提及的所有參考文獻作為參考。Besag,J.,and Kempton,R.(1986),Statistical analysis of field experiments usingneighboring plots,Biometrics,42,231-251.Bhatti,A.U.,Mulla,D.J.,and Frazier,B. E.(1991),"Estimation of soil propertiesand wheat yields on complex eroded hills using geostatistics and thematicmapper images,"Remote Sens. Environ.,37,81-191.Bhatti,A. U.,Mulla,D. J.,Koehler,F.E.,and Gurmani,A.H.(1991),"Identifying andremoving spatial correlation from yield experiments,"Soil Sci.Soc.Am.J.,55,1523-1528.Bradley,J. P.,Knittle,K.H.,and Troyer,A.F.(1988),"Statistical methods in seedcorn product selection,"J.Prod. Agric.,1,34-38.Brownie,C.,Bowman,D. T.,and Burton,J.W.(1993),"Estimating spatial variationin analysis of data from yield trialsa comparison of methods,"Agron.J.,85,1244-1253.Cook,R.D.(1977),“Detection of influential observation in linear regression,”Technometrics,19,15-18.Cressie,N.(1990),"The origins of kriging,"Mathematical Geology,22,239-252.Cressie,N.(1993),Statistics for spatial data,Revised Edition,John Wiley andSons,New York,USA.Cressie,N.,and Lahiri,S.N.(1996),“Asymptotics for REML estimation of spatialcovariance parameters,”Journal of Statistical Planning and Inference,50,327-341.Christensen,R.,Pearson,L.M.,and Johnson,W.(1992),“Case-deletion diagnosticsfor mixed models,”Technometrics,34,38-45.Cullis,B.,Gogel,B.,Verbyla,A.,and Thompson,R.(1998),“Spatial analysis ofmulti-environment early generation variety trials,”Biometrics,54,1-18.Gandin,L.S.(1959),“The problem of optimal interpolation,”Trudy GGO,99,67-75.Gandin,L.S.(1963),Objective Analysis of Meteorological Fields,Gidrometeorologicheskoe Izdatel′stvo(GIMIZ),Leningrad (translated by IsraelProgram for Scientific Translations,Jerusalem,1965).Godambe,V.P.(editor)(1991),Estimating functions,Clarendon Press,Oxford.Goldberger,A.S.(1962),“Best linear unbiased prediction in the generalized linearregression model,”Journal of the American Statistical Association,57,369-375.Goovaerts,P.(1997),Geostatistics for Natural Resources Evaluation,OxfordUniversity Press,New York.Gotway,C.A.,and Hartford,A.H.(1996),“Geostatistical methods for incorporatingauxiliary information in the prediction of spatial variables,”Journal ofAgricultural Biological and Environmental Statistics,1,1,17-39.Gotway,C.A.,and Stroup,W.W.(1997),“A generalized linear model approach tospatial data analysis and prediction,”Journal of Agricultural Biological andEnvironmental Statistics,2,2,157-178.Harville,D.A.(1977),“Maximum likelihood approaches to variance componentestimation and to related problems,”Journal of the American StatisticalAssociation,72,320-340.Heyde,C.C.(1994),“A quasi-likelihood approach to the REML estimatingequations,”Statistics and Probability Letters,21,381-384.Hoaglin,D.C.,and Welsch,R.E.(1978),“The hat matrix in regression andANOVA,”American Statistician,32,17-22.Kolmogorov,A.N.(1941),“Interpolation and extrapolation of stationary randomsequences,”Isvestiia Akademii Nauk SSSR,Seriia Matematicheskiia,5,3-14.
(Translation,1926,Memo RM-3090-PR,Rand Corp.,Santa Monica,CA).Littell,R.C.,Milliken,G.A.,Stroup,W.W.,and Wolfinger,R.D.(1996),SASSystem for Linear Models,SAS Institute Inc.,Cary,North Carolina,USA.Mardia,K.V.,and R.J.Marshall.(1984),“Maximum likelihood estimation of modelsfor residual covarianee in spatial regression,”Biometrika,71,135-146.Matheron,G.(1962),Traitéde Géostatistique Appliquée,Tome I. Mémoires duBureau de Recherches Géologiques et Minières,14,Editions Technip,Paris.Matheron,G.(1969),“Le Krigeage Universel,”Cahiers du Centre de MorphologieMathematique,1,Fontainebleau,France.Montgomery,D.C.,and Peck,E.A.(1992),Introduction to Linear RegressionAnalysis,Second Edition,Wiley,New York.Ovalles,F.A.,and Collins,M.E.(1988).“Evaluation of soil variability in northwestFlorida using geostafistics,”Soil Sci.Soc.Am.J.,52,1702-1708.Patterson,H.D.,and R.Thompson.(1971),“Recovery of interblock informationwhen block sizes are unequal,”Biometrika,58,545-554.Patterson,H.D.,and R.Thompson.(1974),“Maximum likelihood estimation ofcomponents of variance,”Pages 197-207 in Proceedings of the 8thInternational Biometric Conference,Biometric Society,Washington,DC,USA.Searle,S.t.(1971),Linear Models,Wiley,New York.Searle,S.R,Casella,G.,and McCulloch,C.E.(1992),Variance Components,Wiley,New York.Sprague,G.F.,and Eberhart,S.A.(1976),“Corn breeding,”Corn and CornImprovement,J.A.Dudley and G.F.Sprague(eds),Iowa State Univ. Press.Stroup,W.W.,Baenziger,P.S.,and Mulitze,D.K.,(1994)“Removing spatialvariation from wheat yield trialsa comparison of methods,”Crop Science,86,62-66.Ver Hoef,J.M.,and Barry,R.D.(1998),“Modeling crossvariograms for cokrigingand multivariable spatial prediction,”Journal of Statistical Planning andInference,69,275-294.Ver Hoef,J.M.,and Cressie,N.(1993),“Multivariable spatial prediction,”Mathematical Geology,25,219-240.Weiner,N.(1949),Extrapolation,Interpolation,and Smoothing of Stationary Time
Series,MIT Press,Cambridge,MA.Wold,H.(1938),A Study in the Analysis of Stationary Time Series,Almqvist andWiksells,Uppsala.Wolfinger,R.D.,Tobias,R.D.,and Sall,J.(1994),“Computing Gaussian likelihoodsand their derivatives for general linear mixed models,”SIAM Journal onScientific Computing,15,1294-1310.Wu,T.,Mather,D.E.,and Dutilleul,P.(1998),“Application of geostatistical andneighbor analyses to data from plant breeding trials,”Crop Science,38,1545-1553.Yost,R.S.,Uehara,G.,and Fox,R.L.(1982a),“Geostatistical analysis of soilchemical properties of large land areas. I. Semi-variograms,”Soil Sci. Soc. Am.
J.,46,1028-1032.Yost,R.S.,Uehara,G.,and Fox,R.L.(1982b),“Geostatistieal analysis of soilchemical properties of large land areas.II. Kriging,”Soil Sci. Soc. Am.J.,46,1033-1037.Zimmerman,D.L.(1989),“Computationally efficient restricted maximum likelihoodestimation of generalized covariance functions,”Mathematical Geology,21,655-672.Zimmerman,D.L.,and Harville,D.A.(1991),“A random field approach to theanalysis of field-plot experiments and other spatial experiments,”Biometrics,47,1,223-239.
權利要求
1.使用線性混合模型評價農作物品種廣闊區域性能的一種方法,該模型內含地質統計部分并包括固定效應、隨機效應和協方差的參數,本方法包括構建一個廣闊區域數據庫,包括一個或多個農作物品種的試驗點空間坐標和一個或多個農作物品種的性能特征值;利用廣闊區域數據庫中的數據擬合線性混合模型,估計固定效應、隨機效應和協方差的參數;以及使用參數估計結果,對于一個或多個指定空間位置中的每一個,估計農作物品種的長期期望性能。
2.根據權利要求1的方法,其特征在于,數據庫進一步包括共變量數據。
3.根據權利要求1的方法,其特征在于,估計長期期望性能包括估計該農作物品種和另一個農作物品種之間的長期期望性能差異。
4.根據權利要求1的方法,其特征在于,估計參數包括一種約束最大似然方法。
5.根據權利要求1的方法,其特征在于,估計長期期望性能包括一種廣義最小平方方法。
6.根據權利要求1的方法,進一步包括在估計參數之前,使用一種倍率或Cook距離的方法,從數據庫中排除數據。
7.根據權利要求1的方法,其特征在于,數據庫進一步包括試驗點所在的地理區域。
8.根據權利要求7的方法,進一步包括使用參數估計結果,對于一個或多個規定的地理區域中的每一個,估計農作物品種的長期期望性能。
9.根據權利要求1的方法,進一步包括計算與長期期望性能相關聯的一種標準差。
10.根據權利要求1的方法,進一步包括形成長期期望性能的一種輸出。
11.根據權利要求10的方法,其特征在于,該輸出包括文本輸出。
12.根據權利要求10的方法,其特征在于,該輸出包括圖形輸出。
13.根據權利要求12的方法,其特征在于,該輸出包括一個等值線圖,表示長期期望性能的連續曲面。
14.根據權利要求1的方法,進一步包括使用參數估計結果,對于一個或多個指定空間位置中的每一個,以及對于一個或多個指定期間的每一個,預測農作物品種的平均性能。
15.使用線性混合模型評價農作物品種廣闊區域性能的一種方法,該模型內含地質統計部分并包括固定效應、隨機效應和協方差的參數,本方法包括構建一個廣闊區域數據庫,包括一個或多個農作物品種的試驗點空間坐標、試驗點所在的地理區域和一個或多個農作物品種的性能特征值;利用廣闊區域數據庫中的數據擬合線性混合模型,估計固定效應、隨機效應和協方差的參數;以及使用參數估計結果,對于一個或多個指定地理區域中的每一個,估計農作物品種的長期期望性能。
16.根據權利要求15的方法,其特征在于,數據庫進一步包括共變量數據。
17.根據權利要求15的方法,其特征在于,估計長期期望性能包括估計該農作物品種和另一個農作物品種之間的長期期望性能差異。
18.根據權利要求15的方法,其特征在于,估計參數包括一種約束最大似然方法。
19.根據權利要求15的方法,其特征在于,估計長期期望性能包括一種廣義最小平方方法。
20.根據權利要求15的方法,進一步包括在估計參數之前,使用一種倍率或Cook距離的方法,從數據庫中排除數據。
21.根據權利要求15的方法,進一步包括計算與長期期望性能相關聯的一種標準差。
22.根據權利要求15的方法,進一步包括形成長期期望性能的一種輸出。
23.根據權利要求22的方法,其特征在于,該輸出包括文本輸出。
24.根據權利要求22的方法,其特征在于,該輸出包括圖形輸出。
25.根據權利要求24的方法,其特征在于,該輸出包括一個等值區域圖,對于一個或多個指定地理區域中的每一個,表示長期期望性能。
26.根據權利要求15的方法,進一步包括使用參數估計結果,對于一個或多個指定地理區域中的每一個,以及對于一個或多個指定期間的每一個,預測農作物品種的平均性能。
27.使用線性混合模型評價農作物品種廣闊區域性能的一種方法,該模型內含地質統計組件并包括固定效應、隨機效應和協方差的參數,本方法包括構建一個廣闊區域數據庫,包括一個或多個農作物品種的試驗點空間坐標和一個或多個農作物品種的性能特征值;利用廣闊區域數據庫中的數據擬合線性混合模型,估計固定效應、隨機效應和協方差的參數;以及使用參數估計結果,對于一個或多個指定空間位置中的每一個,以及對于一個或多個指定期間的每一個,預測農作物品種的平均性能。
28.根據權利要求27的方法,其特征在于,估計參數包括一種約束最大似然方法。
29.根據權利要求27的方法,其特征在于,預測平均性能包括通用克里金方法。
30.根據權利要求27的方法,其特征在于,數據庫進一步包括共變量數據。
31.根據權利要求30的方法,其特征在于,共變量數據包括一種固定效應。
32.根據權利要求31的方法,其特征在于,預測平均性能包括通用克里金方法。
33.根據權利要求30的方法,其特征在于,共變量數據包括一種反應變量。
34.根據權利要求33的方法,其特征在于,預測平均性能包括通用協克里金方法。
35.根據權利要求31的方法,其特征在于,共變量數據進一步包括一種反應變量。
36.根據權利要求35的方法,其特征在于,預測平均性能包括通用協克里金方法。
37.根據權利要求27的方法,其特征在于,預測平均性能包括預測該農作物品種和另一個農作物品種之間的平均性能差異。
38.根據權利要求27的方法,進一步包括在估計參數之前,使用一種倍率或Cook距離的方法,從數據庫中排除數據。
39.根據權利要求27的方法,其特征在于,數據庫進一步包括試驗點所在的地理區域。
40.根據權利要求39的方法,進一步包括使用參數估計結果,對于一個或多個指定地理區域中的每一個,以及對于一個或多個指定期間的每一個,預測農作物品種的平均性能。
41.根據權利要求27的方法,進一步包括計算與預測平均性能相關聯的一種標準差。
42.根據權利要求27的方法,進一步包括形成預測平均性能的一種輸出。
43.根據權利要求42的方法,其特征在于,該輸出包括文本輸出。
44.根據權利要求42的方法,其特征在于,該輸出包括圖形輸出。
45.根據權利要求44的方法,其特征在于,該輸出包括一個等值線圖,表示預測平均性能的連續曲面。
46.使用線性混合模型評價農作物品種廣闊區域性能的一種方法,該模型內含地質統計部分并包括固定效應、隨機效應和協方差的參數,本方法包括構建一個廣闊區域數據庫,包括一個或多個農作物品種的試驗點空間坐標、試驗點所在的地理區域和一個或多個農作物品種的性能特征值;利用廣闊區域數據庫中的數據擬合線性混合模型,估計固定效應、隨機效應和協方差的參數;以及使用參數估計結果,對于一個或多個指定地理區域中的每一個,以及對于一個或多個指定期間的每一個,預測農作物品種的平均性能。
47.根據權利要求46的方法,其特征在于,估計參數包括一種約束最大似然方法。
48.根據權利要求46的方法,其特征在于,預測平均性能包括通用區塊克里金方法。
49.根據權利要求46的方法,其特征在于,數據庫進一步包括共變量數據。
50.根據權利要求49的方法,其特征在于,共變量數據包括一種固定效應。
51.根據權利要求50的方法,其特征在于,預測平均性能包括通用區塊克里金方法。
52.根據權利要求49的方法,其特征在于,共變量數據包括一種反應變量。
53.根據權利要求52的方法,其特征在于,預測平均性能包括通用區塊協克里金方法。
54.根據權利要求50的方法,其特征在于,共變量數據進一步包括一種反應變量。
55.根據權利要求54的方法,其特征在于,預測平均性能包括通用區塊協克里金方法。
56.根據權利要求46的方法,其特征在于,預測平均性能包括預測該農作物品種和另一個農作物品種之間的平均性能差異。
57.根據權利要求46的方法,進一步包括在估計參數之前,使用一種倍率或Cook距離的方法,從數據庫中排除數據。
58.根據權利要求46的方法,進一步包括計算與預測平均性能相關聯的一種標準差。
59.根據權利要求46的方法,進一步包括形成預測平均性能的一種輸出。
60.根據權利要求59的方法,其特征在于,該輸出包括文本輸出。
61.根據權利要求59的方法,其特征在于,該輸出包括圖形輸出。
62.根據權利要求61的方法,其特征在于,該輸出包括一個等值區域圖,對于一個或多個指定地理區域中的每一個,表示預測平均性能。
63.使用線性混合模型評價農作物品種廣闊區域性能的一種方法,該模型內含地質統計部分并包括固定效應、隨機效應和協方差的參數,本方法包括構建一個廣闊區域數據庫,包括共變量數據、一個或多個農作物品種的試驗點空間坐標、試驗點所在的地理區域和一個或多個農作物品種的性能特征值;利用廣闊區域數據庫中的數據擬合線性混合模型,估計固定效應、隨機效應和協方差的參數;使用參數估計結果,對于一個或多個指定空間位置中的每一個,估計農作物品種的長期期望性能;使用參數估計結果,對于一個或多個指定地理區域中的每一個,估計農作物品種的長期期望性能;使用參數估計結果,對于一個或多個指定空間位置中的每一個,以及對于一個或多個指定期間的每一個,預測農作物品種的平均性能;以及使用參數估計結果,對于一個或多個指定地理區域中的每一個,以及對于一個或多個指定期間的每一個,預測農作物品種的平均性能。
64.根據權利要求63的方法,其特征在于,預測平均性能包括預測該農作物品種和另一個農作物品種之間的平均性能差異。
65.根據權利要求63的方法,其特征在于,估計長期期望性能包括一種廣義最小平方方法。
66.根據權利要求63的方法,其特征在于,預測平均性能包括通用協克里金或通用區塊協克里金方法。
67.根據權利要求63的方法,其特征在于,估計參數包括一種約束最大似然方法。
68.根據權利要求63的方法,進一步包括在估計參數之前,使用一種倍率或Cook距離的方法,從數據庫中排除數據。
69.根據權利要求63的方法,其特征在于,估計長期期望性能包括估計該農作物品種和另一個農作物品種之間的長期期望性能差異。
70.根據權利要求63的方法,進一步包括計算與估計長期期望性能或預測平均性能相關聯的一種標準差。
71.根據權利要求63的方法,進一步包括形成估計長期期望性能或預測平均性能的一種輸出。
72.根據權利要求71的方法,其特征在于,該輸出包括文本輸出。
73.根據權利要求71的方法,其特征在于,該輸出包括圖形輸出。
74.雜交品種培育的一種方法,包括(a)培育一個雜交品種;(b)獲得該雜交品種和一個對照雜交品種的性能數據;(c)使用廣義最小平方方法,對每個雜交品種擬合一個三次多項式曲面到性能數據,使用一個球面變量圖模擬殘差方差;以及(d)對比新的和對照雜交品種的性能。
75.根據權利要求74的方法,其特征在于,在球面變量圖模擬中使用一種熔核效應,其特征在于,使用一種矩量過程的方法估計該熔核效應。
76.雜交品種試驗的一種方法,包括(a)使用廣義最小平方方法,對至少兩個雜交品種擬合一個三次多項式曲面到性能數據,使用一個球面變量圖模擬殘差方差;以及(b)對比這兩個雜交品種的性能。
77.根據權利要求76的方法,其特征在于,在球面變量圖模擬中使用一種熔核效應,其特征在于,使用一種矩量過程的方法估計該熔核效應。
78.一種系統,包括一臺計算機;一個程序,在該計算機上運行并包括以下程序代碼使用廣義最小平方方法,對每個雜交品種擬合一個三次多項式曲面到性能數據;使用一個球面變量圖模擬殘差方差;以及對比新的和對照雜交品種的性能。
79.根據權利要求78的系統,包括一個計算機程序,在該計算機上運行并計算一種估計的熔核效應。
全文摘要
使用一個線性混合模型評價一個農作物品種廣闊區域性能的方法和系統,該模型內含地質統計部分并包括固定效應、隨機效應和協方差的參數(J20)。構建一個廣闊區域數據庫,包括共變量數據、一個或多個農作物品種的試驗點空間坐標、試驗點所在的地理區域和一個或多個農作物品種的性能特征值(H25,H30)。利用廣闊區域數據庫中的數據擬合線性混合模型,估計固定效應、隨機效應和協方差的參數。使用參數估計結果,可以估計農作物品種的長期期望性能。使用參數估計結果,對于一個指定的時間階段,可以預測農作物品種的平均性能。
文檔編號G01V9/00GK1390307SQ00815689
公開日2003年1月8日 申請日期2000年10月14日 優先權日1999年10月15日
發明者斯蒂芬·B·斯塔克, 拉達·G·莫漢提, 杰·麥克爾·維霍夫 申請人:德卡爾博遺傳學公司