分散共分散行列の逆行列

質問

以下の行列を考えます．
行列X（M行N列，成分の値は分散１のガウスノイズ）
行列Xの分散共分散行列S（M行M列），
行列Sの逆行列Y（M行M列）

N＝M+1のときはN≠M+1のときに比べて行列Yの各成分の大きさがかなり大きくなる，という結果が得られました．
（計算にはMATLABの関数covとpinvを用いました．）
これは数学的に正しいことなのでしょうか．
また，正しいとしたらどうしてこのようなことが起こるのでしょうか．

行列Yの計算結果の例を以下に示します．（行列Yを500回算出し，それらの平均値を示しています．）
http://wisteria.orz.ne.jp/download/pinvcovX.jpg
縦・横方向の軸は行列Yのインデックスを，高さ方向の軸は行列Yの成分の値を示しています．
上述した現象の原因についてご教授いただけると幸いです．

回答No.2

noname#227064 · ベストアンサー

Rで試してみたら同じような結果になりました。
その結果が添付画像です。
それぞれのz軸の上限は「z:上限」で表しています。

繰り返し回数を増やしたら各成分が小さくなると思っていたのですが、逆に大きくなりました。

Nが大きくなるにつれ、単位行列に近づいていくようです。
これ自体は分散共分散行列が単位行列に近づくことから当然ですか。

試したコード
func1 <- function(x, M)
{
y <- matrix(x, ncol = M)
z <- cov(y)
c(as.vector(z), as.vector(solve(z)))
}

func2 <- function(M = 20, N = 21, R = 500)
{
x <- matrix(rnorm(M*N*R), nrow = R)
y <- apply(x, 1, function(t) func1(t, M))
z <- rowMeans(y)
list(cov = matrix(z[1:(M*M)], ncol = M), covinv = matrix(z[(M*M+1):(2*M*M)], ncol = M))
}

x <- func2(20, 21, 10000)
y <- func2(20, 21, 500)
z <- func2(20, 22, 500)
w <- func2(20, 100, 500)

par(mfrow = c(2, 2))
persp(x$covinv, zlim = c(-0.5, 10^10), xlab = "x", main = "M = 20, N = 21, 10000回の平均, z:10^10")
persp(y$covinv, zlim = c(-0.5, 10^5), xlab = "x", main = "M = 20, N = 21, 500回の平均, z:10^5")
persp(z$covinv, zlim = c(-0.5, 10^3), xlab = "x", main = "M = 20, N = 22, 500回の平均, z:10^3")
persp(w$covinv, zlim = c(-0.5, 10), xlab = "x", main = "M = 20, N = 100, 500回の平均, z:10")
par(mfrow = c(1, 1))

お礼

検証していただきありがとうございます。
そうですか、計算結果は同じようになりましたか。

> 繰り返し回数を増やしたら各成分が小さくなると思っていたのですが、逆に大きくなりました。

有限長の小数で演算を繰り返すほどに、誤差が膨らんでいくという可能性が考えられますね。

私が用いた言語と異なる言語で検証していただき、たいへん参考になりました。ありがとうございました。

回答No.3

ICE_FALCON

M列N行だと思いますが、
NがMに対して少ないので、
独立したガウスノイズになっているか疑問です。
なので逆行列が正確に求められるかも疑問です。
・・・がinvでも似たような結果になるのでは・・・・。

お礼

> M列N行だと思いますが、
> NがMに対して少ないので、
> 独立したガウスノイズになっているか疑問です。

N行M列でした。失礼しました。独立したガウスノイズになっていることも確認しました。
ご指摘ありがとうございます。

回答No.1

rabbit_cat

＞行列X（M行N列，成分の値は分散１のガウスノイズ）
これは、Ｎ行Ｍ列の間違いですか？
で、発生させたガウスノイズ自体は独立（共分散＝０）ということなんですかね。
ということは、求めた分散共分散行列は、ほとんど単位行列に近いものと思ってよいですか？

多分、その現象は、数値計算上の不安定性（誤差）によるものだと思います。Ｓ行列は、基本的には逆行列が存在するんだと思いますが、pinvではなくて、普通にinvを使えばどうでしょう。

お礼

ご回答ありがとうございます。

> Ｎ行Ｍ列の間違いですか？

はい、間違いでした。N行M列ですね。失礼しました。

> 発生させたガウスノイズ自体は独立（共分散＝０）ということなんですかね。
> ということは、求めた分散共分散行列は、ほとんど単位行列に近いものと思ってよいですか？

はい、仰るとおりです。

> 多分、その現象は、数値計算上の不安定性（誤差）によるものだと思います。Ｓ行列は、基本的には逆行列が存在す> るんだと思いますが、pinvではなくて、普通にinvを使えばどうでしょう。

pinvの場合とinvの場合で計算結果はほぼ同じでした。
ただ、現象の原因としては、rabbit_catさんの仰るように、数値の誤差が膨らんだことは十分に考えられます。