1. 程式人生 > >機器學習筆記--python之scipy

機器學習筆記--python之scipy

學習了numpy和matplotlib,基本上線性代數,概率論的很多計算啊之類的都可以很容易的實現了。此外再學習下scipy這個科學函式庫吧。
scipy包包含致力於科學計算中常見問題的各個工具箱。它的不同子模組相應於不同的應用。像插值,積分,優化,影象處理,特殊函式等等。

1 模組

scipy 由一些特定功能的子模組組成,它們全依賴numpy,但是每個之間基本獨立

模組 功能
scipy.cluster 向量量化 / K-均值
scipy.constants 物理和數學常數
scipy.fftpack 傅立葉變換
scipy.integrate 積分程式
scipy.interpolate 插值
scipy.io 資料輸入輸出
scipy.linalg 線性代數程式
scipy.ndimage n維影象包
scipy.odr 正交距離迴歸
scipy.optimize 優化
scipy.signal 訊號處理
scipy.sparse 稀疏矩陣
scipy.spatial 空間資料結構和演算法
scipy.special 任何特殊數學函式
scipy.stats 統計

之後就學習下部分幾個模組吧。

2 插值 scipy.interpolate

scipy.interpolate對從實驗資料擬合函式來求值沒有測量值存在的點非常有用
正弦函式的實驗:

import numpy as np
from scipy.interpolate import interp1d
import pylab as plt 

measured_time = np.linspace(0, 1, 10) 
noise = (np.random.random(10)*2 - 1) * 1e-1
measures = np.sin(2 * np.pi * measured_time) + noise

#scipy.interpolate.interp1d類會構建線性插值函式:
linear_interp = interp1d(measured_time, measures) #然後scipy.interpolate.linear_interp例項需要被用來求得感興趣時間點的值: computed_time = np.linspace(0, 1, 50) linear_results = linear_interp(computed_time) #三次插值也能通過提供可選關鍵字引數kind來選擇 cubic_interp = interp1d(measured_time, measures, kind='cubic') cubic_results = cubic_interp(computed_time) #Matplotlib影象中顯示如下 plt.plot(measured_time, measures, 'o', ms=6, label='measures') plt.plot(computed_time, linear_results, label='linear interp') plt.plot(computed_time, cubic_results, label='cubic interp') plt.legend() plt.show()

執行結果如下:

3 訊號處理 scipy. signal

3.1 移除訊號的線性趨勢

import numpy as np
import pylab as plt 
from scipy import signal

t = np.linspace(0, 3, 200)
x = t + np.random.normal(size=200)

plt.plot(t, x, linewidth=3)
#移除訊號的線性趨勢
plt.plot(t, signal.detrend(x), linewidth=3)

plt.show()

執行結果:

3.2FFT重取樣

import numpy as np
import pylab as plt 
from scipy import signal

t = np.linspace(0, 5, 100)
x = np.sin(t)

plt.plot(t, x, linewidth=3)
#使用FFT重取樣n個點
plt.plot(t[::2], signal.resample(x, 50), 'ko')
plt.show()

執行結果:

4 統計 scipy.stats

scipy.stats包括統計工具和隨機過程的概率過程。各個隨機過程的隨機數生成器可以從numpy.random中找到。

模組 功能
rvs 隨機變數(就是從這個分佈中抽一些樣本)
pdf 概率密度函式
cdf 累計分佈函式
sf 殘存函式(1-CDF)
ppf 分位點函式(CDF的逆)
isf 逆殘存函式(sf的逆)
stats 返回均值,方差,(費舍爾)偏態,(費舍爾)峰度
moment 分佈的非中心矩

4.1 二項分佈

拋擲10次硬幣,假設在該試驗中正面朝上的概率為0.3。使用stats.binom.pmf計算每次觀測的概率質量函式。

import numpy as np
import matplotlib.pyplot as plt 
from scipy import stats

plt.subplot(121)
n = 10
p = 0.3 
k = np.arange(0, 30) 
binomial = stats.binom.pmf(k, n, p)
plt.plot(k, binomial, 'o-')
#使用rvs函式模擬一個二項隨機變數,其中引數size指定你要進行模擬的次數,這裡為10000次。
plt.subplot(122)
binom_sim = data = stats.binom.rvs(n=10, p=0.3, size=10000)
print "Mean: %g" % np.mean(binom_sim)
print "Sd: %g" % np.std(binom_sim, ddof=1)
plt.hist(binom_sim, bins=10, normed=True)

