NumPy

技术栈
工具链
python科学计算数组线性代数矩阵数值计算

概览

NumPy 技术栈概览

NumPy(Numerical Python)是 Python 科学计算的核心基础库,由 Travis Oliphant 于 2006 年创建。它提供高性能的多维数组对象(ndarray)和丰富的数学运算工具,是 Pandas、SciPy、Matplotlib、scikit-learn、PyTorch 等几乎整个 Python 数据科学生态的基石。

核心特性:

  • 🧱 ndarray — 同构多维数组,是 NumPy 的核心数据结构,内存连续、向量化运算
  • 向量化运算 — 无需 Python 循环即可对整个数组执行批量运算,速度提升 10-100x
  • 🧮 广播机制 — 不同形状的数组之间自动对齐,无需显式复制
  • 📐 线性代数 — 完整的矩阵运算、特征值分解、SVD 等
  • 🎲 随机数生成 — 丰富的概率分布采样
  • 📊 统计函数 — 均值、方差、协方差、直方图等
  • 🔌 C 扩展接口 — 通过 Cython/Numba 无缝衔接 C/C++/Fortran

适用场景: 数值计算、数据预处理、图像处理、信号处理、科学模拟、机器学习数据管道。

安装

1. 环境准备

  • 操作系统: Windows 10+ / macOS 11+ / Linux
  • Python 版本: Python 3.9 - 3.12(推荐 3.10+)
  • 依赖项: pip 或 conda

2. 安装命令

# === pip 安装(推荐)===
pip install numpy

# === conda 安装(自动包含 MKL 加速)===
conda install numpy

# === 验证安装 ===
python -c "import numpy as np; print(f'NumPy {np.__version__}'); print(f'多维数组: {np.array([1,2,3])}')"

# === 安装 NumPy + 常用科学计算栈 ===
pip install numpy scipy matplotlib pandas

# === 查看 NumPy 的 BLAS 后端(性能关键)===
python -c "import numpy; numpy.show_config()"

3. 常见安装问题

问题 1:pip 安装速度慢

pip install numpy -i https://pypi.tuna.tsinghua.edu.cn/simple

问题 2:版本冲突

# 安装特定版本
pip install numpy==1.26.4

# 升级
pip install --upgrade numpy

问题 3:使用 MKL 加速(Intel CPU)

# conda 默认使用 MKL,性能比 pip 的 OpenBLAS 略好
conda install numpy

# 或通过 pip 安装带 MKL 的 wheel(来自 Intel)
pip install intel-numpy

问题 4:Apple Silicon (M1/M2/M3) 加速

NumPy 1.21+ 已原生支持 Apple Silicon 的 Accelerate 框架:

pip install numpy
# 自动使用 vecLib/Accelerate,无需额外配置

问题 5:服务器无网络安装

# 在有网络的环境下载 wheel
pip download numpy -d ./offline_packages
# 复制到服务器后离线安装
pip install numpy --no-index --find-links ./offline_packages

示例

NumPy 数组基础 —— 创建、索引、运算

目标

  • 掌握 ndarray 的创建、属性、索引与切片
  • 理解向量化运算与 Python 循环的性能差异
  • 掌握广播机制

完整代码

import numpy as np

# ============================================================
# 1. 数组创建
# ============================================================

# 从列表创建
a = np.array([1, 2, 3, 4, 5])
print(f"1D 数组: {a}")

b = np.array([[1, 2, 3], [4, 5, 6]])
print(f"2D 数组:\n{b}")

# 特殊数组
zeros = np.zeros((3, 4))           # 全零
ones = np.ones((2, 3))             # 全一
eye = np.eye(3)                     # 单位矩阵
full = np.full((2, 2), 7)          # 填充指定值
print(f"单位矩阵:\n{eye}")

# 序列数组
arange = np.arange(0, 10, 2)      # 步长为2
linspace = np.linspace(0, 1, 5)   # 等分5份
print(f"arange: {arange}")
print(f"linspace: {linspace}")

