view artifacts/src/main/java/org/dive4elements/river/artifacts/model/fixings/fitting/LinearRegressionSqrtErrorFunction.java @ 9646:0380717105ba

Implemented alternative fitting strategy for Log-Linear function.
author Gernot Belger <g.belger@bjoernsen.de>
date Mon, 02 Dec 2019 17:56:15 +0100
parents
children
line wrap: on
line source
/** Copyright (C) 2017 by Bundesanstalt für Gewässerkunde
 * Software engineering by
 *  Björnsen Beratende Ingenieure GmbH
 *  Dr. Schumacher Ingenieurbüro für Wasser und Umwelt
 *
 * This file is Free Software under the GNU AGPL (>=v3)
 * and comes with ABSOLUTELY NO WARRANTY! Check out the
 * documentation coming with Dive4Elements River for details.
 */
package org.dive4elements.river.artifacts.model.fixings.fitting;

import org.apache.commons.math3.analysis.MultivariateFunction;
import org.apache.commons.math3.analysis.UnivariateFunction;
import org.apache.commons.math3.stat.descriptive.DescriptiveStatistics;
import org.apache.commons.math3.stat.regression.SimpleRegression;

/**
 * @author Gernot Belger
 *
 */
final class LinearRegressionSqrtErrorFunction implements MultivariateFunction, UnivariateFunction {

    private final double[] obsDischarges;

    private final double[] obsWaterlevels;

    private final double waterlevelMean;

    private final double waterlevelDeviation;

    private final double[] standardizedWaterlevel;

    private final SQRTFunction sqrtFunction;

    private final TargetFunction targetFunction;

    public LinearRegressionSqrtErrorFunction(final double[] obsDischarges, final double[] obsWaterlevels, final TargetFunction targetFunction) {

        this.obsDischarges = obsDischarges;
        this.obsWaterlevels = obsWaterlevels;

        this.sqrtFunction = new SQRTFunction(obsWaterlevels);
        this.targetFunction = targetFunction;

        final DescriptiveStatistics stats = new DescriptiveStatistics(obsWaterlevels);

        // Compute mean and standard deviation
        this.waterlevelMean = stats.getMean();
        this.waterlevelDeviation = stats.getStandardDeviation();
        // -mittelwert durch standardabweichung /z-transformation)

        // initialize the standardizedSample, which has the same length as the sample
        this.standardizedWaterlevel = new double[obsWaterlevels.length];

        for (int i = 0; i < obsWaterlevels.length; i++)
            this.standardizedWaterlevel[i] = (obsWaterlevels[i] - this.waterlevelMean) / this.waterlevelDeviation;
    }

    @Override
    public double value(final double[] point) {

        final double a = point[0];

        return value(a);
    }

    @Override
    public double value(final double a) {

        final double error = calc_sqrt_trans(a);
        if (Double.isNaN(error))
            return Double.POSITIVE_INFINITY;

        return error;
    }

    private double calc_sqrt_trans(final double a) {

        final double[] waterlevels = calc_stage_trans(a);

        return this.sqrtFunction.calc_sqrt(waterlevels);
    }

    private double[] calc_stage_trans(final double a) {

        final double[] params = calc_parameters_trans(a);

        final double b = params[0];
        final double m = params[1];

        return this.targetFunction.calc_stages(a, b, m);
    }

    public double[] calc_parameters_trans(final double a) {

        // stages_obs_trans = trans_stages(stages_obs, parameter_a);
        final double[] transWaterlevels = trans_stages(this.obsWaterlevels, a);
        final double[] discharges = this.obsDischarges;

        // FIXME: z-nromalizing?

        // printArray(discharges);
        // printArray(transWaterlevels);

        // final DescriptiveStatistics stats = new DescriptiveStatistics(transWaterlevels);
        //
        // // Compute mean and standard deviation
        // final double mean = stats.getMean();
        // final double deviation = stats.getStandardDeviation();
        //
        // final double[] normWaterlevel = new double[transWaterlevels.length];
        // for (int i = 0; i < transWaterlevels.length; i++)
        // normWaterlevel[i] = (transWaterlevels[i] - mean) / deviation;

        final SimpleRegression regression = new SimpleRegression(true);
        for (int i = 0; i < discharges.length; i++)
            regression.addData(discharges[i], transWaterlevels[i]);

        final double normSlope = regression.getSlope();
        final double normIntercept = regression.getIntercept();

        // unnormalize
        // final double slope = normSlope * deviation;
        // final double intercept = normIntercept + mean;
        final double slope = normSlope;
        final double intercept = normIntercept;

        return new double[] { intercept, slope };
    }

    // private void printArray(final double[] discharges) {
    // for (final double d : discharges) {
    // System.out.print(d);
    // System.out.print(", ");
    // }
    // System.out.println(" ");
    // }

    private double[] trans_stages(final double[] waterlevels, final double a) {

        final double[] transWaterlevels = new double[waterlevels.length];
        for (int i = 0; i < transWaterlevels.length; i++)
            transWaterlevels[i] = Math.exp(waterlevels[i] / a);

        return transWaterlevels;
    }
}

http://dive4elements.wald.intevation.org