在R中从相关向量创建相关矩阵

8
我想根据相关系数向量创建相关矩阵,该向量是相关矩阵的上(或下)三角矩阵。目标是将这个向量转换为对角线上为1的相关矩阵。您知道是否有一种方法可以创建一个矩阵,给定三角形在对角线上方,并将对角线设置为1吗?
4个回答

9
您可以说服R,让它将您的向量视为距离对象,然后使用as.matrix进行转换:
> myvec <- c(-.55, -.48, .66, .47, -.38, -.46)
> class(myvec) <- 'dist'
> attr(myvec,'Size') <- 4
> as.matrix(myvec)
      1     2     3     4
1  0.00 -0.55 -0.48  0.66
2 -0.55  0.00  0.47 -0.38
3 -0.48  0.47  0.00 -0.46
4  0.66 -0.38 -0.46  0.00

以下是对@AnandaMahto答案的变体(与上面使用的内部机制类似):

> myvec <- c(-.55, -.48, .66, .47, -.38, -.46)
> mycor <- matrix(0,4,4)
> mycor[ col(mycor) < row(mycor) ] <- myvec
> mycor <- mycor + t(mycor)
> diag(mycor) <- 1
> mycor
      [,1]  [,2]  [,3]  [,4]
[1,]  1.00 -0.55 -0.48  0.66
[2,] -0.55  1.00  0.47 -0.38
[3,] -0.48  0.47  1.00 -0.46
[4,]  0.66 -0.38 -0.46  1.00

1
必须将 diag(1,4) 添加到 as.matrix(myvec) - Ferdinand.kraft

7

我不知道是否有自动化的方法来实现这个,但是可以参考我的评论进行扩展:

myvec <- c(-.55, -.48, .66, .47, -.38, -.46)
mempty <- matrix(0, nrow = 4, ncol = 4)
mindex <- matrix(1:16, nrow = 4, ncol = 4)
mempty[mindex[upper.tri(mindex)]] <- myvec
mempty[lower.tri(mempty)] <- t(mempty)[lower.tri(t(mempty))]
diag(mempty) <- 1
mempty
#       [,1]  [,2]  [,3]  [,4]
# [1,]  1.00 -0.55 -0.48  0.47
# [2,] -0.55  1.00  0.66 -0.38
# [3,] -0.48  0.66  1.00 -0.46
# [4,]  0.47 -0.38 -0.46  1.00

这是一个快速拼凑的函数。我希望我的数学步骤都是正确的!

vec2symmat <- function(invec, diag = 1, byrow = TRUE) {
  Nrow <- ceiling(sqrt(2*length(invec)))

  if (!sqrt(length(invec)*2 + Nrow) %% 1 == 0) {
    stop("invec is wrong length to create a square symmetrical matrix")
  }

  mempty <- matrix(0, nrow = Nrow, ncol = Nrow)
  mindex <- matrix(sequence(Nrow^2), nrow = Nrow, ncol = Nrow, byrow = byrow)
  if (isTRUE(byrow)) {
    mempty[mindex[lower.tri(mindex)]] <- invec
    mempty[lower.tri(mempty)] <- t(mempty)[lower.tri(t(mempty))]
  } else {
    mempty[mindex[upper.tri(mindex)]] <- invec
    mempty[lower.tri(mempty)] <- t(mempty)[lower.tri(t(mempty))]
  }

  diag(mempty) <- diag
  mempty
}

这里是对角线取不同值时的结果。
vec2symmat(1:3, diag = NA)
#      [,1] [,2] [,3]
# [1,]   NA    1    2
# [2,]    1   NA    3
# [3,]    2    3   NA

如果您提供的数据无法创建一个方阵,则会出现以下错误信息。

vec2symmat(1:4)
# Error in vec2symmat(1:4) : 
#   invec is wrong length to create a square symmetrical matrix

而且,使用默认设置。
vec2symmat(1:10)
#      [,1] [,2] [,3] [,4] [,5]
# [1,]    1    1    2    3    4
# [2,]    1    1    5    6    7
# [3,]    2    5    1    8    9
# [4,]    3    6    8    1   10
# [5,]    4    7    9   10    1

还有一个小问题:您的函数是按列填充矩阵的,我该如何按行填充矩阵,就像上面的示例一样?非常感谢! - jeffrey
@jeffrey,为矩阵添加了一个“byrow”参数。不过请查看Greg Snow的回答。 - A5C1D2H2I1M1N2O1R2T1
非常好用,正是我所需要的。你能解释一下你的步骤吗?我对R不太熟悉,只能跟到这里。谢谢。 - AfterWorkGuinness

4

下面是一些辅助函数的答案,这些函数在其他问题中可能会有用:

`lower.tri<-` <- function(x,value){
    x[lower.tri(x)] <- value
    x
}

`upper.tri<-` <- function(x,value){
    y <- t(x)
    lower.tri(y) <- value
    t(y)
}

vec2mat <- function(r){
    n <- (1+sqrt(1+8*length(r)))/2
    x <- diag(1,n)
    lower.tri(x) <- upper.tri(x) <- r
    x
}

编辑:请注意,upper.tri<-并不是简单地在lower.tri<-中将"lower"替换为"upper"而得到的。这会使结果不对称。

结果:

vec2mat(c(-0.55, -0.48, 0.66, 0.47, -0.38, -0.46))

      [,1]  [,2]  [,3]  [,4]
[1,]  1.00 -0.55 -0.48  0.66
[2,] -0.55  1.00  0.47 -0.38
[3,] -0.48  0.47  1.00 -0.46
[4,]  0.66 -0.38 -0.46  1.00

0

在任何编程语言中,我认为你只需要像这样使用几个嵌套的for循环:

给定:

Vector R; // 我将使用圆括号R(3)表示基于1的第三个元素。 // (在具有零基数组的语言中存储在内存位置R[2]中) int N;

Matrix M=Matrix(N,N); // 你的矩阵对象的一个新实例,或者你可以使用数组。

int i,j,k;

k=1;
for(i=1;i<N;i++)
{
   M(i,i)=1;
   for(j=i+1,j<=N;j++)
   { 
      M(i,j)=M(j,i)=R[k];
      k=k+1;
   }
}

这里我假设你知道N是什么,以及你有基本对象,比如向量和矩阵可用。(如果不知道,它们是编写你的第一个“对象”很好的样本问题)复杂的数据结构,如向量、矩阵、复数和直方图都是理想的对象。在科学工作中,正确的面向对象编程方式是使用对象来教会编译器理解你想在实际工作中使用的高级数据类型...对象用于创建适合你工作类型的定制编程语言。任何通用有用的东西都应该放入对象中,因为那些对象将成长和演变为可重用的代码库。

顶层代码可以是非常强大、易读和整洁的应用程序(因为许多细节工作都在对象中完成)。或者,对于快速而粗略的编码,顶层代码是放置所有脆弱的黑客代码的地方,因为它不打算重用。

一旦你调试出这样的东西,你只需要创建一个以相关矢量和N为参数,并为你初始化矩阵的矩阵构造函数。

当然,如果您正在使用一些高级图形数学程序,该程序对您可以使用矩阵和向量的内容有强烈的意见,那么您将不得不通过N个矩阵乘以向量来生成最终矩阵的每个列向量。(否则请阅读手册)
至少,您需要告诉我们这个数学程序叫什么... :)

看一下问题标签,语言已经明确说明了。 - Roland

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接