利用SAS軟件快速方差分析多組資料

宋 釗, 張白鴿, 曹 健

(廣東省農(nóng)業(yè)科學(xué)院蔬菜研究所 廣東省蔬菜新技術(shù)研究重點實驗室，廣州 510640)

SAS的英文全名是Statistical Analysis System，即統(tǒng)計分析系統(tǒng)，從20世紀(jì)60年代末期由美國北卡羅來納州立大學(xué)兩位教授開始發(fā)展的，是當(dāng)今國際上最著名的數(shù)據(jù)分析軟件之一，目前已經(jīng)與SPSS(Statistical product and service solutions)和BMDP(Bio Medical Data Processing)并駕齊驅(qū)[1]，在國際上被譽為數(shù)據(jù)處理和統(tǒng)計分析的標(biāo)準(zhǔn)軟件系統(tǒng)，在國內(nèi)廣泛用于農(nóng)業(yè)試驗統(tǒng)計分析[2]。

對于SAS，SPSS和EXCEL等幾個數(shù)據(jù)統(tǒng)計軟件的使用已經(jīng)有多篇報道從各方面進行了較多研究[3-5]，比較各個軟件在實際使用過程的優(yōu)缺點。由于SAS系統(tǒng)為全英文界面，對于英語基礎(chǔ)比較薄弱的使用者存在一定困難，并且程序模塊化，對于沒有計算機程序語言編程基礎(chǔ)的人員來說不易掌握[5]，所以很長一段時間內(nèi)國內(nèi)研究者在進行單因素或者多因素方差分析時大多使用在農(nóng)業(yè)生物統(tǒng)計中得到廣泛應(yīng)用的Excel[6-7]。但是Excel只能給出方差分析表，不能進行平均數(shù)的多重比較和字母標(biāo)記差異顯著性結(jié)果[5]，不能同時對多組觀察值一次性全部進行方差分析，如果樣本數(shù)較多將費時費力并且容易出錯。本文介紹利用SAS統(tǒng)計軟件同時對20組辣椒品種在3個不同脅迫條件下辣椒葉片葉綠素含量數(shù)據(jù)進行方差分析，檢測同一個品種在不同澇漬脅迫條件葉片葉綠素含量差異顯著性，同時對20個品種間的葉綠素含量差異進行比較。

1 材料與方法

1.1 實驗材料

總計20份辣椒材料，為辣椒屬一年生種(CapsicumannuumL.)，按照果實形狀可以分為尖椒、甜椒、美人椒和線椒4個類型，均為在華南地區(qū)栽培較廣，具有一定代表性的栽培種。

1.2 實驗方法

試驗在廣州市天河區(qū)五山路廣東省農(nóng)業(yè)科學(xué)院蔬菜研究所旁邊科研基地大棚中進行，2015年4月15號播種，5月12間苗，5月15號開始脅迫處理。試驗設(shè)置1個對照和澇漬脅迫3、6和9 d等3個處理，然后恢復(fù)3 d，每個處理3個重復(fù)，每個重復(fù)每個品種10株，穴盤1孔放1株，品種隨機排序。對照進行正常田間管理，缺水后即澆水。澇漬脅迫方法參照尹冬梅等[8]的方法加以改進，水層高于辣椒植株根部2～3 cm，用帶體積刻度的量杯倒入等體積水量開始脅迫處理。

1.3 葉綠素測定

使用日本柯尼卡美能達(Konica Minolta)SPAD-502 PLUS葉綠素?zé)晒鈨x對葉片的葉綠素含量進行測定，該儀器通過測量葉子對兩個波長段里的吸收率，來評估當(dāng)前葉子中的葉綠素的相對含量，用SPAD值表示目前葉子中葉綠素含量相對應(yīng)的參數(shù)。

1.4 數(shù)據(jù)分析軟件和方法

The SAS System for Windows Version 8.01 TS Level 01M0；SAS Institute，Cary，NC，USA(美國北卡羅來納州SAS研究所出品)。本試驗為平衡試驗設(shè)計，采用過程步PROC的ANOVA進行方差分析，同時進行Duncan′s multiple-range test多重比較。

2 結(jié)果與分析

2.1 辣椒葉片葉綠素測量結(jié)果

表1為澇漬脅迫處理后對照和3個處理的植株葉片葉綠素SPAD值，每個品種每個處理測量15個數(shù)據(jù)，20個品種4組數(shù)據(jù)總計1200個數(shù)據(jù)。限于篇幅表1每個品種只列出了5個重復(fù)數(shù)據(jù)，省去了小數(shù)點后面數(shù)值。