# 随机数组
rand = np.random.rand(3, 3)        # [0,1) 均匀分布
randn = np.random.randn(3, 3)      # 标准正态分布
randint = np.random.randint(1, 100, size=(3, 3))  # 随机整数

# 设置随机种子(可复现)
rng = np.random.default_rng(seed=42)
print(f"可复现随机数: {rng.random(5)}")

# ============================================================
# 2. 数组属性
# ============================================================
arr = np.random.randn(4, 3, 28, 28)

print(f"\n形状:   {arr.shape}")      # (4, 3, 28, 28)
print(f"维度数: {arr.ndim}")         # 4
print(f"元素总数: {arr.size}")       # 4*3*28*28 = 9408
print(f"数据类型: {arr.dtype}")     # float64
print(f"元素大小: {arr.itemsize} bytes")
print(f"总内存:  {arr.nbytes} bytes ({arr.nbytes / 1024 / 1024:.2f} MB)")

# 数据类型转换
arr_f32 = arr.astype(np.float32)
print(f"float32 内存减半: {arr_f32.nbytes} bytes")

# ============================================================
# 3. 索引与切片
# ============================================================
arr = np.arange(1, 13).reshape(3, 4)
print(f"\n原始数组:\n{arr}")

# 基本切片
print(f"第2行: {arr[1]}")
print(f"第1列: {arr[:, 0]}")
print(f"子矩阵:\n{arr[0:2, 1:3]}")

# 布尔索引
mask = arr > 6
print(f"大于6的元素: {arr[mask]}")

# 花式索引
print(f"取(0,2)行和(1,3)列:\n{arr[[0, 2], :][:, [1, 3]]}")

# np.ix_ —— 笛卡尔积索引
print(f"np.ix_ 笛卡尔积:\n{arr[np.ix_([0, 2], [1, 3])]}")

# ============================================================
# 4. 形状变换
# ============================================================
a = np.arange(12)

print(f"\nreshape 3×4:\n{a.reshape(3, 4)}")
print(f"reshape 2×2×3:\n{a.reshape(2, 2, 3)}")

# -1 表示自动推断
print(f"自动推断行数:\n{a.reshape(-1, 4)}")

# 转置
m = np.arange(6).reshape(2, 3)
print(f"转置:\n{m.T}")

# 添加/删除维度
x = np.array([1, 2, 3])
print(f"np.newaxis 增加列维度: {x[:, np.newaxis].shape}")  # (3, 1)
print(f"expand_dims: {np.expand_dims(x, axis=0).shape}")     # (1, 3)
print(f"squeeze: {np.ones((1,3,1,5)).squeeze().shape}")     # (3, 5)

# 拼接
a = np.array([[1, 2], [3, 4]])
b = np.array([[5, 6]])
print(f"\n垂直拼接:\n{np.vstack([a, b])}")
print(f"水平拼接:\n{np.hstack([a, b.T])}")

# ============================================================
# 5. 向量化运算(性能对比)
# ============================================================
import time

N = 10_000_000
arr = np.random.randn(N)

# 向量化
start = time.time()
result_np = arr * 2 + 1
print(f"\nNumPy 向量化: {time.time() - start:.4f} 秒")

# Python 循环(不要这样做!)
start = time.time()
result_py = [x * 2 + 1 for x in arr]
print(f"Python 列表推导: {time.time() - start:.4f} 秒")

# 逐元素数学运算
a = np.array([1, 2, 3, 4])
print(f"\n加法: {a + 2}")
print(f"乘法: {a * 3}")
print(f"平方: {a ** 2}")
print(f"指数: {np.exp(a)}")
print(f"对数: {np.log(a)}")

# ============================================================
# 6. 广播机制
# ============================================================

# 标量广播
arr = np.arange(6).reshape(2, 3)
print(f"\n标量广播:\n{arr + 10}")

# 向量广播
row = np.array([10, 20, 30])  # shape (3,)
print(f"向量广播(按行):\n{arr + row}")

col = np.array([[10], [20]])   # shape (2, 1)
print(f"向量广播(按列):\n{arr + col}")

