Come ottenere il pieno distribuzione marginale di un parametro in stan

Question

quando si inizia un esempio standard dal stanwebpage come la seguente:

schools_code <- ' 
    data { 
    int<lower=0> J; // number of schools 
    real y[J]; // estimated treatment effects 
    real<lower=0> sigma[J]; // s.e. of effect estimates 
} 
parameters { 
    real theta[J]; 
    real mu; 
    real<lower=0> tau; 
} 
model { 
    theta ~ normal(mu, tau); 
    y ~ normal(theta, sigma); 
} 
' 

    schools_dat <- list(J = 8, 
       y = c(28, 8, -3, 7, -1, 1, 18, 12), 
       sigma = c(15, 10, 16, 11, 9, 11, 10, 18)) 

fit <- stan(model_code = schools_code, data = schools_dat, 
      iter = 1000, n_chains = 4) 

(questo è stato ottenuto da here)

Tuttavia questo solo fornirmi i quantili del posteriore dei parametri. quindi la mia domanda è: come ottenere altri percentili? I indovinare dovrebbe essere simile a bug

osservazione: ho cercato di introdurre il tag stan però, ho troppo poco reputazione;) spiacente per quello

Answer 1

Ecco il mio tentativo di speranza questo è (?) corretto:

supporre fit è un oggetto ottenuto da stan(...). poi la parte posteriore per qualsiasi percentile derivano:

quantile(fit@sim$sample[[1]]$beta, probs=c((1:100)/100)) 

dove il numero tra parentesi quadre è la catena immagino. nel caso questo non è stato chiaro: i uso rstan

Answer 2

Dal rstan v1.0.3 (non ancora rilasciato), sarà in grado di utilizzare la funzione cavallo di battaglia apply() direttamente su un oggetto di stanfit class che viene prodotto dal stan() function. Se la misura è un oggetto ottenuto da stan(), poi per esempio,

apply(fit, MARGIN = "parameters", FUN = quantile, probs = (1:100)/100) 

o

apply(as.matrix(fit), MARGIN = 2, FUN = quantile, probs = (1:100)/100) 

Il primo riguarda divertente per ciascun parametro in ciascuna catena, mentre il secondo unisce le catene prima di applicare FUN a ciascun parametro. Se ti interessava solo un parametro, allora qualcosa come

beta <- extract(fit, pars = "beta", inc_warmup = FALSE, permuted = TRUE)[[1]] 
quantile(beta, probs = (1:100)/100) 

è un'opzione.