第二章：Scipy 核心模块详解

第二章：Scipy 核心模块详解

第二章：Scipy 核心模块详解

2.1 scipy.integrate：数值积分

scipy.integrate 模块提供了多种数值积分方法，用于计算函数的定积分。

2.1.1 quad 函数：单重积分

quad 函数是最常用的积分函数，用于计算单变量函数的定积分。

import scipy.integrate as integrate
import numpy as np
# 定义被积函数
def f(x):
  return x**2
# 计算积分
result, error = integrate.quad(f, 0, 3)  # 对 x^2 在 [0, 3] 区间积分
print("积分结果:", result)
print("估计误差:", error)
​

代码解释：

• integrate.quad(f, a, b)：计算函数 f 在区间 [a, b] 上的定积分。

• 返回值：返回积分结果和估计误差。

2.1.2 dblquad 函数：二重积分

dblquad 函数用于计算二重积分。

import scipy.integrate as integrate
import numpy as np
# 定义被积函数
def f(x, y):
  return x*y
# 定义 y 的积分上下限（作为 x 的函数）
def y_lower(x):
  return x**2
def y_upper(x):
  return 2*x
# 计算积分
result, error = integrate.dblquad(f, 0, 2, y_lower, y_upper)  # 对 x*y 在特定区域积分
print("积分结果:", result)
print("估计误差:", error)
​

代码解释：

• integrate.dblquad(f, a, b, gfun, hfun)：计算函数 f(x, y) 在区域上的二重积分，其中 x 的积分范围是 [a, b]，y 的积分范围是 [gfun(x), hfun(x)]。

• gfun 和 hfun 是定义 y 积分上下限的函数。

2.1.3 nquad 函数：多重积分

nquad 函数提供了一种通用的方式来计算多重积分。

import scipy.integrate as integrate
import numpy as np
# 定义被积函数
def f(x, y, z):
  return x*y*z
# 定义积分上下限
ranges = [[1, 2], [0, 1], [-1, 1]]
# 计算积分
result, error = integrate.nquad(f, ranges) # 计算三重积分
print("积分结果:", result)
print("估计误差:", error)
​

代码解释：

• integrate.nquad(f, ranges)：计算函数 f 在指定范围内的多重积分。

• ranges 是一个列表，其中每个元素都是一个包含积分上下限的列表。

2.1.4 积分方法选择

scipy.integrate 还提供了其他积分方法，如 fixed_quad, romberg, trapz 等。 选择哪种方法取决于被积函数的性质和所需的精度。

2.2 scipy.optimize：优化

scipy.optimize 模块提供了各种优化算法，用于寻找函数的最小值、最大值或方程的根。

2.2.1 minimize 函数：无约束优化

minimize 函数用于寻找无约束函数的最小值。

import scipy.optimize as optimize
import numpy as np
# 定义目标函数
def f(x):
  return x[0]**2 + (x[1]-1)**2
# 初始猜测值
x0 = np.array([2, 0])
# 优化
result = optimize.minimize(f, x0)
print("最小值:", result.fun)
print("最优解:", result.x)
print("优化是否成功:", result.success)
print("迭代次数:", result.nit)
​

代码解释：

• optimize.minimize(f, x0)：寻找函数 f 的最小值，x0 是初始猜测值。

• 返回值：返回一个包含优化结果的对象，包括最小值、最优解、优化是否成功等信息。

2.2.2 minimize 函数：约束优化

minimize 函数也可以用于约束优化问题。

import scipy.optimize as optimize
import numpy as np
# 定义目标函数
def f(x):
  return x[0]**2 + (x[1]-1)**2
# 定义约束条件
constraints = ({'type': 'eq', 'fun': lambda x: x[0] + x[1] - 1}) # x[0] + x[1] = 1
# 初始猜测值
x0 = np.array([2, 0])
# 优化
result = optimize.minimize(f, x0, constraints=constraints)
print("最小值:", result.fun)
print("最优解:", result.x)
print("优化是否成功:", result.success)
print("迭代次数:", result.nit)
​

代码解释：

• constraints 是一个列表，其中每个元素都是一个字典，用于定义一个约束条件。

• 'type' 指定约束类型，可以是 'eq' (等于), 'ineq' (大于等于)。

• 'fun' 是一个函数，用于计算约束条件的值。

2.2.3 root 函数：求解方程的根

root 函数用于求解方程的根。

import scipy.optimize as optimize
import numpy as np
# 定义方程
def f(x):
  return x**2 - 2
# 初始猜测值
x0 = 1
# 求解方程的根
result = optimize.root(f, x0)
print("根:", result.x)
print("求解是否成功:", result.success)
​

代码解释：

• optimize.root(f, x0)：寻找函数 f 的根，x0 是初始猜测值。

• 返回值：返回一个包含求解结果的对象，包括根、求解是否成功等信息。

2.2.4 优化方法选择

scipy.optimize 提供了多种优化算法，如 BFGS, CG, Newton-CG, L-BFGS-B, SLSQP 等。 选择哪种方法取决于目标函数的性质、约束条件和所需的精度。

2.3 scipy.interpolate：插值

scipy.interpolate 模块提供了各种插值方法，用于根据已知数据点估计未知数据点的值。

2.3.1 interp1d 函数：一维插值

interp1d 函数用于创建一维插值函数。

