使用 R 从列表构造逆矩阵

Construct and inverse matrix from list using R

我有一个从 GCTA 生成的关系矩阵,我可以使用以下函数将其导入 R

ReadGRMBin=function(prefix, AllN=F, size=4){
  sum_i=function(i){
    return(sum(1:i))
  }
  BinFileName=paste(prefix,".grm.bin",sep="")
  NFileName=paste(prefix,".grm.N.bin",sep="")
  IDFileName=paste(prefix,".grm.id",sep="")
  id = read.table(IDFileName)
  n=dim(id)[1]
  BinFile=file(BinFileName, "rb");
  grm=readBin(BinFile, n=n*(n+1)/2, what=numeric(0), size=size)
  NFile=file(NFileName, "rb");
  if(AllN==T){
    N=readBin(NFile, n=n*(n+1)/2, what=numeric(0), size=size)
  }
  else N=readBin(NFile, n=1, what=numeric(0), size=size)
  i=sapply(1:n, sum_i)
  return(list(diag=grm[i], off=grm[-i], id=id, N=N))
}

然后列出非对角线和对角线。

$ diag: num [1:850] 0.878 0.815 1.11 1.161 1.062 ...  
$ off : num [1:360825] 0.0181 -0.0304 -0.0663 -0.0211 -0.0583 ...  
$ n   : int 850

我希望创建一个 grm,我可以从中反转,理想情况下在输出行、列、值中

我试过下面的代码,但它没有以正确的格式读取对角线

  m <- matrix(NA, ncol = length(grm$diag), nrow = length(grm$diag))
    m[lower.tri(m)] <- grm$off
    m[upper.tri(m)] <- t(m)[upper.tri(t(m))]
    diag(m) <- grm$diag
    m
    want=cbind(which(!is.na(m),arr.ind = TRUE),na.omit(as.vector(m)))

而不是将对角线值读取为

2 1、3 1、3 2、4 1、4 2 等

它正在读取对角线的长度方向为

2 1、3 1、4 1、5 1、6 1 等

所以生成的矩阵(缩短)最终是这样的

      [,1]         [,2]        [,3]         [,4]        [,5]
[1,]  0.87798703  0.018129893 -0.03044302 -0.066282429 -0.02106927
[2,]  0.01812989  0.814602911  0.07577287 -0.004078172 -0.03182918
[3,] -0.03044302  0.075772874  1.10976517 -0.055698857 -0.03960679
[4,] -0.06628243 -0.004078172 -0.05569886  1.160611629 -0.01021352
[5,] -0.02106927 -0.031829182 -0.03960679 -0.010213521  1.06245303

如果偏好是这个

      [,1]         [,2]        [,3]         [,4]        [,5]
[1,]  0.87798703   0.018129893  -0.03044302  -0.02106927  -0.04011643
[2,]  0.01812989   0.814602911  -0.06628243  -0.00582625  -0.06237402
[3,] -0.03044302  -0.06628243    1.10976517   0.1315616   -0.1601102
[4,] -0.02106927  -0.00582625    0.1315616    1.160611629 -0.1388046
[5,] -0.04011643  -0.06237402   -0.1601102   -0.1388046    1.06245303

如果您知道如何修改上述代码以提供所需的格式,将不胜感激。 如果可能的话,最终期望的输出将是长格式矩阵的逆矩阵。谢谢

1 1 12456
1 2 78910
1 3 34568
1 4 68942

一个简单的解决方案是调整您的代码以在填充下三角之前先填充上三角(因为应该按列顺序填充上三角):

grm = list(
  diag = 1:5 * 11,
  off  = 0:9)

m <- diag(grm$diag)
m[upper.tri(m)] <- grm$off
m[lower.tri(m)] <- t(m)[lower.tri(t(m))]
#      [,1] [,2] [,3] [,4] [,5]
# [1,]   11    0    1    3    6
# [2,]    0   22    2    4    7
# [3,]    1    2   33    5    8
# [4,]    3    4    5   44    9
# [5,]    6    7    8    9   55