多主成分玉米區(qū)域試驗(yàn)數(shù)據(jù)的GGE 雙標(biāo)圖分析

王大為，王延波，常程，張書萍，孫甲

（遼寧省農(nóng)業(yè)科學(xué)院 玉米研究所，遼寧 沈陽(yáng)110161）

育種者和農(nóng)戶從不同角度關(guān)心同一個(gè)問題“某地種植什么品種更合適”，這實(shí)際上是一個(gè)關(guān)于品種一般適應(yīng)性和特殊適應(yīng)性的典型問題。GGE 雙標(biāo)圖采用“內(nèi)積原理”解釋主成分分析兩向數(shù)據(jù)表的結(jié)果，在一組圖中展示品種的一般適應(yīng)性和特殊適應(yīng)性，給出品種相似性、環(huán)境相關(guān)性和品種環(huán)境互作效應(yīng)等3 方面建議［1］。前人利用Windows 平臺(tái)的GGEbiplot 軟件或R 版本的GGEBiplotGUI 包 分 析 玉 米［2］、水 稻［3］、小 麥［4］、棉 花［5］、馬鈴薯［6］、小豆［7］、烤煙［8］等產(chǎn)量或品質(zhì)相關(guān)數(shù)據(jù)，但是多數(shù)研究只采用了2 個(gè)主成分。由于區(qū)域試驗(yàn)數(shù)據(jù)影響因素復(fù)雜，往往出現(xiàn)多個(gè)主成分或只有1 個(gè)主成分，僅僅采用2 個(gè)主成分分析會(huì)出現(xiàn)擬合不足或過度擬合的情況。當(dāng)2 個(gè)主成分不能充分體現(xiàn)數(shù)據(jù)規(guī)律的時(shí)候，直接制作雙標(biāo)圖可能產(chǎn)生部分信息誤導(dǎo)，需要把數(shù)據(jù)進(jìn)行預(yù)處理再研究［9］。本研究通過采用以GGEBiplotGUI 為主的R語(yǔ)言軟件包對(duì)2018 年和2019 年遼寧玉米中晚熟區(qū)試產(chǎn)量數(shù)據(jù)進(jìn)行分析，探討多于2 個(gè)主成分的區(qū)試數(shù)據(jù)適宜的GGE 雙標(biāo)圖分析方法，旨在為建立包含多個(gè)主成分試驗(yàn)數(shù)據(jù)的可靠分析方法提供理論依據(jù)。

1 材料和方法

1.1 試驗(yàn)材料

試驗(yàn)于2018 年和2019 年在遼寧玉米區(qū)域試驗(yàn)中晚熟組區(qū)域進(jìn)行，每年參試14 個(gè)品種，設(shè)置在11 個(gè)試驗(yàn)地點(diǎn)，包括阜新、昌圖、鐵嶺、于洪、朝陽(yáng)、丹東、法庫(kù)、鳳城、錦州、東陵、營(yíng)口。2018 年只有8個(gè)試驗(yàn)地點(diǎn)統(tǒng)計(jì)數(shù)據(jù)。為了在雙標(biāo)圖中顯示清晰，品種名采用代號(hào)表示（表1）。

1.2 試驗(yàn)方法

試驗(yàn)均采用隨機(jī)區(qū)組排列，3 次重復(fù)，5 行區(qū)，小區(qū)面積20 m2，種植密度是67 500 株·hm-2。對(duì)照品種是鄭單958。施肥水平、田間管理與大田生產(chǎn)水平一致，每次重復(fù)實(shí)收中間3 行風(fēng)干脫粒，產(chǎn)量按14%含水量計(jì)算。

1.3 統(tǒng)計(jì)分析

數(shù)據(jù)處理和統(tǒng)計(jì)分析使用R 語(yǔ)言及其程序包實(shí)現(xiàn)。其中：方差分析、主成分分析、提取及信息比計(jì)算、誤差同質(zhì)性檢測(cè)、相關(guān)性檢測(cè)、熱圖、GGE雙標(biāo)圖分析分別使用了R 語(yǔ)言3.6.3 版及程序包agricolae、psych、reshape2、PerformanceAnalytics、GGEBiplotGUI。

2 結(jié)果與分析

2.1 中晚熟品種在不同環(huán)境下產(chǎn)量表現(xiàn)

