多元正态分布

转自个人微信公众号【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 {MVN}

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实现方法。