当数据分为4部分时,如何将Mardia的偏度系数保存在向量中

问题描述 投票:0回答:1

最小工作示例:

##############################################################################
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 个偏度向量。

r
1个回答
0
投票

稍微改变你的方法,使用上面的参数和矩阵成员

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"

这可能更接近您正在寻找的东西。

© www.soinside.com 2019 - 2024. All rights reserved.