由表2 可 見，2018 年14 個(gè)參 試品種在8 個(gè) 試驗(yàn)點(diǎn)的平均產(chǎn)量變幅為8.69～13.01 t·hm-2。平均產(chǎn)量最高的品種是g2，最低的是g14，而g14 是中晚熟組對(duì)照品種。在不同試驗(yàn)點(diǎn)之間比較，鐵嶺試驗(yàn)點(diǎn)的品種平均產(chǎn)量最高，為11.95 t·hm-2，朝陽(yáng)試驗(yàn)點(diǎn)的品種平均產(chǎn)量最低，為9.56 t·hm-2。用LSD 法進(jìn)行多重比較，品種之間或環(huán)境之間差異均達(dá)到極顯著水平。

表1 遼寧玉米區(qū)域試驗(yàn)中晚熟組品種與代號(hào)Table 1 Codes of varieties of mid-late group in 2018 and 2019 Liaoning maize MET

表2 2018 年參試品種與環(huán)境平均產(chǎn)量與位次Table 2 Yield and rank of varieties and environments in 2018 單位：t·hm-2

由表3 可見，2019 年14 個(gè)參試品種在11 個(gè)試驗(yàn)點(diǎn)的平均產(chǎn)量變幅為8.02～14.17 t·hm-2。平均產(chǎn)量最高的品種是g22，最低的是中晚熟組對(duì)照品種g14。在不同試驗(yàn)點(diǎn)之間比較，昌圖試驗(yàn)點(diǎn)的品種平均產(chǎn)量最高，為13.28 t·hm-2，錦州試驗(yàn)點(diǎn)的品種平均產(chǎn)量最低，為8.73 t·hm-2。品種之間或環(huán)境之間差異均達(dá)到極顯著水平。

表3 2019 年參試品種與環(huán)境平均產(chǎn)量與位次Table 3 Yield and rank of varieties and environments in 2019

2.2 主成分個(gè)數(shù)分析、提取及信息比

采用碎石圖獲取主成分個(gè)數(shù)信息（圖1）。2018年和2019 年都有3 個(gè)特征值大于1，即需要3 個(gè)主成分才能保留樣本中的大部分信息。

擬合度較低的時(shí)候，需要利用信息比（Information Ratio，IR 值）對(duì)特征值排序，其中IR＞1 的特征值含有規(guī)律性信息。2 年試驗(yàn)均是前3 個(gè)特征值IR＞1，需要3 個(gè)主成分來最大擬合原始數(shù)據(jù)，2018年擬合度為81.13%，2019 年為73.85%（表4）。

表4 2018 年和2019 年中晚熟區(qū)試樣本特征值與IR 值Table 4 Eigenvalue and Information Ratio of 2018 and 2019 mid-late MET Samples

2.3 環(huán)境的誤差同質(zhì)性和相關(guān)性檢測(cè)

2.3.1 環(huán)境的誤差同質(zhì)性檢測(cè)

在做雙標(biāo)圖之前，一般需要先檢測(cè)各環(huán)境誤差同質(zhì)性，用來衡量環(huán)境區(qū)分能力。用Bartlett 檢驗(yàn)方法對(duì)各環(huán)境誤差進(jìn)行同質(zhì)性檢測(cè)。結(jié)果表明，2018年，P=0.249 9＞0.05，方差具有齊性，無需定標(biāo)即可分析。2019 年，P=0.013 9＜0.05，初步認(rèn)為方差不具有齊性。又經(jīng)Levene、Fligner-Killeen 檢驗(yàn)，判斷2019 年區(qū)試數(shù)據(jù)呈偏態(tài)分布，需要定標(biāo)后再進(jìn)行分析。

2.3.2 環(huán)境相關(guān)性檢測(cè)

在特征值計(jì)算過程中，同時(shí)得到成分載荷，即觀測(cè)變量與主成分的相關(guān)系數(shù)。成分載荷可幫助解釋主成分的含義。2018 年第一主成分與昌圖、朝陽(yáng)、營(yíng)口、阜新弱相關(guān)（R≥0.70），第二主成分與鐵嶺、丹東、東陵強(qiáng)相關(guān)（R＞0.80），第三主成分與于洪強(qiáng)相關(guān)（R＞0.95），與阜新存 在弱相關(guān)（R＞0.55）。2019 年第一主成分與阜新、昌圖、朝陽(yáng)、法庫(kù)、錦州、營(yíng)口弱相關(guān)（R＞0.60），第二主成分與丹東、鳳城、東陵強(qiáng)相關(guān)（R＞0.80），第三主成分與鐵嶺、于洪弱相關(guān)（R＞0.70）。

2.3.3 環(huán)境聚類分析熱圖可視化

