二维直方图的拟合优度

机器算法验证 拟合优度 直方图
2022-02-15 07:46:53

我有两组代表恒星参数的数据：观察到的一组和建模的一组。通过这些集合，我创建了所谓的双色图 (TCD)。可以在这里看到一个示例：

直方图

A是观察到的数据，B是从模型中提取的数据（别管黑线，点代表数据）我只有一个A图，但可以根据需要生成任意多个不同的B图，我需要的是保留最适合A的那个。

所以我需要一种可靠的方法来检查图B（模型）与图A（观察到的）的拟合优度。

现在我要做的是通过对两个轴进行分箱（每个轴 100 个箱）为每个图表创建一个 2D 直方图或网格（这就是我所说的，也许它有一个更合适的名称）然后我遍历网格的每个单元格我发现该特定单元格的 A和B之间的计数绝对差异。在遍历所有单元格后，我将每个单元格的值相加，因此我最终得到一个表示拟合优度的正参数（gf) 在A和B之间。越接近零，拟合越好。基本上，这就是该参数的样子：

gf=∑ij|aij−bij|; 在哪里aij是图表A中该特定单元格的星数（由ij） 和bij是B的数字。

这就是那些(aij−bij)每个单元格中的计数差异看起来像我创建的网格（请注意，我没有使用(aij−bij)在这张图片中，但我在计算gf范围）：

赫斯

问题是我被告知这可能不是一个好的估计器，主要是因为除了说这个拟合比另一个拟合更好，因为参数较低，我真的不能再说什么了。

重要：

（感谢@PeterEllis 提出这个问题）

1- B中的点与A中的点不是一一对应的。这是在寻找最佳拟合时要记住的重要一点： A和B中的点数不一定相同，拟合优度测试也应该考虑这种差异并尽量减少它。

2-我尝试拟合到A的每个B数据集（模型输出）中的点数不固定。

我见过在某些情况下使用卡方检验：

∑i(Oi−Ei)2/Ei; 在哪里Oi是观察到的频率（模型）和Ei是预期频率（观察）。

但问题是：如果Ei是零？如上图所示，如果我在该范围内创建这些图表的网格，则会有很多单元格Ei为零。

另外，我读过一些人建议在涉及直方图的情况下应用对数似然泊松检验。如果这是正确的，如果有人可以指导我如何使用该测试来解决这个特殊情况，我将非常感激（请记住，我的统计知识很差，所以请尽可能简单:)

1个回答

好的，我已经广泛修改了这个答案。我认为，与其对您的数据进行分箱并比较每个箱中的计数，不如在我原来的答案中提出的建议，即拟合二维内核密度估计并比较它们是一个更好的主意。更好的是，在 Tarn Duong 的R ks 包中有一个函数 kde.test()可以轻松完成这项工作。

查看 kde.test 的文档以获取更多详细信息以及您可以调整的参数。但基本上它几乎完全符合您的要求。它返回的 p 值是在零假设下生成您正在比较的两组数据的概率，即它们是从同一分布中生成的。所以 p 值越高，A 和 B 之间的拟合越好。参见下面的示例，这很容易发现 B1 和 A 不同，但 B2 和 A 似乎是相同的（这就是它们的生成方式） .

# generate some data that at least looks a bit similar
generate <- function(n, displ=1, perturb=1){
    BV <- rnorm(n, 1*displ, 0.4*perturb)
    UB <- -2*displ + BV + exp(rnorm(n,0,.3*perturb))
    data.frame(BV, UB)
}
set.seed(100)
A <- generate(300)
B1 <- generate(500, 0.9, 1.2)
B2 <- generate(100, 1, 1)
AandB <- rbind(A,B1, B2)
AandB$type <- rep(c("A", "B1", "B2"), c(300,500,100))

# plot
p <- ggplot(AandB, aes(x=BV, y=UB)) + facet_grid(~type) + 
    geom_smooth() +     scale_y_reverse() + theme_grey(9)
win.graph(7,3)
p +geom_point(size=.7)

在此处输入图像描述

> library(ks)
> kde.test(x1=as.matrix(A), x2=as.matrix(B1))$pvalue
[1] 2.213532e-05
> kde.test(x1=as.matrix(A), x2=as.matrix(B2))$pvalue
[1] 0.5769637

