如何估計相關觀測值的方差？

假設我們有 n 個觀察值 xi （i 從 1 到 n），每個都來自均值為 0 和一些方差分量的正態分佈： Xi∼N(0,σ2) . 隨機變量 Xi s 有一些（假設已知）相關結構。IE： Corr(Xi,Xj)=ρij .

我們應該如何估計 σ2 ?

相關結構將如何影響估計量，例如 ˆσ2=∑(x2i)n ? 我們能否設計一個更好的估計器（假設相關結構已知）？

讓 R=(ρij) 是相關矩陣，使得協方差矩陣是 Σ=σ2R. 考慮 x=(x1,…,xn)′ 成為一個單一的觀察 n - 均值為零的變量正態分佈和 Σ 協方差。因為對數似然 N≥1 獨立的這種觀察是，直到一個附加常數（僅取決於 N, n, 和 R ) 由

Λ(σ)=N∑i=1(−nlogσ)−12σ2N∑i=1x′i,R−1,xi,

它有臨界點 σ→0, σ→∞, 並在任何解決方案

0=ddσΛ(σ)=−nNσ+1σ3N∑i=1x′i,R−1,xi.

除非表格 ∑Ni=1x′i,R−1,xi 為零，有一個唯一的全局最大值

ˆσ2=1nNN∑i=1x′i,R−1,xi.

這是最大似然估計。 即使在 N=1 （這是問題中提出的情況）。

直觀地說，這將優於任何忽略假設的相關性的估計 R. 為了檢查，我們可以計算這個估計值的 Fisher 信息矩陣——但我會把它留給你，部分原因是結果不應該令人信服 N （最大似然漸近結果可能不適用）。

**為了說明實際發生的情況，**這裡有一千個估計 σ2 從實驗 n=4 和 N=1. 的價值 σ 被設置為 1 始終。在這些實驗中，相關性總是

R=(10.30555690.55133770.5100989 0.305556910.12401510.09634469 0.55133770.124015111−0.4209064 0.51009890.09634469−0.42090641)

（一開始隨機生成）。在圖中，最左邊的面板是上述最大似然估計的直方圖；中間面板是使用通常（無偏）方差估計量的估計直方圖；右側面板是兩組估計值的 QQ 圖。斜線是平等線。您可以看到通常的方差估計器往往會產生更多的極值。 它也是有偏差的（由於忽略了相關性）：MLE 的平均值是 0.986 - 驚人地接近真實值 σ2=12=1 而通常估計的平均值僅為 0.791。（我寫“令人驚訝”是因為它是眾所周知的通常的最大似然估計量 σ2, 在不涉及相關性的情況下，具有順序偏差 1/(nN), 在這種情況下相當大。）

您可以通過修改、、、、和隨機數種子R的值來試驗生成這些數字的代碼。n``sigma``N``n.sim``Rho``17

f <- function(x, Rho) { # The MLE of sigma^2 given data `x` and correlation `Rho`
 S <- solve(Rho)
 sum(apply(x, 1, function(x) x %*% S %*% x)) / length(c(x))
}
n <- 4
sigma <- 1
N <- 1
n.sim <- 1e3
set.seed(17)
#
# Generate a random correlation matrix. Larger values of `d` yield more
# spherical matrices in general.
#
d <- 1
Rho <- cor(matrix(rnorm(n*(n+d)), ncol=n))
(ev <- eigen(Rho, only.values=TRUE)$values) 
#
# Run the experiments.
#
library(MASS)
sim <- replicate(n.sim, {
 x <- matrix(mvrnorm(N, rep(0,n), sigma^2 * Rho), N)
 c(f(x, Rho),  var(c(x)))
})
(rowMeans(sim))
#
# Plot the results.
#
par(mfrow=c(1,3))
hist(sim[1,], col=gray(.93), xlab="Estimate", 
    main=expression(paste("Histogram of Estimates of ", sigma^2)))
abline(v = sigma^2, col="Red", lwd=2)

hist(sim[2,],col=gray(.93), xlab="Estimate",
    main=expression(paste("Histogram of Independent Estimates of ", sigma^2)))
abline(v = sigma^2, col="Red", lwd=2)

y1 <- sort(sim[1,])
y2 <- sort(sim[2,])
plot(y1, y2,asp=1, xlab="Correlation-based estimate", ylab="Independent estimate")
abline(0:1, col="Red", lwd=2)
par(mfrow=c(1,1))

引用自：https://stats.stackexchange.com/questions/535727