# 广播规则:从后往前对齐维度,不匹配的维度必须为 1
# (2, 3) + (1, 3) → (2, 3)  ✅
# (2, 3) + (2, 1) → (2, 3)  ✅
# (2, 3) + (2,  ) → (2, 3)  ✅(自动在前面补 1)

# ============================================================
# 7. 聚合操作
# ============================================================
arr = np.random.randn(4, 5)

print(f"\n全局求和: {arr.sum():.4f}")
print(f"全局均值: {arr.mean():.4f}")
print(f"全局最大值: {arr.max():.4f}")
print(f"全局最小值: {arr.min():.4f}")
print(f"标准差: {arr.std():.4f}")

# 按轴聚合
print(f"\n按列求和 (axis=0): {arr.sum(axis=0)}")  # 对每列求和
print(f"按行求和 (axis=1): {arr.sum(axis=1)}")     # 对每行求和
print(f"每列均值: {arr.mean(axis=0)}")
print(f"每行最大值: {arr.max(axis=1)}")

# argmax / argmin
print(f"每列最大值索引: {arr.argmax(axis=0)}")

运行输出示例

NumPy 向量化: 0.0123 秒
Python 列表推导: 3.4567 秒
# 向量化比列表推导快 ~280 倍!

关键要点

概念 说明
ndarray 同构多维数组,所有元素类型相同
shape 元组,描述各维度大小
dtype 数据类型,如 float64、int32
向量化 一次操作作用于整个数组,底层 C 实现
广播 不同形状数组的自动对齐运算
axis=0 沿行方向(对列操作),结果减少第0维
np.newaxis 增加一个大小为1的维度

NumPy 线性代数与数据应用

目标

  • 掌握矩阵运算、特征值分解、SVD
  • 求解线性方程组
  • 使用 NumPy 进行数据清洗与特征工程
  • 多项式拟合与信号处理入门

完整代码

1. 矩阵运算

import numpy as np

# ============================================================
# 矩阵乘法
# ============================================================
A = np.array([[1, 2], [3, 4]])
B = np.array([[5, 6], [7, 8]])

# 三种等价写法
C1 = A @ B              # Python 3.5+ 推荐
C2 = np.matmul(A, B)
C3 = np.dot(A, B)
print(f"矩阵乘法:\n{C1}")

# 逐元素乘法
print(f"逐元素乘法:\n{A * B}")

# 向量点积
a = np.array([1, 2, 3])
b = np.array([4, 5, 6])
print(f"点积: {np.dot(a, b)}")          # 32
print(f"内积: {np.inner(a, b)}")        # 32
print(f"外积:\n{np.outer(a, b)}")       # 3×3 矩阵

# ============================================================
# 矩阵属性
# ============================================================
print(f"\n行列式: {np.linalg.det(A):.2f}")
print(f"逆矩阵:\n{np.linalg.inv(A)}")
print(f"转置: {(A @ np.linalg.inv(A)).round(10)}")  # 应等于单位矩阵

# ============================================================
# 求解线性方程组 Ax = b
# ============================================================
# 3x + y = 9
# x + 2y = 8
A = np.array([[3, 1], [1, 2]])
b = np.array([9, 8])

x = np.linalg.solve(A, b)
print(f"\n方程组的解: x={x[0]:.1f}, y={x[1]:.1f}")  # x=2, y=3
print(f"验证 Ax-b = {A @ x - b}")                    # 应接近 [0,0]

2. 特征值与特征向量

# ============================================================
# 特征值分解
# ============================================================
A = np.array([[4, 2], [1, 3]])
eigenvalues, eigenvectors = np.linalg.eig(A)

print(f"特征值: {eigenvalues}")
print(f"特征向量:\n{eigenvectors}")

# 验证: A @ v = λ @ v
for i in range(len(eigenvalues)):
    v = eigenvectors[:, i]
    λ = eigenvalues[i]
    print(f"验证 A·v{i}: {A @ v}")
    print(f"验证 λ{i}·v: {λ * v}")
    print(f"误差: {np.linalg.norm(A @ v - λ * v):.2e}\n")

