格式化分类预测器以在STAN中使用的正确方法是什么？我似乎无法将分类预测器输入为正常因子变量,那么转换正常分类变量以便STAN可以接受它的最快方法是什么？

例如,假设我有一个继续预测器和类似这样的分类预测器:

     income      country
1  62085.59      England
2  60806.33      England
3  60527.27      England
4  67112.64          USA
5  57675.92          USA
6  58128.44          USA
7  60822.47 South Africa
8  55805.80 South Africa
9  63982.99 South Africa
10 64555.45      Belgium

我如何准备将其输入rstan？

我想将权重结合到WINBUGS模型brms用权重做的可能性中.

通常的BUGS方法来实现这一目标dnorm并且dpois不起作用dbin.

正如@ paul.buerkner 在这里所说的那样,这是用这样的Stan代码完成的:

vector[N] weights; \\ model weights 
target += weights[n] * neg_binommial_2_log_lpmf(Y[n] | mu[n], shape);

当我在我的BUGS模型中实现这种方法时,我得到下面的详细错误(参见编辑).

以下是数据和模型:

library(R2WinBUGS)
dat <- data.frame(
  A = c(1, 1, 0, 0), B = c(1, 0, 1, 0),
  Pass = c(278, 100, 153, 79), Fail = c(743, 581, 1232, 1731), Weights= c(3, 1, 12, 3))
N <- length(dat$Pass)
case <- dat$Pass
nn <- dat$Fail+dat$Pass
A <- dat$A
B <- dat$B
data …

我在Stan使用最大似然优化,但不幸的是该optimizing()函数没有报告标准错误:

> MLb4c <- optimizing(get_stanmodel(fitb4c), data = win.data, init = inits)  
STAN OPTIMIZATION COMMAND (LBFGS)
init = user
save_iterations = 1
init_alpha = 0.001
tol_obj = 1e-012
tol_grad = 1e-008
tol_param = 1e-008
tol_rel_obj = 10000
tol_rel_grad = 1e+007
history_size = 5
seed = 292156286
initial log joint probability = -4038.66
    Iter      log prob        ||dx||      ||grad||       alpha      alpha0  # evals  Notes 
      13      -2772.49  9.21091e-005     0.0135987     0.07606      0.9845       15   
Optimization terminated normally: 
  Convergence detected: relative gradient magnitude is below …

我想rstan在rnotebook 中获得模型的估计系数

我有以下stan块:

