双求和代码的奇怪输出

Question

我正在尝试编写代码来计算以下两个方程式：

到目前为止，我有以下代码：

set.seed(123)
a = rnorm(100, 0, 0.02)
b = rnorm(100, 0, 0.03)
c = rnorm(100, 0, 0.04)
dt = cbind(a, b, c)
wg = runif(ncol(dt))
wg = wg/sum(wg)
wg = round(wg, digits = 3)

options(scipen = 999)
cv = cov(dt)
VARs = diag(cv)
cr = cor(dt)
length(wg)

sum_1 = 0
for (i in 1 : length(wg) ) {
  for (j in min(i+1, length(wg)) ) {
    sum_1 = sum_1 + 2 * wg[i] * wg[j] * cr[i, j]
  }}
rho_1 = sum_1/(1 - sum(wg^2))

sum_2 = 0
for (i in 1 : length(wg) ) {
  for (j in min(i+1, length(wg)) ) {
    sum_2 = sum_2 + 2 * wg[i] * wg[j] * sqrt(VARs[i]) * sqrt(VARs[j])
  }}

PV = t(wg) %*% cv %*% wg
ss = sum(wg^2 * VARs)
rho_2 = (PV - ss)/sum_2

rhos = c(rho_1, rho_2)
rhos

[1]  1.388473356 -0.004104948

我希望这两个彼此接近。我认为我的代码中可能有错误。如果有人可以验证代码，我将不胜感激。

谢谢。

感谢李哲源 Zheyuan Li 的代码，我对这两个相关性度量模拟了 100 个观察结果，它们看起来确实很接近：

Answer 1

复杂求和的矢量化能力。

第一个等式：

# double summation in numerator (including multiplier 2)
diag(cr) <- 0
double_sum.1 <- c(crossprod(wg, cr %*% wg))

# single summation in denominator
single_sum.1 <- c(crossprod(wg))

# result
double_sum.1 / (1 - single_sum.1)

第二个方程：

dcv <- sqrt(diag(cv))

# single summation in numerator
single_sum.2 <- c(crossprod(wg * dcv))

# double summation in denominator (including multiplier 2)
double_sum.2 <- sum(tcrossprod(wg * dcv)) - single_sum.2

# result
PV <- c(crossprod(wg, cv %*% wg))
(PV - single_sum.2) / double_sum.2

修复原始循环：

i in 1:(length(wg)-1)
j in (i+1):length(wg)

备注

在我看来，这两个方程式应该不会给出接近的结果。

为了解这一点，我们现在在您的示例中使用真协方差矩阵和真相关矩阵。当您生成 3 组独立的正态样本时，它们将具有 对角线 协方差矩阵和 同一性 相关矩阵：

dcv <- c(0.02, 0.03, 0.04)
cr <- diag(3)

这立即意味着第一个等式会给出 0。但是第二个是（做一个快速的文书工作）

sum_i (1 - w_i) ^ 2 * sigma_i^2

这怎么可能是零？

双求和代码的奇怪输出

Weird output from double summation code

r

calculus

备注