# ============================================================
# SVD 奇异值分解
# ============================================================
M = np.array([[1, 2, 3], [4, 5, 6]])
U, S, Vt = np.linalg.svd(M, full_matrices=False)

print(f"U (左奇异向量):\n{U}")
print(f"Σ (奇异值): {S}")
print(f"Vt (右奇异向量转置):\n{Vt}")

# 重构
M_reconstructed = U @ np.diag(S) @ Vt
print(f"重构矩阵:\n{M_reconstructed}")
print(f"重构误差: {np.linalg.norm(M - M_reconstructed):.2e}")

3. 数据清洗与统计

# ============================================================
# 处理缺失值
# ============================================================
data = np.array([1.0, 2.0, np.nan, 4.0, np.nan, 6.0])

# 检测 NaN
nan_mask = np.isnan(data)
print(f"\nNaN 位置: {nan_mask}")
print(f"NaN 数量: {np.sum(nan_mask)}")

# 填充 NaN(用均值)
clean_mean = np.where(nan_mask, np.nanmean(data), data)
print(f"均值填充: {clean_mean}")

# 填充 NaN(用前向填充)
from numpy.lib.stride_tricks import sliding_window_view
# 简化:用 np.nan_to_num
clean_zero = np.nan_to_num(data, nan=0)
print(f"零填充: {clean_zero}")

# ============================================================
# 描述性统计
# ============================================================
scores = np.array([78, 85, 92, 88, 76, 95, 89, 83, 91, 87])

print(f"\n均值: {np.mean(scores):.1f}")
print(f"中位数: {np.median(scores):.1f}")
print(f"标准差: {np.std(scores):.1f}")
print(f"方差: {np.var(scores):.1f}")

# 分位数
print(f"25%: {np.percentile(scores, 25):.1f}")
print(f"50%: {np.percentile(scores, 50):.1f}")
print(f"75%: {np.percentile(scores, 75):.1f}")

# 检测异常值(Z-score 方法)
z_scores = (scores - scores.mean()) / scores.std()
outliers = np.abs(z_scores) > 2
print(f"异常值 (|z| > 2): {scores[outliers]}")

# ============================================================
# 相关系数
# ============================================================
x = np.array([1, 2, 3, 4, 5])
y = np.array([2, 4, 5, 4, 6])

corr_matrix = np.corrcoef(x, y)
print(f"\n相关系数矩阵:\n{corr_matrix}")
print(f"皮尔逊相关系数: {corr_matrix[0, 1]:.4f}")

# ============================================================
# 直方图与分箱
# ============================================================
data = np.random.randn(1000)
hist, bin_edges = np.histogram(data, bins=20)
print(f"\n直方图计数: {hist}")
print(f"分箱边界: {bin_edges}")

# 数字分箱
ages = np.array([15, 22, 35, 42, 55, 68, 73])
bins = [0, 18, 35, 60, 100]
labels = ["未成年", "青年", "中年", "老年"]
categories = np.digitize(ages, bins) - 1
print(f"年龄分组索引: {categories}")
print(f"年龄分组: {np.array(labels)[categories]}")

4. 多项式拟合

# ============================================================
# 多项式回归
# ============================================================
x = np.array([0, 1, 2, 3, 4, 5])
y = np.array([1.1, 3.5, 8.2, 15.9, 26.2, 39.1])  # 大致 y = x² + 2x + 1 + noise

# 拟合二次多项式
coeffs = np.polyfit(x, y, deg=2)  # 返回 [a, b, c] 对应 ax² + bx + c
print(f"\n拟合系数: a={coeffs[0]:.4f}, b={coeffs[1]:.4f}, c={coeffs[2]:.4f}")

# 预测
x_new = np.linspace(0, 6, 50)
y_pred = np.polyval(coeffs, x_new)

# 计算 R²
y_fit = np.polyval(coeffs, x)
ss_res = np.sum((y - y_fit) ** 2)
ss_tot = np.sum((y - np.mean(y)) ** 2)
r2 = 1 - ss_res / ss_tot
print(f"R² = {r2:.6f}")

