-
-
Notifications
You must be signed in to change notification settings - Fork 269
RStan Getting Started (Português)
(Translation: Luiz Max Carvalho; revision: Marco Inácio, Daniel Marcelino)
RStan é a interface para o Stan no R. Para mais informações sobre o Stan e a sua linguagem de modelagem, visite o website em http://mc-stan.org/
Quase todas as instruções de instalação abaixo são para a versão do RStan citada acima, que necessita que você tenha o R na versão 3.4.0 ou mais recente. Se necessário, você pode baixar a versão mais recente do R aqui.
Além disso, recomendamos fortemente que você use o RStudio, versão 1.2.x ou mais recente dado que este tem ótimo suporte ao Stan.
Por precaução, muitas vezes é necessário remover qualquer Rstan existente no sistema executando o comando a seguir:
remove.packages("rstan")
if (file.exists(".RData")) file.remove(".RData")
Depois disso, reabra o R.
Na maioria dos casos você pode simplesmente digitar (exatamente dessa forma)
install.packages("rstan", repos = "https://cloud.r-project.org/", dependencies = TRUE)
No entanto, se você usa o Linux, ou se getOption("pkgType")
está configurado para "source"
, ou ainda, se o R perguntar se você deseja instalar a versão mais recente do Rstan a partir do código-fonte ("from source"), vá até a página correspondente para cada plataforma: Windows, Mac ou Linux.
No RStudio (preferencialmente) ou no R, execute uma vez
pkgbuild::has_build_tools(debug = TRUE)
para se certificar que o seu conjunto de ferramentas C++ (toolchain) está usando o pacote pkgbuild que é instalado junto com o Rstan.
Se essa linha de código retorna TRUE
, então sua toolchain está instalada corretamente e você pode pular para a próxima seção.
Caso contrário,
- Se você utiliza Windows e RStudio (recomendado), uma janela pop-up vai aparecer perguntando se você deseja instalar o Rtools. Clique em Yes/Sim e espere a instalação terminar.
- Se você utiliza Windows mas não o RStudio, uma mensagem aparecerá no console do R dizendo-lhe para instalar o Rtools. Mais informações sobre baixar e instalar rtools podem ser úteis, mas você normalmente NÃO precisa ir para a seção intitulada "Instalando o RStan a partir da fonte".
- Se você utiliza Mac, um link irá aparecer, mas não clique nele. Em vez disso, vá aqui
- Se você utiliza Linux (incluindo o subsistema Windows para Linux), então vá aqui.
Se você seguir as instruções acima mas não obtiver sucesso, pode obter ajuda no Forum do Stan no Discourse mas por favor se certifique de mencionar qual o seu sistema operacional e também de informar se você está utilizando ou não o RStudio e qual o output quando você tenta as instruções acima.
Esse passo é opcional, mas pode resultar em programas Stan compilados que rodam muito mais rápido do que rodariam de outra forma. Simplesmente cole o código abaixo no R
dotR <- file.path(Sys.getenv("HOME"), ".R")
if (!file.exists(dotR)) dir.create(dotR)
M <- file.path(dotR, ifelse(.Platform$OS.type == "windows", "Makevars.win", "Makevars"))
if (!file.exists(M)) file.create(M)
cat("\nCXX14FLAGS=-O3 -march=native -mtune=native",
if( grepl("^darwin", R.version$os)) "CXX14FLAGS += -arch x86_64 -ftemplate-depth-256" else
if (.Platform$OS.type == "windows") "CXX11FLAGS=-O3 -march=native -mtune=native" else
"CXX14FLAGS += -fPIC",
file = M, sep = "\n", append = TRUE)
No entanto, esteja ciente que mudar o nível de otimização para O3
pode causar problemas para outros pacotes além do RStan e que, em casos raros, especificar -march=native -mtune=native
pode fazer com que programas Stan não funcionem. Caso você futuramente venha a mudar as configurações da sua toolchain C++, execute
M <- file.path(Sys.getenv("HOME"), ".R", ifelse(.Platform$OS.type == "windows", "Makevars.win", "Makevars"))
file.edit(M)
O restante deste documento assume que você já instalou o Rstan seguindo as instruções acima.
O nome do pacote é rstan (minúsculas), então começaremos executando
library("rstan") # observe as mensagens de inicialização
Como diz a mensagem de inicialização, se você utiliza o rstan localmente numa máquina com múltiplos processadores e tem bastante memória RAM para estimar o seu modelo em paralelo, execute
options(mc.cores = parallel::detectCores())
Além disso, você deve seguir a segunda mensagem de inicialização que pede que você execute
rstan_options(auto_write = TRUE)
que permite que você salve uma versão compilada do seu programa Stan no disco rígido, de forma que este não precise ser compilado novamente (a menos que você o mude). Por último, se você utiliza Windows, haverá uma terceira mensagem de inicialização pedindo que execute
Sys.setenv(LOCAL_CPPFLAGS = '-march=native')
que não é necessário se você seguiu os conselhos sobre a toolchain C++ na seção anterior.
Esse é um exemplo da seção 5.5 de Gelman et al (2003), que estudam os efeitos de programas educacionais sobre a performance dos estudantes no Scholastic Aptitude Test (SAT) em oito escolas. Chamaremos este exemplo de "eight schools".
Começaremos escrevendo um programa Stan com o modelo num arquivo de texto. Se você está utilizando o RStudio na versão 1.2.x ou mais recente, clique em Arquivo -> Novo arquivo -> Arquivo Stan.
Se não, abra seu editor de texto favorito.
De qualquer forma, cole o código abaixo e salve seu trabalho em um arquivo chamado 8schools.stan
no diretório de trabalho (nota: você pode obter qual o diretório de trabalho atual rodando getwd()
).
// saved as 8schools.stan
data {
int<lower=0> J; // número de escolas
real y[J]; // estimativas dos efeitos dos tratamentos
real<lower=0> sigma[J]; // erro padrão das estimativas
}
parameters {
real mu; // efeito populacional do tratamento
real<lower=0> tau; // desvio padrão dos efeitos do tratamento
vector[J] eta; // desvio de mu por escola (não-normalizado)
}
transformed parameters {
vector[J] theta = mu + tau * eta; // efeitos do tratamento por escola
}
model {
target += normal_lpdf(eta | 0, 1); // log-densidade a priori
target += normal_lpdf(y | theta, sigma); // log-verossimilhança
}
Certifique-se de que seu programa Stan termina com uma linha em branco, sem nenhum caractere, incluindo espaços e comentários.
Nesse programa Stan, nós vamos usar theta
como a "transformada" de mu
, eta
e tau
(no bloco transformed parameters
) em vez de declarar theta
diretamente no bloco parameters
, o que permite que o amostrador rode mais eficientemente (veja aqui a explicação detalhada).
Podemos preparar os dados -- que geralmente são uma lista com nomes (named list) -- no R usando:
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))
E nós podemos obter um ajuste com o seguinte comando R. Note que o argumento file =
deve apontar para onde o arquivo 8schools.stan
está no seu sistema, a menos que você o tenha posto no diretório onde está o R, caso em que o comando abaixo deverá funcionar.
fit <- stan(file = '8schools.stan', data = schools_dat)
O objeto fit
, retornado pela função stan
é um objeto S4 de classe stanfit
.
Métodos como print
, plot
, e pairs
estão associados com o resultado ajustado de modo que podemos utilizar o código que se segue para avaliar os resultados em fit
.
O método print
oferece um sumário do parâmetro do modelo e também da densidade log-posteriori sob o nome de lp__
(ver o output de exemplo que se segue).
Para mais métodos e detalhes da classe stanfit
, veja a ajuda da classe stanfit
.
Em particular, podemos usar a função extract
em objetos stanfit
para obter as amostras.
O método extract
extrai as amostras de um objeto stanfit
na forma de uma lista de arrays para os parâmetros de interesse, ou apenas uma array.
Além disso, as funções S3 as.array
, as.matrix
, e as.data.frame
estão definidas para objetos stanfit
(utilize help("as.array.stanfit")
para olhar a documentação de ajuda no R).
print(fit)
plot(fit)
pairs(fit, pars = c("mu", "tau", "lp__"))
la <- extract(fit, permuted = TRUE) # retorna uma lista de arrays
mu <- la$mu
### retorna um array de três dimensões: iterações, cadeias, parâmetros (iterations, chains, parameters)
a <- extract(fit, permuted = FALSE)
### use funções S3 em objetos
a2 <- as.array(fit)
m <- as.matrix(fit)
d <- as.data.frame(fit)
Outro exemplo popular é o "Rats". Podemos encontrar, por exemplo, a versão OpenBUGS neste link, que é originalmente de Gelfand et al (1990).
Os dados são sobre o crescimento semanal de 30 ratos por cinco semanas.
Listamos os dados na tabela abaixo, em que que utilizamos x
para denotar as datas em que os dados foram coletados. Podemos testar esse exemplo usando os dados em rats.txt e o código Stan para modelo em rats.stan.
Rat | x=8 | x=15 | x=22 | x=29 | x=36 | Rat | x=8 | x=15 | x=22 | x=29 | x=36 | |
---|---|---|---|---|---|---|---|---|---|---|---|---|
1 | 151 | 199 | 246 | 283 | 320 | 16 | 160 | 207 | 248 | 288 | 324 | |
2 | 145 | 199 | 249 | 293 | 354 | 17 | 142 | 187 | 234 | 280 | 316 | |
3 | 147 | 214 | 263 | 312 | 328 | 18 | 156 | 203 | 243 | 283 | 317 | |
4 | 155 | 200 | 237 | 272 | 297 | 19 | 157 | 212 | 259 | 307 | 336 | |
5 | 135 | 188 | 230 | 280 | 323 | 20 | 152 | 203 | 246 | 286 | 321 | |
6 | 159 | 210 | 252 | 298 | 331 | 21 | 154 | 205 | 253 | 298 | 334 | |
7 | 141 | 189 | 231 | 275 | 305 | 22 | 139 | 190 | 225 | 267 | 302 | |
8 | 159 | 201 | 248 | 297 | 338 | 23 | 146 | 191 | 229 | 272 | 302 | |
9 | 177 | 236 | 285 | 350 | 376 | 24 | 157 | 211 | 250 | 285 | 323 | |
10 | 134 | 182 | 220 | 260 | 296 | 25 | 132 | 185 | 237 | 286 | 331 | |
11 | 160 | 208 | 261 | 313 | 352 | 26 | 160 | 207 | 257 | 303 | 345 | |
12 | 143 | 188 | 220 | 273 | 314 | 27 | 169 | 216 | 261 | 295 | 333 | |
13 | 154 | 200 | 244 | 289 | 325 | 28 | 157 | 205 | 248 | 289 | 316 | |
14 | 171 | 221 | 270 | 326 | 358 | 29 | 137 | 180 | 219 | 258 | 291 | |
15 | 163 | 216 | 242 | 281 | 312 | 30 | 153 | 200 | 244 | 286 | 324 |
y <- as.matrix(read.table('https://raw.github.com/wiki/stan-dev/rstan/rats.txt', header = TRUE))
x <- c(8, 15, 22, 29, 36)
xbar <- mean(x)
N <- nrow(y)
T <- ncol(y)
rats_fit <- stan('https://raw.githubusercontent.com/stan-dev/example-models/master/bugs_examples/vol1/rats/rats.stan')
Você pode rodar muitos dos exemplos do BUGS - e outros que nós criamos - no Stan executando
model <- stan_demo()
e escolhendo um modelo de exemplo a partir da lista que aparece. A primeira vez que você chamar stan_demo()
, o programa irá lhe perguntar se você deseja baixar esses exemplos. Você deve escolher a opção 1 para colocar esses exemplos no diretório em que o rstan foi instalado de modo que eles possam ser utilizados no futuro sem que você precise baixá-los novamente.
O objeto model
acima é uma instância da classe stanfit
, então você pode chamar print
, plot
, pairs
, extract
, etc no mesmo posteriormente.
Mais detalhes sobre o Rstan podem ser encontrados na documentação, incluindo a vignette do pacote rstan.
Por exemplo, você pode usar help(stan)
e help("stanfit-class")
para ver a ajuda para a função stan
e a classe S4 stanfit
.
Veja também Stan's modeling language manual (em inglês) para mais detalhes sobre as rotinas de amostragem, otimização e a linguagem de modelagem Stan.
Além disso, a Lista de email dos usuários de Stan pode ser utilizada para discutir o uso do Stan, postar exemplos ou perguntar sobre o Stan ou o Rstan. Quando precisar de ajuda, é importante que você ofereça informação suficiente, por exemplo:
- sintaxe formatada apropriadamente na linguagem de modelagem Stan (utilize "```" para delimitar o seu código)
- dados
- código R necessário
- Um dump da mensagem de erro utilizando
verbose=TRUE
ecores=1
quando chamar as funçõesstan
ousampling
- informação sobre o R que você utiliza, com o output da função
sessionInfo()
no R
- Gelman, A., Carlin, J. B., Stern, H. S., and Rubin, D. B. (2003). Bayesian Data Analysis, CRC Press, London, 2nd Edition.
- Stan Development Team. Stan Modeling Language User's Guide and Reference Manual.
- Gelfand, A. E., Hills S. E., Racine-Poon, A., and Smith A. F. M. (1990). "Illustration of Bayesian Inference in Normal Data Models Using Gibbs Sampling", Journal of the American Statistical Association, 85, 972-985.
- Stan
- R
- BUGS
- OpenBUGS
- JAGS
- Rcpp