轉自個人微信公眾號【Memo_Cleon】的統計學習筆記:多元正態分布檢驗的R實現方法。
多元正態分布也稱多元高斯分布。如同正態分布在單變量分析中的地位類似,在對多個因變量(多元)同時進行分析時,常常假設因變量組合成的向量服從一個多元正態分布。比如重複測量數據將重複的測量結果(比如各個時間點上的測量結果)視為不同的因變量,可以採用多元方差分析,此時就要求各個因變量的組合向量服從多元正態分布。對多元正態分布的判斷通常採用的邊際分布來判斷,即每個因變量的分布呈正態或近似正態。但實際上單因變量正態是多因變量多元正態的必要非充分條件:所有因變量的組合服從多元正態分布,每個因變量的分布(邊際分布)必然呈正態;每個因變量的分布呈正態分布,所有因變量的組合未必呈正態分布。只要有一個因變量不服用正態分布則組合分布肯定不服從多元正態分布。
雖然我們可以通過邊際分布來大體判斷,有沒有統計方法直接判斷是否滿足多元正態分布呢?當然,R有無所不能的「包」!比如程序包mvnormtest【函數mshapiro.test】、程序包mvShapiroTest【mvShapiro.Test】、程序包MVN【函數mvn】、程序包mvnTest【函數AD.test、CM.test、DH.test、HZ.test、R.test、S2.test等】、程序包mvtnorm【函數pmvnorm】、程序包energy【函數mvnorm.e、mvnorm.test、mvnorm.etest】……
我們以mshapiro.test {mvnormtest}和mvn {MVN}為例演示多元正態分布檢驗的R實現方法。
mshapiro.test {mvnormtest}mshapiro.test(U):Performs the Shapiro-Wilk test for multivariate normality。U要求是數值型矩陣,且a matrix with number of columns (sample size) between 3 and 5000。因為在一般在錄入時行表示觀測/記錄,列表示變量/欄位,因此使用此函數需要先對行和列進行轉置,轉置函數為t(「矩陣或數據框」)。
示例:兩組新生兒出生時的體重與身長數據如下。預對兩組的體重和身高同時進行檢驗,可以考慮多元方差分析,應考察多元正態性。
文件導入:
library(openxlsx) #調用程序包openxlsx
mn<-read.xlsx(「D:/Temp/multivnorm.xlsx」,1) #從名稱為multivnorm的excel文件中導入第1個sheet的數據到數據框mn中
**將數據調整成mshapiro.test的分析格式:**對數據框mn的行列進行轉置
A<-t(mn[1:8,3:4])
B<-t(mn[9:16,3:4])
多元正態分析:
library(mvnormtest) #調用程序包mvnormtest
mshapiro.test(A)
mshapiro.test(B)
結果顯示A組多元正態分布SW檢驗W=0.912,P
=0.366>0.05,呈二元正態分布;B組多元正態分布SW檢驗W=0.861,P=0.122>0.05,呈二元正態分布。
mvn(data, subset = NULL, mvnTest = c(「mardia」, 「hz」, 「royston」, 「dh」, 「energy」), covariance = TRUE, tol = 1e-25, alpha = 0.5, scale = FALSE, desc = TRUE, transform = 「none」, R = 1000, univariateTest = c(「SW」, 「CVM」, 「Lillie」, 「SF」, 「AD」), univariatePlot = 「none」, multivariatePlot = 「none」, multivariateOutlierMethod = 「none」, bc = FALSE, bcType = 「rounded」, showOutliers = FALSE, showNewData = FALSE)。Data為矩陣或數據框。
示例:10名肥胖患者在醫生指導下服用藥物減肥,按統一標準記錄服藥前和服藥後1-4周的體重。
這是一個典型的無對照的重複測量數據,如將各個時點的體重視為不同的因變量,因變量(結果變量)就不止一個,此時可考慮多元方差分析(MANOVA)對它們同時進行分析。多元方差分析要求多元正態性、組間方差-協方差矩陣同質性、各因變量間有一定的相關性(個體內不獨立,但個體間獨立)。
library(openxlsx) #調用程序包openxlsx
mn2<-read.xlsx(「D:/Temp/multivnorm.xlsx」,2) #從名稱為multivnorm的excel文件中導入第2個sheet的數據到數據框mn2中
U<-mn2[1:10,2:6]
library(MVN)
mvn(U, mvnTest = c(「dh」), multivariatePlot = 「qq」)
註:本例10行6列,列名為id、W0、W1、W2、W3、W4。分析時僅篩選列,行不進行篩選,U<-multivnorm[1:10,2:6]可以直接U<-multivnorm[2:6]或者U<-multivnorm[-1]都可以。
結果:本例採用Doornik-Hansen了多變量正態性檢驗,當然也可以在語句中修改命令,換成Marida, Royston, Henze-Zirkler’s, E-Statistics等方法。結果顯示E=6.79,P=0.75>0.05,數據滿足多元正態分布。同時結果還給出了每個變量的Shapiro-Wilk檢驗結果,結果顯示W0、W1、W2、W3、W4均滿足正態分布。另外本例同時給出了多元正態性的QQ圖。
轉自個人微信公眾號【Memo_Cleon】的統計學習筆記:多元正態分布檢驗的R實現方法。