RStan – это реализация языка программирования Stan, используемого для построения статистических моделей на основе байесовского подхода, в R. Подробнее обо всех этих страшных словах можно почитать здесь, здесь и здесь (осторожно, возможен английский).

При составлении данной инструкции использовались материалы англоязычной инструкции по установке с сайта проекта mc-stan.org, а также материалы семинаров Б. Б. Демешева по эконометрике.

Как установить RStan

Если у Вас Windows

Есть мнение, что проприетарные системы (читай – Windows и Mac OS) вредят открытости научного знания и разобщают исследователей, однако если они установлены на большинстве компьютеров в мире, то с ними придётся считаться. Данный раздел посвящён процессу установки среды RStan на компьютер с ОС Windows (начиная с Windows 7; более старые версии могут поддерживаться, но работа на них не гарантируется).

Требования:

Сам R. Если у Вас его нет, то ~~что вы тут делаете? :)~~ идите и качайте последнюю версию с https://www.r-project.org.
RStudio – предполагается, что он у вас тоже уже есть; в противном случае качаем (хотя его наличие/отсутствие не критично).
Набор инструментов Rtools. Этот набор содержит компилятор C++, необходимый для создания исполняемого файла на основе модели на языке Stan (и не только для этого).

Скачиваем набор отсюда: https://cran.r-project.org/bin/windows/Rtools/. Нужно выбрать наиболее свежую версию, поддерживающую Вашу версию R (например, Rtools34.exe, если у Вас R 3.3.0 и выше).
Устанавливаем. Как уже было сказано, Rtools предназначен не только для работы со Stan, поэтому установочный комплект содержит массу всего нужного и не совсем нужного. Дабы не захлямлять свой жёсткий диск лишними гига(?)байтами, рекомендую ставить следующие параметры установки:

Для 64-битной версии Windows

Для 32-битной версии Windows

На следующем этапе (в окне “Выберите дополнительные задачи”) отметьте галочкой опцию “Edit the system PATH” (нужно, чтобы впоследствии R “увидел” компилятор C++):

Именно так

Дальше ничего не меняем и ждем завершения процесса установки.
Проверяем правильность установки.

Перезапускаем R/RStudio и вбиваем следующие команды:

Sys.getenv('PATH')
system('g++ -v')
system('where make')

В итоге должно получиться следующее:

Должно быть так (выделено существенное)

Если что-то идёт не так (например, появляются сообщения об ошибках и пр.) – смело переустанавливайте Rtools.

Установка RStan (наконец-то!)

Откройте/перезапустите R/RStudio.
Установите Rstan из репозитория CRAN, запустив эту команду (в точности!):

install.packages("rstan", repos = "https://cloud.r-project.org", dependencies = TRUE)

В случае облома попробуйте другой вариант:

Sys.setenv(MAKEFLAGS = "-j4") # вместо 4 поставьте количество ядер Вашего процессора
install.packages("rstan", type = "source")

Снова перезапустите R/RStudio.
Наконец, проверьте правильность установки RStan, запустив ещё одну команду:

fx <- inline::cxxfunction( signature(x = "integer", y = "numeric" ) , '
    return ScalarReal( INTEGER(x)[0] * REAL(y)[0] ) ;
' )
fx( 2L, 5 ) # должно получиться 10

Готово!

Если у Вас Mac OS или Linux

Данный раздел посвящён процессу установки среды RStan на компьютеры от Apple, а также на истинно свободные компьютеры с истинно свободной ОС Linux :)

MacOS

Огромное спасибо Никите Герману за предоставленный MacBook :)

Требования:

R.
RStudio (наличие/отсутствие не критично).
Набор инструментов Xcode. Этот набор для Mac OS содержит компилятор C++, необходимый для создания исполняемого файла на основе модели на языке Stan. Установить его можно следующим образом:

Сначала нужно проверить используемую версию Mac OS. Для этого нужно нажать на логотип Apple в левом верхнем углу и в выпадающем меню выбрать пункт “Об этом Mac” (“About this Mac”).

Пример

