几何中位数

在正常中位数中，我们找到一个距离总和最小的点。类似的概念适用于二维空间。
给定二维空间中的N个点，任务是找出一个点(x, y)，从该点到输入点的距离总和最小（也称为最小距离的中心）。
例子：

Input: (1, 1), (3, 3) 
Output: Geometric Median = (2, 2) with minimum distance = 2.82843
Input: (0, 0), (0, 0), (0, 12) 
Output: Geometric Median = (0, 0) with minimum distance = 12 
 

方法：
乍一看，这个问题似乎要求我们找到给定输入点的几何中心点（换句话说，质心）的中点。由于它是输入的“中心”点，从中心到所有给定输入点的距离总和应该自动最小化。这个过程类似于找到重心$N$离散质量粒子。第一个示例测试用例甚至给出了正确的答案。但是当我们将相同的逻辑应用于第二个示例时会发生什么？
我们可以清楚地看到几何中心或质心$(0, 0), (0, 0), (0, 12)$在$(0, 4)$ .所以根据欧几里德距离公式，从质心到所有 3 个输入点的总距离为$4+4+8 = 16$但最佳点应该是$(0, 0)$ ，给我们一个总距离$12$那么，我们错在哪里呢？
直观地，您可以认为输入点的质心为我们提供了输入点的算术平均值。但是我们需要的是输入点的中心趋势，以便达到该中心趋势（或换句话说，欧几里得距离）的成本最小化。这被称为一组点的几何中值。从概念上讲，中值与给定输入的均值有很大不同。
没有任何定义的正确算法来查找几何中位数。我们为解决这类问题所做的是逼近一个解决方案并确定我们的解决方案是否确实是几何中值。
算法
有两个重要的变量：

• current_point – 存储可能是几何中值的点的 x 和 y 坐标。
• minimum_distance – 存储从 current_point 到所有输入点的欧几里得距离的总和。

在每次逼近之后，如果我们找到一个距离总和更低的新点，那么我们更新当前点的值和到新点和新距离的最小距离的值。
首先，我们找到给定点的质心，将其作为当前点（或中位数）并将距离之和存储在最小距离中。然后，我们遍历给定的输入点，依次假设每个输入点是中位数，然后计算到其他点的距离。如果这个距离小于 minimum_distance，那么我们将 current_point 和 minimum_distance 的旧值更新为新值。否则，旧值保持不变。
然后我们进入一个while循环。在该循环内，我们总共从 current_point 移动 test_distance 的距离（我们假设本例中的 test_distance 为 1000） $4$方向（左、上、右、下）。因此我们得到$4$新点。然后我们计算这些新点到给定输入点的距离。如果这个距离总和低于之前的 minimum_distance ，那么我们将 current_point 和 minimum_distance 的旧值更新为新值并重复 while 循环。否则，我们将 test_distance 除以$2$然后重复while循环。
while 循环的终止条件是一个称为“lower_limit”的特定值。值越低，我们的近似精度越高。当lower_limit 超过test_distance 时循环终止。
下面是上述方法的实现：

CPP

// C++ implementation of the approach
#include 
using namespace std;
 
// To store a point in 2-D space
struct Point {
    double x, y;
};
 
// Test points. These points are the left,
// up, right and down relative neighbours
// (arranged circularly) to the
// current_point at a distance of
// test_distance from current_point
Point test_point[] = { { -1.0, 0.0 },
                       { 0.0, 1.0 },
                       { 1.0, 0.0 },
                       { 0.0, -1.0 } };
 
// Lowest Limit till which we are going
// to run the main while loop
// Lower the Limit higher the accuracy
double lower_limit = 0.01;
 
// Function to return the sum of Euclidean
// Distances
double distSum(Point p,
                        Point arr[], int n)
{
    double sum = 0;
    for (int i = 0; i < n; i++) {
        double distx = abs(arr[i].x - p.x);
        double disty = abs(arr[i].y - p.y);
        sum += sqrt((distx * distx) + (disty * disty));
    }
 
    // Return the sum of Euclidean Distances
    return sum;
}
 
// Function to calculate the required
// geometric median
void geometricMedian(Point arr[], int n)
{
 
    // Current x coordinate and y coordinate
    Point current_point;
 
    for (int i = 0; i < n; i++) {
        current_point.x += arr[i].x;
        current_point.y += arr[i].y;
    }
 
    // Here current_point becomes the
    // Geographic MidPoint
    // Or Center of Gravity of equal
    // discrete mass distributions
    current_point.x /= n;
    current_point.y /= n;
 
    // minimum_distance becomes sum of
    // all distances from MidPoint to
    // all given points
    double minimum_distance =
       distSum(current_point, arr, n);
 
    int k = 0;
    while (k < n) {
        for (int i = 0; i < n, i != k; i++) {
            Point newpoint;
            newpoint.x = arr[i].x;
            newpoint.y = arr[i].y;
            double newd =
                   distSum(newpoint, arr, n);
            if (newd < minimum_distance) {
                minimum_distance = newd;
                current_point.x = newpoint.x;
                current_point.y = newpoint.y;
            }
        }
        k++;
    }
 
    // Assume test_distance to be 1000
    double test_distance = 1000;
    int flag = 0;
 
    // Test loop for approximation starts here
    while (test_distance > lower_limit) {
 
        flag = 0;
 
        // Loop for iterating over all 4 neighbours
        for (int i = 0; i < 4; i++) {
 
            // Finding Neighbours done
            Point newpoint;
            newpoint.x = current_point.x
                 + (double)test_distance * test_point[i].x;
            newpoint.y = current_point.y
                 + (double)test_distance * test_point[i].y;
 
            // New sum of Euclidean distances
            // from the neighbor to the given
            // data points
            double newd = distSum(newpoint, arr, n);
 
            if (newd < minimum_distance) {
 
                // Approximating and changing
                // current_point
                minimum_distance = newd;
                current_point.x = newpoint.x;
                current_point.y = newpoint.y;
                flag = 1;
                break;
            }
        }
 
        // This means none of the 4 neighbours
        // has the new minimum distance, hence
        // we divide by 2 and reiterate while
        // loop for better approximation
        if (flag == 0)
            test_distance /= 2;
    }
 
    cout << "Geometric Median = ("
         << current_point.x << ", "
         << current_point.y << ")";
    cout << " with minimum distance = "
         << minimum_distance;
}
 
// Driver code
int main()
{
 
    int n = 2;
    Point arr[n];
    arr[0].x = 1;
    arr[0].y = 1;
    arr[1].x = 3;
    arr[1].y = 3;
    geometricMedian(arr, n);
 
    return 0;
}

输出：

Geometric Median = (2, 2) with minimum distance = 2.82843

参考资料：几何中值，最小距离中心