Python scipy.optimize.leastsq to Java org.apache.commons.mat

Question

我尝试模仿在Python中开发的该算法，该算法基于所见WiFi站点的位置来计算地理位置。

该算法首先使用numpy函数来计算观察到的纬度和纵向的基本加权平均值。为了最大程度地减少可能的WiFi位置错误的影响，它还使用" scipy.optimize.leastsq "方法，以便以统计方式计算，并在可能的情况下进行计算。

我想在Java Android平台上实现相同的行为。

对于所有其他计算，我成功依赖于org.apache.commons.math3。因此，对于最小二乘问题，我从逻辑上尝试依靠https://commons.apache.org/proper/commons-math/userguide/leastsquares.html。

我的问题，如果我很明白，Scipy为我管理雅各布功能定义的复杂性和我差的数学技能不允许我正确地定义了最少的Quarparesproblem模型。我尝试了一些基于此示例的实验，这似乎对我需要的东西很固定，但是结果并不好，因为我不知道如何处理"雅各布"部分。

作为某人为这篇文章所做的，有人可以为我做同样的事情，并尝试以简单的方式解释它？

有关Python部分的工作方式的更多详细信息：

所使用的" scipy.optimize.leastsq"语句是：

(lat, lon), cov_x, info, mesg, ier = 
scipy.optimize.leastsq(func, initial, args=data, full_output=True)

数据是：毫秒/信号强度的纬度/经度/年龄，例如：data = numpy.array([(43.48932915, 1.66561772, 1000, -20), (43.48849093, 1.6648176, 2000, -10), (43.48818612, 1.66615113, 3000, -50)])

初始计算的加权平均纬度/经度，在此示例中：initial = 43.48864654, 1.66550075