Через Spotlight проверяем, установлен ли Xcode на Вашем Mac:
- Нажимаем клавишу control и пробел.
- В появившуюся строку вбиваем Xcode.

Если появляется что-то подобное изображённому на картинке ниже, то эта программа у Вас есть и Вам нужно просто её обновить через App Store. Если же нет – переходим к следующему шагу.

Пример

Для OS X 10.8 “Mountain Lion” и выше:
- Откройте App Store и найдите через поиск программу Xcode.
- Выберите программу и выполните все шаги по её установке (можно оставить все параметры по умолчанию).
- Запустите Xcode и примите лицензионное соглашение.
Для OS X 10.7 “Lion” или Mac OS X 10.6 “Snow Leopard”

Этот раздел перекочевал сюда из оригинальной инструкции.
- Зарегистрируйтесь на сайте для разработчиков, если Вы ещё не сделали этого.
- Перейдите на сайт Downloads Page for Developers; если нужно – войдите с помощью Apple ID.
- Скачайте и установите наиболее свежую версию Xcode.
- Откройте R/RStudio и наберите команду:

system('clang++ -v')

Если компилятор C++ установлен, то появится номер версии; если же нет, то вернитесь на страницу Downloads Page for Developers и скачайте последнюю версию Command Line Tools.

Переходите к установке RStan.

Установка RStan

Откройте/перезапустите R/RStudio.
Установите Rstan из репозитория CRAN, запустив эту команду (в точности!):

install.packages("rstan", repos = "https://cloud.r-project.org", dependencies = TRUE)

В случае облома попробуйте другой вариант:

install.packages("rstan", type = "source")

Снова перезапустите R/RStudio.
Наконец, проверьте правильность установки RStan, запустив ещё одну команду:

fx <- inline::cxxfunction( signature(x = "integer", y = "numeric" ) , '
    return ScalarReal( INTEGER(x)[0] * REAL(y)[0] ) ;
' )
fx( 2L, 5 ) # должно получиться 10

Готово!

Linux

Эх, всё через одно место – консоль… или менеджер пакетов.

Требования:

R.
RStudio (наличие/отсутствие не критично).
С помощью менеджера пакетов ставим build-essential, libssl-dev (строго до версии 1.0), а также свежую версию g++ или clang++ (на выбор).

Установка RStan

Откройте/перезапустите R/RStudio.
В R/RStudio задайте количество ядер процессора, которое Вы хотите выделить для программы:

Sys.setenv(MAKEFLAGS = "-j4") ## поставьте нужное количество вместо 4

Установите Rstan из репозитория CRAN, запустив эту команду (в точности!):

install.packages("rstan", repos = "https://cloud.r-project.org", dependencies = TRUE)

В случае облома попробуйте другой вариант:

install.packages("rstan", type = "source")

Снова перезапустите R/RStudio.
Наконец, проверьте правильность установки RStan, запустив ещё одну команду:

fx <- inline::cxxfunction( signature(x = "integer", y = "numeric" ) , '
    return ScalarReal( INTEGER(x)[0] * REAL(y)[0] ) ;
' )
fx( 2L, 5 ) # должно получиться 10

Готово!

Как пользоваться RStan

Запускаем R/RStudio.
Необходимый пакет называется rstan. Подгружаем его с помощью следующей команды:

library('rstan')

Появляющееся после выполнения команды сообщение предлагает установить режим параллельных вычислений для ускорения процесса (рекомендуется для компьютеров с многоядерными процессорами и большим объёмом оперативной памяти (от 6 Гбайт)). Эти команды также активируют режим записи скомпилированных программ на жёсткий диск (полезно, если не хотите каждый раз компилировать программу заново для разных данных в одной модели).

rstan_options(auto_write = TRUE)
options(mc.cores = parallel::detectCores())

В дальнейшем можете прописывать эти команды непосредственно в файле скрипта R, с помощью которого выполняется анализ данных.

Создаём текстовый файл и сохраняем его с расширением .stan, затем создаём стандартный скрипт R.

Структура файла .stan включает в себя несколько основных блоков:

data{содержимое} – определяет внешние данные, которые мы “скармливаем” модели;
parameters{содержимое} – параметры модели (случайные величины), подчиняющиеся некоему заранее определённому (априорному) распределению;
model{содержимое} – сама модель;
generated quantities{содержимое} – блок для построения прогнозных значений.

В файле .R должны содержаться команды для компиляции программы, передачи ей данных и анализа полученных результатов.

Для примера решим в STAN задачу Сэра Томаса Байеса.

Монетку с неизвестной вероятностью выпадения орла \(p\) подбросили три раза. И все три раза она выпала орлом. Какова вероятность того, что она выпадет орлом в следующий раз? Какова вероятность того, что \(p\) больше половины с учётом имеющихся наблюдений?

Модель для данных проста:

\[ y_n|p \sim Bernoulli(p) \]

Априорное мнение о \(p\) возьмём равномерное на \([0;1]\):

\[ p \sim Uniform[0;1] \]

Создаём файл bayes_coins.stan:

data {
  int N; // число наблюдений
  int y[N]; // отдельные наблюдения
}
parameters {
  real<lower=0, upper=1> p; // вероятность орла
}
model {
  // априорно:
  p ~ uniform(0, 1);
  // модель: как наблюдения связаны с параметром
  for (n in 1:N) {
    y[n] ~ bernoulli(p);
  }
}
generated quantities {
  int y_new; // результат следующего броска
  y_new = bernoulli_rng(p); // генерация будущего игрека
}

Благородные доны и доньи не забывают перенос строки в последней строке :)

И теперь из R сделаем всё нужное:

library(rstan)
library(bayesplot)

model <- stan_model(file = "coins.stan")

## In file included from file2284732993c4.cpp:8:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/src/stan/model/model_header.hpp:4:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math.hpp:4:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/rev/mat.hpp:4:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/rev/core.hpp:12:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/rev/core/gevv_vvv_vari.hpp:5:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/rev/core/var.hpp:7:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/BH/include/boost/math/tools/config.hpp:13:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/BH/include/boost/config.hpp:39:
## /Library/Frameworks/R.framework/Versions/3.3/Resources/library/BH/include/boost/config/compiler/clang.hpp:196:11: warning: 'BOOST_NO_CXX11_RVALUE_REFERENCES' macro redefined [-Wmacro-redefined]
## #  define BOOST_NO_CXX11_RVALUE_REFERENCES
##           ^
## <command line>:6:9: note: previous definition is here
## #define BOOST_NO_CXX11_RVALUE_REFERENCES 1
##         ^
## In file included from file2284732993c4.cpp:8:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/src/stan/model/model_header.hpp:4:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math.hpp:4:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/rev/mat.hpp:4:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/rev/core.hpp:42:
## /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/rev/core/set_zero_all_adjoints.hpp:14:17: warning: unused function 'set_zero_all_adjoints' [-Wunused-function]
##     static void set_zero_all_adjoints() {
##                 ^
## In file included from file2284732993c4.cpp:8:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/src/stan/model/model_header.hpp:4:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math.hpp:4:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/rev/mat.hpp:4:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/rev/core.hpp:43:
## /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/rev/core/set_zero_all_adjoints_nested.hpp:17:17: warning: 'static' function 'set_zero_all_adjoints_nested' declared in header file should be declared 'static inline' [-Wunneeded-internal-declaration]
##     static void set_zero_all_adjoints_nested() {
##                 ^
## In file included from file2284732993c4.cpp:8:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/src/stan/model/model_header.hpp:4:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math.hpp:4:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/rev/mat.hpp:11:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/prim/mat.hpp:60:
## /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/prim/mat/fun/autocorrelation.hpp:17:14: warning: function 'fft_next_good_size' is not needed and will not be emitted [-Wunneeded-internal-declaration]
##       size_t fft_next_good_size(size_t N) {
##              ^
## In file included from file2284732993c4.cpp:8:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/src/stan/model/model_header.hpp:4:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math.hpp:4:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/rev/mat.hpp:11:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/prim/mat.hpp:299:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/prim/arr.hpp:36:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/StanHeaders/include/stan/math/prim/arr/functor/integrate_ode_rk45.hpp:13:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/BH/include/boost/numeric/odeint.hpp:61:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/BH/include/boost/numeric/odeint/util/multi_array_adaption.hpp:29:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/BH/include/boost/multi_array.hpp:21:
## In file included from /Library/Frameworks/R.framework/Versions/3.3/Resources/library/BH/include/boost/multi_array/base.hpp:28:
## /Library/Frameworks/R.framework/Versions/3.3/Resources/library/BH/include/boost/multi_array/concept_checks.hpp:42:43: warning: unused typedef 'index_range' [-Wunused-local-typedef]
##       typedef typename Array::index_range index_range;
##                                           ^
## /Library/Frameworks/R.framework/Versions/3.3/Resources/library/BH/include/boost/multi_array/concept_checks.hpp:43:37: warning: unused typedef 'index' [-Wunused-local-typedef]
##       typedef typename Array::index index;
##                                     ^
## /Library/Frameworks/R.framework/Versions/3.3/Resources/library/BH/include/boost/multi_array/concept_checks.hpp:53:43: warning: unused typedef 'index_range' [-Wunused-local-typedef]
##       typedef typename Array::index_range index_range;
##                                           ^
## /Library/Frameworks/R.framework/Versions/3.3/Resources/library/BH/include/boost/multi_array/concept_checks.hpp:54:37: warning: unused typedef 'index' [-Wunused-local-typedef]
##       typedef typename Array::index index;
##                                     ^
## 8 warnings generated.

df <- list(N = 2, y = c(1, 1))
fit <- sampling(model, data = df)

## 
## SAMPLING FOR MODEL 'coins' NOW (CHAIN 1).
## 
## Chain 1, Iteration:    1 / 2000 [  0%]  (Warmup)
## Chain 1, Iteration:  200 / 2000 [ 10%]  (Warmup)
## Chain 1, Iteration:  400 / 2000 [ 20%]  (Warmup)
## Chain 1, Iteration:  600 / 2000 [ 30%]  (Warmup)
## Chain 1, Iteration:  800 / 2000 [ 40%]  (Warmup)
## Chain 1, Iteration: 1000 / 2000 [ 50%]  (Warmup)
## Chain 1, Iteration: 1001 / 2000 [ 50%]  (Sampling)
## Chain 1, Iteration: 1200 / 2000 [ 60%]  (Sampling)
## Chain 1, Iteration: 1400 / 2000 [ 70%]  (Sampling)
## Chain 1, Iteration: 1600 / 2000 [ 80%]  (Sampling)
## Chain 1, Iteration: 1800 / 2000 [ 90%]  (Sampling)
## Chain 1, Iteration: 2000 / 2000 [100%]  (Sampling)
##  Elapsed Time: 0.010419 seconds (Warm-up)
##                0.010141 seconds (Sampling)
##                0.02056 seconds (Total)
## 
## 
## SAMPLING FOR MODEL 'coins' NOW (CHAIN 2).
## 
## Chain 2, Iteration:    1 / 2000 [  0%]  (Warmup)
## Chain 2, Iteration:  200 / 2000 [ 10%]  (Warmup)
## Chain 2, Iteration:  400 / 2000 [ 20%]  (Warmup)
## Chain 2, Iteration:  600 / 2000 [ 30%]  (Warmup)
## Chain 2, Iteration:  800 / 2000 [ 40%]  (Warmup)
## Chain 2, Iteration: 1000 / 2000 [ 50%]  (Warmup)
## Chain 2, Iteration: 1001 / 2000 [ 50%]  (Sampling)
## Chain 2, Iteration: 1200 / 2000 [ 60%]  (Sampling)
## Chain 2, Iteration: 1400 / 2000 [ 70%]  (Sampling)
## Chain 2, Iteration: 1600 / 2000 [ 80%]  (Sampling)
## Chain 2, Iteration: 1800 / 2000 [ 90%]  (Sampling)
## Chain 2, Iteration: 2000 / 2000 [100%]  (Sampling)
##  Elapsed Time: 0.010196 seconds (Warm-up)
##                0.009739 seconds (Sampling)
##                0.019935 seconds (Total)
## 
## 
## SAMPLING FOR MODEL 'coins' NOW (CHAIN 3).
## 
## Chain 3, Iteration:    1 / 2000 [  0%]  (Warmup)
## Chain 3, Iteration:  200 / 2000 [ 10%]  (Warmup)
## Chain 3, Iteration:  400 / 2000 [ 20%]  (Warmup)
## Chain 3, Iteration:  600 / 2000 [ 30%]  (Warmup)
## Chain 3, Iteration:  800 / 2000 [ 40%]  (Warmup)
## Chain 3, Iteration: 1000 / 2000 [ 50%]  (Warmup)
## Chain 3, Iteration: 1001 / 2000 [ 50%]  (Sampling)
## Chain 3, Iteration: 1200 / 2000 [ 60%]  (Sampling)
## Chain 3, Iteration: 1400 / 2000 [ 70%]  (Sampling)
## Chain 3, Iteration: 1600 / 2000 [ 80%]  (Sampling)
## Chain 3, Iteration: 1800 / 2000 [ 90%]  (Sampling)
## Chain 3, Iteration: 2000 / 2000 [100%]  (Sampling)
##  Elapsed Time: 0.010261 seconds (Warm-up)
##                0.010913 seconds (Sampling)
##                0.021174 seconds (Total)
## 
## 
## SAMPLING FOR MODEL 'coins' NOW (CHAIN 4).
## 
## Chain 4, Iteration:    1 / 2000 [  0%]  (Warmup)
## Chain 4, Iteration:  200 / 2000 [ 10%]  (Warmup)
## Chain 4, Iteration:  400 / 2000 [ 20%]  (Warmup)
## Chain 4, Iteration:  600 / 2000 [ 30%]  (Warmup)
## Chain 4, Iteration:  800 / 2000 [ 40%]  (Warmup)
## Chain 4, Iteration: 1000 / 2000 [ 50%]  (Warmup)
## Chain 4, Iteration: 1001 / 2000 [ 50%]  (Sampling)
## Chain 4, Iteration: 1200 / 2000 [ 60%]  (Sampling)
## Chain 4, Iteration: 1400 / 2000 [ 70%]  (Sampling)
## Chain 4, Iteration: 1600 / 2000 [ 80%]  (Sampling)
## Chain 4, Iteration: 1800 / 2000 [ 90%]  (Sampling)
## Chain 4, Iteration: 2000 / 2000 [100%]  (Sampling)
##  Elapsed Time: 0.010533 seconds (Warm-up)
##                0.010804 seconds (Sampling)
##                0.021337 seconds (Total)

fit

## Inference for Stan model: coins.
## 4 chains, each with iter=2000; warmup=1000; thin=1; 
## post-warmup draws per chain=1000, total post-warmup draws=4000.
## 
##        mean se_mean   sd  2.5%   25%   50%   75% 97.5% n_eff Rhat
## p      0.75    0.00 0.19  0.29  0.64  0.80  0.91  0.99  1833    1
## y_new  0.75    0.01 0.43  0.00  1.00  1.00  1.00  1.00  3610    1
## lp__  -2.86    0.02 0.82 -5.26 -3.07 -2.53 -2.31 -2.25  1481    1
## 
## Samples were drawn using NUTS(diag_e) at Sat Feb 11 23:54:54 2017.
## For each parameter, n_eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor on split chains (at 
## convergence, Rhat=1).

fit_array <- as.array(fit)
mcmc_hist(fit_array)

Программируем!

…что значит “как”?

Полный список всех возможных команд для языка Stan ищите на mc-stan.org :)

RStan: инструкция по применению

Андрей Горюшев

Как установить RStan

Если у Вас Windows

Требования:

Установка RStan (наконец-то!)

Если у Вас Mac OS или Linux

MacOS

Требования:

Установка RStan

Linux

Требования:

Установка RStan

Как пользоваться RStan