```{stan output.var="rats"}
data {
  int<lower=0> N;
  int<lower=0> T;
  real x[T];
  real y[N,T];
  real xbar;
}
parameters {
  real alpha[N];
  real beta[N];

  real mu_alpha;
  real mu_beta;          // beta.c in original bugs model

  real<lower=0> sigmasq_y;
  real<lower=0> sigmasq_alpha;
  real<lower=0> sigmasq_beta;
}
transformed parameters {
  real<lower=0> sigma_y;       // sigma in original bugs model
  real<lower=0> sigma_alpha;
  real<lower=0> sigma_beta;

  sigma_y = sqrt(sigmasq_y);
  sigma_alpha = sqrt(sigmasq_alpha);
  sigma_beta =  sqrt(sigmasq_beta);
}
model {
  mu_alpha ~ normal(0, 100);
  mu_beta ~ normal(0, 100); …

我正在使用很棒的绘图库bayesplot来可视化我正在使用的模型的后验概率区间rstanarm。我想通过将系数的后验区间放到同一图上来以图形方式比较来自不同模型的绘制。

例如，想象一下，对于beta1, beta2, beta3两个不同模型的三个参数，我从后验中抽取了 1000 次：

# load the plotting library
library(bayesplot)
#> This is bayesplot version 1.6.0
#> - Online documentation and vignettes at mc-stan.org/bayesplot
#> - bayesplot theme set to bayesplot::theme_default()
#>    * Does _not_ affect other ggplot2 plots
#>    * See ?bayesplot_theme_set for details on theme setting
library(ggplot2)

# generate fake posterior draws from model1
fdata <- matrix(rnorm(1000 * 3), ncol = 3)
colnames(fdata) <- c('beta1', 'beta2', 'beta3')

# fake posterior …

我正在尝试通过 rstan 来学习 Stan （因为我熟悉 R）。我尝试运行一个简单的混合帕累托和正态模型。它编译得很好（据我所知），但它无法采样，给我错误：

“(-2, 2) 之间的初始化在 100 次尝试后失败。尝试指定初始值、减小约束值的范围或重新参数化模型。

调用采样器时发生错误；采样未完成”

可以说我已经尝试了各种方法来参数化事物，并尝试设置初始值，但都无济于事。

我的 R+rstan 代码如下：

library(rstan)
rpareto = function(n, location, shape){location/runif(n)^(1/shape)}
sdvec=runif(1e3,0.1,1)
HMFtest=list(x=rpareto(1e3,10,2)+rnorm(1e3,0,sdvec), sdev=sdvec, N=1e3)

HMF.stan <- "
data {
  int<lower=0> N;
  real x[N];
  real sdev[N];
}
parameters {
  real<lower=0,upper=20> y_min;
  real<lower=0,upper=4> alpha;
  real xtrue[N];
}
model {
  y_min ~ lognormal(1, 1);
  alpha ~ lognormal(1, 1);
  xtrue ~ pareto(y_min, alpha);
  for(i in 1:N){
    x[i] ~ normal(xtrue[i], sdev[i]);
  }
}
"

stan.test <- stan(model_code=HMF.stan, data=HMFtest, pars=c('y_min','alpha'), chains=1, iter=30000, warmup=10000)

这个例子与 …

安装rstan时出现以下错误：

Error in .shlib_internal(args) :
C++14 standard requested but CXX14 is not defined

从研究中了解到 C++14 编译器应该可用。如何在配置 R 时安装相同的。使用以下命令配置 R

./configure --with-readline=no --with-x=no

并安装

yum install -y devtoolset-6

但仍然无法更新 C++14 并且 rstan 给出错误

Default C++ compiler:      g++   -g -O2
C++98 compiler:            g++  -g -O2
C++11 compiler:            g++ -std=gnu++11 -g -O2
C++14 compiler:            g++   -g -O2  
C++17 compiler:              
Fortran 90/95 compiler:    gfortran -g -O2
Obj-C compiler: 

设置文件

 yum -y update
 yum install -y centos-release-scl
 yum install -y devtoolset-6
 yum install -y devtoolset-6-gcc-gfortran …

问题总结

我正在将brms::brm_multiple()模型拟合到大型数据集，其中使用该包估算了缺失的数据mice。数据集的大小使得并行处理的使用非常可取。但是，我不清楚如何最好地使用计算资源，因为我不清楚如何brms在核心之间划分估算数据集的采样。

如何选择以下选项以最大限度地有效利用计算资源？

• 插补数 ( m)
• 链数 ( chains)
• 核心数 ( cores)

概念示例

假设我天真地（或者为了举例而故意愚蠢地）选择m = 5, chains = 10, cores = 24。因此，需要在 HPC 上保留的 24 个核心之间分配 5 x 10 = 50 个链。如果没有并行处理，这将需要约 50 个时间单位（不包括编译时间）。

我可以想象 的三种并行化策略brms_multiple()，但可能还有其他策略：

场景 1：并行估算数据集，串行关联链

这里，5 个插补中的每一个都分配给它自己的处理器，该处理器串行运行 10 个链。处理时间为 10 个单位（与非并行处理相比，速度提高了 5 倍），但糟糕的规划浪费了 19 个核心 x 10 个时间单位 = 190 个核心时间单位（ctu；= 80% 的预留计算资源）。有效的解决方案是设置cores= m。

场景 2：串行估算数据集，并行关联链

在这里，采样首先获取第一个估算数据集，并在 …

我想根据这些参数的先验来绘制来自 stan 模型的参数估计的直方图。我曾尝试通过在 stan 中运行模型，用 ggplot2 绘制它，然后使用 R 的随机生成器函数（例如rnorm()，rbinom()）覆盖先验分布的近似值来做到这一点，但我遇到了许多缩放问题，使图形无法查看对。

我在想一个更好的方法是直接从先验分布中采样，然后根据参数估计绘制这些样本，但是运行一个完整的单独模型只是为了从先验中采样似乎非常耗时。我想知道是否有一种方法可以在现有模型内或与之并行执行此操作。

这是一个示例脚本。

# simulate linear model
a <- 3 # intercept
b <- 2 # slope

# data
x <- rnorm(28, 0, 1)
eps <- rnorm(28, 0, 2)
y <- a + b*x + eps

# put data into list
data_reg <- list(N = 28, x = x, y = y)

# create the model string

ms <- "
    data {
    int<lower=0> N;
    vector[N] x; …

我正在尝试在对象bayesplot上实现包中的函数INLA，并且有点不确定如何从后验预测分布中提取。我想我几乎已经拿到了，但是rstan 抽签比抽签更加多变INLA。

在 中，使用小插图rstan中的简化示例我可以：bayesplot

library(bayesplot)
library(ggplot2)
library(rstanarm)
library(ggpubr)
library(tidyverse)


#rstan model set up
roaches$roach100 <- roaches$roach1 / 100 # pre-treatment number of roaches (in 100s)
fit_poisson <- stan_glm(y ~ roach100 + treatment + senior, offset = log(exposure2), family = poisson(link = "log"), data = roaches,  seed = 1111,  refresh = 0) 


#In order to use the PPC functions from the bayesplot package we need a vector y of outcome …