R Academy · Урок

Написание моделей Stan в R

Определяйте блоки данных, параметры и блок модели в синтаксисе Stan

Урок 2 из 413 шагов

«Написание моделей Stan в R» — бесплатный урок R Academy на CoddyKit. Это урок 2 из 4. Ты можешь прочитать весь урок бесплатно ниже — а потом практиковать его прямо в браузере с встроенным редактором кода и ИИ-репетитором 24/7. Это часть пути обучения R Academy, и твой прогресс синхронизируется между веб-версией и приложением CoddyKit. Курс R Academy содержит 4 уроков всего.

Что такое Stan

Stan — это язык вероятностного программирования для байесовского статистического моделирования. RStan — интерфейс R для Stan. Вы записываете модель на языке Stan, похожем на C++, а Stan компилирует её в эффективный код C++, который выполняет выборку методом гамильтонова Монте-Карло (HMC) для приближения апостериорного распределения.

# Stan installation check
library(rstan)

# Check version
cat('RStan version:', as.character(packageVersion('rstan')), '\n')

# Enable parallel chains
options(mc.cores = parallel::detectCores())

# Reuse compiled models across sessions
rstan_options(auto_write = TRUE)

cat('Stan ready. Cores:', parallel::detectCores(), '\n')

Структура модели Stan

Модель Stan может содержать до шести именованных блоков: data, transformed data, parameters, transformed parameters, model и generated quantities. Три обязательных блока — это data, parameters и model.

library(rstan)

# Stan model as a character string in R
stan_code <- '
data {
  int<lower=0> N;       // number of observations
  vector[N] x;          // predictor
  vector[N] y;          // response
}
parameters {
  real alpha;           // intercept
  real beta;            // slope
  real<lower=0> sigma;  // noise (must be positive)
}
model {
  // Priors
  alpha ~ normal(0, 10);
  beta  ~ normal(0, 10);
  sigma ~ exponential(1);

  // Likelihood
  y ~ normal(alpha + beta * x, sigma);
}
'

cat('Stan model defined as a string in R\n')
cat('Blocks: data, parameters, model\n')

data{}: объявление входных данных

Блок data объявляет все внешние данные, которые модель получает из R. Доступные типы включают int, real, vector[N], matrix[M,N] и array. Ограничения, например <lower=0>, проверяются во время выполнения.

library(rstan)

# Comprehensive data block examples
data_block_examples <- '
data {
  // Scalars
  int<lower=1> N;               // at least 1 observation
  int<lower=2> K;               // at least 2 groups

  // Constrained scalars
  real<lower=0, upper=1> rate;  // probability

  // Vectors
  vector[N] y;                  // continuous response
  array[N] int<lower=0, upper=1> z;  // binary outcomes

  // Matrix
  matrix[N, K] X;               // design matrix

  // Integer array
  array[N] int group;           // group membership 1..K
}
'

cat(data_block_examples)

parameters{}: Типы параметров

Блок parameters объявляет неизвестные величины, которые Stan будет выбирать методом выборки. Ограничения в этом блоке определяют пространство параметров: <lower=0> — для положительных величин, <lower=0, upper=1> — для вероятностей. Stan автоматически применяет якобиановы поправки к логарифму вероятности для параметров с ограничениями.

library(rstan)

# Common parameter declarations
param_examples <- '
parameters {
  // Unconstrained
  real mu;                        // mean
  vector[K] beta;                 // regression coefficients

  // Positive (sigma, lambda, variance)
  real<lower=0> sigma;
  real<lower=0> lambda;

  // Probability
  real<lower=0, upper=1> theta;

  // Simplex (sums to 1, for mixture weights)
  simplex[K] pi;

  // Correlation matrix
  corr_matrix[K] Omega;

  // Cholesky factor of covariance
  cholesky_factor_cov[K] L_Sigma;
}
'
cat(param_examples)

model{}: Априорные распределения и правдоподобие

Блок model накапливает логарифм апостериорного распределения с помощью синтаксиса ~ (тильда). Запись y ~ normal(mu, sigma) является сокращением для добавления логарифма нормальной PDF к target. Явные приращения логарифма вероятности можно записывать с помощью target += normal_lpdf(y | mu, sigma).

library(rstan)

model_block_example <- '
model {
  // --- Priors ---
  mu    ~ normal(0, 10);          // weakly informative
  sigma ~ cauchy(0, 2.5);         // half-Cauchy for scale
  beta  ~ normal(0, 1);           // standardised coefficients

  // --- Likelihood ---
  // Tilde notation (most common)
  y ~ normal(mu + X * beta, sigma);

  // Equivalent explicit notation:
  // target += normal_lpdf(y | mu + X * beta, sigma);

  // For loop (less common but valid)
  // for (i in 1:N)
  //   target += normal_lpdf(y[i] | mu, sigma);
}
'
cat(model_block_example)

Подгонка простой нормальной модели