表1 20份辣椒材料澇漬脅迫下葉片葉綠素?zé)晒釹PAD值Table 1 Chlorophyll SPAD values of leaves of 20 pepper varieties under waterlogging stress

2.2 比較每個品種處理間的差異顯著性

2.2.1 建立SAS數(shù)據(jù)集

如表1數(shù)據(jù)在SAS系統(tǒng)中是不能直接使用的，需要轉(zhuǎn)換成SAS的格式，也就是建立數(shù)據(jù)集，我們通過SAS語言的數(shù)據(jù)步DATA中的INPUT函數(shù)和CARDS選項建立數(shù)據(jù)集。

核心代碼及說明如下：

DATA Work.yelvsu; /*輸入20個辣椒品種葉綠素數(shù)據(jù)*/

DO trt=1 TO 4; /*4個處理數(shù)序號，保留*/

DO t=1 TO 15;

INPUT brd1-brd20@@;

OUTPUT;

END;

END;

DROP t; /*15個重復(fù)數(shù)序號，丟棄*/

CARDS;

38 32 32 42 37 34 29 32 32 38 37 37 35 44 40 35 28 30 37 30

40 35 33 38 40 36 37 29 39 38 36 36 38 41 42 39 31 33 37 34

40 34 35 39 32 36 27 27 43 34 39 35 37 38 35 39 38 33 33 33

39 39 29 35 37 39 31 29 36 33 42 36 34 36 38 39 39 37 34 31

41 36 38 37 32 33 33 29 31 25 38 37 39 37 37 40 28 28 36 34

……(省略部分類似數(shù)據(jù))

;

RUN;

“DATA Work.yelvsu;”語句建立臨時庫Work的數(shù)據(jù)集yelvsu；“DO trt=1 TO 4;”語句表示在最外層進行4次循環(huán)，因為有4組數(shù)據(jù)；“DO t=1 TO 15;”語句表示連續(xù)寫入15個重復(fù)值；“DROP t;”語句則表示重復(fù)次數(shù)1到15不需要寫入數(shù)據(jù)集，只寫入相對應(yīng)的SPAD值；“INPUT brd1-brd20@@;”語句表示一行寫入20個品種的SPAD值；“OUTPUT;”語句表示將過程步中的數(shù)據(jù)輸出到數(shù)據(jù)集中保存起來；“CARDS;”語句則列出INPUT所有變量對應(yīng)的數(shù)據(jù)，以“;”表示數(shù)據(jù)的結(jié)束。執(zhí)行之后自動建立trt和brd1,brd2,brd3…,brd20等21列數(shù)據(jù)。只有當(dāng)原始數(shù)據(jù)排列格式如表1時，使用上述SAS程序語句就可以直接建立適合SAS進行統(tǒng)計分析的數(shù)據(jù)集。

2.2.2 對20個辣椒品種同時進行方差分析

在Excel中只能對20組辣椒澇漬脅迫實驗數(shù)據(jù)逐個進行方差分析，在SAS中可以一次性對20組數(shù)據(jù)同時進行分析。使用SAS的過程步PROC對數(shù)據(jù)進行分析，同時進行DUNCAN多重比較，顯著性水平0.05。

核心代碼及注釋如下：

PROC ANOVA data=WORK.Yelvsu; /*調(diào)用上一步建立的數(shù)據(jù)集Yelvsu*/

class TRT;

model BRD1 BRD2 BRD3 BRD4 BRD5 BRD6 BRD7 BRD8 BRD9 BRD10 BRD11 BRD12 BRD13

BRD14 BRD15 BRD16 BRD17 BRD18 BRD19 BRD20=TRT; /*對20個品種進行方差分析*/

means TRT / HOVTEST=LEVENE; /*進行方差齊性檢驗*/

means trt / DUNCAN alpha=0.05; /*進行DUNCAN多重比較*/

RUN;