5. 信号处理(简易)

# ============================================================
# FFT 快速傅里叶变换
# ============================================================
sr = 1000  # 采样率 1000Hz
t = np.linspace(0, 1, sr, endpoint=False)
# 50Hz + 120Hz 的两个正弦波叠加
signal = np.sin(2 * np.pi * 50 * t) + 0.5 * np.sin(2 * np.pi * 120 * t)

# FFT
fft_vals = np.fft.fft(signal)
freqs = np.fft.fftfreq(len(signal), 1 / sr)

# 取正频率部分
positive_mask = freqs > 0
freqs_pos = freqs[positive_mask]
magnitude = np.abs(fft_vals[positive_mask])

# 找峰值频率
peak_idx = np.argsort(magnitude)[-5:]  # 前5个最大幅值
print(f"\n主要频率成分 (Hz): {freqs_pos[peak_idx]}")
print(f"对应幅值: {magnitude[peak_idx].round(2)}")
# 应能找到 50Hz 和 120Hz

关键要点

函数 说明
np.linalg.solve(A, b) 解 Ax = b
np.linalg.eig(A) 特征值与特征向量
np.linalg.svd(M) 奇异值分解
np.polyfit(x, y, deg) 多项式拟合
np.fft.fft(signal) 快速傅里叶变换
np.corrcoef(x, y) 相关系数矩阵
np.isnan(x) 检测 NaN

教程

NumPy 入门教程 —— 从 Python 列表到向量化思维

本章目标

  • 理解 NumPy 在 Python 生态系统中的定位
  • 掌握 ndarray 的内存模型与性能原理
  • 学会用向量化思维替代显式循环
  • 掌握常见数据操作技巧

1. NumPy 的定位

Python 科学计算栈层级:

应用层:    scikit-learn | TensorFlow/PyTorch | SciPy
             ↓              ↓                ↓
中间层:    Pandas      |  NumPy (ndarray)  |  Matplotlib
             ↓              ↓                ↓
底层:      NumPy (C/Fortran 实现 + Python API)
             ↓
         BLAS/LAPACK (MKL / OpenBLAS / Accelerate)

一句话:NumPy 是 Python 数据科学的事实标准"通用语言"。

2. ndarray 内存模型

为什么 ndarray 这么快?

# Python 列表:每个元素是 PyObject 指针
py_list = [1, 2, 3, 4]
# 内存布局: [obj_ptr] → [int_obj 1]
#             [obj_ptr] → [int_obj 2]  # 内存碎片化
#             ...

# NumPy ndarray:连续内存块
np_arr = np.array([1, 2, 3, 4], dtype=np.int64)
# 内存布局: [1][2][3][4]  # 连续 32 bytes

三大性能来源:

  1. 连续内存 — CPU 缓存友好,SIMD 向量指令
  2. C 循环 — 运算代码在 C 层执行,绕过 Python 解释器
  3. BLAS — 线性代数调用高度优化的 BLAS 库(MKL/OpenBLAS)

数据排布:C-order vs Fortran-order

arr = np.arange(12).reshape(3, 4)

# C-order(默认,按行存储)
print(arr.flags['C_CONTIGUOUS'])  # True
# 内存: [0 1 2 3 | 4 5 6 7 | 8 9 10 11]

# Fortran-order(按列存储)
arr_f = np.asfortranarray(arr)
# 内存: [0 4 8 | 1 5 9 | 2 6 10 | 3 7 11]

# 性能影响:沿连续维度的操作更快

3. 向量化思维

从"循环思维"到"向量化思维"

# 问题:计算 1 到 100 万每个数的 sin 值

# ❌ Python 循环思维
import math
result = [math.sin(i) for i in range(1, 1_000_001)]

# ✅ NumPy 向量化思维
x = np.arange(1, 1_000_001)
result = np.sin(x)

常见向量化模式

# ① 条件赋值
# ❌ 循环
for i in range(len(data)):
    if data[i] > 0:
        data[i] = 1
    else:
        data[i] = -1