Ниже приведён полный минимальный рабочий процесс Stan: определите строку с моделью, подготовьте список данных, вызовите stan() и изучите результаты с помощью print(). При первом запуске Stan компилирует модель, а затем кэширует её для последующих запусков.

library(rstan)

# Stan model: estimate mean and SD of a normal distribution
normal_model <- '
data {
  int<lower=0> N;
  vector[N] y;
}
parameters {
  real mu;
  real<lower=0> sigma;
}
model {
  mu    ~ normal(0, 10);
  sigma ~ exponential(0.1);
  y     ~ normal(mu, sigma);
}
'

# Simulate data
set.seed(42)
y_data <- rnorm(50, mean = 5, sd = 2)

# Fit the model
fit <- stan(
  model_code = normal_model,
  data       = list(N = length(y_data), y = y_data),
  chains     = 2,
  iter       = 1000,
  warmup     = 500,
  refresh    = 0  # suppress iteration output
)

print(fit, pars = c('mu', 'sigma'))

Линейная регрессия в Stan

Stan упрощает байесовскую линейную регрессию: задайте нормальные априорные распределения для коэффициентов и экспоненциальное или полунормальное распределение Коши для sigma. Апостериорное распределение задаёт полное распределение неопределённости для каждого коэффициента, а не только точечные оценки.

library(rstan)

lin_reg_model <- '
data {
  int<lower=0> N;
  vector[N] x;
  vector[N] y;
}
parameters {
  real alpha;
  real beta;
  real<lower=0> sigma;
}
model {
  alpha ~ normal(0, 10);
  beta  ~ normal(0, 10);
  sigma ~ exponential(1);
  y ~ normal(alpha + beta * x, sigma);
}
generated quantities {
  vector[N] y_rep;  // posterior predictive
  for (i in 1:N)
    y_rep[i] = normal_rng(alpha + beta * x[i], sigma);
}
'

set.seed(7)
n <- 80
x <- rnorm(n); y <- 2 + 3 * x + rnorm(n, 0, 1)

fit <- stan(model_code = lin_reg_model,
            data = list(N = n, x = x, y = y),
            chains = 2, iter = 1000, refresh = 0)
print(fit, pars = c('alpha', 'beta', 'sigma'))

Блок transformed parameters{}

Блок transformed parameters вычисляет производные величины на основе параметров, полученных методом выборки. Это удобно для задания моделей численно устойчивым способом (например, с помощью разложений Холецкого и логарифмических преобразований), сохраняя при этом понятный смысл основных параметров.

library(rstan)

# Model with transformed parameters
trans_param_example <- '
data {
  int<lower=0> N;
  array[N] int<lower=0> y;  // counts
}
parameters {
  real log_lambda;           // work in log space for stability
}
transformed parameters {
  real<lower=0> lambda;
  lambda = exp(log_lambda);  // transform back to original scale
}
model {
  log_lambda ~ normal(1, 2);  // prior on log scale
  y ~ poisson(lambda);
}
'

set.seed(42)
y_counts <- rpois(30, lambda = 5)

fit_poisson <- stan(
  model_code = trans_param_example,
  data       = list(N = length(y_counts), y = y_counts),
  chains     = 2, iter = 1000, refresh = 0
)
print(fit_poisson, pars = c('log_lambda', 'lambda'))

Блок generated quantities{}

Блок generated quantities выполняется после выборки и вычисляет дополнительные величины: апостериорные прогнозные выборки (y_rep), логарифм правдоподобия для LOO-CV или преобразованные параметры для сводных показателей. Значения в этом блоке получают выборкой из апостериорного прогнозного распределения.

library(rstan)

# Use generated quantities for posterior predictive checks
model_with_gq <- '
data {
  int<lower=0> N;
  vector[N] y;
}
parameters {
  real mu;
  real<lower=0> sigma;
}
model {
  mu    ~ normal(0, 10);
  sigma ~ exponential(0.5);
  y     ~ normal(mu, sigma);
}
generated quantities {
  vector[N] y_rep;           // replicated datasets
  real mean_y_rep;           // mean of replicated data
  for (i in 1:N)
    y_rep[i] = normal_rng(mu, sigma);
  mean_y_rep = mean(y_rep);
}
'

set.seed(1)
y_obs <- rnorm(40, 3, 1.5)
fit <- stan(model_code = model_with_gq,
            data = list(N = length(y_obs), y = y_obs),
            chains = 2, iter = 1000, refresh = 0)
print(fit, pars = c('mu', 'sigma', 'mean_y_rep'))

Передача данных из R в Stan

Аргумент data функции stan() представляет собой именованный список R. Имена должны в точности совпадать с именами переменных, объявленных в блоке Stan data{}. Векторы преобразуются в vector[N] Stan; целые числа — в int; матрицы R — в matrix[M,N] Stan.

library(rstan)

# Data preparation: names must match Stan data block exactly
set.seed(42)
n <- 60
x1 <- rnorm(n)
x2 <- rnorm(n)
y  <- 1.5 + 2 * x1 - 0.8 * x2 + rnorm(n, 0, 0.5)

