[Octave]Использование fminunc не всегда дает последовательное решение

Я пытаюсь найти коэффициенты в уравнении для моделирования переходной характеристики двигателя, которая имеет форму 1-e^x. Уравнение, которое я использую для моделирования, имеет вид

a(1)*t^2 + a(2)*t^3 + a(3)*t^3 + ...

(Это получено в исследовательской статье, используемой для определения параметров двигателя)

Иногда использование fminunc для нахождения коэффициентов работает нормально, и я получаю хороший результат, и он достаточно хорошо соответствует обучающим данным. В других случаях возвращаемые коэффициенты ужасны (будут чрезвычайно выше, чем должен быть вывод, и на порядки). Это особенно происходит, когда я начал использовать термины более высокого порядка: использование любой модели, использующей x^8 или выше (x^9, x^10, x^11 и т. д.), всегда приводит к плохим результатам.

Поскольку иногда это работает, я не могу понять, почему моя реализация была бы неправильной. Я пробовал fminunc, предоставляя градиенты, а также не предоставляя градиенты, но разницы нет. Я рассматривал возможность использования других функций для расчета коэффициентов, таких как polyfit, но в этом случае у него должны быть члены, которые повышаются от 1 до члена высшего порядка, но модель, которую я использую, имеет наименьшую мощность при 2 .

Вот основной код:

clear;

%Overall Constants
max_power = 7;

%Loads in data
%data = load('TestData.txt');
load testdata.mat

%Sets data into variables
indep_x = data(:,1); Y = data(:,2);

%number of data points
m = length(Y);

%X is a matrix with the independant variable
exps = [2:max_power];
X_prime = repmat(indep_x, 1, max_power-1); %Repeats columns of the indep var
X = bsxfun(@power, X_prime, exps);

%Initializes theta to rand vals
init_theta = rand(max_power-1,1);

%Sets up options for fminunc
options = optimset( 'MaxIter', 400, 'Algorithm', 'quasi-newton');

%fminunc minimizes the output of the cost function by changing the theta paramaeters
[theta, cost] = fminunc(@(t)(costFunction(t, X, Y)), init_theta, options)

%
Y_line = X * theta;

figure;
hold on; plot(indep_x, Y, 'or');
hold on; plot(indep_x, Y_line, 'bx');

А вот и costFunction:

function [J, Grad] = costFunction (theta, X, Y)
   %# of training examples

   m = length(Y);

    %Initialize Cost and Grad-Vector
    J = 0;
    Grad = zeros(size(theta));

    %Poduces an output based off the current values of theta
    model_output = X * theta;

    %Computes the squared error for each example then adds them to get the total error
    squared_error = (model_output - Y).^2;
    J = (1/(2*m)) * sum(squared_error);

    %Computes the gradients for each theta t
    for t = 1:size(theta, 1)
        Grad(t) = (1/m) * sum((model_output-Y) .* X(:, t));
    end

endfunction

Любая помощь или совет будут оценены.