經(jīng)過聚類分析并制作heatmap 熱圖（圖2），2018年，于洪、鐵嶺、丹東、東陵聚為一類，朝陽(yáng)、阜新、營(yíng)口、昌圖聚為另一類。2019 年，東陵、丹東、鳳城聚為一類，鐵嶺、于洪、朝陽(yáng)、阜新、昌圖、法庫(kù)、錦州、營(yíng)口聚為另一類。

圖1 利用碎石圖判斷主成分個(gè)數(shù)Fig.1 Judging principal component numbers using scree plot

圖2 環(huán)境聚類分析熱圖Fig.2 Clustering analysis heatmap of environments

2.4 2 種分析方法的雙標(biāo)圖比較

2.4.1 分組分析雙標(biāo)圖

利用前2 個(gè)主成分先做出初始雙標(biāo)圖，參照各區(qū)垂線劃分顯示的主要模式，結(jié)合前述聚類分析熱圖和環(huán)境相關(guān)系數(shù)，把試驗(yàn)點(diǎn)分組進(jìn)行研究。

由圖3A、圖3B 可見，2018 年環(huán)境分為2 組：于洪、鐵嶺、丹東、東陵為一組，簡(jiǎn)稱：2018 遼東組；朝陽(yáng)、阜新、營(yíng)口、昌圖為一組，簡(jiǎn)稱：2018 遼西組。由圖3C、圖3D 可見，2019 年環(huán)境 也分為2 組：東陵、丹東、鳳城為一組，簡(jiǎn)稱2019 遼東組；鐵嶺、于洪、朝陽(yáng)、阜新、昌圖、法庫(kù)、錦州、營(yíng)口為一組，簡(jiǎn)稱2019 遼西組。分組后2018 年2 組雙標(biāo)圖擬合度是87.39%和83.07%，2019 年2 組雙標(biāo)圖擬合度是91.02%和68.62%。

圖3 GGE 雙標(biāo)圖的“哪個(gè)贏在哪里”功能形態(tài)Fig.3 Which-won-where view of the GGE biplot

在2018 遼東組中，丹東、鐵嶺被劃為一小組，最高產(chǎn)品種為g1；東陵為一小組，最高產(chǎn)品種為g8；于洪為一小組，最高產(chǎn)品種為g13。在2018 遼西組中，昌圖、朝陽(yáng)、阜新、營(yíng)口被劃為同一小組內(nèi)，最高產(chǎn)品種g2。

在2019 遼東組中，環(huán)境因素幾乎掩蓋了品種因素，丹東和鳳城小組最高產(chǎn)品種是g15，東陵最高產(chǎn)品種是g22。在2019 遼西組中，昌圖、錦州、法庫(kù)、營(yíng)口、鐵嶺被劃為一小組，g12 是最高產(chǎn)品種。g21 是阜新和朝陽(yáng)小組最高產(chǎn)品種，g18 是于洪最高產(chǎn)品種。

由圖4A、圖4B 可見，在2018 遼東組中，東陵的代表性較好，于洪與其它地點(diǎn)相關(guān)性較小，丹東的區(qū)分能力和代表性均低。在2018 遼西組中，昌圖、阜新的區(qū)分力好，阜新與其它地點(diǎn)相關(guān)性較小，昌圖代表性高。由圖4C、圖4D 可見，在2019遼東組中，鳳城與丹東相關(guān)性強(qiáng)，東陵與鳳城、丹東相關(guān)性弱，鳳城的代表性較強(qiáng)。在2019 遼西組中，法庫(kù)、錦州、營(yíng)口的代表性和區(qū)分力均較強(qiáng)，阜新與朝陽(yáng)相關(guān)性強(qiáng)，于洪與鐵嶺相關(guān)性強(qiáng)，阜新、朝陽(yáng)與于洪、鐵嶺之間出現(xiàn)負(fù)相關(guān)。

由圖5A、圖5B 可見，在2018 遼東組中，g9 產(chǎn)量表現(xiàn)最好，最不穩(wěn)定品種是g1，比較穩(wěn)定的是g8。在2018 遼西組中，g2 平均產(chǎn)量最高，后面依次是g12、g3，最不穩(wěn)定的是g13；g3 是比較高產(chǎn)穩(wěn)產(chǎn)的品種。由圖5C、圖5D 可見，在2019 遼東組中，最高產(chǎn)品種是g15，后面依次是g22、g24；各品種均勻分布在環(huán)境平均值附近，穩(wěn)定性相似度較高。在2019 遼西組中，g12 平均產(chǎn)量最高且穩(wěn)定，后面依次是g11、g22；最不穩(wěn)定品種是g24。

圖4 GGE 雙標(biāo)圖的“區(qū)分力與代表性”功能形態(tài)Fig.4 Discrimination and representativeness view of the GGE biplot

