"近似的"最大公约数

例如,假设您有一个浮点数列表,它大约是常见数量的倍数

2.468,3.700,6.1699

这大约是1.234的倍数.您如何描述这种"近似gcd",您将如何进行计算或估算呢？

与我对这个问题的回答密切相关.

您可以运行Euclid的gcd算法,其中任何小于0.01(或您选择的少数)都是伪0.使用您的数字:

3.700 = 1 * 2.468 + 1.232,
2.468 = 2 * 1.232 + 0.004. 

所以前两个数字的伪gcd是1.232.现在你用你的最后一个数字取gcd:

6.1699 = 5 * 1.232 + 0.0099.

所以1.232是伪gcd,而mutiples是2,3,5.要改善此结果,您可以对数据点进行线性回归:

(2,2.468), (3,3.7), (5,6.1699).

斜率是改进的伪gcd.

警告:第一部分是算法在数值上不稳定 - 如果从非常脏的数据开始,你就遇到了麻烦.

将您的测量值表示为最低值的倍数.因此,您的列表变为1.00000,1.49919,2.49996.这些值的小数部分将非常接近1/Nth,因为N的某个值取决于您的最低值与基频的接近程度.我建议循环增加N,直到找到足够精确的匹配.在这种情况下,对于N = 1(即,假设X = 2.468是您的基频),您会发现标准偏差为0.3333(三个值中的两个是X*1的0.5.),这是不可接受的高.对于N = 2(即,假设2.468/2是您的基频),您会发现标准差几乎为零(所有三个值都在.X/2倍数的.001内),因此2.468/2是您的近似值GCD.

我的计划中的主要缺陷是,当最低测量值最准确时,它最有效,但情况可能并非如此.这可以通过多次执行整个操作来减轻,每次丢弃测量列表上的最低值,然后使用每个过程的结果列表来确定更精确的结果.改善结果的另一种方法是调整GCD以最小化GCD的整数倍与测量值之间的标准偏差.

这让我想起了找到实数的理性数近似值的问题.标准技术是持续分数扩展:

def rationalizations(x):
    assert 0 <= x
    ix = int(x)
    yield ix, 1
    if x == ix: return
    for numer, denom in rationalizations(1.0/(x-ix)):
        yield denom + ix * numer, numer

我们可以直接将它应用于Jonathan Leffler和Sparr的方法:

>>> a, b, c = 2.468, 3.700, 6.1699
>>> b/a, c/a
(1.4991896272285252, 2.4999594813614263)
>>> list(itertools.islice(rationalizations(b/a), 3))
[(1, 1), (3, 2), (925, 617)]
>>> list(itertools.islice(rationalizations(c/a), 3))
[(2, 1), (5, 2), (30847, 12339)]

从每个序列中挑选出第一个足够好的近似值.(3/2和5/2在这里.)或代替直接比较3.0/2.0〜1.499189 ......,你会发现比六百十七分之九百二十五使用更大的整数比3/2,使3/2的绝佳场所停止.

你所分割的数字中的哪一个并不重要.(使用A/B和C/B你2/3和5/3,例如).一旦你有整数比,你可以细化基本使用shsmurfy的线性回归的隐含估值.大家都赢了!

我假设你的所有数字都是整数值的倍数.对于我的其余解释,A将表示您尝试查找的"根"频率,B将是您必须开始的数字数组.

你要做的是表面上类似于线性回归.您正在尝试找到线性模型y = mx + b,它可以最小化线性模型和一组数据之间的平均距离.在您的情况下,b = 0,m是根频率,y代表给定值.最大的问题是没有明确给出自变量X. 我们对X的唯一了解是它的所有成员都必须是整数.

您的第一个任务是尝试确定这些自变量.我现在能想到的最好的方法是假定给定频率具有几乎连续的索引(x_1=x_0+n).所以B_0/B_1=(x_0)/(x_0+n)给出一个(希望)小整数n.然后,您可以利用以下事实:x_0 = n/(B_1-B_0)从n = 1开始,并保持棘轮直到k-rnd(k)在某个阈值内.拥有x_0(初始索引)后,您可以近似根频率(A = B_0/x_0).然后,您可以通过查找来估计其他索引x_n = rnd(B_n/A).此方法不是很健壮,如果数据中的错误很大,可能会失败.

如果您想要更好地近似根频率A,则可以使用线性回归来最小化线性模型的误差,因为您具有相应的因变量.最简单的方法是使用最小二乘拟合.Wolfram的Mathworld对这个问题进行了深入的数学处理,但是通过一些谷歌搜索可以找到一个相当简单的解释.