import scipy.interpolate as interpolate
import numpy as np
import matplotlib.pyplot as plt
# 已知数据点
x = np.array([0, 1, 2, 3, 4])
y = np.array([0, 1, 4, 9, 16])
# 创建插值函数
f = interpolate.interp1d(x, y, kind='quadratic') # 线性，二次，三次插值
# 估计未知数据点的值
x_new = np.linspace(0, 4, 100)
y_new = f(x_new)
# 绘图
plt.plot(x, y, 'o', label='原始数据')
plt.plot(x_new, y_new, '-', label='插值结果')
plt.legend()
plt.show()
​

代码解释：

• interpolate.interp1d(x, y, kind='linear')：创建一个线性插值函数，x 和 y 是已知数据点的坐标。

• kind 参数指定插值类型，可以是 'linear', 'quadratic', 'cubic' 等。

2.3.2 griddata 函数：多维插值

griddata 函数用于创建多维插值函数。

import scipy.interpolate as interpolate
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
# 已知数据点
points = np.random.rand(10, 2)
values = np.random.rand(10)
# 创建网格
grid_x, grid_y = np.mgrid[0:1:100j, 0:1:100j]
# 插值
grid_z = interpolate.griddata(points, values, (grid_x, grid_y), method='cubic')
# 绘图
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_surface(grid_x, grid_y, grid_z, cmap='viridis')
plt.show()
​

代码解释：

• interpolate.griddata(points, values, (grid_x, grid_y), method='linear')：创建一个多维插值函数，points 是已知数据点的坐标，values 是已知数据点的值，(grid_x, grid_y) 是要估计的网格坐标。

• method 参数指定插值方法，可以是 'linear', 'cubic', 'nearest' 等。

2.3.3 插值方法选择

选择哪种插值方法取决于数据的性质和所需的精度。 线性插值简单快速，但精度较低； 立方插值精度较高，但计算量较大。

2.4 scipy.fft：快速傅里叶变换

scipy.fft 模块提供了快速傅里叶变换 (FFT) 的功能，用于将信号从时域转换到频域。

import scipy.fft as fft
import numpy as np
import matplotlib.pyplot as plt
# 生成信号
N = 600
T = 1.0 / 800.0
x = np.linspace(0.0, N*T, N, endpoint=False)
y = np.sin(50.0 * 2.0*np.pi*x) + 0.5*np.sin(80.0 * 2.0*np.pi*x)
# 计算 FFT
yf = fft.fft(y)
xf = fft.fftfreq(N, T)[:N//2]
# 绘图
plt.plot(xf, 2.0/N * np.abs(yf[0:N//2]))
plt.grid()
plt.show()
​

代码解释：

• fft.fft(y)：计算信号 y 的 FFT。

• fft.fftfreq(N, T)：计算 FFT 的频率。

• N 是信号的长度，T 是采样间隔。

2.5 scipy.signal：信号处理

scipy.signal 模块提供了各种信号处理功能，如滤波、窗函数、谱分析等。

import scipy.signal as signal
import numpy as np
import matplotlib.pyplot as plt
# 生成信号
t = np.linspace(0, 1, 1000, endpoint=False)
sig  = np.sin(2*np.pi*10*t) + np.sin(2*np.pi*20*t)
# 设计滤波器
sos = signal.butter(10, 15, 'low', fs=1000, output='sos')
# 应用滤波器
filtered = signal.sosfilt(sos, sig)
# 绘图
plt.plot(t, sig, label="原始信号")
plt.plot(t, filtered, label="滤波后信号")
plt.legend()
plt.show()
​

代码解释：

• signal.butter(N, Wn, btype, fs, output)：设计一个 Butterworth 滤波器。

• N 是滤波器的阶数，Wn 是截止频率，btype 是滤波器类型 ('low', 'high', 'bandpass', 'bandstop')，fs 是采样频率，output 指定输出类型 ('sos' 是二阶节)。

• signal.sosfilt(sos, sig)：应用滤波器到信号 sig。

2.6 scipy.stats：统计

scipy.stats 模块提供了各种统计函数，如概率分布、统计检验、描述性统计等。

import scipy.stats as stats
import numpy as np
import matplotlib.pyplot as plt
# 生成随机数
data = np.random.normal(size=1000)
# 计算统计量
mean = np.mean(data)
std = np.std(data)
# 拟合正态分布
params = stats.norm.fit(data)
print("拟合参数 (mean, std):", params)
# 绘制直方图和拟合曲线
plt.hist(data, bins=30, density=True, alpha=0.6, label="直方图")
x = np.linspace(-4, 4, 100)
p = stats.norm.pdf(x, *params)
plt.plot(x, p, 'k', linewidth=2, label="拟合曲线")
plt.legend()
plt.show()
​

代码解释：

• stats.norm.fit(data)：将正态分布拟合到数据 data。

• stats.norm.pdf(x, *params)：计算正态分布的概率密度函数。

2.7 总结

Scipy 提供了许多强大的核心模块，可以用于解决各种科学计算问题。 本章介绍了 scipy.integrate, scipy.optimize, scipy.interpolate, scipy.fft, scipy.signal 和 scipy.stats 等核心模块，并通过代码示例和图表来阐述其功能和用法。 掌握这些模块将使您能够更有效地使用 Scipy 进行科学计算。