R дисперсия на многочленно разпределение

В кода по-долу p_hat съдържа MLE на вероятностите за X1, X2 и X3 в дадената извадка от данни. Съгласно страницата многочленно разпределение в Уикипедия ковариационната матрица за приблизителните вероятности се изчислява по следния начин:

set.seed(102)
X <- rmultinom(n=1, size=100, prob =c(0.1,0.3,0.6))
p_hat <- X/sum(X)

# print covariance matrix
cov_matrix <- matrix(0, nrow=length(p_hat), ncol=length(p_hat))
rownames(cov_matrix) <- c("X1","X2","X3"); colnames(cov_matrix) <- c("X1","X2","X3");
for (r in 1: length(p_hat)){
  for (c in 1: length(p_hat)){
    if(r==c){cov_matrix[r,c] <- p_hat[r] * (1-p_hat[r])}
    else{cov_matrix[r,c] <- -p_hat[r] *p_hat[c]}
  }
}

Правилно ли е това изпълнение?

Има ли R функция, която може да произведе тази ковариационна матрица, дадена prob =c(0.1,0.3,0.6) за мултиномиално разпределение?


person Zhubarb    schedule 13.11.2013    source източник


Отговори (1)


Можете дори да използвате outer и diag, за да получите същия резултат

> p <- drop(p_hat)
> variance         <-  p*(1-p)
> covariance       <- -outer(p, p)
> diag(covariance) <-  variance
> covariance
       [,1]    [,2]    [,3]
[1,]  0.090 -0.0290 -0.0610
[2,] -0.029  0.2059 -0.1769
[3,] -0.061 -0.1769  0.2379
person Jilber Urbina    schedule 13.11.2013