R - 标准化矩阵以具有单位对角线

Posted 2023-02-14

【中文标题】R - 标准化矩阵以具有单位对角线

【英文标题】：R - Standardize matrix to have unit diagonals

【发布时间】：2022-01-10 14:30:05

【问题描述】：

我正在寻求生成以下矩阵：

Θ = B + δIp ∈ Rp×p，其中 Ip 是单位矩阵，每个非对角项 在 B（对称矩阵）中独立生成，概率等于 0.5 0.1 或 0，概率为 0.9。选择参数 δ > 0 使得 Θ 是正定的。该矩阵被标准化为具有单位对角线（从协方差矩阵转换为相关矩阵）。

我认为我拥有大部分代码，但我不确定如何标准化矩阵以在 R 中在语法上具有单位对角线（理论上，为什么这是矩阵的一个有用特性）。

# set number of cols/rows
p <- 5
set.seed(123)
# generate matrix B with values of 0.5 given probabilities
B <- matrix(sample(c(0,0.5), p^2, replace=TRUE, prob=c(0.9,0.1)), p)
# call the matrix lower triangle, need a symmetric matrix
i <- lower.tri(B)
B[i] <- t(B)[i]
diag(B) <- rep(0, p)
# finding parameter delta, such that Θ is positive definite.
(delta <- -min(eigen(B, symmetric=TRUE, only.values=TRUE)$values))
# set theta (delta is 2.8802)
theta <- B + 2.89*(diag(p))
# now to standardize the matrix to have unit diagonals ?

【问题讨论】：

你只想要theta <- B/delta + diag(p) 吗？

【参考方案1】：

有很多方法可以做到这一点，但以下在计时实验中非常快：

v <- 1/sqrt(diag(theta))
B <- theta * outer(v, v)

这将所有行和列除以其标准差，即对角线元素的平方根。

只要任何对角线为零（或负数），它就会失败：但在这种情况下，这种标准化是不可能的。首先计算平方根及其倒数可以让您以最少的计算尽快了解该过程是否会成功。

顺便说一句，在代码的第一部分（对角线为零）中计算 B 的直接方法是

B <- as.matrix(structure(sample(c(0,1/2), p*(p-1)/2, replace=TRUE, prob=c(.9,.1),
               Size=p, Diag=TRUE, class="dist"))

这消除了多余的采样。

【讨论】：

一般来说这很好，但我们不知道对角线元素是rep(delta, p)吗？

@Brian True，但缩小该情况的解决方案是没有意义的，因为与应用 rep 相比，提取矩阵的对角线没有比较成本。如果稍后将相同的方法应用于一般协方差矩阵，这也可以防止出现可能的错误。在这种特殊情况下，p 的全局划分可能会略有加速。但是，我将这个问题解释为一般性地询问如何将协方差矩阵转换为其相关的相关矩阵。