随机化对称矩阵的非对角线元素

Randomize non-diagonal elements of symmetric matrix

我有一个对称矩阵,我想在保持对角线元素不变的情况下随机打乱它。所有行总和为 1,洗牌后总和仍应为 1。

下面的玩具示例:

A <- rbind(c(0.6,0.1,0.3),c(0.1,0.6,0.3),c(0.1,0.3,0.6))
A
#      [,1] [,2] [,3]
# [1,]  0.6  0.1  0.3
# [2,]  0.1  0.6  0.3
# [3,]  0.1  0.3  0.6

我想要一个矩阵 B,其对角线元素与 A 相同且仍然对称,但元素随机打乱以生成类似于

的矩阵
B <- rbind(c(0.6,0.3,0.1), c(0.3,0.6,0.1), c(0.3,0.1,0.6))
B
#      [,1] [,2] [,3]
# [1,]  0.6  0.3  0.1
# [2,]  0.3  0.6  0.1
# [3,]  0.3  0.1  0.6

我的目标是在 24 * 24 矩阵上执行此操作,因此代码可能会很混乱,并且不需要具有低计算成本的东西。到目前为止,我已经尝试过循环,但代码很快变得过于复杂,我想知道是否有更直接的方法来做到这一点。

一个选项可以是:

set.seed(123)
t(mapply(function(x, y) {
    ind <- which(seq_along(x) != y)
    `[<-`(x, ind, sample(x[ind]))
    },
    x = asplit(A, 1),
    y = 1:nrow(A)))

     [,1] [,2] [,3]
[1,]  0.6  0.1  0.3
[2,]  0.1  0.6  0.3
[3,]  0.1  0.3  0.6

由于您希望将行向总和保持为 1,因此您只能打乱每行中不包括对角线元素的每行元素。

set.seed(2021)

t(sapply(seq(nrow(A)), function(x) {
  tmp <- A[x, ]
  tmp[-x] <- sample(tmp[-x])
  tmp
}))

#     [,1] [,2] [,3]
#[1,]  0.6  0.1  0.3
#[2,]  0.3  0.6  0.1
#[3,]  0.1  0.3  0.6

试试下面的代码

t(mapply(
  function(x, k) replace(x, k, sample(x[k])),
  asplit(A, 1),
  asplit(row(A) != col(A), 1)
))

这给出了

     [,1] [,2] [,3]
[1,]  0.6  0.1  0.3
[2,]  0.3  0.6  0.1
[3,]  0.1  0.3  0.6

基础 R 解决方案:

n <- 3

t(apply(cbind(A, 1:n), 1, 
  function(x) {x[-c(x[n+1], n+1)] <- sample(x[-c(x[n+1], n+1)]); x[1:3]}))

另一个解决方案,这次基于tidyverse/purrr

library(tidyverse)

A <- rbind(c(0.6,0.1,0.3),c(0.1,0.6,0.3),c(0.1,0.3,0.6))
n <- 3

set.seed(23)

t(A) %>% as.data.frame %>% 
  map2_dfr(1:n, ~ {.x[-.y] <- sample(.x[-.y], n-1); .x}) %>%
  unname %>% as.matrix %>% t

#>      [,1] [,2] [,3]
#> [1,]  0.6  0.1  0.3
#> [2,]  0.3  0.6  0.1
#> [3,]  0.3  0.1  0.6

获取非对角元素的索引。子集值和行索引。在每一行中,打乱值并重新分配。

i = row(A) != col(A)
A[i] = ave(A[i], row(A)[i], FUN = sample)
A
#      [,1] [,2] [,3]
# [1,]  0.6  0.1  0.3
# [2,]  0.3  0.6  0.1
# [3,]  0.3  0.1  0.6

如果您不想覆盖原始矩阵,请分配给一个副本。

A = rbind(c(0.6,0.1,0.3), c(0.1,0.6,0.3), c(0.1,0.3,0.6))
i = row(A) != col(A)
A2 = A

set.seed(1)
A2[i] = ave(A[i], row(A)[i], FUN = sample)
A2
#      [,1] [,2] [,3]
# [1,]  0.6  0.1  0.3
# [2,]  0.1  0.6  0.3
# [3,]  0.3  0.1  0.6

set.seed(12)
A2[i] = ave(A[i], row(A)[i], FUN = sample)
A2
#      [,1] [,2] [,3]
# [1,]  0.6  0.3  0.1
# [2,]  0.3  0.6  0.1
# [3,]  0.1  0.3  0.6