19

次のようなプロットがいくつかあります。

ここに画像の説明を入力

x 軸の約 5.5 から 8 の間の勾配を見つけるには、どのような方法があるのでしょうか。このようなプロットがいくつかある場合、傾きの値を自動的に見つける方法があるかどうかもっと疑問に思っています。

助言がありますか?

ployfit()、または線形回帰を考えています。問題は、値を自動的に見つける方法がわからないことです。

4

4 に答える 4

38

A generic way to find linear parts in data sets is to calculate the second derivative of the function, and see where it is (close to) zero. There are several things to consider on the way to the solution:

  • How to calculate the second derivative of noisy data? One fast and simple method, that can easily be adapted to different noise levels, data set sizes and expected lengths of the linear patch, is to convolve the data with a convolution kernel equal to the second derivative of a Gaussian. The adjustable part is the width of the kernel.

  • What does "close to zero" mean in your context? To answer this question, you have to experiment with your data.

  • The results of this method could be used as an input to the chi^2-method described above, to identify candidate regions in the data set.

Here some source code that will get you started:

from matplotlib import pyplot as plt

import numpy as np

# create theoretical data
x_a = np.linspace(-8,0, 60)
y_a = np.sin(x_a)
x_b = np.linspace(0,4,30)[1:]
y_b = x_b[:]
x_c = np.linspace(4,6,15)[1:]
y_c = np.sin((x_c - 4)/4*np.pi)/np.pi*4. + 4
x_d = np.linspace(6,14,120)[1:]
y_d = np.zeros(len(x_d)) + 4 + (4/np.pi)

x = np.concatenate((x_a, x_b, x_c, x_d))
y = np.concatenate((y_a, y_b, y_c, y_d))


# make noisy data from theoretical data
y_n = y + np.random.normal(0, 0.27, len(x))

# create convolution kernel for calculating
# the smoothed second order derivative
smooth_width = 59
x1 = np.linspace(-3,3,smooth_width)
norm = np.sum(np.exp(-x1**2)) * (x1[1]-x1[0]) # ad hoc normalization
y1 = (4*x1**2 - 2) * np.exp(-x1**2) / smooth_width *8#norm*(x1[1]-x1[0])



# calculate second order deriv.
y_conv = np.convolve(y_n, y1, mode="same")

# plot data
plt.plot(x,y_conv, label = "second deriv")
plt.plot(x, y_n,"o", label = "noisy data")
plt.plot(x, y, label="theory")
plt.plot(x, x, "0.3", label = "linear data")
plt.hlines([0],-10, 20)
plt.axvspan(0,4, color="y", alpha=0.2)
plt.axvspan(6,14, color="y", alpha=0.2)
plt.axhspan(-1,1, color="b", alpha=0.2)
plt.vlines([0, 4, 6],-10, 10)
plt.xlim(-2.5,12)
plt.ylim(-2.5,6)
plt.legend(loc=0)
plt.show()

This is the result: enter image description here

smooth_width is the width of the convolution kernel. In order to adjust the amount of noise, change the value 0.27 in random.normal to different values. And please note, that this method does not work well close to the border of the data space.

As you can see, the "close-to-zero" requirement for the second derivative (blue line) holds quite well for the yellow parts, where the data is linear.

于 2012-12-05T16:40:30.743 に答える
5

Ramer Douglas Peucker アルゴリズムを使用して、データをより小さな線分のセットに単純化できます。このアルゴリズムを使用すると、すべてのデータ ポイントが特定の線分からepsilon離れないように指定できます。epsilon線分の傾きから、曲線の傾きを概算できます。

ここには、RDP アルゴリズムの Python 実装があります。

于 2012-12-03T21:42:13.097 に答える
1

これは可能な解決策にすぎません。プリセットされた最小値よりも長い最小カイ ^ 2 値を持つポイントの直線セグメントが見つかります。

from matplotlib.pyplot import figure, show
from numpy import pi, sin, linspace, exp, polyfit
from matplotlib.mlab import stineman_interp

x = linspace(0,2*pi,20);
y = x + sin(x) + exp(-0.5*(x-2)**2);

num_points = len(x)

min_fit_length = 5

chi = 0

chi_min = 10000

i_best = 0
j_best = 0

for i in range(len(x) - min_fit_length):
    for j in range(i+min_fit_length, len(x)):

        coefs = polyfit(x[i:j],y[i:j],1)
        y_linear = x * coefs[0] + coefs[1]
        chi = 0
        for k in range(i,j):
            chi += ( y_linear[k] - y[k])**2

        if chi < chi_min:
            i_best = i
            j_best = j
            chi_min = chi
            print chi_min

coefs = polyfit(x[i_best:j_best],y[i_best:j_best],1)
y_linear = x[i_best:j_best] * coefs[0] + coefs[1]


fig = figure()
ax = fig.add_subplot(111)
ax.plot(x,y,'ro')
ax.plot(x[i_best:j_best],y_linear,'b-')


show()

ここに画像の説明を入力

ただし、より大きなデータセットでは問題が発生することがわかります...

于 2012-12-03T22:00:16.800 に答える
0

データの「モデル」がほとんど直線に適合するデータで構成されており、最後にいくつかの外れ値または波状のビットがある場合は、RANSACアルゴリズムを試すことができます。

(非常に冗長で申し訳ありませんが)疑似コードは次のようになります。

choose a small threshold distance D

for N iterations:
    pick two random points from your data, a and b
    fit a straight line, L, to a and b
    count the inliers: data points within a distance D of the line L
    save the parameters of the line with the most inliers so far

estimate the final line using ALL the inliers of the best line
于 2012-12-03T21:30:33.803 に答える