對本試驗而言，辣椒葉片葉綠素SPAD值是因變量，進行澇漬脅迫使得辣椒葉片葉綠素含量在處理間出現(xiàn)差異，那么脅迫處理為“因素”，加上對照，此“因素”具有4個“水平”，在單因素方差分析中該“因素”為自變量，自變量是唯一的。在SAS語句中，通過“class TRT;”語句將TRT設(shè)置為水平，對4個水平進行方差分析和比較。按照蓋鈞鎰[9]對實驗統(tǒng)計方法的描述，本試驗中對1個辣椒品種進行4個水平方差分析為“組內(nèi)觀察值數(shù)目相等的單向分組資料的方差分析”，那么上述SAS語句“model BRD1 BRD2 BRD3 BRD4 BRD5 BRD6 BRD7 BRD8 BRD9 BRD10 BRD11 BRD12 BRD13 BRD14 BRD15 BRD16 BRD17 BRD18 BRD19 BRD20=TRT; ”可以一次性完成20組單向分組資料的方差分析，這樣就一次性對20個品種在澇漬脅迫3、6和9 d后葉片的葉綠素含量差異進行了顯著性分析，而不用重復(fù)20次相同的工作?！癿eans TRT / HOVTEST=LEVENE;”表示在方差分析之前先進行方差的齊性檢驗。“means trt / DUNCAN alpha=0.05;”表示執(zhí)行唐肯氏多范圍檢定(Duncan’s Multiple-Range Test)，當(dāng)進行顯著性水平為0.05的多重比較時，“alpha=0.05”參數(shù)可以省略。

2.3 比較辣椒品種間的差異顯著性

2.3.1 建立數(shù)據(jù)集

比較20個辣椒品種在3個處理與對照之間的葉片葉綠素SPAD含量差異顯著性，需要重新建立數(shù)據(jù)集。

核心語句及注釋如下：

DATA Work.yelvsu1; /*輸入所有品種葉綠素數(shù)據(jù)*/

DO trt=1 TO 4; /*4個處理數(shù)，保留*/

DO t=1 TO 15; /*15個重復(fù)數(shù)，不寫入數(shù)據(jù)庫*/

DO brd=1 TO 20; /*20個品種數(shù)，保留*/

INPUT spad@@;

OUTPUT;

END;

END;

END;

DROP t;

CARDS;

38 32 32 42 37 34 29 32 32 38 37 37 35 44 40 35 28 30 37 30

40 35 33 38 40 36 37 29 39 38 36 36 38 41 42 39 31 33 37 34

40 34 35 39 32 36 27 27 43 34 39 35 37 38 35 39 38 33 33 33

39 39 29 35 37 39 31 29 36 33 42 36 34 36 38 39 39 37 34 31

41 36 38 37 32 33 33 29 31 25 38 37 39 37 37 40 28 28 36 34

……(省略部分類似數(shù)據(jù))

;

RUN;

由于是比較品種間在不同澇漬脅迫條件下葉片葉綠素含量的差異，那么這里可以把20個品種看成一個因素的20個水平，所以需要將所有20個品種葉片葉綠素SPAD值放到SAS數(shù)據(jù)集單獨一列中，以便進行處理間的方差比較，在SAS語句中通過添加一個循環(huán)來實現(xiàn)，具體語句為“DO brd=1 TO 20;”。執(zhí)行之后數(shù)據(jù)集yelvsu1會產(chǎn)生trt，brd和spad共3列數(shù)據(jù)。

2.3.2進行方差分析

比較澇漬脅迫下辣椒品種間葉片葉綠素含量差異，并進行DUNCAN多重比較。

核心代碼和注釋如下：

PROC anova data=Work.yelvsu1; /*調(diào)用上一步建立的數(shù)據(jù)集Yelvsu1*/

class BRD;

model spad=BRD;

means BRD / HOVTEST=LEVENE;

means brd / DUNCAN; /*默認(rèn)顯著性水平為0.05*/

RUN;

“class BRD;”語句是把所有20個辣椒品種看作一個因素的20個水平，比較20個辣椒品種間在澇漬脅迫下葉綠素含量是否存在差異顯著性。如果把所有的辣椒品種看成一個整體，比較4個處理之間在澇漬脅迫下葉片葉綠素含量是否有差異，在程序中把class的值修改為TRT即可。在進行多重比較時，不寫alpha參數(shù)值則默認(rèn)顯著性水平為0.05。如果想同時比較不同品種間和不同處理間辣椒葉片SPAD值差異，按照蓋鈞鎰[9]對實驗統(tǒng)計方法的描述，此類分析屬于“組內(nèi)有重復(fù)觀察值的兩向分組資料的方差分析”[9]，即多因素方差分析。與楊孔雀等[10]使用的方法不同，SAS系統(tǒng)PROC步的anova分析主要研究單個因素對因變量的影響[11]，多因素方差分析建議使用PROC步的GLM模型進行參數(shù)估計。

3 討論

