R coda "Il principale minore di ordine 3 non è definito positivo"

Question

Qualcuno può spiegarmi cosa significa questo messaggio di errore?

ho un metodo di campionamento MCMC chiamato hitandrun (in pacchetto non finito https://github.com/davidkane9/kmatching) che mi dà un elenco di matrici aventi colonne come uscita multivariata dei singoli campioni, qui è un esempio della sua azione:

> A = matrix(1, ncol = 3) 
> b = 1 
> ## gives me solutions of Ax = b (a.k.a x + y + z = 1) 
> h = hitandrun(A,b, n=10, chains = 2) 
> h 
[[1]] 
      [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] 
[1,] 0.1804431 0.3340590 0.4195820 0.2061222 0.3591085 0.09984353 0.6707110 0.3926639 0.1283919 
[2,] 0.6135745 0.4256909 0.3619727 0.2918238 0.5057426 0.81919629 0.2368842 0.1178713 0.2666737 
[3,] 0.2059824 0.2402501 0.2184453 0.5020541 0.1351489 0.08096018 0.0924048 0.4894647 0.6049344 
     [,10] 
[1,] 0.1322112 
[2,] 0.4736057 
[3,] 0.3941831 

[[2]] 
      [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] 
[1,] 0.32883534 0.1284182 0.1735151 0.2005726 0.94511422 0.61653717 0.5130324 0.33228224 0.2088865 
[2,] 0.65868549 0.3066952 0.5182009 0.3065610 0.01214334 0.07007548 0.1191157 0.01137002 0.3311197 
[3,] 0.01247917 0.5648866 0.3082840 0.4928664 0.04274244 0.31338735 0.3678519 0.65634774 0.4599938 
      [,10] 
[1,] 0.61412223 
[2,] 0.32289039 
[3,] 0.06298738 

Volevo vedere la diagnostica Gelman-Rubin di questi dati, per vedere quanto ho bisogno di diluirlo, ma quando lo metto in funzione, ottengo un errore oscuro, e non so cosa significhi. Ecco che:

> mclist = lapply(h, function(x) mcmc(t(x), thin = 5)) 
> gelman.diag(mclist) 
Error in chol.default(W) : 
    the leading minor of order 1 is not positive definite 

(immagino ora è ordine 1, ma prima che fosse ordine 3) Ci sono esperti di Coda su SO? Ho provato a eseguire il debug, ma mi porta a una funzione interna La_chol e non so cosa fare da lì.

Answer 1

Sembra esserci qualcosa di sbagliato nell'ottenere la stima multivariata del diagnostico Gelman-Rubin. L'impostazione multivariate = FALSE risolve il problema e genera una stima di una singola variabile per ogni variabile. Tuttavia, la maggior parte delle mie variabili sono correlate a causa della natura del problema che sto cercando di risolvere, quindi mi chiedo (e spero) che questo mi stia dando una sopravvalutazione della diagnostica.

Answer 2

Recentemente ho avuto questo problema con il mio modello - ho trovato la mia matrice non era simmetrica, in quanto era una matrice GRM asimmetrica. Ho usato questa funzione per renderlo simmetrico e il problema è stato risolto (tratto da Most Efficient way to create a symmetric matrix):

symmetrise <- function(mat){ rowscols <- which(lower.tri(mat), arr.ind=TRUE) sparseMatrix(i=rowscols[,1], 
       j=rowscols[,2], 
       x=mat[lower.tri(mat)], 
       symmetric=TRUE, 
       dims=c(nrow(mat),ncol(mat))) } 

Questa replicato la parte inferiore sinistra alla parte superiore destra.