将此向量转换为对称矩阵的快速方法是什么？

Question

大家！假设我有一个长度为 n(n+1)/2 的向量：

a = (a_11, a_12, a_22, ...., a_nn) 

现在我想把它变成一个对称矩阵，意思是

我可以一一分配值，但我想知道是否有更快的方法来创建这个矩阵？非常感谢！！

Answer 1

您可以编写一个函数来执行此操作：

如果你的向量为

a11,a12,a13..a1n,a22,a23..a2n, a33,..a3n,..ann

你可以做的：

vec2mat <- function(x){
  p <- sqrt(1 + 8 * length(x))/ 2 - 0.5
  m <- matrix(0, p, p)
  m[lower.tri(m, diag = TRUE)] <- x
  m[upper.tri(m)] <- (t(m))[upper.tri(m)]
  m
}

现在：

vec2mat(1:6)
     [,1] [,2] [,3]
[1,]    1    2    3
[2,]    2    4    5
[3,]    3    5    6
vec2mat(1:10)
     [,1] [,2] [,3] [,4]
[1,]    1    2    3    4
[2,]    2    5    6    7
[3,]    3    6    8    9
[4,]    4    7    9   10

如果你有

a11, a12,a22,a31, a32, a33...

vec2mat <- function(x){
  p <- sqrt(1 + 8 * length(x))/ 2 - 0.5
  m <- matrix(0, p, p)
  m[upper.tri(m, diag = TRUE)] <- x
  m[lower.tri(m)] <- t(m)[lower.tri(m)]
  m
}
 vec2mat(1:10)
     [,1] [,2] [,3] [,4]
[1,]    1    2    4    7
[2,]    2    3    5    8
[3,]    4    5    6    9
[4,]    7    8    9   10

Answer 2

你可以尝试一下：

#define n
n <- 2

#Create n(n+1)/2 objects in global environment
#OP already has that
a_11 <- 5
a_12 <- 9
a_22 <- 8

#Create n X n matrix with NA
mat <- matrix(nrow = n, ncol = n)
#Get all the individual objects in one vector
vec <- unlist(mget(ls(pattern = 'a_')), use.names = FALSE)
#Replace upper (or lower) triangular elements with it
mat[upper.tri(mat, diag = TRUE)] <- vec
#Copy the elements to other half.
mat[lower.tri(mat)] <- mat[upper.tri(mat)]

#      [,1] [,2]
#[1,]    5    9
#[2,]    9    8

Answer 3

您可以像这样创建一个稀疏矩阵：

a <- 1:6
n <- as.integer(-0.5 + sqrt(0.25 + 2 * length(a)))

library(Matrix)
sparseMatrix(x = a, dims = c(n, n), symmetric = TRUE, 
             i = sequence(1:n), j = rep(1:n, 1:n))
#3 x 3 sparse Matrix of class "dsCMatrix"
#          
#[1,] 1 2 4
#[2,] 2 3 5
#[3,] 4 5 6

如果需要密集矩阵，请在结果上使用

as.matrix

。如果向量的顺序与您显示的不同（例如，如其他一些答案所假设的行优先），则需要稍微调整

i

和

j

的计算。

Answer 4

a <- 1:25
sq.m <- matrix( a, ncol=sqrt(length(a) ) )

> sq.m
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    6   11   16   21
[2,]    2    7   12   17   22
[3,]    3    8   13   18   23
[4,]    4    9   14   19   24
[5,]    5   10   15   20   25

显然这不是一个对称矩阵，但如果初始向量的阶数为 1，那么它就成功了。

如果您想强制非方阵，使其上三角元素与下三角元素相同，则可以在赋值的两侧进行索引：

 sq.m[ upper.tri(sq.m) ] <- sq.m[lower.tri(sq.m)]
> sq.m
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    2    3    5   10
[2,]    2    7    4    8   14
[3,]    3    8   13    9   15
[4,]    4    9   14   19   20
[5,]    5   10   15   20   25

Answer 5

如果你有大向量并且速度很重要，那么你可以使用这个使用 Rfast 包的函数（取自@Onyambu）。

vec2mat_rfast <- function(x){
    p <- sqrt(1 + 8 * length(x))/ 2 - 0.5
    m <- matrix(0, p, p)
    m[Rfast::upper_tri(m, diag = TRUE)] <- x
    m[Rfast::lower_tri(m)] <- Rfast::transpose(m)[Rfast::lower_tri(m)]
    m
}

vec2mat <- function(x){
    p <- sqrt(1 + 8 * length(x))/ 2 - 0.5
    m <- matrix(0, p, p)
    m[upper.tri(m, diag = TRUE)] <- x
    m[lower.tri(m)] <- t(m)[lower.tri(m)]
    m
}

y=numeric(500500)


> microbenchmark::microbenchmark(a<-vec2mat(y),b<-vec2mat_rfast(y),times = 10)
Unit: milliseconds
                    expr     min       lq      mean    median       uq      max neval
         a <- vec2mat(y) 88.9461 101.5294 114.96752 108.86680 112.9854 180.9679    10
   b <- vec2mat_rfast(y) 33.1351  34.5294  49.61295  45.26955  62.8214  80.5069    10

> all.equal(a,b)
[1] TRUE

Answer 6

将速度与上述方法进行比较会很有趣：

# For i,j element in pxp symmetric matrix y compute subscript in compact vector x (columns moving fastest)
ssub <- function(i, j, p) {
  ii <- pmin(i, j)
  jj <- pmax(i, j)
  (ii - 1) * p + jj - (ii - 1) - (ii - 1) * (ii - 2) / 2
}

p <- 5
y <- matrix(NA, p, p)
nv <- p * (p + 1) / 2
x <- 1 : nv
y[] <- x[ssub(row(y), col(y), p)]
y