2017-04-11 4 views
0

J'utilise rjags en tant qu'échantillonneur. Le modèle a 3 matrices définies. La fonction coda.samples renvoie une liste d'échantillons. Si je prends la première liste d'exemples les noms de colonnes ressemblent à ceci:Reconstitution de variables à partir d'objets mcmc

> colnames(output[[1]]) 
"A[1,1]" "A[2,1]" "A[1,2]" "A[2,2]" ... 
"B[1,1]" "B[2,1]" "B[3,1]" "B[4,1]" ... 
"C[1,1]" "C[2,1]" 

De toute évidence, A, B et C sont des matrices dans mon modèle. Je veux les reconstruire en fonction de la moyenne de ces échantillons. Je peux facilement obtenir les moyens avec colMeans(output[[1]]) mais je n'ai aucune idée de comment reconstruire facilement les matrices à partir de ce vecteur.

Un bon moyen de reconstruction serait la fonction relist(). Donc, si j'avais les matrices A, B et C dans une liste L = list(A=A,B=B,C=C) alors je pourrais transformer cette liste en un vecteur avec unlist() et reconvertir avec relist(). Je suis à la recherche de quelque chose de similaire/readymade pour les objets mcmc, mais sans succès jusqu'à présent - je ne peux pas croire que je suis le premier à en avoir besoin. Évidemment, relist(colMeans(output[[1]])) ne fonctionne pas.

Quelqu'un peut m'aider à reconstruire?

Modifier: notez également que la fonction relist() n'a besoin que d'un squelette, donc l'extraction du squelette de colnames(output[[1]]) ferait également l'affaire. Ou est-ce que je me complique?

Répondre

0

Je ne pense pas que relist() fera l'affaire ...

Je suppose que votre objet output est un objet de classe mcmc.list, tel que défini dans le package R coda et output[[1]] est un objet de classe mcmc représentant l'échantillon pour la première chaîne MCMC.

Je suis assez sûr coda ne comprend pas que, par exemple. "A[1,1]" est une matrice JAGS, elle le gère simplement comme un nom de variable. Vous devrez donc parcourir les variables pertinentes et imposer la structure vous-même.

Idéalement, vous souhaitez envelopper cela dans une fonction comme ce qui suit:

getMatrix <- function(output, varname, rows, cols) { 
    unname(
    sapply(1:cols, function(j) 
     sapply(1:rows, function(i) 
     summary(output[,sprintf("%s[%s,%s]", varname, i, j)])[[1]][1] 
    ) 
    ) 
) 
} 

Ainsi, par exemple, si votre matrice B stockée dans output[[1]] a 3 lignes et 4 colonnes écrirait:

getMatrix(output[[1]], "B", 3, 4) 

pour obtenir les moyennes comme objet de matrice dans R.