功能是

  def func(initial, data):
        return numpy.array([
            geographic distance((float(point[latitude]), float(point[longitude])), (initial[latitude], initial[longitude])).meters * min(math.sqrt(2000.0 / float(point[age])), 1.0) / math.pow(float(point[signal strength]), 2)

结果是：43.4885401095, 1.6648660983

我在Java中进行的实验，我替换了数据值，并且已经更改了计算"模型"的方式。我简化了信号强度和年龄值。但这是一个事实，结果表明这还不够。

double modelI = calculateVincentyDistance(o.getY(), o.getX(), center.getY(), center.getX())* Math.min(Math.sqrt(2000.0/1000.0), 1.0) / Math.pow(-10, 2);

我还将尝试https://github.com/odinsbane/least-squares-in-java，但是我不确定是否正确使用它，因为我不掌握它的工作方式。

fyi，我使用vincenty距离计算，例如，可以用haversine或Euclidean代替。

感谢您的任何帮助！

Answer 1

该代码不容易端口，因为Scipy提供了更通用的最小二乘最小化接口，而Apache Commons Math Math提供了曲线拟合。仍然可以作为曲线拟合来重述许多优化问题。在Python代码中，您最小化

F(current_point) = Sum{ (distance(known_point[i], current_point) * weight[i])^2 } -> min

Java曲线拟合问题有些不同：

F(current_point) = Sum{ (target_value[i] - model[i](current_point))^2  } -> min

可以通过将所有target_value s分配给0并使model[i]计算从current_point到known_point[i]的加权距离来创建等效拟合问题。

在一般情况下，这种问题没有使用公式的精确解决方案，并且使用了一些数值优化方法。这是另一个区别：Java实现明确要求您为优化器提供手段来计算优化函数的导数。如果未提供Dfun，则Python代码似乎使用某种差异。您可以手动地在Java中执行类似的操作，也可以使用FinitedIfferencesDifferentediator，但是对于简单的公式，使用derivativativertucture

明确地对其进行编码可能更容易。

static class PositionInfo {
    public final double latitude;
    public final double longitude;
    public final int ageMs;
    public final int strength;
    public PositionInfo(double latitude, double longitude, int ageMs, int strength) {
        this.latitude = latitude;
        this.longitude = longitude;
        this.ageMs = ageMs;
        this.strength = strength;
    }
    public double getWeight() {
        return Math.min(1.0, Math.sqrt(2000.0 / ageMs)) / (strength * strength);
    }
}

static DerivativeStructure getWeightedEuclideanDistance(double tgtLat, double tgtLong, PositionInfo knownPos) {
    DerivativeStructure varLat = new DerivativeStructure(2, 1, 0, tgtLat); // latitude is 0-th variable of 2 for derivatives up to 1
    DerivativeStructure varLong = new DerivativeStructure(2, 1, 1, tgtLong); // longitude is 1-st variable of 2 for derivatives up to 1
    DerivativeStructure latDif = varLat.subtract(knownPos.latitude);
    DerivativeStructure longDif = varLong.subtract(knownPos.longitude);
    DerivativeStructure latDif2 = latDif.pow(2);
    DerivativeStructure longDif2 = longDif.pow(2);
    DerivativeStructure dist2 = latDif2.add(longDif2);
    DerivativeStructure dist = dist2.sqrt();
    return dist.multiply(knownPos.getWeight());
}
// as in https://en.wikipedia.org/wiki/Haversine_formula
static DerivativeStructure getWeightedHaversineDistance(double tgtLat, double tgtLong, PositionInfo knownPos) {
    DerivativeStructure varLat = new DerivativeStructure(2, 1, 0, tgtLat);
    DerivativeStructure varLong = new DerivativeStructure(2, 1, 1, tgtLong);
    DerivativeStructure varLatRad = varLat.toRadians();
    DerivativeStructure varLongRad = varLong.toRadians();
    DerivativeStructure latDifRad2 = varLat.subtract(knownPos.latitude).toRadians().divide(2);
    DerivativeStructure longDifRad2 = varLong.subtract(knownPos.longitude).toRadians().divide(2);
    DerivativeStructure sinLat2 = latDifRad2.sin().pow(2);
    DerivativeStructure sinLong2 = longDifRad2.sin().pow(2);
    DerivativeStructure summand2 = varLatRad.cos().multiply(varLongRad.cos()).multiply(sinLong2);
    DerivativeStructure sum = sinLat2.add(summand2);
    DerivativeStructure dist = sum.sqrt().asin();
    return dist.multiply(knownPos.getWeight());
}

使用这种准备工作，您可以做这样的事情：

public static void main(String[] args) {
    // latitude/longitude/age in milliseconds/signal strength
    final PositionInfo[] data = new PositionInfo[]{
            new PositionInfo(43.48932915, 1.66561772, 1000, -20),
            new PositionInfo(43.48849093, 1.6648176, 2000, -10),
            new PositionInfo(43.48818612, 1.66615113, 3000, -50)
    };

    double[] target = new double[data.length];
    Arrays.fill(target, 0.0);
    double[] start = new double[2];
    for (PositionInfo row : data) {
        start[0] += row.latitude;
        start[1] += row.longitude;
    }
    start[0] /= data.length;
    start[1] /= data.length;
    MultivariateJacobianFunction distancesModel = new MultivariateJacobianFunction() {
        @Override
        public Pair<RealVector, RealMatrix> value(final RealVector point) {
            double tgtLat = point.getEntry(0);
            double tgtLong = point.getEntry(1);
            RealVector value = new ArrayRealVector(data.length);
            RealMatrix jacobian = new Array2DRowRealMatrix(data.length, 2);
            for (int i = 0; i < data.length; i++) {
                DerivativeStructure distance = getWeightedEuclideanDistance(tgtLat, tgtLong, data[i]);
                //DerivativeStructure distance = getWeightedHaversineDistance(tgtLat, tgtLong, data[i]);
                value.setEntry(i, distance.getValue());
                jacobian.setEntry(i, 0, distance.getPartialDerivative(1, 0));
                jacobian.setEntry(i, 1, distance.getPartialDerivative(0, 1));
            }
            return new Pair<RealVector, RealMatrix>(value, jacobian);
        }
    };

    LeastSquaresProblem problem = new LeastSquaresBuilder()
            .start(start)
            .model(distancesModel)
            .target(target)
            .lazyEvaluation(false)
            .maxEvaluations(1000)
            .maxIterations(1000)
            .build();
    LeastSquaresOptimizer optimizer = new LevenbergMarquardtOptimizer().
            withCostRelativeTolerance(1.0e-12).
            withParameterRelativeTolerance(1.0e-12);
    LeastSquaresOptimizer.Optimum optimum = optimizer.optimize(problem);
    RealVector point = optimum.getPoint();
    System.out.println("Start = " + Arrays.toString(start));
    System.out.println("Solve = " + point);
}

P.S。体重的逻辑对我来说似乎是可疑的。在您引用的问题中，OP具有半径的一些估计，然后是一个明显的重量。使用信号强度的反向平方，该信号强度在对数 dbm中都对我来说似乎很奇怪。