在SAS數(shù)據(jù)集的建立方法上，除了本文所述采用DATA步的input函數(shù)方法外，還可以使用從外部導(dǎo)入數(shù)據(jù)文件的方法，但是在導(dǎo)入數(shù)據(jù)之前需要將數(shù)據(jù)整理成SAS所需要的格式。在農(nóng)業(yè)生物統(tǒng)計中常常需要對大量的數(shù)據(jù)進行統(tǒng)計分析，在本篇文章中有1500個葉綠素SPAD值數(shù)據(jù)，在將該數(shù)據(jù)在Excel中整理成SAS所需格式時，由于數(shù)據(jù)超長，1500個數(shù)據(jù)就需要占用1500行，超過一個屏幕的高度，這樣在整理數(shù)據(jù)容易造成數(shù)據(jù)的丟失、重疊等問題。當(dāng)使用SAS分析較大數(shù)量的數(shù)據(jù)時建議使用DATA步來建立數(shù)據(jù)集，盡量不要采用人工整理數(shù)據(jù)然后導(dǎo)入SAS的方式，避免人為錯誤導(dǎo)致分析結(jié)果的失真。

4 小結(jié)

使用SAS同時對20個辣椒品種進行差異顯著性分析，明顯比用Excel逐個進行方差分析更為方便快捷，并且能用字母標(biāo)注法顯示各個品種在各處理間的顯著性差異。在SAS系統(tǒng)中Work數(shù)據(jù)庫是一個臨時庫，當(dāng)SAS系統(tǒng)關(guān)閉時，臨時數(shù)據(jù)庫中的內(nèi)容會被自動清空，所以在實際使用過程中建議單獨建立一個永久庫，這樣避免每次進行數(shù)據(jù)分析時都要重新建立數(shù)據(jù)集的麻煩。DATA步INPUT語句中的“@@”符號表示順序依次連續(xù)讀入數(shù)據(jù)，如果沒有此符號則按照行讀入數(shù)據(jù)。在實際SAS分析數(shù)據(jù)時我們往往遇到組內(nèi)觀察值數(shù)據(jù)不完整的情況，這時我們可以將空格數(shù)據(jù)用英文狀態(tài)的“.”來代替，但是過程步PROC這里不能使用ANOVA，ANOVA主要對平衡實驗設(shè)計的數(shù)據(jù)執(zhí)行變異數(shù)分析，對于不平衡試驗設(shè)計數(shù)據(jù)則建議用一般線性模型GLM進行變異數(shù)分析。

[1]彭紹英.世界統(tǒng)計與分析全才 SAS 系統(tǒng)應(yīng)用指南[M].北京: 希望電子出版社, 2000: 1-868.

[2]顏亭玉, 杜曉林.SAS聚類分析過程在農(nóng)業(yè)試驗統(tǒng)計分析中的應(yīng)用[J].北京農(nóng)學(xué)院學(xué)報, 2009,24(1): 73-75.

[3]周 倩, 張晉昕.含缺失值的重復(fù)測量資料分析在SPSS和SAS中的實現(xiàn)[J].循證醫(yī)學(xué), 2013,13(2): 120-123.

[4]林 潔, 孫志明.SAS、PASS、Stata三種常用軟件樣本量計算方法及結(jié)果差異的比較[J].中國醫(yī)藥導(dǎo)報, 2015,12(18): 133-137.

[5]詹秋文.Excel和SAS在生物統(tǒng)計學(xué)的應(yīng)用比較[J].生物學(xué)雜志, 2009,26(1): 74-76.

[6]霍志軍, 李菊艷, 潘曉琳.Excel在農(nóng)業(yè)生物統(tǒng)計分析中的應(yīng)用[J].現(xiàn)代化農(nóng)業(yè), 2003(9): 28-30.

[7]范 平, 崔黨群, 詹克慧, 等.Excel軟件在生物統(tǒng)計實驗教學(xué)中的綜合開發(fā)應(yīng)用[J].實驗技術(shù)與管理, 2003,20(2): 65-69.

[8]尹冬梅, 管志勇, 陳素梅, 等.菊花及其近緣種屬植物耐澇評價體系建立及耐澇性鑒定[J].植物遺傳資源學(xué)報, 2009, 10(3): 399-404.

[9]蓋鈞鎰.試驗統(tǒng)計方法[M].北京: 中國農(nóng)業(yè)出版社, 2000: 100-127.

[10]楊孔雀, 孫占育, 蒙惠軍.用SAS軟件進行組內(nèi)有重復(fù)的兩向分組資料的方差分析[J].河北農(nóng)業(yè)科學(xué), 2010,14(10): 165-166.

[11]阮 敬.SAS統(tǒng)計分析從入門到精通[M].北京: 人民郵電出版社, 2009: 99-122.