薛亞宏,王 嘉,程喜林
(1.甘肅工業職業技術學院 經管學院,甘肅 天水 7410252;2.蘭州大學 信息科學與工程學院,蘭州 7300003;3.西北大學 城市與環境學院,西安 741069)
方差分析(Analysis of variance,ANOVA)是英國統計學家兼遺傳學家Gitbert提出的一種分析方法,在林業遺傳、科學試驗、醫學研究等眾多領域有極其廣泛的應用[1]。方差是統計量分析中的一類,是假設檢驗與區間估計的推廣和延伸,由于試驗樣本的成因方式不同,一般采用的分法方法也有所差異,在統計學界,通常用單因子(one-way)、復因子(double-way)或多因子(N-way)方法。
主成分分析是根據已有數據推斷假設數據受某主要因素影響程度的一種分析方法,其本質仍為假設檢驗[2]。通常需要提供理論樣本結果和實測樣本結果,兩者通過矩陣排列得到回歸模型,最終得到影響總體數據分布的主要因素及數值,一般要采用二維或三維曲線進行二次以上模擬,在誤差允許范圍內滿足達到精度即終止計算。
復因子分析即有兩個影響因子。通常應用于醫學療效的驗證,其一般的分析步驟為:首先將病人隨機(一般、等距、整群)分為k組,每組有x人,將每位病人的療效監測指標記為tm,n,其中下標m、n分別表示第m組,(i=1,2,…,k),n表示某組內病人的診療編號,(n=1,2,…,w),則第m組的第n個病人的監測指標即為tm,n。按照MATLAB特有的表達式記號,記第k組所有病人的所有監測指標為tm,或各組的第w位病人的監測標為tn,則這樣計算出向量,這涉及算術平均數的范籌,如此則構造出標準方差分析表,并根據所給出的監測數據找出效果分析數據。
以上所采用藥物作為分組的依據,稱為復因子(Complex factor),它們的差異均值稱為復水平。其中m值與p值較為關鍵,直接影響概率值p<a的置信度及拒絕假設目標H0,否則假設不成立,療效分析驗證為假。
在MATLAB中,anoval0、anova2可分別進行單因子、復因子分析,并有效地給出較為精確的結果,其基本格式為[p,Tab,Stats]=anova1(Q)[3],其中,Q為需要分析的數據,該數據是一個k×w矩陣,其行對應于分組號,運算結果會返回檢驗值s,以及檢測數據表Tab;該函數還將打開兩個主程序窗口viewer和power,分別以表式、盒式結構呈現。
案例1:以非嗎啡類中樞型鎮痛藥物鹽酸曲馬多(Tramadol)為例[4],現將40個病人(醫學低于30為小樣本)樣本分為6組,每組5人,患者(patient)使用同一藥物(假定無其他輔助藥物),記錄從用藥到痊愈時間(h),觀測所用藥物的療效是否存在顯著差異,觀測數據如表1所示。

表1 痊愈時間觀測數據表
算法設計:
基于以上監測數據,現構造出一個5×6型矩陣,命名為矩陣Q,對各組數據采用復因子方差分析,得出以下分析結果:

在程序運行過程中,anoval()會自動呈現兩個窗口,分別是盒式圖、分析表,同時顯示概率值p,其中a=0.03或0.04,表示置信水平,顯然從結果來看應拒絕假設,即藥物對痊愈時間有顯著影響。
案例2:以巖松、油松、赤松3種松樹樹種在甘肅省天水市小隴山林區黨川、利橋、草川、草灘4地(林場)的生長情況為例,每地每類樹種選擇6株,測量其胸徑,并進行雙因子方差分析,觀測數據如表2所示:

表2 甘肅省小隴山林區巖松等3類松樹生長觀測數據
Q算法設計:
基于以上監測數據,現調用anoval2()函數,命名為矩陣H,對各組數據采用雙因子方差分析,anoval2()命令及矩陣列排列如下:

從結果來看,由于PA=0.01393,所以應該拒絕H1假設。可以初步推斷,列數據對監測結果有顯著影響,即小隴山林區下轄黨川等4地3類松樹樹種對其胸徑有顯著影響。
以下計算均值,以反映不同樹種在同一林場生長差異:

根據結果分析,赤松胸徑明顯大于巖松和油松。PH與PHQ的差距較大,從而判斷假設為真,故接受假設[5]。即黨川等4地各自對3類松樹樹種的胸徑有輕微影響,不同區域(林場)對不同松樹樹種胸徑成長無顯著影響。
主成分分析是一種常見的多因素分析方法,在信息模擬、疾病預防、地理信息采集、工程造價測算、農作物產量分析等領域有著廣泛應用[6]。通常采用SPSS、R等平臺進行分析,但由于原數據類型的多樣性,輸出圖形特征的兩極分化(異端非同步)現象較為普遍,經與實際監測比對出現較大偏差,結論不穩定,因此不具有代表性。在這種情況下,利用MATLAB在數據降維處理方面的精度、效度以及在圖象表現上的多維仿真優勢,通過調用Corr()函數,建立協方差矩陣及特征向量、主成分貢獻率及累計貢獻率、建立變量指標方程組對數據進行降維處理,實現源數據分析從多維到低維的轉化。
首先建立協方差矩陣

再根據R計算出特征向量ei及特征值,做反序處理;
最后,通過轉換變量(降維),構建平面坐標方程zi=azi+bzi+Ri,獲取主體變量影響因素間的關系構成,即主成分分析的基本表達式。
案例3:以甘肅省天水市李子園鉛鋅礦第四紀淺層地貌特征分析為例,某測量點三維坐標參數分別為 x=ωcos2ω, y=ωsin2ω, z=0.887x+3.463y,現通過MATLAB生成一維數組,并輸出以2個測量數位為基本單位矢量模擬表達式。
算法設計如下:

顯然,基于對降維矢量輸出原理的分析,進一步利用空間坐標變換,對三維原數據做放樣投影,得到二維數組[7]。
執行以上命令,輸出結果為:

值得注意的是,結果中的σ向量、e向量非測量高程測序排列,要通過fliplr()函數和real()函數執行反序和翻轉輸出,目地是使特征值按常規測序呈現,為RNSS測繪系統數據導入做必要的前期配置。
具體語句如下:

轉換后的3×3矩陣提供二維數據(z列為0)如下:

故新坐標系可表示為:

該坐標方程實現了對三維高程測量數組的降維(縱向投影),通過坐標轉化使RNSS源數據壓縮于二維平面上,一方面能準確表現該區域第四紀地貌分布特征,另一方面在同類型礦區主要作業區域地形圖繪制(表層、淺層)中提供了滿足繪制精度要求的一種新的計算途徑,其誤差范圍與多基點均勻采樣在同一水平[8],但其在數據生成原理、仿真形式以及中間變量轉換等多個方面集成了ArcGIS、C++的優勢,有效降低了測圖成本。
基于MATLAB的復合因子方差分析與主成分分析計算原理清晰,算法邏輯性強,語法調用靈活。在實踐中,以應用統計學基本計算理論為基礎,結合矩陣運算、坐標變換等數學手段,最終通過MATLAB實現對樣本數據的終端處理,有效彌補了SPSS、R等分析工具在圖像擬合優度與坐標維度無法兼顧的不足,特別是三維數字測圖、工程概預算、造價分析等領域內能大大地降低數據交叉,有效降低項目成本,有較強的實用意義。