r - MASS 包中 MVRNorm 背后的代数
问题描述
这是代码:
function (n = 1, mu, Sigma, tol = 1e-06, empirical = FALSE, EISPACK = FALSE)
{
p <- length(mu)
if (!all(dim(Sigma) == c(p, p)))
stop("incompatible arguments")
if (EISPACK)
stop("'EISPACK' is no longer supported by R", domain = NA)
eS <- eigen(Sigma, symmetric = TRUE)
ev <- eS$values
if (!all(ev >= -tol * abs(ev[1L])))
stop("'Sigma' is not positive definite")
X <- matrix(rnorm(p * n), n)
if (empirical) {
X <- scale(X, TRUE, FALSE)
X <- X %*% svd(X, nu = 0)$v
X <- scale(X, FALSE, TRUE)
}
X <- drop(mu) + eS$vectors %*% diag(sqrt(pmax(ev, 0)), p) %*%
t(X)
nm <- names(mu)
if (is.null(nm) && !is.null(dn <- dimnames(Sigma)))
nm <- dn[[1L]]
dimnames(X) <- list(nm, NULL)
if (n == 1)
drop(X)
else t(X)
}
我很好奇的问题是:
x <- eS$vectors %*% diag(sqrt(ev)) %*% t(x) # ignoring drop(mu)
...
t(x)
为什么会这样
X^T = UVZ^T,其中 Z 是标准化的 MVN?
我原以为这将是 X = UVZ,其中 X ~ MVN(0, UV(I)(UV)^T) = MVN(0, Sigma)?
回应 Siong Thye Goh 的回答:
我可以看到代数,并且仅通过考虑维度就可以做到这一点,但是考虑到多元正态的属性,转置一切的整个行为似乎很奇怪。即 X = UVZ
我做了一些审查,发现这实际上是一个Matrix Normal,并且那里的仿射变换以类似的方式工作。也就是说,X = Z (UV)^T。
我不确定我在理解这一点时是否缺少一些愚蠢的东西,或者我是否完全错过了关于为什么所有内容都被转置的图片,例如MVN 的 Wikipedias 仿射变换
解决方案
U 是 Sigma 的特征向量。即 Sigma = UV^2 U^T,其中 V 是对角矩阵。
让我们计算协方差矩阵 E[X^TX],看看它是否等于 Sigma,其中 X=UVZ^T 和 Z^T 满足 E[Z^TZ]=I,即单位矩阵。
我们有
E[X^TX]=E[UVZ^TZVU^T]=UVE[Z^TZ]VU^T=UV^2U^T=Sigma
推荐阅读
- continuous-integration - CircleCI 保存工作流中相关作业的输出
- webpack - 热模块更换 HMR 找不到更新。需要完全重新加载!创建反应应用程序和电子
- amazon-web-services - 共享的公共映像或社区映像是否可以像 AWS 中的私有映像那样是特定于区域的?
- angular - 如何让 uiSignal 处理嵌套的异步调用?
- python - Kolmogorov Smirnov 检验 python 中的拟合优度
- cuda - 如何在远程服务器上调试 CUDA 代码?
- django - 您是否忘记注册或加载此标签?
- php - 在 laravel-5.6 中找不到类“应用程序/帖子”
- azure - Microsoft Azure - 面向公众的负载均衡器后面的 2 个 VM
- asp.net-core - 通过托管 IdentityServer 的单独微服务使用 OATH 保护 WebApi 的集成测试