2.4.2 部分高產(chǎn)品種分析雙標(biāo)圖

由圖6A 可見，2018 年 刪 除5 個(gè)低產(chǎn) 品 種g6、g7、g8、g13、g14，擬合度下降為68.78%。8 個(gè)地點(diǎn)被分為2 小組，昌圖、阜新、于洪、朝陽(yáng)分為一個(gè)小組，營(yíng)口、東陵、鐵嶺、丹東分為一個(gè)小組。各組名義上高產(chǎn)品種分別為g2、g12 和g1。

由圖6B 可見，2019 年刪除5 個(gè)低產(chǎn)品種g16、g19、g23、g13、g14，擬合度下降為66.56%。11 個(gè)地點(diǎn)被分為3 組，營(yíng)口、錦州、阜新、昌圖、法庫(kù)、朝陽(yáng)、鐵嶺分為一個(gè)小組，東陵、鳳城、丹東分為一個(gè)小組，于洪單獨(dú)為一個(gè)小組。各組名義上最高產(chǎn)品種分別為g12、g22、g15。

2.4.3 同時(shí)參加兩年區(qū)試的4 個(gè)品種雙標(biāo)圖

采用同時(shí)參加了兩年區(qū)試的g11、g12、g13、g14 等4 個(gè)品種和8 個(gè)重復(fù)性環(huán)境進(jìn)行了年際間重復(fù)GGE 雙標(biāo)圖分析。由圖7 可見，g12 和g11 在兩年區(qū)試大部分地點(diǎn)表現(xiàn)高產(chǎn)。8 個(gè)試驗(yàn)環(huán)境中，營(yíng)口、阜新、鐵嶺等3 個(gè)地點(diǎn)年際間向量夾角較小，重復(fù)性較強(qiáng)。

圖5 GGE 雙標(biāo)圖的“高產(chǎn)性與穩(wěn)定性”功能形態(tài)Fig.5 Mean vs.stability view of the GGE biplot

圖6 只包含較高產(chǎn)品種的“哪個(gè)贏在哪里”功能形態(tài)Fig.6 Which-won-where view of GGE biplot excluding some low yielding genotypes

圖7 4 個(gè)品種2 年區(qū)試GGE 雙標(biāo)圖功能形態(tài)Fig.7 GGE biplot for 4 genotypes in 2 years MET

3 討論與結(jié)論

3.1 品種適宜區(qū)域推薦

GGE 雙標(biāo)圖中“哪個(gè)贏在哪里”圖按照品種-環(huán)境的相互關(guān)系以發(fā)自環(huán)境中心點(diǎn)的矢量垂線來對(duì)區(qū)試點(diǎn)分組，位于各組內(nèi)多邊形頂點(diǎn)即為該組最高產(chǎn)品種?！案弋a(chǎn)性與穩(wěn)定性”圖中，平均環(huán)境軸所指的方向是品種在所有環(huán)境下的近似平均產(chǎn)量走向，品種與平均環(huán)境軸之間的垂線表示品種的穩(wěn)定性。綜合產(chǎn)量和穩(wěn)定性指標(biāo)，2018 年既高產(chǎn)又穩(wěn)產(chǎn)的品種是遼東組的g9、g8，遼西組的g2、g12。連續(xù)2 年在遼西組表現(xiàn)較高產(chǎn)且穩(wěn)定的是g12。去除數(shù)據(jù)中的低產(chǎn)品種之后，可以對(duì)相對(duì)高產(chǎn)品種的適宜環(huán)境給出建議。g12 在部分高產(chǎn)品種的分析中表現(xiàn)突出，尤其在阜新、朝陽(yáng)、昌圖試驗(yàn)點(diǎn)表現(xiàn)出較高的G+GE 效應(yīng)。g12（東1775）于2020 年初通過了省級(jí)品種審定，適宜在阜新、朝陽(yáng)、昌圖等類似生態(tài)區(qū)種植

3.2 試驗(yàn)點(diǎn)綜合評(píng)價(jià)

GGE 雙標(biāo)圖“區(qū)分力與代表性”圖中，區(qū)試點(diǎn)線段和平均環(huán)境軸的夾角度數(shù)是其對(duì)目標(biāo)環(huán)境代表性的度量，矢量長(zhǎng)度則是環(huán)境區(qū)分力的度量。2018 年，阜新、于洪區(qū)分力強(qiáng)但代表性弱，可以發(fā)揮淘汰品種的功能；昌圖是既有區(qū)分力又有代表性的區(qū)試點(diǎn)。2019 年，東陵區(qū)分力較強(qiáng)但代表性弱，可以發(fā)揮淘汰品種的功能；法庫(kù)、錦州、營(yíng)口是既有區(qū)分力又有代表性的區(qū)試點(diǎn)。

