5

ここに画像の説明を入力

上記のような単純な対数曲線があります。セグメント化された線でこの曲線に適合し、これらの線分の開始点と終了点を表示できる Matlab の関数はありますか? matlab でカーブ フィッティング ツールボックスを確認しました。彼らは、1行またはいくつかの関数のいずれかでカーブフィッティングを行っているようです. 1行だけカーブフィッティングしたくありません。

直接的な機能がない場合、同じ目標を達成するための代替手段は何でも構いません。私の目標は、セグメント化された線で曲線を当てはめ、これらのセグメントの終点の位置を取得することです。

4

3 に答える 3

7

まず第一に、あなたの問題はカーブフィッティングとは呼ばれていません。カーブ フィッティングとは、データがあり、ある意味でそれを説明する最適な関数を見つけることです。一方、関数の区分線形近似を作成したいとします。

次の戦略を提案します。

  1. 手動でセクションに分割します。セクション サイズは導関数に依存する必要があります。大きな導関数 -> 小さなセクション
  2. セクション間のノードで関数をサンプリングします
  3. 上記のポイントを通過する線形補間を見つけます。

これを行うコードの例を次に示します。セクションの量が少ないにもかかわらず、赤い線 (補間) が元の関数に非常に近いことがわかります。これは、アダプティブ セクション サイズが原因で発生します。

ここに画像の説明を入力

function fitLogLog()
   x = 2:1000;
   y = log(log(x));

   %# Find section sizes, by using an inverse of the approximation of the derivative
   numOfSections = 20;
   indexes = round(linspace(1,numel(y),numOfSections));
   derivativeApprox = diff(y(indexes));
   inverseDerivative = 1./derivativeApprox;
   weightOfSection =  inverseDerivative/sum(inverseDerivative);   
   totalRange = max(x(:))-min(x(:));   
   sectionSize = weightOfSection.* totalRange;

   %# The relevant nodes
   xNodes = x(1) + [ 0 cumsum(sectionSize)];
   yNodes = log(log(xNodes));

   figure;plot(x,y);   
   hold on;
   plot (xNodes,yNodes,'r');
   scatter (xNodes,yNodes,'r');
   legend('log(log(x))','adaptive linear interpolation');
end
于 2012-09-23T22:29:59.580 に答える
5

アンドレイの適応ソリューションは、より正確な全体的な適合を提供します。ただし、必要なのが固定長のセグメントである場合は、すべての近似値の完全なセットも返すメソッドを使用して、機能するはずの何かがあります。速度が必要な場合はベクトル化できます。

Nsamp = 1000;     %number of data samples on x-axis
x = [1:Nsamp];    %this is your x-axis
Nlines = 5;       %number of lines to fit

fx = exp(-10*x/Nsamp);  %generate something like your current data, f(x)
gx = NaN(size(fx));     %this will hold your fitted lines, g(x)

joins = round(linspace(1, Nsamp, Nlines+1));  %define equally spaced breaks along the x-axis

dx = diff(x(joins));   %x-change
df = diff(fx(joins));  %f(x)-change

m = df./dx;   %gradient for each section

for i = 1:Nlines
   x1 = joins(i);   %start point
   x2 = joins(i+1); %end point
   gx(x1:x2) = fx(x1) + m(i)*(0:dx(i));   %compute line segment
end

subplot(2,1,1)
h(1,:) = plot(x, fx, 'b', x, gx, 'k', joins, gx(joins), 'ro');
title('Normal Plot')

subplot(2,1,2)
h(2,:) = loglog(x, fx, 'b', x, gx, 'k', joins, gx(joins), 'ro');
title('Log Log Plot')

for ip = 1:2
    subplot(2,1,ip)
    set(h(ip,:), 'LineWidth', 2)
    legend('Data', 'Piecewise Linear', 'Location', 'NorthEastOutside')
    legend boxoff
end

MATLABプロット出力

于 2012-09-23T23:33:32.597 に答える
0

これはこの質問に対する正確な回答ではありませんが、検索に基づいてここにたどり着いたので、平均 (または散布図の間隔データの中央値、またはその他の関数)。

まず、回帰を使用した関連しているがより洗練された代替手段は、明らかにウィキペディアのページにリストされている MATLAB コードがあり、多変量適応回帰スプラインです。

ここでの解決策は、重複する間隔の平均を計算してポイントを取得することです

function [x, y] = intervalAggregate(Xdata, Ydata, aggFun, intStep, intOverlap)
% intOverlap in [0, 1); 0 for no overlap of intervals, etc.
% intStep    this is the size of the interval being aggregated.

minX = min(Xdata);
maxX = max(Xdata);

minY = min(Ydata);
maxY = max(Ydata);

intInc = intOverlap*intStep; %How far we advance each iteraction.
if intOverlap <= 0
   intInc = intStep;
end
nInt = ceil((maxX-minX)/intInc); %Number of aggregations

parfor i = 1:nInt
    xStart = minX + (i-1)*intInc;
    xEnd   = xStart + intStep;
    intervalIndices = find((Xdata >= xStart) & (Xdata <= xEnd));
    x(i) = aggFun(Xdata(intervalIndices));
    y(i) = aggFun(Ydata(intervalIndices));
end

たとえば、ペアになった X と Y のデータの平均を計算するために、長さ 0.1 の間隔が互いに約 1/3 重なっていると便利でした (散布図を参照)。

Xdat と Ydat の散布図の例

[x,y] = intervalAggregate(Xdat, Ydat, @平均, 0.1, 0.333)

×＝

列 1 ～ 8

0.0552    0.0868    0.1170    0.1475    0.1844    0.2173    0.2498    0.2834

9列目から15列目

0.3182    0.3561    0.3875    0.4178    0.4494    0.4671    0.4822

y =

列 1 ～ 8

0.9992    0.9983    0.9971    0.9955    0.9927    0.9905    0.9876    0.9846

9列目から15列目

0.9803    0.9750    0.9707    0.9653    0.9598    0.9560    0.9537

x が増加すると、y がわずかに減少する傾向があることがわかります。そこから、線分を描画したり、他の種類のスムージングを実行したりするのは簡単です。

(このソリューションをベクトル化しようとしたわけではないことに注意してください。Xdata がソートされている場合は、はるかに高速なバージョンが想定される可能性があります。)

于 2014-01-22T00:01:12.703 に答える