最小工作示例:
##############################################################################
set.seed(10)
##############################################################################
library(Matrix)
library(psych)
##############################################################################
N0 <- 1
n0 <- 5
p0 <- 2
q0 <- 4
n <- n0*q0
##############################################################################
m2 <- matrix(c(0,0),p0,1)
s2 <- matrix(c(1,0,0,1),p0,p0)
##############################################################################
Dat <- array(data=NA,dim=c(n0,p0,q0,N0))
m1 <- array(data=NA,dim=c(q0,1,1))
#MNA<-array(data=NA,dim=c(q0,1,1))
for (i in 1:N0){
for (j in 1:q0){
Dat[,,j,i] <- mvrnorm(n=n0,m2,s2)
MNA[,,j,i]<- mardia(Dat[,,j,i],plot=FALSE)
m1<-MNA$b1p
}
}
I
想问如何获得一个名为 MNA 的向量,其中将显示我的所有 4 个偏度系数。 因为当我在 R 中编写 MNA 时,它只给出最后一个数据片段的偏度系数。
MNA<- mardia(x = Dat[, , j, i], plot = FALSE)
Mardia tests of multivariate skew and kurtosis
Use describe(x) the to get univariate tests
n.obs = 5 num.vars = 2
b1p = 1.12 skew = 0.93 with probability <= 0.92
small sample skew = 2.23 with probability <= 0.69
b2p = 3.35 kurtosis = -1.3 with probability <= 0.19>
I
仅获取最后一条数据的结果。 我期望得到一个名为 MNA 的向量,其中呈现所有 4 个偏度向量。
稍微改变你的方法,使用上面的参数和矩阵成员
library(Matrix)
library(psych)
library(MASS) # for mvrnorm()
set.seed(10)
mardia_lst <- list() # my suggested change as receiver
for (i in 1:q0) {
for (j in 1:N0) {
Dat[,,j,i] <- MASS::mvrnorm(n = n0, m2, s2)
mardia_lst[[i]] <- mardia(Dat[,,j,i], plot = FALSE)
}
}
str(mardia_lst)
List of 4
$ :List of 12
..$ n.obs : int 5
..$ n.var : int 2
..$ b1p : num 0.593
..$ b2p : num 3.35
..$ skew : num 0.494
..$ small.skew: num 1.19
..$ p.skew : num 0.974
..$ p.small : num 0.88
..$ kurtosis : num -1.3
..$ p.kurt : num 0.194
..$ d : num [1:5] 1.556 1.244 0.417 1.102 1.625
..$ Call : language mardia(x = Dat[, , j, i], plot = FALSE)
..- attr(*, "class")= chr [1:2] "psych" "mardia"
$ :List of 12
..$ n.obs : int 5
..$ n.var : int 2
..$ b1p : num 1.76
..$ b2p : num 4.13
..$ skew : num 1.47
..$ small.skew: num 3.52
..$ p.skew : num 0.833
..$ p.small : num 0.475
..$ kurtosis : num -1.08
..$ p.kurt : num 0.279
..$ d : num [1:5] 0.393 1.717 1.699 1.375 0.35
..$ Call : language mardia(x = Dat[, , j, i], plot = FALSE)
..- attr(*, "class")= chr [1:2] "psych" "mardia"
$ :List of 12
..$ n.obs : int 5
..$ n.var : int 2
..$ b1p : num 0.79
..$ b2p : num 2.93
..$ skew : num 0.659
..$ small.skew: num 1.58
..$ p.skew : num 0.956
..$ p.small : num 0.812
..$ kurtosis : num -1.42
..$ p.kurt : num 0.156
..$ d : num [1:5] 1.12 1.12 1.26 1.67 1.06
..$ Call : language mardia(x = Dat[, , j, i], plot = FALSE)
..- attr(*, "class")= chr [1:2] "psych" "mardia"
$ :List of 12
..$ n.obs : int 5
..$ n.var : int 2
..$ b1p : num 1.35
..$ b2p : num 3.2
..$ skew : num 1.13
..$ small.skew: num 2.7
..$ p.skew : num 0.89
..$ p.small : num 0.609
..$ kurtosis : num -1.34
..$ p.kurt : num 0.18
..$ d : num [1:5] 0.943 0.92 1.75 1.263 1.268
..$ Call : language mardia(x = Dat[, , j, i], plot = FALSE)
..- attr(*, "class")= chr [1:2] "psych" "mardia"
这可能更接近您正在寻找的东西。