标签: bayesian

我正在寻找一个R可用于根据计数数据训练狄利克雷先验的包。我正在询问一位正在使用 的同事R，而我自己并没有使用它，所以我不太清楚如何寻找软件包。搜索起来有点困难，因为“R”是一个非特定的搜索字符串。CRAN上好像没有什么，但是还有其他地方可以看吗？

我正在使用 OpenBUGS 和 R 包研究二项式混合模型R2OpenBUGS。我已经成功构建了更简单的模型，但是一旦我添加另一个级别来进行不完善的检测，我就会不断收到错误variable X is not defined in model or in data set。我尝试了许多不同的方法，包括更改数据结构以及将数据直接输入 OpenBUGS。我发布此内容是希望其他人有此错误的经验，并且也许知道为什么 OpenBUGS 无法识别变量 X，尽管据我所知它已明确定义。

我也收到了错误expected the collection operator c error pos 8- 这不是我以前收到的错误，但我也同样感到困惑。

模型和数据模拟函数均来自 Kery 的《生态学家 WinBUGS 简介》(2010)。我要注意的是，这里的数据集是代替我自己的数据的，这是相似的。

我包括构建数据集和模型的函数。抱歉长度。

# Simulate data: 200 sites, 3 sampling rounds, 3 factors of the level 'trt', 
# and continuous covariate 'X'

