matlab - 係数にいくつかの制約がある多項式近似 matlab

Question

次の種類の関数で補間する必要があるデータがあります。

f(x) = ax4 + bx2 + c

と。a > 0_ b ≤ 0残念ながら、MATLAB ではpolyfit多項式の係数に対する制約は許可されていません。これを行うMATLAB関数があるかどうかは誰にもわかりませんか? それ以外の場合、どのように実装できますか?

事前にどうもありがとうございました、

エリザベッタ

Answer 1

を使用して、目的関数を手動で定義するfminsearchことができます。fminunc

または、問題を少し異なる方法で定義することもできます。

f(x) = a2x4 - b2x2 + c

これで、新しいaとを制約なしbで最適化できますが、探している最終のとが正 (負の場合) であることを確認できます。ab

Answer 2

制約がなければ、問題は単純な線形システムとして記述および解決できます。

% Your design matrix ([4 2 0] are the powers of the polynomial)
A = bsxfun(@power, your_X_data(:), [4 2 0]);

% Best estimate for the coefficients, [a b c], found by 
% solving A*[a b c]' = y in a least-squares sense
abc = A\your_Y_data(:)

これらの制約は、制約されたモデルが実際にデータの基礎となる場合、もちろん自動的に満たされます。例えば、

% some example factors
a = +23.9;
b = -15.75;
c = 4;

% Your model
f = @(x, F) F(1)*x.^4 + F(2)*x.^2 + F(3);

% generate some noisy XY data
x = -1:0.01:1;
y = f(x, [a b c]) + randn(size(x));

% Best unconstrained estimate a, b and c from the data
A = bsxfun(@power, x(:), [4 2 0]);
abc = A\y(:);

% Plot results
plot(x,y, 'b'), hold on
plot(x, f(x, abc), 'r')
xlabel('x (nodes)'), ylabel('y (data)')

ただし、その制約付きモデルによって正確に記述されていないデータに制約を課すと、問題が発生する可能性があります。

% Note: same data, but flipped signs 
a = -23.9;
b = +15.75;
c = 4;

f = @(x, F) F(1)*x.^4 + F(2)*x.^2 + F(3);

% generate some noisy XY data
x = -1:0.01:1;
y = f(x, [a b c]) + randn(size(x));

% Estimate a, b and c from the data, Forcing a>0 and b<0
abc = fmincon(@(Y) sum((f(x,Y)-y).^2), [0 0 0], [-1 0 0; 0 +1 0; 0 0 0], zeros(3,1));

% Plot results
plot(x,y, 'b'), hold on
plot(x, f(x, abc), 'r')
xlabel('x (nodes)'), ylabel('y (data)')

(このソリューションには がありa == 0、モデルの選択が正しくないことを示しています)。

の正確な等価性がa == 0問題になる場合: もちろん、 を設定しても違いはありませんa == eps(0)。数値的には、これは実際のデータでは目立ちませんが、それでもゼロではありません。

とにかく、モデルが適切に選択されておらず、制約がすべてを機能させるための「修正」であるか、適合させる前にデータを実際にバイアス解除/再スケーリングする必要があるか、または同様の前提条件が適用されるのではないかという疑いがあります(私は人々がこの種のことをしているのをよく見てきましたので、はい、私はこの点で少し偏っています:)。

では、これらの制約の背後にある本当の理由は何ですか?

Answer 3

カーブ フィッティング ツールボックスがある場合、fitでは、'upper' および 'lower' オプションを使用して制約を設定できます。あなたは次のようなものが欲しいでしょう。

M=fit(x, f, 'poly4', 'upper', [-inf, 0, -inf, 0, -inf], 'lower', [0, 0, 0, 0, -inf]);

-inf を使用して、特定の係数を制約なしに設定することに注意してください。

これにより、関連する係数を持つ cfit オブジェクトが得られます。たとえば、x^4 項の M.p1 を使用してこれらにアクセスできます。または、 fevalを使用して、任意のポイントで関数を評価できます。

最適化ツールボックスでも lsqcurvefit を使用して同様のことができると思います。