3.3 分組分析的優(yōu)點(diǎn)與注意事項(xiàng)

在區(qū)試數(shù)據(jù)存在多個(gè)主成分的情況下，如果不分組，使用只有2 個(gè)主成分的初始GGE 雙標(biāo)圖分析，所展示的主要模式會(huì)被隨機(jī)誤差所掩蓋或扭曲。就本試驗(yàn)數(shù)據(jù)而言，會(huì)出現(xiàn)3 種弊端：第一，在可信度方面，擬合度較低。2018 年為62.94%，2019 年為62.89%。第二，在品種推薦方面，只能顯示出g2 和g9 等并不突出的品種。第三，在環(huán)境評(píng)價(jià)方面，2018 年于洪與東陵矢量線幾乎重疊，會(huì)得出可以去掉其中一個(gè)試驗(yàn)點(diǎn)的錯(cuò)誤結(jié)論。分組后擬合度提高，圖中所顯示的規(guī)律更接近真實(shí)情況。分組分析可以更多地保留數(shù)據(jù)原始信息，是解決數(shù)據(jù)較大、關(guān)系較復(fù)雜的有效途徑，是針對(duì)多主成分?jǐn)?shù)據(jù)的適宜分析方法。但是由于區(qū)域試驗(yàn)點(diǎn)往往數(shù)量有限，也不宜過度細(xì)分組別，否則得到的GGE 雙標(biāo)圖內(nèi)容含義不清。建議參照聚類分析和熱圖結(jié)果，選取盡量少的組別進(jìn)行分組分析。

3.4 多年重復(fù)試驗(yàn)分析的意義

鑒于區(qū)域試驗(yàn)的參試規(guī)則特殊性，難以獲得多個(gè)品種連續(xù)多年的重復(fù)試驗(yàn)數(shù)據(jù)。能夠?qū)ι贁?shù)重復(fù)參加2 年區(qū)試的品種進(jìn)行年際間GGE 雙標(biāo)圖分析很有意義。在本試驗(yàn)中，g12 和g11 表現(xiàn)年際間穩(wěn)定高產(chǎn)，營(yíng)口、阜新、鐵嶺這3 個(gè)地點(diǎn)在2 年區(qū)試中表現(xiàn)出較大的環(huán)境穩(wěn)定性。品種與環(huán)境間的GE 互作在年際間具有一定的重復(fù)性，這為遼寧玉米劃分生態(tài)區(qū)、進(jìn)而針對(duì)性地選擇區(qū)域特殊適應(yīng)性品種提供了參考。

3.5 數(shù)據(jù)聚集的特殊情況

2019 年遼東組各品種靠近環(huán)境中心位置集中分布，難以區(qū)分。究其原因，雙標(biāo)圖是一種降維分析方法，要么該組數(shù)據(jù)變異接近于環(huán)境平均值，要么該組數(shù)據(jù)變異位于其它維度。從本試驗(yàn)現(xiàn)狀分析，可能是前者，即遼東組各試驗(yàn)點(diǎn)環(huán)境與參試品種的GE 互作效應(yīng)較小。

GGE 雙標(biāo)圖在主效可加互作可乘模型（Additive main multiplicative interaction，AMMI）考 察G×E 效應(yīng)的基礎(chǔ)上又更多考慮了基因G 的作用［10，11］。GGE 雙標(biāo)圖主要用于直觀分析基因在不同環(huán)境下表達(dá)的相對(duì)水平而不是絕對(duì)水平，原則上應(yīng)當(dāng)把其當(dāng)作一個(gè)建立假說的工具而不是做出最終結(jié)論的工具［1，12］。本試驗(yàn)是對(duì)2018 年和2019年遼寧玉米中晚熟組區(qū)試產(chǎn)量性狀的分析，僅代表1 種熟期類型品種產(chǎn)量對(duì)2 年光溫反應(yīng)的差異。品種與區(qū)試點(diǎn)推薦結(jié)果只具有限意義，不宜擴(kuò)展到其它熟期或者多年應(yīng)用。今后需要繼續(xù)深入對(duì)不同熟期類型品種產(chǎn)量、品質(zhì)、抗性等多種性狀和多年區(qū)試結(jié)果進(jìn)行分析，從中挖掘出更多具有規(guī)律性的信息，供育種者和品種推廣者參考。