# ✅ 向量化
data = np.where(data > 0, 1, -1)

# ② 分段函数
conditions = [data < 0, data == 0, data > 0]
choices = [-1, 0, 1]
result = np.select(conditions, choices)

# ③ 分组聚合(没有 groupby 但有替代方案)
labels = np.array(["a", "b", "a", "c", "b"])
values = np.array([10, 20, 30, 40, 50])

# 用 unique + 循环(分组数少时仍很快)
for label in np.unique(labels):
    mask = labels == label
    print(f"{label}: sum={values[mask].sum()}, mean={values[mask].mean()}")

4. 高级索引技巧

# ============================================================
# 布尔索引 —— 条件筛选
# ============================================================
data = np.random.randn(1000, 5)

# 筛选第3列 > 0 的所有行
filtered = data[data[:, 2] > 0]

# 多条件
mask = (data[:, 0] > 0) & (data[:, 1] < 0)  # 用 & 不是 and
filtered = data[mask]

# ============================================================
# 花式索引 —— 任意顺序选择
# ============================================================
arr = np.array(["a", "b", "c", "d", "e"])
indices = [0, 2, 4, 1, 3]
print(arr[indices])  # ['a' 'c' 'e' 'b' 'd']

# 行列同时花式索引
matrix = np.arange(25).reshape(5, 5)
rows = [0, 2, 4]
cols = [1, 3]
# 需要 ix_ 取笛卡尔积
print(matrix[np.ix_(rows, cols)])

# ============================================================
# where 寻址
# ============================================================
arr = np.array([5, 2, 8, 1, 9, 3])
indices = np.where(arr > 4)
print(indices)  # (array([0, 2, 4]),)
print(arr[indices])  # [5 8 9]

5. 性能优化技巧

5.1 避免不必要的复制

# ❌ 切片后再赋值会创建临时数组
a = np.random.randn(1000000)
b = a[::2].copy()  # copy() 强制复制

# ✅ 直接使用视图
b = a[::2]  # 视图,不复制数据

# ❌ 增量改变形状
for i in range(100):
    a = a.reshape(100, 100)  # 每次都复制

# ✅ 一次性 reshape
a = a.reshape(100, 100)

5.2 out 参数避免分配新内存

# ❌ 创建3个临时数组
result = np.sin(a) + np.cos(a) + np.tan(a)

# ✅ 使用 out 参数复用内存
result = np.empty_like(a)
np.sin(a, out=result)
np.cos(a, out=result)   # 覆盖 result
# 这个例子不太合适——如需累加:
result = np.sin(a)
np.add(result, np.cos(a), out=result)
np.add(result, np.tan(a), out=result)

5.3 选择正确的数据类型

# float64 (8 bytes) → float32 (4 bytes) 内存减半
a = np.random.randn(10_000_000).astype(np.float32)

# float64 → int8 内存减至 1/8
labels = np.random.randint(0, 10, 10_000_000, dtype=np.int8)

6. NumPy vs Pandas:何时用哪个?

场景 推荐
纯数值矩阵运算 NumPy
图像/信号处理 NumPy
混合类型表格数据 Pandas
数据清洗/探索 Pandas
深度学习张量 PyTorch/TensorFlow
简单统计聚合 都可以

思考题

  1. 为什么 np.sin(arr)[math.sin(x) for x in arr] 快得多?
  2. ndarray 的 viewcopy 有什么区别?什么时候它们是危险的?
  3. 广播机制在什么情况下会失败?如何 debug 形状不匹配?
  4. 什么时候不应该用 NumPy(即纯 Python 反而更好)?

参考资料

  1. [1] Wes McKinney. Python for Data Analysis (3rd Edition). 2022.
  2. [2] NumPy Community. NumPy 官方文档. 2024.
  3. [3] Nicolas P. Rougier. From Python to NumPy. 2017.
  4. [4] Juan Nunez-Iglesias, Stéfan van der Walt, Harriet Dashnow. Elegant SciPy: The Art of Scientific Python. 2017.