最新のリリースでは、このページがまだ翻訳されていません。 このページの最新版は英語でご覧になれます。

曲線近似と分布近似

この例では、曲線近似と分布近似を実行する方法を示し、どのような場合にこれらの手法が適しているかについて説明します。

曲線近似と分布近似の選択

曲線近似と分布近似は、異なるタイプのデータ解析です。

  • 曲線近似は、応答変数を予測子変数の関数としてモデル化する場合に使用します。

  • 分布近似は、単一変数の確率分布をモデル化する場合に使用します。

曲線近似

以下の実験データでは、予測子変数は、薬物の摂取後の時間 time です。応答変数は、血流中の薬物濃度 conc です。応答データ conc のみが実験誤差の影響を受けると仮定します。

time = [ 0.1   0.1   0.3   0.3   1.3   1.7   2.1   2.6   3.9   3.9 ...
         5.1   5.6   6.2   6.4   7.7   8.1   8.2   8.9   9.0   9.5 ...
         9.6  10.2  10.3  10.8  11.2  11.2  11.2  11.7  12.1  12.3 ...
        12.3  13.1  13.2  13.4  13.7  14.0  14.3  15.4  16.1  16.1 ...
        16.4  16.4  16.7  16.7  17.5  17.6  18.1  18.5  19.3  19.7]';
conc = [0.01  0.08  0.13  0.16  0.55  0.90  1.11  1.62  1.79  1.59 ...
        1.83  1.68  2.09  2.17  2.66  2.08  2.26  1.65  1.70  2.39 ...
        2.08  2.02  1.65  1.96  1.91  1.30  1.62  1.57  1.32  1.56 ...
        1.36  1.05  1.29  1.32  1.20  1.10  0.88  0.63  0.69  0.69 ...
        0.49  0.53  0.42  0.48  0.41  0.27  0.36  0.33  0.17  0.20]';

血中濃度を時間の関数としてモデル化したいとします。conctime に対してプロットします。

plot(time,conc,'o');
xlabel('Time');
ylabel('Blood Concentration');

conctime の関数として 2 パラメーターのワイブル曲線に従うと仮定します。ワイブル曲線には、次の形状とパラメーターがあります。

y=c(x/a)(b-1)e-(x/a)b,

ここで、a は水平方向のスケーリング、b は形状パラメーター、c は垂直方向のスケーリングです。

非線形最小二乗法を使用して、ワイブル モデルをあてはめます。

modelFun =  @(p,x) p(3) .* (x./p(1)).^(p(2)-1) .* exp(-(x./p(1)).^p(2));
startingVals = [10 2 5];
nlModel = fitnlm(time,conc,modelFun,startingVals);

ワイブル曲線をデータ上にプロットします。

xgrid = linspace(0,20,100)';
line(xgrid,predict(nlModel,xgrid),'Color','r');

あてはめたワイブル モデルには問題があります。fitnlm は、実験誤差が加法的であり、一定している分散をもつ対称的な分布に由来すると仮定します。しかし、散布図には、誤差分散が曲線の高さに比例することが示されています。さらに、誤差が加法的かつ対称的であると、負の血中濃度測定値が可能になります。

乗法的誤差が対数スケールで対称的であると仮定する方が現実的です。この仮定に基づき、両辺の対数をとって、ワイブル曲線をデータにあてはめます。非線形最小二乗法を使用して曲線をあてはめます。

log(y)=log(c)+(b-1)log(x/a)-(x/a)b.

nlModel2 = fitnlm(time,log(conc),@(p,x) log(modelFun(p,x)),startingVals);

新しい曲線を既存のプロットに追加します。

line(xgrid,exp(predict(nlModel2,xgrid)),'Color',[0 .5 0],'LineStyle','--');
legend({'Raw Data','Additive Errors Model','Multiplicative Errors Model'});

モデル オブジェクト nlModel2 には、精度の推定値が含まれています。モデルの適合度をチェックすることをお勧めします。たとえば、対数スケール上に残差をプロットして、乗法的誤差の分散が一定であるという仮定をチェックします。

この例では、乗法的誤差モデルを使用しても、モデルの予測にはほとんど影響がありません。モデルのタイプによる影響が大きい例については、線形性への変換による非線形モデルの近似の注意点を参照してください。

曲線近似のための関数

  • Statistics and Machine Learning Toolbox™ には、モデルをあてはめるための関数として、非線形最小二乗モデル用のfitnlm、一般化線形モデル用のfitglm、ガウス過程回帰モデル用のfitrgp、サポート ベクター マシン回帰モデル用のfitrsvmが含まれています。

  • Curve Fitting Toolbox™ には、曲線近似のタスクを簡素化するコマンド ライン ツールおよびグラフィカル ツールが用意されています。たとえば、さまざまなモデルに対する開始係数値の自動選択や、堅牢でノンパラメトリックな近似法が用意されています。

  • Optimization Toolbox™ には、係数に制約があるモデルの解析など、複雑なタイプの曲線近似解析を実行するための関数があります。

  • MATLAB® 関数polyfitは多項式モデルをあてはめます。MATLAB 関数fminsearchは他の種類の曲線近似に役立ちます。

