Почему я получаю System.ArithmeticException в этом примере, выполняя нелинейную регрессию к гауссову

Я переношу некоторый код, который я прототипировал в Matlab, на C #, и мне нужно выполнить нелинейную регрессию гауссова к моим демонстрационным данным. В настоящее время я пробую Accord.Net и придумал приведенный ниже код, следуя примеру здесь.

К сожалению, я получаю эти исключения и не понимаю, почему:

System.ArithmeticException: при вычислении ошибки получено не конечное число. Пожалуйста, убедитесь, что во входных данных нет постоянных столбцов."

с сообщением об ошибке:

Расчет ошибки дал не конечное число. Пожалуйста, убедитесь, что во входных данных нет постоянных столбцов.

Это код, который я придумал:

    public static double[] FitGaussian()
    {
        var nls = CreateGaussianFitObject(new[] { 1.43e+04, 0.05093, 0.8098 }, 100, 0);

        double[,] data =
        {
            { -0.2,    12973.3071 },
            { -0.1,    13846.1569 },
            { -0,      14243.9094 },
            {  0.1,    14215.6044 },
            {  0.2,    13840.9077 },
        };

        // Extract inputs and outputs
        double[][] inputs = data.GetColumn(0).ToJagged();
        double[] outputs = data.GetColumn(1);

        var regression = nls.Learn(inputs, outputs);

        var result = regression.Coefficients;
        return result;
    }

    public static NonlinearLeastSquares CreateGaussianFitObject(double[] StartValues, int maxIterations, double tolerance)
    {
        Func<double, double, double, double> expFnc = (mu, sig, x) => m.Exp(-m.Pow((x - mu) / sig, 2));

        var nls = new NonlinearLeastSquares()
        {
            NumberOfParameters = 3,

            // Reference Gaussian Distributioon used in Matlab: f(x) =  a1*exp(-((x-b1)/c1)^2)
            //Function = (w, x) => w[0] * System.Math.Exp(-System.Math.Pow((x[0] - w[1]) / w[2], 2)),
            Function = (w, x) => w[0] * expFnc(w[1],w[2],x[0]),

            // Derivative in respect to the weights:
            Gradient = (w, x, r) =>
            {
                r[0] = expFnc(w[1], w[2], x[0]); // e^(-(x - b1)^2/c1^2)   <=>   diff a1*exp(-((x-b1)/c1)^2) w.r.t. a1
                r[1] = (2 * w[0] * (x[0] - w[1]) * expFnc(w[1], w[2], x[0])) / m.Pow(w[2],2); // (2 a1 (x - b1) e^(-(x - b1)^2/c1^2))/c1^2   <=>   diff a1*exp(-((x-b1)/c1)^2) w.r.t. b1
                r[2] = (2 * w[0] * m.Pow((x[0] - w[1]),2) * expFnc(w[1], w[2], x[0])) / m.Pow(w[2],3); // (2 a1 (x - b1)^2 e^(-(x - b1)^2/c1^2))/c1^3   <=>   diff a1*exp(-((x-b1)/c1)^2) w.r.t. c1
            },

            Algorithm = new LevenbergMarquardt()
            {
                MaxIterations = maxIterations,
                Tolerance = tolerance
            }
        };

        return nls;
    }

Я дважды проверил функцию градиента (которую я подозреваю как виновника), и все в порядке.

Соответствующий код Matlab таков:

x = [-0.2,-0.1,0,0.1,0.2];
y = [12973.3071,13846.1569,14243.9094,14215.6044,13840.9077];
f = fit(x.',y.','gauss1');
fitObject = f;

который дает этот вывод:

f = 

     General model Gauss1:
     f(x) =  a1*exp(-((x-b1)/c1)^2)
     Coefficients (with 95% confidence bounds):
       a1 =    1.43e+04  (1.419e+04, 1.441e+04)
       b1 =     0.05093  (0.03501, 0.06686)
       c1 =      0.8098  (0.7263, 0.8932)

Как видите, я даже пытался дать ожидаемые коэффициенты из Matlab в качестве начального значения для реализации на C#.

Кроме того, меня смущает приведенный выше пример, поскольку он кажется мне поменяли местами коэффициенты и входы(?).


person packoman    schedule 29.11.2017    source источник


Ответы (1)


Я забыл использовать аргумент StartValues в конструкторе NonlinearLeastSquares(), что вызвало исключение. Следует читать:

    var nls = new NonlinearLeastSquares()
    {
        NumberOfParameters = 3,
        StartValues = StartValues,
        // Reference Gaussian Distributioon used in Matlab: f(x) =  a1*exp(-((x-b1)/c1)^2)
        //Function = (w, x) => w[0] * System.Math.Exp(-System.Math.Pow((x[0] - w[1]) / w[2], 2)),
        Function = (w, x) => w[0] * expFnc(w[1],w[2],x[0]),

        // Derivative in respect to the weights:
        Gradient = (w, x, r) =>
        {
            r[0] = expFnc(w[1], w[2], x[0]); // e^(-(x - b1)^2/c1^2)   <=>   diff a1*exp(-((x-b1)/c1)^2) w.r.t. a1
            r[1] = (2 * w[0] * (x[0] - w[1]) * expFnc(w[1], w[2], x[0])) / m.Pow(w[2],2); // (2 a1 (x - b1) e^(-(x - b1)^2/c1^2))/c1^2   <=>   diff a1*exp(-((x-b1)/c1)^2) w.r.t. b1
            r[2] = (2 * w[0] * m.Pow((x[0] - w[1]),2) * expFnc(w[1], w[2], x[0])) / m.Pow(w[2],3); // (2 a1 (x - b1)^2 e^(-(x - b1)^2/c1^2))/c1^3   <=>   diff a1*exp(-((x-b1)/c1)^2) w.r.t. c1
        },

        Algorithm = new LevenbergMarquardt()
        {
            MaxIterations = maxIterations,
            Tolerance = tolerance
        }
    };

По-видимому, StartValues установлен в массив нулей, если он не установлен явно. Это вызвало деление на ноль и, следовательно, исключение. Я узнал это, добавив

Debug.WriteLine($"Weights: w[0]: {w[0]}, w[1]: {w[1]}, w[2]: {w[2]}");
Debug.WriteLine($"Result: r[0]: {r[0]}, r[1]: {r[1]}, r[2]: {r[2]}");

к лямбда-выражению для свойства Gradient.

person packoman    schedule 02.12.2017