algorithm - 将浮点向量舍入到整数向量

Question

作为一个更大的模拟程序的一部分(我很乐意分享其背景，但与问题无关)，我遇到了以下问题，正在寻找一个好的算法。

问题:给定一个长度为 n 的 float 组 f(包含元素 f_1, ..., f_n)，指定 n 维空间中的一个点。可以公平地假设 f 的总和为 0.0(取决于浮点精度)。

寻求:求长度为n的整数数组i(有元素i_1, ..., i_n)，在n维空间指定一个格点，使得i和d(f,i)之和正好为0，合适的测度f 和 i 之间的距离最小化。

至于合适的度量，我认为最好的方法是最小化相对误差(即对 (i_j/f_j-1)^2 的 j 求和)，但最小化普通欧几里德距离(即求和(i_j-f_j)^2) 中的 j 以上也可能有效。

我很容易想到的最好的算法是在第 i 个网格(和为 0)上猜测一个合适的点，然后重复切换到具有最小距离的相邻网格点(和为 0)，直到所有邻居都有一个更大的距离。给定距离函数的凹性，这应该收敛于解。

但该算法似乎很笨拙。谁能做得更好？

How to round floats to integers while preserving their sum? 有相关讨论但它没有达到我正在寻找的最优点。

上下文附录:如果您对此感兴趣(也因为我觉得它很酷)，让我具体说明问题出现的背景。

该问题作为交易模拟的一部分出现(这是更大模拟的一部分)。在每个地点，代理商提供交易多种商品的服务。由于每个地点和商品都是单独处理的，我们可以专注于一个地点和商品并对其进行顺序处理。

每个代理人 j 都有一定数量的货币 c_j 和一定数量的商品 q_j，它们是并且必须保持完整。每个代理人还指定一个实值的、连续的、非负的、非单调递减的需求函数 d_j(p)，它实质上表示代理人在任何给定价格下想要拥有多少单位的商品。

交易按以下步骤执行。

1. 为每个代理人 j 计算预算约束需求函数 b_j(p) = min (d_j(p), q_j + c_j/p)
2. 定义总量 q_tot(q_j 的 j 求和)和总预算约束需求函数 b_tot(p)(b_j(p) 的 j 求和)。
3. 使用牛顿法求 b_tot(p_eq) = q_tot 的均衡价格 p_eq。如果没有均衡价格，返回。
4. 将每个代理人的交易量 f_j 定义为 b_j(p_eq)-q_j(净购买量为正，净销售量为负)。
5. 从计算中删除 f_j 非常小(例如，小于 1/10)的代理
6. 每个q_j应该调整为q_j+f_j，每个c_j应该调整为c_j-p_eq*f_j
7. 问题就在此时出现。 f_j 和 p_eq 是 float ，但 q_j 和 c_j 需要保持整数。为避免创建或销毁货币或商品，数组 f_j 和 (-p_eq*f_j) 需要一致地舍入为整数

Answer 1

当 n = 2

这很简单。定义要使用以下算法计算的 i_k。 i_k 的总和将始终为 0，欧氏距离将始终最小。

if (f[k] < 0)
   i[k] = int(f[k]-0.5);
else
   i[k] = int(f[k]+0.5);

int() 是一个返回 float 的整数部分的函数。它会截断 float 。

说明

让我们定义 e_k 使得 f_k = i_k + e_k。

对于 n = 2，f₀ = -f₁。两个f_k大小相同但符号相反。通过向 0 舍入，这两个误差也具有相同的大小但符号相反。因为 sum(f) = sum(i) + sum(e) 并且 sum(e) 和 sum(f) 等于 0，所以 sum(i) = 0。

除了平衡两个误差的大小外，通过四舍五入到最接近的整数，我们将误差最小化。 sum(e²) 将是最小的。

当 n = 3

我们像上面那样计算 i_k。然后 sum(i) 可能取值 -1、0 或 1。

当 sum(i) = -1 时，我们必须增加 i_k 之一。选择具有最大 e_k 的 i_k(所有 e_k 都是正值)。

