python中的scipy库_SciPy库学习
Python庫之SciPy
工具:Pycharm2019.01,Python3.5
關(guān)于ScipyScipy包中有許多工具用來解決科學(xué)計算的一些共同問題。不同的子模塊負(fù)責(zé)不同的應(yīng)用。比如插值、擬合、最優(yōu)化、圖像處理、數(shù)值計算以及特殊函數(shù)等等。
Scipy中許多模塊都是基于numpy的,scipy命名空間中許多主要的函數(shù)實際上就是numpy函數(shù)比如scipy.cos就是np.cos,故通常一起導(dǎo)入numpy庫和scipy庫
File input/output: scipy.io
import numpy as np
from scipy import io as siosio.savemat('file_path',{'key',value}) : 將一個用命名和數(shù)組的字典保存到matlab樣式的.mat文件中。
sio.loadmat('file_path'):導(dǎo)入.mat格式的文件1
2
3
4
5
6
7
8
9
10
11import numpy as np
import matplotlib.pyplot as plt
from scipy import io as sio
a = np.ones((3, 3))
sio.savemat('file.mat',{'a':a}) # 將命名為‘a(chǎn)’的字典保存為.mat格式的文件,默認(rèn)保存到當(dāng)前目錄下,(.mat)后綴可以不加
data = sio.loadmat('file.mat')
print(data['a'])
- - - - - - - - -
[[1. 1. 1.]
[1. 1. 1.]
[1. 1. 1.]]
Linear algebra operations: scipy.linalg
線形代數(shù)運(yùn)算
import numpy as np
from scipy import linalg1
2
3
4
5
6
7
8
9a=np.array([[1,2],[3,4]])
print(linalg.det(a)) ?# 求方陣的行列式,要求必須是方陣
print(linalg.inv(a)) # 求逆,必須可逆,否則報錯提示奇異(singular不可逆)矩陣
print(np.allclose(np.dot(a,linalg.inv(a)),np.eye(2))) # 判斷是否為單位矩陣
- - - - - - - -
-2.0
[[-2. 1. ]
[ 1.5 -0.5]]
True
奇異值分解:
關(guān)于奇異值分解:SVD簡介1
2
3
4
5import numpy as np
from scipy import linalg
arr = np.arange(9).reshape((3, 3)) + np.diag([1, 0, 1])
print(arr)
uarr, spec, vharr = linalg.svd(arr) ?# uarr應(yīng)該是偽逆
Optimization and fit: scipy.optimize
scipy.optimize提供一種算法來最小化函數(shù)值(規(guī)模或多維度),曲線擬合和根查找1
2
3import numpy as np
import matplotlib.pyplot as plt
import scipy.optimize as opt
Curve fitting
曲線擬合,這里以sin()曲線為例1
2
3
4
5
6
7
8
9
10def test_func(x, a, b): ?# 定義函數(shù)
return a * np.sin(b * x)
x_data = np.linspace(-5, 5, num=50)
y_data = 2.9 * np.sin(1.5 * x_data) + np.random.normal(size=50) # sin函數(shù)存在偏差
params, params_covariance =opt.curve_fit(test_func, x_data, y_data, p0=[2, 2]) # 求上面函數(shù)參數(shù)a,b
plt.scatter(x_data,y_data,label='Data')
plt.legend(loc='upper right')
plt.plot(x_data,test_func(x_data,params[0],params[1]),label='fitted function') # 擬合
plt.legend(loc='upper right')
plt.show() # 顯示效果如下
Finding the minimum of a scalar function
尋找標(biāo)量函數(shù)最小值。
例如:
下面函數(shù)全局最優(yōu)(小)值為-1.3,而局部最優(yōu)值大概3.8。
用opt.minimize()函數(shù),需要定義函數(shù)f(),以及傳入起點坐標(biāo)x0,返回找到的最小值信息(包括函數(shù)值,坐標(biāo)等)。
為什么說是局部最優(yōu)呢,當(dāng)x0=5時,返回的函數(shù)值是8.3,位于x=3.83處。1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18def f(x):
return x**2 + 10*np.sin(x)
x = np.arange(-10, 10, 0.1)
plt.plot(x, f(x))
plt.show()
result=opt.minimize(f,x0=0) # 尋找局部最小值的x坐標(biāo)
print(result) # 打印返回結(jié)果信息,
- - - - - - - - - - -
fun: -7.945823375615215 ?# 局部最小值
hess_inv: array([[0.08589237]])
jac: array([-1.1920929e-06])
message: 'Optimization terminated successfully.'
nfev: 18
nit: 5
njev: 6
status: 0
success: True
x: array([-1.30644012]) ?# 局部最小值的坐標(biāo)
如果我們不知道全局最優(yōu)點附近信息也就無法合理選擇起點X0了,而opt.basinhopping()函數(shù)可以找到全局最優(yōu)解。將上面opt.minimize(f,x0=0)替換成opt.basinhopping(f,x0)即可,x0還是需要賦初值,經(jīng)測試當(dāng)x0=7時找到的是全局最優(yōu)(-1.3,-7.9),而當(dāng)x0=5時找到的是局部最優(yōu)值(3.8,8.3),故為了找到全局最優(yōu)保險的方法還是需要多試幾個參數(shù)。
試試opt.minimize_scalar()方法:能找到全局最優(yōu)解
opt.minimize_scalar()是只有一個變量時的最小函數(shù)值,而opt.minimize()是變量作為一個參數(shù)向量傳進(jìn)去的。1
2
3
4
5
6
7def f(x):
return x**2 + 10*np.sin(x)
x = np.arange(-10, 10, 0.1)
result=opt.minimize_scalar(f)
print(result.x, result.fun)
- - - - - - - - - - - - - - - - ?- - - - -
-1.3064400120612139 -7.945823375615284
Finding the roots of a scalar function
需求:想找到x=?時,f(x)=0。
方法:opt.root(f,x0),x0需要賦初值1
2
3
4
5
6
7
8
9
10def f(x):
return x**2 + 10*np.sin(x)
x = np.arange(-10, 10, 0.1)
root=opt.root(f,x0=-3) # 設(shè)置起點,只能找到一個,要找到另外一個需要調(diào)參使得x0=1
print(root.x,root.fun)
root=opt.root(f,x0=1)
print(root.x,root.fun)
- - - - - - - - - - - - - - - - -
[-2.47948183] [-1.77635684e-15]
[0.] [0.]
將以上三個函數(shù)繪制到一張圖中:1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43def f(x): # 定義函數(shù)
return x**2 + 10*np.sin(x)
x = np.arange(-10, 10, 0.1)
# Global optimization
grid = (-10, 10, 0.1)
xmin_global = opt.brute(f, (grid, )) # 暴力找出參數(shù)列表中所得解的最優(yōu)值
print("Global minima found %s" % xmin_global)
# Constrain optimization
xmin_local = opt.fminbound(f, 0, 10) # 在取區(qū)間(0,10)找局部最優(yōu)解
print("Local minimum found %s" % xmin_local)
root = opt.root(f, 1) # our initial guess is 1
print("First root found %s" % root.x)
root2 = opt.root(f, -2.5)
print("Second root found %s" % root2.x)
fig = plt.figure(figsize=(6, 4)) # 新建畫布對象指定size
ax = fig.add_subplot(111) # 將畫布分成一行一列指定第一塊區(qū)域作為對象賦給ax
# Plot the function
ax.plot(x, f(x), 'b-', label="f(x)")
# Plot the minima
xmins = np.array([xmin_global[0], xmin_local])
ax.plot(xmins, f(xmins), 'go', label="Minima")
# Plot the roots
roots = np.array([root.x, root2.x])
ax.plot(roots, f(roots), 'kv', label="Roots")
# Decorate the figure
ax.legend(loc='best') # 自適應(yīng)選擇位置
ax.set_xlabel('x')
ax.set_ylabel('f(x)')
ax.axhline(0, color='gray') # 設(shè)置y=0畫橫線
plt.show()
- - - - - - - - - - - - - - - - - -
Global minima found [-1.30641113]
Local minimum found 3.8374671194983834
First root found [0.]
Second root found [-2.47948183]
未完待續(xù)…..
總結(jié)
以上是生活随笔為你收集整理的python中的scipy库_SciPy库学习的全部內(nèi)容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: QIIME 2教程. 22命令行界面q2
- 下一篇: 药店零售管理php系统,医药POS零售管