plt.show()

執行結果:
Mean: 2.9956
Sd: 1.44187

4.2 泊松分佈

import numpy as np
import matplotlib.pyplot as plt 
from scipy import stats

rate =2
n = np.arange(0, 10) 
y = stats.poisson.pmf(n, rate)

plt.subplot(121)
plt.plot(n, y, 'o-')
#模擬1000個服從泊松分佈的隨機變數
plt.subplot(122)
data = stats.poisson.rvs(mu=2, loc=0, size=1000)
print "Mean: %g" % np.mean(data)
print "Sd: %g" % np.std(data, ddof=1)
plt.hist(data, bins=9, normed=True)

plt.show()

執行結果:
Mean: 2.105
Sd: 1.4677

4.3 正態分佈

import numpy as np
from scipy import stats
import matplotlib.pyplot as plt 

plt.subplot(121)
mu = 0 
sigma = 1 
x = np.arange(-5, 5, 0.1)
y = stats.norm.pdf(x, 0, 1)
plt.plot(x, y)

plt.subplot(122)

data = stats.norm.rvs(0, 1, size=1000)
plt.hist(data, bins=10, normed=True)

plt.show()

4.4 指數分佈

import numpy as np
import matplotlib.pyplot as plt 
from scipy import stats

plt.subplot(121)
lambd = 0.5 
x = np.arange(0, 15, 0.1)
y = lambd*np.exp(-lambd*x)
plt.plot(x, y)

#下模擬1000個隨機變數
plt.subplot(122)
data = stats.expon.rvs(scale=2, size=1000)
print "Mean: %g" % np.mean(data)
print "Sd: %g" % np.std(data, ddof=1)
plt.hist(data, bins=20, normed=True)

plt.show()

執行結果如下:
Mean: 2.01166
Sd: 2.02959

5 優化和擬合 scipy.optimize

優化是找到最小值或等式的數值解的問題

import numpy as np
import matplotlib.pyplot as plt
from scipy import optimize
#定義以下函式x^2+2sin(x)
def f(x):
    return x**2+10*np.sin(x)

x = np.arange(-10, 10, 0.1)

#找到這個函式最小值一般而有效的方法是從初始點使用梯度下降法,這裡用擬牛頓法
print optimize.fmin_bfgs(f, 0)

#如果函式有區域性最小值,演算法會因初始點不同找到這些區域性最小而不是全域性最小
print optimize.fmin_bfgs(f, 3, disp=0)

#為了找到全域性最小點,最簡單的演算法是蠻力演算法
grid = (-10, 10, 0.1)
xmin_global = optimize.brute(f, (grid,))
print xmin_global

#為了找到區域性最小,我們把變數限制在(0, 10)之間
xmin_local =optimize.fminbound(f, 0, 10) 
print xmin_local

#找到標量函式的根
root1 = optimize.fsolve(f, 1)
print root1

root2 = optimize.fsolve(f, -2) 
print root2

#曲線擬合
xdata = np.linspace(-10, 10, num=20)
ydata = f(xdata) + np.random.random(xdata.size)

def f2(x, a, b): 
    return a*x**2+b*np.sin(x)

guess = [2, 2]
#通過最小二乘擬合擬合來找到幅度
params, params_convariance = optimize.curve_fit(f2, xdata, ydata, guess)
print params

#畫圖顯示所有的資訊
plt.plot(x, f(x), c='r', label="f(x)")
plt.plot(xmin_global, f(xmin_global), '^', label='global minima')
plt.plot(xmin_local, f(xmin_local), '^', label='local minima')
plt.plot(root1, f(root1), 'o', label='root1')
plt.plot(root2, f(root2), 'o', label='root2')
plt.plot(x, f2(x, params[0], params[1]), '--', c='b', label="curve fit result")
plt.legend()
plt.show()

執行結果:
Optimization terminated successfully.
Current function value: -7.945823
Iterations: 5
Function evaluations: 24
Gradient evaluations: 8
[-1.30644003]
[ 3.83746663]
[-1.30641113]
3.8374671195
[ 0.]
[-2.47948183]
[ 1.0089413 9.92366871]


基本上把python用到的數學庫簡單的過了一遍,相信之後遇到一些函式庫也不會那麼陌生了。接下去就開始真正的學習機器學習演算法了。