在合理的假设下,一个好的判别式是唯一的数据值,它导致 B 的概率密度在分割点左侧的面积等于 A 在右侧的面积(反之亦然,给出相同的点) .
一个简单的方法是计算两个经验累积分布函数 (CDF) 作为shown here 并搜索它们以提供分割点。这是两个 CDF 之和为 1 的点。
简单地说,构建经验 CDF 只是对每个数据集进行排序并将数据用作 x 轴值。从左到右绘制曲线,从 y=0 开始,在每个 x 值处向上走 1/n 步。这样的曲线从 x 1 的 0 逐渐上升到 x >= data[n] 的 y = CDF(x) = 1。有一种稍微复杂的方法可以给出连续的逐步线性曲线而不是阶梯,在某些假设下,它是对真实 CDF 的更好估计。在
请注意,上面的讨论只是为了提供直觉。 CDF 由已排序的数据数组完美表示。不需要新的数据结构;即 x[i], i=1,2,...,n 是曲线到达 y = i/n 处的 x 值。
使用两个 CDF,R(x) 和 B(x),根据您的图表,您想要找到唯一点 x,使得 |1 - R(x) - B(x)|被最小化(使用分段线性 CDF,您将始终能够将其设为零)。这可以通过二分搜索很容易地完成。
这种方法的一个好处是您可以通过将两个 CDF 维护在有序集合(平衡二叉搜索树)中来使其动态化。随着点的增加,新的分割点很容易找到。
有序集需要“订单统计”。 Here is a reference。我的意思是您需要能够查询排序集以检索 CDF 中任何存储的 x 值的序数。这可以通过跳过列表和树来完成。
我编写了这个算法的一个变体。它使用分段 CDF 近似,但也允许在重复数据点处进行“垂直步骤”。这使算法有些复杂,但还不错。然后我使用二分法(而不是组合二分搜索)来找到分割点。正常的二分算法需要修改以适应 CDF 中的垂直“步骤”。我认为这一切都是正确的,但经过了轻微的测试。
不处理的一个极端情况是数据集具有不相交的范围。这将在较低的顶部和较高的底部之间找到 a 点,这是一个完全有效的鉴别器。但是你可能想做一些更有趣的事情,比如返回某种加权平均值。
请注意,如果您对数据可以达到的实际最小值和最大值有一个很好的概念,并且它们不会出现在数据中,则应考虑添加它们,以免 CDF不经意间偏了。
在您的示例数据上,代码生成 4184.76,看起来非常接近您在图表中选择的值(略低于最小和最大数据的一半)。
请注意,我没有对数据进行排序,因为它已经是。排序绝对是必要的。
public class SplitData {
// Return: i such that a[i] <= x < a[i+1] if i,i+1 in range
// else -1 if x < a[0]
// else a.length if x >= a[a.length - 1]
static int hi_bracket(double[] a, double x) {
if (x < a[0]) return -1;
if (x >= a[a.length - 1]) return a.length;
int lo = 0, hi = a.length - 1;
while (lo + 1 < hi) {
int mid = (lo + hi) / 2;
if (x < a[mid])
hi = mid;
else
lo = mid;
}
return lo;
}
// Return: i such that a[i-1] < x <= a[i] if i-1,i in range
// else -1 if x <= a[0]
// else a.length if x > a[a.length - 1]
static int lo_bracket(double[] a, double x) {
if (x <= a[0]) return -1;
if (x > a[a.length - 1]) return a.length;
int lo = 0, hi = a.length - 1;
while (lo + 1 < hi) {
int mid = (lo + hi) / 2;
if (x <= a[mid])
hi = mid;
else
lo = mid;
}
return hi;
}
// Interpolate the CDF value for the data a at value x. Returns a range.
static void interpolate_cdf(double[] a, double x, double[] rtn) {
int lo_i1 = lo_bracket(a, x);
if (lo_i1 == -1) {
rtn[0] = rtn[1] = 0;
return;
}
int hi_i0 = hi_bracket(a, x);
if (hi_i0 == a.length) {
rtn[0] = rtn[1] = 1;
return;
}
if (hi_i0 + 1 == lo_i1) { // normal interpolation
rtn[0] = rtn[1]
= (hi_i0 + (x - a[hi_i0]) / (a[lo_i1] - a[hi_i0]))
/ (a.length - 1);
return;
}
// we're on a joint or step; return range answer
rtn[0] = (double)lo_i1 / (a.length - 1);
rtn[1] = (double)hi_i0 / (a.length - 1);
assert rtn[0] <= rtn[1];
}
// Find the data value where the two given data set's empirical CDFs
// sum to 1. This is a good discrimination value for new data.
// This deals with the case where there's a step in either or both CDFs.
static double find_bisector(double[] a, double[] b) {
assert a.length > 0;
assert b.length > 0;
double lo = Math.min(a[0], b[0]);
double hi = Math.max(a[a.length - 1], b[b.length - 1]);
double eps = (hi - lo) * 1e-7;
double[] a_rtn = new double[2], b_rtn = new double[2];
while (hi - lo > eps) {
double mid = 0.5 * (lo + hi);
interpolate_cdf(a, mid, a_rtn);
interpolate_cdf(b, mid, b_rtn);
if (1 < a_rtn[0] + b_rtn[0])
hi = mid;
else if (a_rtn[1] + b_rtn[1] < 1)
lo = mid;
else
return mid; // 1 is included in the interpolated range
}
return 0.5 * (lo + hi);
}
public static void main(String[] args) {
double split = find_bisector(a, b);
System.err.println("Split at x = " + split);
}
static final double[] a = {
385, 515, 975, 1136, 2394, 2436, 4051, 4399, 4484, 4768, 4768, 4849,
4856, 4954, 5020, 5020, 5020, 5020, 5020, 5020, 5020, 5020, 5020, 5052,
5163, 5200, 5271, 5421, 5421, 5442, 5746, 5765, 5903, 5992, 5992, 6046,
6122, 6205, 6208, 6239, 6310, 6360, 6416, 6512, 6536, 6543, 6581, 6609,
6696, 6699, 6752, 6796, 6806, 6855, 6859, 6886, 6906, 6911, 6923, 6953,
7016, 7072, 7086, 7089, 7110, 7232, 7278, 7293, 7304, 7309, 7348, 7367,
7378, 7380, 7419, 7453, 7454, 7492, 7506, 7549, 7563, 7721, 7723, 7731,
7745, 7750, 7751, 7783, 7791, 7813, 7813, 7814, 7818, 7833, 7863, 7875,
7886, 7887, 7902, 7907, 7935, 7942, 7942, 7948, 7973, 7995, 8002, 8013,
8013, 8015, 8024, 8025, 8030, 8038, 8041, 8050, 8056, 8060, 8064, 8071,
8081, 8082, 8085, 8093, 8124, 8139, 8142, 8167, 8179, 8204, 8214, 8223,
8225, 8247, 8248, 8253, 8258, 8264, 8265, 8265, 8269, 8277, 8278, 8289,
8300, 8312, 8314, 8323, 8328, 8334, 8363, 8369, 8390, 8397, 8399, 8399,
8401, 8436, 8442, 8456, 8457, 8471, 8474, 8483, 8503, 8511, 8516, 8533,
8560, 8571, 8575, 8583, 8592, 8593, 8626, 8635, 8635, 8644, 8659, 8685,
8695, 8695, 8702, 8714, 8715, 8717, 8729, 8732, 8740, 8743, 8750, 8756,
8772, 8772, 8778, 8797, 8828, 8840, 8840, 8843, 8856, 8865, 8874, 8876,
8878, 8885, 8887, 8893, 8896, 8905, 8910, 8955, 8970, 8971, 8991, 8995,
9014, 9016, 9042, 9043, 9063, 9069, 9104, 9106, 9107, 9116, 9131, 9157,
9227, 9359, 9471
};
static final double[] b = {
12, 16, 29, 32, 33, 35, 39, 42, 44, 44, 44, 45, 45, 45, 45, 45, 45, 45,
45, 45, 47, 51, 51, 51, 57, 57, 60, 61, 61, 62, 71, 75, 75, 75, 75, 75,
75, 76, 76, 76, 76, 76, 76, 79, 84, 84, 85, 89, 93, 93, 95, 96, 97, 98,
100, 100, 100, 100, 100, 102, 102, 103, 105, 108, 109, 109, 109, 109,
109, 109, 109, 109, 109, 109, 109, 109, 110, 110, 112, 113, 114, 114,
116, 116, 118, 119, 120, 121, 122, 124, 125, 128, 129, 130, 131, 132,
133, 133, 137, 138, 144, 144, 146, 146, 146, 148, 149, 149, 150, 150,
150, 151, 153, 155, 157, 159, 164, 164, 164, 167, 169, 170, 171, 171,
171, 171, 173, 174, 175, 176, 176, 177, 178, 179, 180, 181, 181, 183,
184, 185, 187, 191, 193, 199, 203, 203, 205, 205, 206, 212, 213, 214,
214, 219, 224, 224, 224, 225, 225, 226, 227, 227, 228, 231, 234, 234,
235, 237, 240, 244, 245, 245, 246, 246, 246, 248, 249, 250, 250, 251,
255, 255, 257, 264, 264, 267, 270, 271, 271, 281, 282, 286, 286, 291,
291, 292, 292, 294, 295, 299, 301, 302, 304, 304, 304, 304, 304, 306,
308, 314, 318, 329, 340, 344, 345, 356, 359, 363, 368, 368, 371, 375,
379, 386, 389, 390, 392, 394, 408, 418, 438, 440, 456, 456, 458, 460,
461, 467, 491, 503, 505, 508, 524, 557, 558, 568, 591, 609, 622, 656,
665, 668, 687, 705, 728, 817, 839, 965, 1013, 1093, 1126, 1512, 1935,
2159, 2384, 2424, 2426, 2484, 2738, 2746, 2751, 3006, 3184, 3184, 3184,
3184, 3184, 4023, 5842, 5842, 6502, 7443, 7781, 8132, 8237, 8501
};
}