分布近似

電気部品の寿命の分布をモデル化したいとします。変数 life は、同じ電気部品 50 個が故障するまでの時間の測定値です。

life = [ 6.2 16.1 16.3 19.0 12.2  8.1  8.8  5.9  7.3  8.2 ...
        16.1 12.8  9.8 11.3  5.1 10.8  6.7  1.2  8.3  2.3 ...
         4.3  2.9 14.8  4.6  3.1 13.6 14.5  5.2  5.7  6.5 ...
         5.3  6.4  3.5 11.4  9.3 12.4 18.3 15.9  4.0 10.4 ...
         8.7  3.0 12.1  3.9  6.5  3.4  8.5  0.9  9.9  7.9]';

ヒストグラムでデータを可視化します。

binWidth = 2;
lastVal = ceil(max(life));
binEdges = 0:binWidth:lastVal+1;
h = histogram(life,binEdges);
xlabel('Time to Failure');
ylabel('Frequency');
ylim([0 10]);

多くの場合に寿命データはワイブル分布に従うので、1 つのアプローチとして、前述した曲線近似の例のワイブル曲線をヒストグラムにあてはめることが考えられます。このアプローチを試すには、ヒストグラムを一連の点 (x,y) に変換してから、曲線をこれらの点にあてはめます。ここで、x はビンの中心、y はビンの高さです。

counts = histcounts(life,binEdges);
binCtrs = binEdges(1:end-1) + binWidth/2;
h.FaceColor = [.9 .9 .9];
hold on
plot(binCtrs,counts,'o');
hold off

ただし、曲線をヒストグラムにあてはめることには問題があり、通常は推奨されません。

  1. このプロセスは、最小二乗近似の基本的な仮定に反しています。ビンの個数は非負なので、測定誤差は対称的にはなりません。また、分布の裾と中央ではビンの個数の変動性が異なります。そして、ビンの個数は、合計が一定なので、独立した測定値ではありません。

  2. ヒストグラムは経験的確率密度関数 (pdf) をスケーリングしたバージョンなので、ワイブル曲線をバーの高さにあてはめる場合、曲線を制約しなければなりません。

  3. 連続的なデータの場合、データではなくヒストグラムに曲線をあてはめると、情報が失われます。

  4. ヒストグラムのバーの高さは、ビンのエッジと幅の選択に依存します。

多くのパラメトリック分布の場合、これらの問題が回避されるので、パラメーターの推定には最尤法の方が優れています。ワイブル pdf の形状は、ワイブル曲線とほとんど同じです。

y=(b/a)(x/a)(b-1)e-(x/a)b.

ただし、この関数を積分すると 1 になる必要があるので、スケール パラメーター cb/a に置き換えられています。最尤法を使用してワイブル分布をデータにあてはめるには、fitdist を使用し、分布名として 'Weibull' を指定します。最小二乗法と異なり、最尤法では、pdf とバーの高さとの差分二乗和を最小化せずに、スケーリングされたヒストグラムに最も一致するワイブル pdf を求めます。

pd = fitdist(life,'Weibull');

データのスケーリングされたヒストグラムをプロットし、あてはめた pdf を重ね合わせます。

h = histogram(life,binEdges,'Normalization','pdf','FaceColor',[.9 .9 .9]);
xlabel('Time to Failure');
ylabel('Probability Density');
ylim([0 0.1]);
xgrid = linspace(0,20,100)';
pdfEst = pdf(pd,xgrid);
line(xgrid,pdfEst)

モデルの適合度をチェックすることをお勧めします。

通常は曲線をヒストグラムにあてはめることは推奨されませんが、このプロセスが適しているケースもあります。たとえば、カスタム一変量分布の近似を参照してください。

分布近似のための関数

  • Statistics and Machine Learning Toolbox™ には、確率分布オブジェクトをデータにあてはめる関数fitdistが含まれています。また、最尤法を使用してパラメトリック分布をあてはめるための専用の近似関数 (wblfitなど)、専用の近似関数がないカスタム分布をあてはめるための関数mle、ノンパラメトリック分布モデルをデータにあてはめるための関数ksdensityも含まれています。

  • Statistics and Machine Learning Toolbox にはDistribution Fitterアプリもあります。このツールは、可視化や診断プロットの生成など、分布近似における多くの作業を簡素化します。

  • Optimization Toolbox™ の関数を使用すると、パラメーターに制約がある分布など、複雑な分布をあてはめることができます。

  • MATLAB® 関数fminsearchは、最尤法による分布近似を提供します。

参考

| | | | | | | | |

関連するトピック