data.fn <- function(nsite = 180, nrep = 3, xmin = -1, xmax = 1, alpha.vec = c(0.01,0.2,0.4,1.1,0.01,0.2), beta0 = 1, beta1 …

我正在使用贝叶斯推理做课堂作业。为此，我使用MCMCregress来自 的函数MCMCpack。当我想要获得残差时，问题就出现了，因为该函数不提供残差，所以我必须“手动”计算它们（在 R 中）。

我的模型是：

Model <- MCMCregress(Y~X1+X2+X3+X4+X5, data=DATA)

其中X1和X5是连续的，而X2和X3是X4二分的。模型输出为我提供了每个变量的估计：

(Intercept) = 1.90, 
X1 = -0.02, 
X2 = -0.05, 
X3 = 0.32,
X4 = 0.61,
X5 = -0.003,
sigma2 = 0.54

我想我必须这样做：

1.90 - 0.02*X1 - 0.05*X2 + 0.32*X3 ...

但我知道我在 R 代码中遗漏了一些重要的东西，所以我想知道哪个是正确的 R 代码以获得残差。

这是一个可重现的示例（尽管它与原始数据不对应）：

Y  <- c(0.2,0.8,1,4.3,5,3.5,3.2,1,3.3,1,2,4,3.6,3,5,4.3,3.2,4,1,2)
X1 <- c(17,13,15,NA,12,24,15,NA,12,22,14,12,18,NA,10,13,12,11,26,10)
X2 <- c(0,0,0,1,0,0,1,1,1,1,0,1,0,0,0,0,0,1,NA,NA)
X3 <- c(0,0,1,1,1,0,1,0,0,0,0,1,0,NA,0,NA,NA,0,1,0)
X4 <- c(1,0,1,0,0,1,0,0,NA,0,NA,NA,0,0,1,0,1,1,0,1)
X5 <- c(2.46,4.56,32.1,NA,NA,NA,NA,NA,NA,3.76,5.67,4.56,NA,
        17.32,12.2,4.56,7.2,1.2,NA,9.2)

X2 …

我使用以下行来更新每次试验中的 beta 发行版并给出 arm 推荐（我使用 scipy.stats.beta）：

self.prior = (1.0,1.0)
def get_recommendation(self):
    sampled_theta = []
    for i in range(self.arms):
        #Construct beta distribution for posterior
        dist = beta(self.prior[0]+self.successes[i],
                    self.prior[1]+self.trials[i]-self.successes[i])
        #Draw sample from beta distribution
        sampled_theta += [ dist.rvs() ]
    # Return the index of the sample with the largest value
    return sampled_theta.index( max(sampled_theta) )

但目前，它只适用于奖励是二元的（成功或失败）。我想修改它，使其适用于非二元奖励。（例如奖励：2300、2000，...）。我怎么做？

背景：

通常，R 给出众所周知的分布的分位数。在这些分位数中，较低的 2.5% 到较高的 97.5% 覆盖了这些分布下 95% 的面积。

问题：

假设我有一个 F 分布（df1 = 10，df2 = 90）。在 R 中，如何确定此分布下的 95% 区域，使得该 95% 仅覆盖高密度区域，而不是 R 通常给出的 95%（请参阅下面的 R 代码）？

注意：显然，最高密度是“模式”（下图中的虚线）。所以我想，人们必须从“模式”转向尾部。

这是我的 R 代码：

curve(df(x, 10, 90), 0, 3, ylab = 'Density', xlab = 'F value', lwd = 3)

Mode = ( (10 - 2) / 10 ) * ( 90 / (90 + 2) )

abline(v = Mode, lty = 2)

CI = qf( c(.025, .975), 10, …

我正在努力实现一个模型，其中狄利克雷变量的集中因子依赖于另一个变量。

情况如下：

系统由于组件故障而失败（共有三个组件，每次测试/观察时只有一个组件失败）。

组件发生故障的概率取决于温度。

这是该情况的（已注释的）简短实现：

import numpy as np
import pymc3 as pm
import theano.tensor as tt


# Temperature data : 3 cold temperatures and 3 warm temperatures
T_data = np.array([10, 12, 14, 80, 90, 95])

# Data of failures of 3 components : [0,0,1] means component 3 failed
F_data = np.array([[0, 0, 1], \
       [0, 0, 1], \
       [0, 0, 1], \
       [1, 0, 0], \
       [1, 0, 0], \
       [1, 0, 0]])

n_component = 3

# When temperature …

我发现到处都发布了同样的问题，而且我似乎无法找到任何适合我的数据的解决方案，我想知道我是否正在尝试将我的数据拟合到一个模型太复杂了。

我正在尝试将我的数据拟合到 MCMCglmm 包中的多项逻辑回归模型。我查看了许多不同的文档、教程和 MCMCglmm 手册本身，主要是Florian Jaeger 的教程，它非常出色且全面。然而，我在他选择先验的 G 结构和 R 结构的值时迷失了方向，并且不断收到此错误消息

Error in priorformat(if (NOpriorG) { : \n  V is the wrong dimension for some prior$G/prior$R elements\n

特别是，我不确定n在这两个数据中应该给我的数据赋予什么值，但有些不透明的错误消息表明这是一个问题V

这是我的数据的（相关）子集：

CG_imm       locuteur   enquete       loc_age   loc_sexe      left    liquid   right    articulation_C1    voice_C1   NC_C1    NC_right    voice_right    right2          pos    logfreq      realization\nabordable    44ajs1     Nantes        79         M            bl      l        p        stop               V          Vstop    stop        NV             NVstop          adj    NA           2\nadmettre     91adb1     Brunoy        54         M            tR      R        E …

我正在为我的主项目演示 hyperopt 的 TPE 算法，但似乎无法让算法收敛。根据我从原始论文和 YouTube讲座中了解到的，TPE 算法按以下步骤工作：

（下文中，x=超参数，y=损失）

1. 首先创建 [x,y] 的搜索历史记录，例如 10 个点。
2. 根据超参数的损失对超参数进行排序，并使用一些分位数 γ 将它们分为两组（γ = 0.5 意味着这些组的大小相同）
3. 对较差的超参数组 (g(x)) 和良好的超参数组 (l(x)) 进行核密度估计
4. 良好的估计在 g(x) 中的概率较低，在 l(x) 中的概率较高，因此我们建议在 argmin(g(x)/l(x)) 处评估函数
5. 在建议的点评估 (x,y) 对并重复步骤 2-5。

我已经在 python 中对目标函数 f(x) = x^2 实现了这一点，但算法无法收敛到最小值。

import numpy as np
import scipy as sp
from matplotlib import pyplot as plt
from scipy.stats import gaussian_kde


def objective_func(x):
    return x**2

def measure(x):
    noise = np.random.randn(len(x))*0
    return x**2+noise

def split_meassures(x_obs,y_obs,gamma=1/2):
    #split x and …

我使用贝叶斯优化包（https://github.com/fmfn/BayesianOptimization）进行参数优化。默认情况下，该库迭代浮点数，但我需要迭代整数，我该如何执行？

def black_box_function(x, y):
    return -x ** 2 - (y - 1) ** 2 + 1
from bayes_opt import BayesianOptimization
pbounds = {'x': (2, 4), 'y': (-3, 3)}

optimizer = BayesianOptimization(
    f=black_box_function,
    pbounds=pbounds,
    verbose=2, 
    random_state=1,
)

optimizer.maximize(
    init_points=2,
    n_iter=3,
)
# and as you can see it is iterated not over integers.
|   iter    |  target   |     x     |     y     |
-------------------------------------------------
|  1        | -7.135    |  2.834    |  1.322    |
|  2        | -7.78     |  2.0      | …

我试图拟合一个多项式逻辑回归模型，其rjags结果是一个具有 3 个级别的分类（名义）变量（结果），解释变量是年龄（连续）和组（具有 3 个级别的分类）。在此过程中，我想获得Age和Group的后验均值和基于分位数的 95% 区域。

我不太擅长，for loop我认为这就是我为模型编写的代码无法正常工作的原因。

我的 beta 先验遵循正态分布，\xce\xb2j \xe2\x88\xbc Normal(0,100) for j \xe2\x88\x88 {0, 1, 2}。

可重现的 R 代码

library(rjags)\n\nset.seed(1)\ndata <- data.frame(Age = round(runif(119, min = 1, max = 18)),\n                   Group = c(rep("pink", 20), rep("blue", 18), rep("yellow", 81)), \n                   Outcome = c(rep("A", 45), rep("B", 19), rep("C", 55)))\n\nX <- as.matrix(data[,c("Age", "Group")]) \nJ <- ncol(X)\nN <- nrow(X)\n\n## Step …