当 sum(i) = 1 时，我们必须减少 i_k 之一。选择e_k最小的i_k(所有e_k都是负数)。

说明

当 n=3 时，我们有 3 个错误值 |e_k| < 0.5。结果 |sum(e)|永远不能是 2 或更多。由于sum(i)只能取整数值，sum(e)只能取值-1、0或1，sum(i)也是如此。

|总和(e)| = 1 当所有 e_k 具有相同的符号时。这是因为 |e_k| < 0.5。您总是需要三个相同符号的错误才能达到 1 或 -1。请注意，与符号不同且 sum(i) = 0 的情况相比，这种情况在统计上较少出现。

我们如何决定选择哪个 i_k？

当 sum(i) = 1 时，sum(e) = -1 并且所有 e_k 都为负。我们必须减一 i_k。减少一个 i_k 将通过增加其 e_k 来平衡，因为 sum(i) + sum(e) = 0。因此我们应该选择 e_k 以便增加它，产生最小幅度的错误。这是最接近 -0.5 的 e_k，因此也是最小的 e_k。这确保了 sum(e²) 最小。

当 sum(i) = -1 时，同样的逻辑适用。在这种情况下，sum(e) = 1 并且所有 e_k 都是正数。增加一个 i_k 通过减少它的 e_k 来平衡。通过选择最接近 0.5 的 e_k，我们得到最小的总和 (e²)。

当 n = 4

在这种情况下，sum(i) 仍然限于值 -1、0 和 1。

当sum(i) = -1时，取最大的e_k。

当sum(i) = 1时，取最小的e_k。

说明

|总和(e)|无法达到 2。这就是为什么 sum(i) 被限制为值 -1、0 和 1 的原因。

与 n = 3 的不同之处在于，现在我们可能有 3 或 4 个具有相同符号的 e_k。但是通过应用上述规则，sum(e²) 保持最小值。

泛化

当 n > 4 |sum(e)|可以大于 1。在这种情况下，我们必须修改多个 i_k。

一般算法如下

sum(i) -> m
when m = 0, 
    we are done.
when m < 0, 
    increment the m i<sub>k</sub> with biggest e<sub>k</sub>. 
when m > 0, 
    decrement the m i<sub>k</sub> with smallest e<sub>k</sub>.

这是 python 2.7 代码

def solve(pf):
    n = len(pf)
    
    # construct array pi from pf
    pi = [round(val) for val in pf]
    print "pi~:", pi
    
    # compute the sum of the integers
    m = sum(val for val in pi)
    print "m :", m
    
    # if the sum is zero, we are done
    if m == 0:
        return pi
        
    # compute the errors
    pe = [pf[k]-pi[k] for k in xrange(n)]
    print "pe :", pe

    # correct pi when m is negative    
    while m < 0:
        # locate the pi with biggest error
        biggest = 0
        for k in xrange(1,n):
            if pe[k] > pe[biggest]:
                biggest = k
                
        # adjust this integer i
        pi[biggest] += 1
        pe[biggest] -= 1
        m += 1
    
    # correct pi when m is positive    
    while m > 0:
        # locate the pi with smallest error
        smallest = 0
        for k in xrange(1,n):
            if pe[k] < pe[smallest]:
                smallest = k
        
        # adjust this integer i
        pi[smallest] -= 1
        pe[smallest] += 1
        m -= 1
   
    return pi
    
    
if __name__ == '__main__': 
    print "Example case when m is 0"    
    pf = [1.1, 2.2, -3.3]
    print "pf :", pf
    pi = solve( pf )
    print "pi :", pi
    
    print "Example case when m is 1"    
    pf = [0.6, 0.7, -1.3]
    print "pf :", pf
    pi = solve( pf )
    print "pi :", pi
    
    print "Example case when m is -1"    
    pf = [0.4, 1.4, -1.8]
    print "pf :", pf
    pi = solve( pf )
    print "pi :", pi