我在下面的原始答案仅保留，因为现在有来自其他地方的链接，这些链接没有意义

首先，可能有其他方法可以解决这个问题。

Justel 等人提出了 Kolmogorov-Smirnov 拟合优度检验的多元扩展，我认为可以在您的案例中使用它，以测试每组建模数据与原始数据的拟合程度。我找不到这个的实现（例如在 R 中），但也许我看起来不够努力。

或者，可能有一种方法可以通过将copula拟合到原始数据和每组建模数据，然后比较这些模型来做到这一点。在 R 和其他地方有这种方法的实现，但我对它们不是特别熟悉，所以没有尝试过。

但是要直接解决您的问题，您采取的方法是合理的。有几点建议自己：

  • 除非您的数据集比看起来要大，否则我认为 100 x 100 的网格太多了。直观地说，我可以想象您得出的结论是各种数据集比它们更加不同，这仅仅是因为您的 bin 的精度意味着您有很多 bin 中的点数很少，即使数据密度很高。然而，这最终是一个判断的问题。我当然会用不同的分箱方法检查你的结果。

  • 完成分箱并将数据转换为（实际上）具有两列且行数等于分箱数（在您的情况下为 10,000）的计数列联表后，您就有一个比较两列的标准问题的计数。卡方检验或拟合某种泊松模型都可以，但正如您所说，由于大量的零计数存在尴尬。这些模型中的任何一个通常都是通过最小化差异的平方和来拟合的，由预期计数的倒数加权；当这接近零时，可能会导致问题。

编辑 - 这个答案的其余部分我现在不再认为是一种合适的方法。

我认为 Fisher 的精确检验在这种情况下可能没有用或不合适，因为交叉表中行的边际总数不固定。它会给出一个似是而非的答案，但我发现很难将其使用与实验设计的原始推导相协调。我将原始答案留在这里，因此评论和后续问题是有意义的。另外，可能仍然有一种方法可以回答 OP 对数据进行分箱的所需方法，并通过基于平均绝对值或平方差的一些测试来比较这些分箱。这种方法仍将使用ng×2下面提到的交叉表并测试独立性，即寻找 A 列与 B 列具有相同比例的结果。

我怀疑上述问题的解决方案是使用Fisher 精确检验，将其应用于您的ng×2交叉表，其中ng是 bin 的总数。尽管由于表中的行数，完整的计算不太可能是实用的，但您可以使用蒙特卡罗模拟获得对 p 值的良好估计（Fisher 检验的 R 实现将其作为表的一个选项：大于 2 x 2，我怀疑其他包也是如此）。这些 p 值是第二组数据（来自您的一个模型）在您的 bin 中与原始数据具有相同分布的概率。因此，p 值越高，拟合越好。

我模拟了一些数据看起来有点像你的，发现这种方法非常有效地识别我的“B”组数据中的哪些是从与“A”相同的过程中生成的，哪些略有不同。肯定比肉眼更有效。

  • 用这种方法测试变量的独立性ng×2列联表，A中的点数与B中的点数不同并不重要（尽管注意它是如果您按照最初的建议仅使用绝对差或平方差的总和，则会出现问题）。但是，您的每个版本的 B 具有不同数量的点确实很重要。基本上，较大的 B 数据集将倾向于返回较低的 p 值。我可以想到几个可能的解决方案来解决这个问题。1. 您可以通过从大于该大小的所有 B 组中随机抽取该大小的样本，将所有 B 组数据减少到相同大小（B 组中最小的大小）。2. 您可以首先为每个 B 组拟合一个二维核密度估计，然后从该估计中模拟大小相等的数据。3.您可以使用某种模拟来计算p值与大小的关系并使用它来“纠正” 您从上述过程中获得的 p 值，因此它们具有可比性。可能还有其他选择。你做哪一个将取决于 B 数据是如何生成的，大小有多么不同等。

希望有帮助。

其它你可能感兴趣的问题
上一篇在统计学中的重要性是什么？e-X2e−x2 下一篇哪些 MCMC 算法/技术用于离散参数？