# Build the design matrix
X <- cbind(x1, x2)  # 60 x 2 matrix

# Named list passed to stan(data = ...)
stan_data <- list(
  N = n,          # int
  K = ncol(X),    # int
  X = X,          # matrix[N, K]
  y = y           # vector[N]
)

cat('Stan data list elements:\n')
for (nm in names(stan_data)) {
  cat(' ', nm, ': class =', class(stan_data[[nm]]),
      'dim =', paste(dim(stan_data[[nm]]), collapse = 'x'),
      '\n')
}

Просмотр скомпилированной модели

После подгонки используйте print(fit), чтобы просмотреть апостериорные сводные показатели (среднее, se_mean, sd, квантили, Rhat, n_eff) для всех параметров. Используйте stan_plot(fit) для построения графика интервалов и traceplot(fit) для оценки перемешивания цепей.

library(rstan)

# Reuse the simple normal model from scene 6
normal_model <- '
data { int<lower=0> N; vector[N] y; }
parameters { real mu; real<lower=0> sigma; }
model {
  mu ~ normal(0, 10);
  sigma ~ exponential(0.1);
  y ~ normal(mu, sigma);
}
'

set.seed(5)
fit <- stan(model_code = normal_model,
            data = list(N = 50, y = rnorm(50, 7, 3)),
            chains = 2, iter = 1000, refresh = 0)

# Detailed summary table
print(fit)

# Extract as data frame
posterior_df <- as.data.frame(fit)
cat('\nPosterior samples shape:', nrow(posterior_df),
    'rows x', ncol(posterior_df), 'cols\n')

Быстрая проверка

В модели Stan необходимо ограничить параметр строго положительными значениями (например, стандартное отклонение). Какое объявление является правильным?

Повторение: написание моделей Stan

Главные выводы:

  • Модель Stan содержит три основных блока: data{}, parameters{}, model{}
  • Ограничения: <lower=0>, <upper=1>, <lower=0, upper=1> для вероятностей
  • Типы: int, real, vector[N], matrix[M,N], simplex[K]
  • Блок модели: используйте синтаксис с тильдой y ~ normal(mu, sigma) для априорных распределений и правдоподобия
  • transformed parameters{} — для производных величин; generated quantities{} — для вычислений после выборки
  • Передавайте данные в виде именованного списка R, имена которого в точности соответствуют именам блоков Stan
  • rstan_options(auto_write = TRUE) кэширует скомпилированные модели
library(rstan)

# Stan model skeleton
model_skeleton <- '
data    { int N; vector[N] y; }
parameters { real mu; real<lower=0> sigma; }
model   { mu ~ normal(0,10); sigma ~ exponential(1); y ~ normal(mu, sigma); }
'
cat(model_skeleton)
cat('\nFit with: stan(model_code = model_skeleton, data = list(N=..., y=...), chains=4)\n')
Можно начать бесплатно

Изучай R с ИИ-репетитором — бесплатно

Пиши и запускай код прямо в браузере, получай мгновенную помощь от ИИ-репетитора 24/7 и продолжи учиться на сайте или в приложении.

Курсы
43
Уроки
159

Часто задаваемые вопросы

Урок «Написание моделей Stan в R» бесплатный?

Да — полный текст урока «Написание моделей Stan в R» бесплатно доступен здесь в веб-версии. Чтобы практиковать его интерактивно (встроенный редактор кода и ИИ-репетитор 24/7) и разблокировать остальной курс R Academy, подпишись на CoddyKit PRO. Курс R Academy содержит 4 уроков всего.

Чему я научусь в уроке «Написание моделей Stan в R»?

Определяйте блоки данных, параметры и блок модели в синтаксисе Stan Ты практикуешь R Academy с помощью реального кода, который запускаешь прямо в браузере, и ИИ-репетитор 24/7 отвечает на твои вопросы во время урока.

Нужен ли мне опыт, чтобы начать R Academy?

Предыдущий опыт не требуется. R Academy на CoddyKit структурирован для всех уровней — от новичков до продвинутых, поэтому ты можешь начать отсюда или с самого начала и учиться в своем темпе. Это урок 2 из 4.

Сколько времени занимает урок «Написание моделей Stan в R»?

Большинство уроков CoddyKit занимают около 5–10 минут. Каждый из них компактный и интерактивный, поэтому ты постоянно делаешь прогресс и продолжаешь с того же места в веб-версии и приложении.

Можно ли писать и запускать код в этом уроке R Academy?

Да. Каждый урок R Academy включает встроенный редактор кода, поэтому ты пишешь и запускаешь реальный код прямо в браузере и получаешь моментальную обратную связь от AI — локальная установка не требуется.

Все уроки этого курса

  1. Введение в байесовское мышление
  2. Написание моделей Stan в R
  3. Выборка MCMC и диагностика
  4. Проверка апостериорных предсказаний
← Назад к R Academy