2017-01-30 58 views
3

請看下面的例子:查找花鍵的平滑使用SciPy的

import numpy as np 
import math 
import matplotlib.pyplot as plt 
from scipy import interpolate 
xs = np.linspace(1,10,500) 
ys = [0.92 * x ** 2.3 + 0.0132 * x ** 4 + 0.0743 * (x - 9) ** 3 - 4 * (x -3) ** 2 + 80 * math.sin(math.sin(x)) + 10 * math.sin(x*5) + 1.2* np.random.normal(-4,4,1) for x in xs] 
ys[200] = ys[200] + 130 
ys[201] = ys[201] + 135 
ys[202] = ys[202] + 129 
ys[203] = ys[203] + 128 
ys[204] = ys[204] + 131 
ys[205] = ys[205] + 130 
ys[206] = ys[206] + 129 
ys[207] = ys[207] + 129 
ys[208] = ys[208] + 128 
ys[209] = ys[209] + 130 

如果我在這一點上繪製xsys,它產生了很好的圖形: a oisy dataset for testing

現在我使用scipy.interpolate.splrep爲這個數據擬合樣條曲線。我使用了兩種不同的樣條曲線,以適應數據的兩個不同的部分:

tck = interpolate.splrep(xs[0:199], ys[0:199], s = 1000) 
ynew2 = interpolate.splev(xs[0:199], tck, der = 0) 

和:

tck = interpolate.splrep(xs[210:500], ys[210:500], s = 9000) 
ynew3 = interpolate.splev(xs[210:500], tck, der = 0) 

然後我們有: Sample spline fit of the same data as above

現在我想以編程方式檢測的質量適合。貼合不應該太直 - 即保留特徵,也不應該將嘈雜的變化過度檢測爲特徵。

我打算使用送入人工神經網絡峯值計數器。

然而,在這一點上,我的問題是:

  • 是否SciPy的/ numpy的有一個內置的功能,在那裏我可以在splrep輸出飼料,它會計算的最低或最高和在任何特定間隔的最大/最小密度?

注:
我知道R**2價值,我希望找到另一項措施來檢測的功能保存。

回答

0

SciPy的不具有查找三次樣條的臨界點的方法。最近我們有sproot它找到一個三次樣條的根。爲了在這裏有用,我們必須擬合階4的樣條,使得導數是三次樣條。這是我做下面

from scipy.interpolate import splrep, splev, splder, sproot 

tck1 = splrep(xs[0:199], ys[0:199], k=4, s=1000) 
tck2 = splrep(xs[210:500], ys[210:500], k=4, s=9000) 
roots1 = sproot(splder(tck1), 1000)  # 1000 is an upper bound for the number of roots 
roots2 = sproot(splder(tck2), 1000) 

x1 = np.linspace(xs[0], xs[198], 1000)  # plot both splines 
plt.plot(x1, splev(x1, tck1)) 
x2 = np.linspace(xs[210], xs[499], 1000) 
plt.plot(x2, splev(x2, tck2))    

plt.plot(roots1, splev(roots1, tck1), 'ro')  # plot their max/min points 
plt.plot(roots2, splev(roots2, tck2), 'ro') 
plt.show() 

critical points

差別是顯而易見的。

您還可以找到在任何特定的時間間隔根數,如[3,4]:

np.where((3 <= roots1) & (roots1 <= 4))[0].size # 29 

或等價地,np.sum((3 <= roots1) & (roots1 <= 4))