多因素一元方差分析
(1)多因素一元方差分析的MATLAB實現
MATLAB工具工具箱中提供了anovan函數,用來根據樣本觀測值向量進行均衡或非均衡試驗的多因素一元方差分析,檢驗多個(N個)因素的主效應或交互式效應是否顯著,調用格式如下:
<1>p=anovan(y,group)
根據樣本觀測值向量y進行均衡或非均衡試驗的多因素一元方差分析,檢驗多個因素的主效應是否顯著。輸入參數group是一個元胞數組,它的每一個元胞對應一個因素,是該因素的水平列表,與y等長,用來标記y中的每個觀測所對應的因素的水平。
anovan函數還生成1個圖形,用來顯示一個标準的多因素一元方差分析表。
<2>p=anovan(y,group,param1,val1,param2,val2,......)
通過指定一個或多個成對出現的參數名與參數值來控制多因素一元方差分析。可用的參數名與參數值如下表
參數名 | 參數值 | 說明 |
‘alpha’ | (0,1)内的标量 | 指定置信水平 |
‘continuous’ | 下标向量 | 用來指明哪些分組變量被作為連續變量,而不是離散的分類變量 |
‘display’ | ‘on’‘off’ | 用來指定是否顯示方差分析表 |
‘model’ | ‘linear’,‘interaction’‘full’ | 用來指定所用模型的類型,‘linear’(默認)隻對N個主效應進行檢驗,不考慮交互效應,‘interaction’對主效應和兩個因素的交互效應減小檢驗‘full’對N個主效應和全部的交互效應進行檢驗。也可以通過0和1的矩陣自定義效應項 |
‘nested’ | 由0和1構成的矩陣M | 指定分組變量之間的嵌套關系 |
‘random’ | 下标向量 | 用來指明哪些分組變量是随機的。 |
‘sstype’ | 1,2或3 | 指定平方和的類型,默認值為3 |
‘varnames’ | 字符矩陣或字符串元胞數組 | 指定分組變量的名稱。當沒有指定時,默認用‘X1’,‘X2’,…’XN’做為名稱 |
<3>[p,table]=anovan(......)
返回元胞數組形式的方差分析表table。
<4>[p,table,stats]=anovan(.....)
返回一個結構體變量stats,用于進行後續的多重比較,。當某因素對實驗指标的影響顯著時,在後續的分析中,可用調用multcompare函數,把stats作為其輸入,進行多重比較。
調用anovan函數進行分析,在不考慮區組因素的情況下,分析氮、磷兩種肥料的施用量對水稻的産量是否有顯著影響,并分析交互作用是否顯著,然後找出在因素A,B的哪種水平組合下水稻的平均産量高,顯著性水平為0.05.
%定義一個矩陣,輸入原始數據
yield=[38 29 36 40
45 42 37 43
58 46 52 51
67 70 65 71
62 64 61 70
58 63 71 69];
yield=yield' %矩陣轉置
%将數據矩陣yield按列拉長成24行1列的向量
yield=yield(:);
%定義因素A(氮)的水平列表向量
A=strcat({'N'},num2str([ones(8,1);2*ones(8,1);3*ones(8,1)]));
%定義因素B(磷)的水平列表向量
B=strcat({'P'},num2str([ones(4,1);2*ones(4,1)]));
B=[B;B;B];
%将因素A、B的水平列表向量與yield向量放在一起構成一個元胞數組,以元胞數組形式顯示出來
[A,B,num2cell(yield)]
%指定因素名稱,A表示氮肥施用量,B表示磷肥施用量
varnames={'A','B'};
%調用anovan函數作雙因素一元方差分析,返回主效應A、B和交互效應AB所對應的p值向量
%還返回方差分析表table,結構體變量stats,标識模型效應項的矩陣term
[p,table,stats,term]=anovan(yield,{A,B},'model','full','varnames',varnames)
ans =
'N1' 'P1' [38]
'N1' 'P1' [29]
'N1' 'P1' [36]
'N1' 'P1' [40]
'N1' 'P2' [45]
'N1' 'P2' [42]
'N1' 'P2' [37]
'N1' 'P2' [43]
'N2' 'P1' [58]
'N2' 'P1' [46]
'N2' 'P1' [52]
'N2' 'P1' [51]
'N2' 'P2' [67]
'N2' 'P2' [70]
'N2' 'P2' [65]
'N2' 'P2' [71]
'N3' 'P1' [62]
'N3' 'P1' [64]
'N3' 'P1' [61]
'N3' 'P1' [70]
'N3' 'P2' [58]
'N3' 'P2' [63]
'N3' 'P2' [71]
'N3' 'P2' [69]
p =
0.0000
0.0004
0.0080
table =
'Source' 'Sum Sq.' 'd.f.' 'Singular?' 'Mean Sq.' 'F' 'Prob>F'
'A' [3.0670e 03] [ 2] [ 0] [1.5335e 03] [78.3064] [1.3145e-09]
'B' [ 368.1667] [ 1] [ 0] [ 368.1667] [18.8000] [3.9813e-04]
'A*B' [ 250.3333] [ 2] [ 0] [ 125.1667] [ 6.3915] [ 0.0080]
'Error' [ 352.5000] [ 18] [ 0] [ 19.5833] [] []
'Total' [ 4038] [ 23] [ 0] [] [] []
stats =
source: 'anovan'
resid: [24x1 double]
coeffs: [12x1 double]
Rtr: [6x6 double]
rowbasis: [6x12 double]
dfe: 18
mse: 19.5833
nullproject: [12x6 double]
terms: [3x2 double]
nlevels: [2x1 double]
continuous: [0 0]
vmeans: [2x1 double]
termcols: [4x1 double]
coeffnames: {12x1 cell}
vars: [12x2 double]
varnames: {2x1 cell}
grpnames: {2x1 cell}
vnested: []
ems: []
denom: []
dfdenom: []
msdenom: []
varest: []
varci: []
txtdenom: []
txtems: []
rtnames: []
term =
1 0
0 1
1 1
返回的向量p是主效應A,B和交互效應AB所對應的檢驗的p值,table是元胞數組的方差分析表,stats是一個結構體變量,可用于後續的分析中(如多重比較),矩陣term的3行分布表示了3個效應項:主效應項A,主效應項B和交互式效應項AB,還生成了一個方差分析表。
返回的結果與調用anova2函數得到的結果一樣,因素A,因素B以及他們的交互作用所對應的p值均小于給定的顯著想水平0.05,所以可以認為氮、磷兩種肥料的施用量對水稻的産量均有非常顯著的影響,并且他們之間的交互作用也是非常顯著的。
(3)多重比較。
調用multcompare函數,把anovan函數返回的結構體變量stats作為它的輸入,對因素A、B的每種水平的組合進行多重比較,找出在因素A、B的哪種水平組合下水稻的平均産量最高。
%多重比較
[c,m,h,gnames]=multcompare(stats,'dimension',[1 2])
%查看各處理的均值
[gnames,num2cell(m)]
c =
1.0000 2.0000 -25.9446 -16.0000 -6.0554 0.0009
1.0000 3.0000 -38.4446 -28.5000 -18.5554 0.0000
1.0000 4.0000 -15.9446 -6.0000 3.9446 0.4236
1.0000 5.0000 -42.4446 -32.5000 -22.5554 0.0000
1.0000 6.0000 -39.4446 -29.5000 -19.5554 0.0000
2.0000 3.0000 -22.4446 -12.5000 -2.5554 0.0093
2.0000 4.0000 0.0554 10.0000 19.9446 0.0483
2.0000 5.0000 -26.4446 -16.5000 -6.5554 0.0006
2.0000 6.0000 -23.4446 -13.5000 -3.5554 0.0047
3.0000 4.0000 12.5554 22.5000 32.4446 0.0000
3.0000 5.0000 -13.9446 -4.0000 5.9446 0.7926
3.0000 6.0000 -10.9446 -1.0000 8.9446 0.9995
4.0000 5.0000 -36.4446 -26.5000 -16.5554 0.0000
4.0000 6.0000 -33.4446 -23.5000 -13.5554 0.0000
5.0000 6.0000 -6.9446 3.0000 12.9446 0.9251
m =
35.7500 2.2127
51.7500 2.2127
64.2500 2.2127
41.7500 2.2127
68.2500 2.2127
65.2500 2.2127
h =
2
gnames =
'A=N1,B=P1'
'A=N2,B=P1'
'A=N3,B=P1'
'A=N1,B=P2'
'A=N2,B=P2'
'A=N3,B=P2'
ans =
'A=N1,B=P1' [35.7500] [2.2127]
'A=N2,B=P1' [51.7500] [2.2127]
'A=N3,B=P1' [64.2500] [2.2127]
'A=N1,B=P2' [41.7500] [2.2127]
'A=N2,B=P2' [68.2500] [2.2127]
'A=N3,B=P2' [65.2500] [2.2127]
返回的矩陣c是6種處理(N1P1,N2P1,N3P1,N1P2,N2P2,N3P2)間多重比較的結果矩陣,每一行的前兩列是進行比較的兩個處理的編号,第4列是兩個處理的均值之差,第3列是兩個處理均值差的95%置信下限,第5列是兩個處理均值差的95%置信下限,當兩個處理均值差的95%置信區間不包含0時,說明在顯著性水平0.05下,這兩個處理均值間差異是顯著的。m矩陣列出了6種處理的平均值,很明顯第5個處理(即N2P2)的平均值最大,由矩陣c或交互式多重比較的圖中可以得到,處理5與處理3,6差異不顯著,所以可以認為第3個和第6個處理也是可以的,所以,綜上,可以在處理3,5,6中做出選擇,即N3P1,N2P2,N3P2。
,更多精彩资讯请关注tft每日頭條,我们将持续为您更新最新资讯!