1

我有一个SampleOpenTURNS，我想在它上面安装一个发行版。为了考虑参数的数量，我想使用贝叶斯信息准则（BIC）。贝叶斯信息准则 (BIC) 根据加权最大似然准则对模型列表进行排名，该准则考虑了样本大小和每个分布的参数数量。BIC 分数越低越好。

我知道FittingTest.BestModelBIC返回最适合数据的模型。然而，我希望看到的不仅仅是最合适的：也许第二好的 BIC 对我来说有更多的物理意义？

如何在 OpenTURNS 中执行此操作？

PS 这里是BestModelBIC的一个例子。

4

1 回答 1

1

我建议使用GetContinuousUniVariateFactories来获取连续分布工厂的列表，并将其与FittingTest.BIC计算分布的 BIC 分数的函数结合起来。最后，可以使用 的sortAccordingToAComponent方法对表格进行排序Sample。

这是此方法的详细示例。我们首先生成一个样本。

import openturns as ot
import tqdm
import pandas as pd
sample = ot.Normal().getSample(100)

GetContinuousUniVariateFactories静态方法返回用于连续分布的所有可用工厂的列表。我们可以使用这个列表而无需进一步处理，但直方图会在排名中排在第一位，因为它是专门为此目的而设计的。因此，我们不将其包含在 BIC 分数的计算中。

factories = ot.DistributionFactory.GetContinuousUniVariateFactories()
marginalFactories = []
for factory in factories:
    if str(factory) != "HistogramFactory":
        print(factory)
        marginalFactories.append(factory)
number_of_factories = len(marginalFactories)
print("Number of selected factories:", number_of_factories)

这打印：

ArcsineFactory
BetaFactory
BurrFactory
[...]
Number of selected factories: 29

在以下脚本中，我们for对之前创建的列表中的所有工厂执行循环。稍后我们将按升序对 BIC 分数进行排序。这就是我们将 BIC 分数和边际指数存储在score_array样本中的原因。对于某些分布，计算可能会很长。因此我们使用该tqdm模块来打印进度条。最后，某些发行版并非基于此特定样本。为了避免破坏 for 循环，我们将对BIC方法的调用包装到一个try/except. 如果分布拟合失败，我们将 BIC 分数设置为浮点数的最大有限值（即MaxScalar），大约等于 $10^{308}$。

score_array = ot.Sample(number_of_factories, 2)
for i in tqdm.tqdm(range(number_of_factories)):
    try:
        factory = marginalFactories[i]
        fitted_dist, bic_score = ot.FittingTest.BIC(sample, factory)
        score_array[i] = [i, bic_score]
    except TypeError:
        print("Cannot build ", factory)
        score_array[i] = [i, ot.SpecFunc.MaxScalar]

关键步骤是对包含 BIC 分数的数组进行排序。

sorted_BIC_scores = score_array.sortAccordingToAComponent(1)

可能有超过 30 个分布可以构建到样本上。在这里，我们将列表限制为 BIC 分数最低的前 10 个分布。我们将使用 Pandas 来很好地打印 BIC 分数。为此，我们创建了一个BIC_data列表，其中包含工厂名称和相应的 BIC 分数。这是使用第一列中的分布索引的地方sorted_BIC_scores。但是，Samplestores floats: 在将其用作索引之前，我们必须将它们转换为整数。

BIC_data = []
rank = list(range(min(number_of_factories, 10)))
for i in rank:
    distribution_index = int(sorted_BIC_scores[i, 0])
    factory = marginalFactories[distribution_index]
    BIC_score = sorted_BIC_scores[i, 1]
    BIC_data.append([factory, BIC_score])

现在是最简单的部分，我们最终使用 Pandas 的DataFrame.

columns = ["BIC", "Factory"]
df = pd.DataFrame(BIC_data, index=rank, columns=columns,)
print(df)

这打印的内容如下所示：

秩 BIC 工厂
0 普通工厂 2.902476
1 威布尔麦克斯工厂 2.933988
2 威布尔MinFactory 2.938930
3 截断法线工厂 2.939140
4 物流工厂 2.945102
5 日志正常工厂 2.948479
6 学生工厂 2.948733
7 三角工厂 2.968305
8 贝塔工厂 3.033244
9 瑞利工厂 3.051117

幸运的是，我们看到NormalFactory符合高斯样本。截断正态工厂排在第 4 位，在两种 Weibull 分布之后。

于 2021-01-30T20:35:15.650 回答