NumPy随机数模块

数据分析专题 · 随机数生成与概率分布

专题：Python数据分析系统学习

关键词：数据分析, NumPy, 随机数, random, 正态分布, 蒙特卡洛, 种子, Generator

一、NumPy随机数模块概述

NumPy的 numpy.random 模块是Python数据分析生态中最常用的随机数生成工具。它为科学计算和数据分析提供了高效、灵活的随机数生成能力，涵盖了从简单的均匀分布到复杂的高维随机抽样等多种功能。无论是模拟实验、蒙特卡洛方法、统计推断，还是机器学习中的数据随机打乱和批量采样，都离不开这个模块。

NumPy随机数模块的核心价值在于：第一，它基于高效的C语言实现，能快速生成大量随机数，远快于Python内置的 random 模块；第二，它支持多维数组的批量生成，与NumPy的数组操作无缝衔接；第三，它提供了丰富的概率分布函数，满足各种统计需求；第四，它通过种子机制确保结果可重复，这对科学研究和调试至关重要。

自NumPy 1.17版本起，随机数模块经历了一次重大架构升级，引入了基于 Generator 对象的新风格API。新风格API提供了更好的性能、更灵活的位生成器选择和更强的可扩展性。虽然旧风格的 RandomState API依然可用，官方已推荐在新代码中使用新风格API。

本笔记将系统讲解NumPy随机数模块的核心功能，包括种子设置、各类分布函数、随机抽样操作、新旧API对比以及蒙特卡洛模拟的实际应用。

环境准备

本笔记所有代码基于NumPy 1.17+版本。建议使用 pip install numpy>=1.17 安装最新版本。导入约定：import numpy as np。

二、随机数种子设置与可重复性

计算机生成的随机数是伪随机数，即通过确定性算法生成的看似随机的数列。伪随机数的特点是：给定相同的初始状态（种子），后续产生的随机数序列完全一致。这一特性在科学计算中非常重要——它允许实验可复现，方便调试和结果验证。

2.1 旧风格API的种子设置

在使用 np.random.RandomState 旧风格API时，通过 seed() 方法或构造函数设置种子。

import numpy as np # 方式一：设置全局种子 np.random.seed(42) print(np.random.rand(3)) # [0.37454012 0.95071431 0.73199394] # 方式二：创建独立的RandomState对象 rng = np.random.RandomState(2024) print(rng.randn(3)) # [-0.85044496 0.52900849 0.56419963] # 相同种子产生相同序列 rng2 = np.random.RandomState(2024) print(rng2.randn(3)) # [-0.85044496 0.52900849 0.56419963] # 全局种子会影响所有np.random.*调用 np.random.seed(42) a = np.random.uniform(0, 1, 3) np.random.seed(42) b = np.random.uniform(0, 1, 3) print(np.allclose(a, b)) # True

2.2 新风格API的种子设置

新风格API使用 default_rng() 函数创建 Generator 对象，通过 BitGenerator 管理随机数生成算法的状态。

from numpy.random import default_rng # 创建默认的Generator（不设种子，每次运行结果不同） rng = default_rng() print(rng.random(3)) # 设置种子确保可重复 rng = default_rng(12345) print(rng.random(3)) # [0.22733602 0.31675834 0.79736546] # 相同种子产生相同序列 rng2 = default_rng(12345) print(rng2.random(3)) # [0.22733602 0.31675834 0.79736546] # 指定位生成器（BitGenerator） from numpy.random import PCG64, MT19937, Philox # PCG64是默认的位生成器 rng_pcg = default_rng(PCG64(12345)) # MT19937是Mersenne Twister（旧风格使用的算法） rng_mt = default_rng(MT19937(12345)) # Philox支持并行计算中的独立流 rng_ph = default_rng(Philox(12345))

2.3 种子设置的最佳实践

最佳实践建议

1. 避免使用全局 np.random.seed()，在多模块项目中可能导致不可预见的交互。2. 推荐创建独立的 Generator 对象管理各自的状态。3. 在并行计算中使用 Philox 或 PCG64 位生成器以保证各进程的随机数流独立。4. 调试和测试时固定种子，生产环境中不设种子以保持随机性。

三、均匀分布随机数生成

均匀分布是最基本的概率分布，指在指定区间内每个值被选中的概率相等。NumPy提供了多种生成均匀分布随机数的函数。

3.1 标准均匀分布 [0, 1)

# 旧风格API print(np.random.rand(5)) # 5个[0,1)均匀分布随机数 print(np.random.rand(2, 3)) # 2x3的[0,1)均匀分布随机数矩阵 # 新风格API rng = default_rng(42) print(rng.random(5)) # 5个[0,1)均匀分布随机数 print(rng.random((2, 3))) # 2x3矩阵 print(rng.random(size=10)) # 使用size参数

3.2 指定区间的均匀分布

# randint: 生成指定范围内的随机整数 [low, high) # 旧风格 print(np.random.randint(0, 10, 5)) # 5个[0,10)的随机整数 print(np.random.randint(5, size=3)) # 3个[0,5)的随机整数 # 新风格 rng = default_rng(42) print(rng.integers(0, 10, 5)) # 5个[0,10)的随机整数 print(rng.integers(5, size=3)) # 3个[0,5)的随机整数 print(rng.integers(1, 100, (3, 4))) # 3x4的随机整数矩阵 # uniform: 生成指定范围的均匀分布浮点数 [low, high) # 旧风格 print(np.random.uniform(0, 1, 5)) # 5个[0,1)浮点数 print(np.random.uniform(-1, 1, (2, 3))) # 2x3的[-1,1)矩阵 # 新风格 rng = default_rng(42) print(rng.uniform(0, 10, 5)) # 5个[0,10)浮点数 print(rng.uniform(-1, 1, (2, 3))) # 2x3的[-1,1)矩阵

3.3 random_integers vs randint

函数对比

random_integers(low, high, size) 是旧风格API中已废弃的函数，其区间为闭区间 [low, high]，包括上下界。而 randint(low, high, size) 是半开区间 [low, high)，不包括上界。在新风格API中统一使用 integers() 方法，并通过 endpoint=True 参数控制是否包含上界。

# random_integers（已废弃，仅旧风格） print(np.random.random_integers(1, 6, 10)) # 模拟掷骰子，范围1-6闭区间 # 新风格API的等效用法 rng = default_rng(42) print(rng.integers(1, 6, 10)) # [1,6) 不包含6 print(rng.integers(1, 6, 10, endpoint=True)) # [1,6] 包含6，模拟骰子
函数API风格区间类型说明
rand(d0,d1,...,dn)旧[0,1)float标准均匀分布
randint(low,high,size)旧[low,high)int随机整数
uniform(low,high,size)旧/新[low,high)float指定范围均匀分布
random_integers(low,high,size)旧(已废弃)[low,high]int闭区间随机整数
random(size)新[0,1)floatGenerator的默认方法
integers(low,high,size)新[low,high)int可指定endpoint参数

四、正态分布随机数生成

正态分布（高斯分布）是自然界中最常见的概率分布之一，在统计学、机器学习和信号处理等领域有广泛应用。正态分布由均值（mu）和标准差（sigma）两个参数唯一确定，概率密度函数呈钟形曲线。

4.1 标准正态分布

标准正态分布的均值为0，标准差为1。

# 旧风格API print(np.random.randn(5)) # 5个标准正态分布随机数 print(np.random.randn(2, 3)) # 2x3的标准正态分布矩阵 # randn的命名含义：random + normal distribution # 新风格API rng = default_rng(42) print(rng.standard_normal(5)) # 5个标准正态分布随机数 print(rng.standard_normal((2,3))) # 2x3的标准正态分布矩阵 # 验证：大量样本的均值和标准差应接近0和1 samples = rng.standard_normal(100000) print(f"均值: {samples.mean():.4f}") # 接近0 print(f"标准差: {samples.std():.4f}") # 接近1 print(f"在[-3,3]内的比例: {((samples > -3) & (samples < 3)).mean():.4f}") # 约99.7%

4.2 指定均值和标准差的正态分布

# 旧风格API # normal(loc=均值, scale=标准差, size=形状) data = np.random.normal(loc=100, scale=15, size=10000) # IQ分数模拟：均值100，标准差15 print(f"均值: {data.mean():.2f}") # 约100 print(f"标准差: {data.std():.2f}") # 约15 print(f"最小值: {data.min():.2f}") print(f"最大值: {data.max():.2f}") # 新风格API rng = default_rng(42) data = rng.normal(loc=50, scale=10, size=(1000, 5)) # 1000行5列 print(f"整体形状: {data.shape}") # (1000, 5) print(f"整体均值: {data.mean():.2f}") # 约50 print(f"每列均值: {data.mean(axis=0)}") # 多维正态分布应用：模拟考试成绩 exam_scores = rng.normal(loc=75, scale=12, size=500) # 截断到有效范围[0, 100] exam_scores = np.clip(exam_scores, 0, 100) print(f"平均分: {exam_scores.mean():.1f}") print(f"最高分: {exam_scores.max():.0f}") print(f"及格率: {(exam_scores >= 60).mean()*100:.1f}%")

4.3 三种正态分布函数的对比

函数API风格含义参数
randn(d0,...,dn)旧标准正态 N(0,1)维度参数
normal(loc,scale,size)旧/新一般正态 N(loc,scale^2)均值、标准差、形状
standard_normal(size)新标准正态 N(0,1)形状参数
# 三种函数的关系：normal(loc,scale) = loc + scale * randn/standard_normal loc, scale = 10, 2 rng = default_rng(42) a = rng.normal(loc, scale, 100000) b = loc + scale * rng.standard_normal(100000) print(f"方法一均值: {a.mean():.4f}, 标准差: {a.std():.4f}") print(f"方法二均值: {b.mean():.4f}, 标准差: {b.std():.4f}")

五、其他常见概率分布

NumPy的随机数模块支持数十种概率分布，涵盖了科学计算和数据分析领域的大部分需求。以下介绍几种最常用的分布。

5.1 泊松分布 (Poisson)

泊松分布用于描述单位时间内随机事件发生的次数，如网站访问量、客服电话数、放射性衰变等。

# 模拟一个书店每小时的顾客到达数（平均每小时20人） rng = default_rng(42) hourly_customers = rng.poisson(lam=20, size=24*7) # 一周7天×24小时 print(f"平均每小时: {hourly_customers.mean():.1f} 人") print(f"最忙时段: {hourly_customers.max()} 人") print(f"空闲时段: {hourly_customers.min()} 人") print(f"超过30人的时段比例: {(hourly_customers > 30).mean()*100:.1f}%") # 核电站故障模拟（罕见事件，平均每年0.5次） faults = rng.poisson(lam=0.5, size=1000) # 模拟1000年 fault_years = np.sum(faults > 0) print(f"1000年中有故障的年份数: {fault_years}") print(f"最大年故障次数: {faults.max()}")

5.2 二项分布 (Binomial)

二项分布描述n次独立伯努利试验中成功次数的分布，每个试验的成功概率为p。

# 模拟投篮命中：罚球10次，命中率70% rng = default_rng(42) shots = rng.binomial(n=10, p=0.7, size=10000) print(f"平均命中数: {shots.mean():.2f}") # 约7次 print(f"命中10次(全中)概率: {(shots==10).mean()*100:.2f}%") print(f"命中≤4次概率: {(shots<=4).mean()*100:.2f}%") # 质量控制：抽样检测，抽20件，次品率5% samples = rng.binomial(n=20, p=0.05, size=100000) print(f"平均次品数: {samples.mean():.2f}") print(f"零次品概率: {(samples==0).mean()*100:.1f}%") print(f"含2件以上次品概率: {(samples>=2).mean()*100:.1f}%")

5.3 Beta分布 (Beta)

Beta分布是定义在[0,1]区间上的连续概率分布，常用于建模概率本身的分布，如点击率、转化率等比率参数。

# Beta分布模拟点击率的不确定性 # a=α, b=β 参数控制分布形状 rng = default_rng(42) # 点击次数多(α=100, β=900) → 集中在0.1附近 ctr_high_data = rng.beta(a=100, b=900, size=10000) print(f"高数据量 - 平均点击率: {ctr_high_data.mean():.4f}") print(f"高数据量 - 标准差: {ctr_high_data.std():.4f}") # 点击次数少(α=2, β=18) → 分布更分散 ctr_low_data = rng.beta(a=2, b=18, size=10000) print(f"低数据量 - 平均点击率: {ctr_low_data.mean():.4f}") print(f"低数据量 - 标准差: {ctr_low_data.std():.4f}") # Beta分布作为共轭先验的典型应用 # 先验: Beta(1,1) 均匀分布 → 观察到5次点击、95次未点击 → 后验: Beta(6,96) posterior = rng.beta(a=6, b=96, size=10000) print(f"后验点击率: {posterior.mean()*100:.2f}%")

5.4 Gamma分布 (Gamma)

Gamma分布是正偏态分布，常用于建模等待时间、降水量、保险赔付金额等正值数据。

# Gamma分布模拟日降水量 rng = default_rng(42) # shape=k, scale=θ rainfall = rng.gamma(shape=2.0, scale=5.0, size=365*3) # 3年的日降水量 print(f"日均降水量: {rainfall.mean():.1f} mm") print(f"日最大降水量: {rainfall.max():.1f} mm") print(f"无雨日比例: {(rainfall < 0.5).mean()*100:.1f}%") print(f"暴雨日(>30mm)比例: {(rainfall > 30).mean()*100:.2f}%")

5.5 指数分布 (Exponential)

指数分布用于建模事件发生的时间间隔，如客户到达间隔、设备故障间隔等。

# 指数分布模拟服务时间：平均服务时间2分钟 rng = default_rng(42) # scale = 1/λ, λ是平均发生率 service_times = rng.exponential(scale=2.0, size=10000) print(f"平均服务时间: {service_times.mean():.2f} 分钟") print(f"中位服务时间: {np.median(service_times):.2f} 分钟") print(f"服务时间>5分钟概率: {(service_times > 5).mean()*100:.1f}%") print(f"服务时间<1分钟概率: {(service_times < 1).mean()*100:.1f}%") # 指数分布的无记忆性 # P(X > t+s | X > s) = P(X > t) # 即已经等待了5分钟，再等3分钟的概率 = 新来者等3分钟的概率
分布函数参数典型应用
泊松poisson(lam, size)lam=平均发生率计数数据、到达次数
二项binomial(n, p, size)n=试验次数, p=成功概率AB测试、质量控制
Betabeta(a, b, size)a,b=形状参数概率建模、贝叶斯统计
Gammagamma(shape, scale, size)shape=形状,scale=尺度降水量、保险精算
指数exponential(scale, size)scale=平均间隔排队论、可靠性分析
卡方chisquare(df, size)df=自由度假设检验、拟合优度
F分布f(dfnum, dfden, size)dfnum/dfden=自由度方差分析(ANOVA)
多项multinomial(n, pvals, size)n=试验数, pvals=概率向量分类抽样

六、随机抽样与排列

在实际数据分析中，经常需要从数据集中随机抽样、打乱数据顺序或生成随机排列。NumPy提供了高效的函数来完成这些操作。

6.1 随机抽样 (choice)

# 从一维数组中随机抽样 rng = default_rng(42) population = np.arange(100) # 有放回抽样 sample_with = rng.choice(population, size=10, replace=True) print(f"有放回抽样: {sample_with}") print(f"是否有重复: {len(np.unique(sample_with)) < 10}") # 无放回抽样 sample_without = rng.choice(population, size=10, replace=False) print(f"无放回抽样: {sample_without}") print(f"是否有重复: {len(np.unique(sample_without)) < 10}") # 设置概率权重 items = np.array(['A', 'B', 'C']) weights = np.array([0.7, 0.2, 0.1]) samples = rng.choice(items, size=1000, p=weights) print(f"A比例: {(samples=='A').mean()*100:.1f}%") print(f"B比例: {(samples=='B').mean()*100:.1f}%") print(f"C比例: {(samples=='C').mean()*100:.1f}%") # 旧风格API print(np.random.choice(10, 5)) # 从[0,10)抽5个 print(np.random.choice(5, 3, replace=False)) # 无放回 print(np.random.choice(['a','b','c'], 2)) # 从列表抽2个

6.2 数组打乱 (shuffle)

# shuffle: 原地打乱数组（仅沿第一个轴） rng = default_rng(42) data = np.arange(10) print(f"打乱前: {data}") rng.shuffle(data) print(f"打乱后: {data}") # 多维数组的打乱（打乱行） matrix = np.arange(20).reshape(5, 4) print(f"打乱前:\n{matrix}") rng.shuffle(matrix) print(f"打乱后:\n{matrix}") # 旧风格API arr = np.arange(5) np.random.shuffle(arr) print(arr) # 注意：shuffle在原地操作，不返回新数组！ # 如需保留原始数组，先复制再打乱 data = np.arange(10) data_copy = data.copy() rng.shuffle(data_copy) print(f"原始: {data}") print(f"打乱副本: {data_copy}")

6.3 随机排列 (permutation)

# permutation: 返回打乱后的新数组或随机排列的索引 rng = default_rng(42) # 生成0-9的随机排列 perm = rng.permutation(10) print(f"随机排列: {perm}") # 打乱数组并返回副本（不修改原数组） arr = np.array([10, 20, 30, 40, 50]) perm_arr = rng.permutation(arr) print(f"原始: {arr}") print(f"排列后: {perm_arr}") # 多维数组：permutation仅打乱第一轴 matrix = np.arange(12).reshape(4, 3) print(f"原始矩阵:\n{matrix}") print(f"排列后:\n{rng.permutation(matrix)}") # 旧风格API print(np.random.permutation(5)) print(np.random.permutation([1, 2, 3, 4]))

6.4 实用示例：数据集划分

# 机器学习中的训练集/测试集划分 rng = default_rng(42) n_samples = 1000 indices = np.arange(n_samples) # 随机打乱索引 rng.shuffle(indices) # 70%训练集, 30%测试集 split = int(0.7 * n_samples) train_idx = indices[:split] test_idx = indices[split:] print(f"训练集大小: {len(train_idx)}") print(f"测试集大小: {len(test_idx)}") # 确保无重叠 print(f"交集是否为空: {len(np.intersect1d(train_idx, test_idx)) == 0}")

七、新风格API：Generator与BitGenerator

自NumPy 1.17起引入的 Generator 对象是新风格随机数API的核心组件。它在架构、性能和可扩展性方面都有了显著提升。

7.1 架构对比

旧风格 RandomState

  • 基于Mersenne Twister (MT19937)算法
  • 全局状态管理
  • 难以并行化
  • 速度相对较慢
  • 已冻结，不再增加新功能

新风格 Generator

  • 可插拔BitGenerator架构
  • 对象化管理状态
  • 支持独立流（并行友好）
  • 速度提升2-5倍
  • 持续更新，支持新分布
# Generator的创建方式 from numpy.random import default_rng, Generator from numpy.random import PCG64, MT19937, Philox, SFC64 # 方式一：使用default_rng（推荐，自动选择最佳算法） rng = default_rng() # 方式二：直接指定BitGenerator rng_pcg64 = Generator(PCG64()) # 默认算法，高质量高速度 rng_mt19937 = Generator(MT19937()) # Mersenne Twister，兼容旧版 rng_philox = Generator(Philox()) # 支持独立流，适合并行 rng_sfc64 = Generator(SFC64()) # 轻量级算法 # PCG64的状态大小 pcg = PCG64(seed=42) print(f"PCG64状态: {pcg.state}") # 生成独立随机数流（并行计算关键特性） child_rngs = [default_rng(PCG64(seed).jumped(i)) for i in range(4)] # 每个child_rng都拥有独立的随机数流

7.2 新旧API函数对照表

功能旧风格 (RandomState)新风格 (Generator)
创建对象np.random.RandomState(seed)np.random.default_rng(seed)
[0,1)均匀分布rand(d0,d1,...)random(size)
指定范围均匀uniform(low,high,size)uniform(low,high,size)
随机整数randint(low,high,size)integers(low,high,size)
标准正态randn(d0,d1,...)standard_normal(size)
一般正态normal(loc,scale,size)normal(loc,scale,size)
随机抽样choice(a,size,...)choice(a,size,...)
打乱shuffle(x)shuffle(x)
排列permutation(x)permutation(x)
设置种子seed(s) 或 RandomState(s)default_rng(s)

7.3 性能对比

# 新旧API性能对比 import time n = 10_000_000 # 旧风格 start = time.time() np.random.uniform(0, 1, n) old_time = time.time() - start # 新风格 rng = default_rng(42) start = time.time() rng.uniform(0, 1, n) new_time = time.time() - start print(f"旧风格 RandomState: {old_time:.3f} 秒") print(f"新风格 Generator: {new_time:.3f} 秒") print(f"加速比: {old_time/new_time:.1f}x") # 生成不同类型的随机数 # uniform, normal, integers, poisson rng = default_rng(42) floats = rng.uniform(0, 1, 1000000) # 100万个均匀分布 normals = rng.normal(0, 1, 1000000) # 100万个正态分布 ints = rng.integers(0, 100, 1000000) # 100万个整数 poissons = rng.poisson(5, 1000000) # 100万个泊松分布 print(f"完成: 生成了4组100万随机数, 形状: {floats.shape}")

八、蒙特卡洛模拟应用

蒙特卡洛方法是一类通过随机采样来解决确定性问题的数值计算方法。NumPy的随机数模块是实现蒙特卡洛模拟的基础工具。以下通过几个经典案例演示其应用。

8.1 计算圆周率π

通过单位正方形内随机投点，统计落在内切圆中的比例来估算π值。

# 蒙特卡洛法估算π值 rng = default_rng(42) def estimate_pi(n_points): x = rng.uniform(-1, 1, n_points) y = rng.uniform(-1, 1, n_points) inside = (x**2 + y**2) <= 1.0 pi_estimate = 4 * inside.sum() / n_points return pi_estimate for n in [1000, 10000, 100000, 1000000]: pi_est = estimate_pi(n) error = abs(pi_est - np.pi) / np.pi * 100 print(f"投点数: {n:>8d} π估算: {pi_est:.6f} 误差: {error:.4f}%") # 多次实验观察收敛性 n_experiments = 1000 estimates = np.array([estimate_pi(10000) for _ in range(n_experiments)]) print(f"1000次实验平均: {estimates.mean():.6f}") print(f"1000次实验标准差: {estimates.std():.6f}") print(f"95%置信区间: [{np.percentile(estimates, 2.5):.6f}, {np.percentile(estimates, 97.5):.6f}]")

8.2 模拟股票价格路径

使用几何布朗运动模型模拟股票价格的随机游走过程。

# 几何布朗运动模拟股票价格 rng = default_rng(42) def simulate_stock_price(S0, mu, sigma, T, n_steps, n_paths=1): """ 模拟股票价格路径 S0: 初始价格 mu: 年化收益率 sigma: 年化波动率 T: 时间跨度（年） n_steps: 时间步数 n_paths: 模拟路径数 """ dt = T / n_steps # 生成标准正态分布随机数 z = rng.standard_normal((n_steps, n_paths)) # 计算价格变化 returns = (mu - 0.5 * sigma**2) * dt + sigma * np.sqrt(dt) * z # 累积乘积得到价格序列 price = S0 * np.exp(np.cumsum(returns, axis=0)) return np.vstack([np.full(n_paths, S0), price]) # 包含初始价格 # 模拟参数 S0 = 100.0 # 初始价格100元 mu = 0.08 # 年化收益率8% sigma = 0.25 # 年化波动率25% T = 1.0 # 模拟1年 n_steps = 252 # 交易日数 n_paths = 1000 # 模拟1000条路径 prices = simulate_stock_price(S0, mu, sigma, T, n_steps, n_paths) final_prices = prices[-1, :] print(f"最终价格均值: {final_prices.mean():.2f}") print(f"最终价格中位数: {np.median(final_prices):.2f}") print(f"最终价格标准差: {final_prices.std():.2f}") print(f"亏损概率(P<100): {(final_prices < S0).mean()*100:.1f}%") print(f"涨幅超20%概率: {(final_prices > S0*1.2).mean()*100:.1f}%") print(f"95%价格区间: [{np.percentile(final_prices, 2.5):.1f}, {np.percentile(final_prices, 97.5):.1f}]")

8.3 投掷硬币与中心极限定理

# 通过蒙特卡洛模拟验证中心极限定理 rng = default_rng(42) def coin_toss_experiment(n_tosses, n_experiments): """进行n_experiments次投掷n_tosses次硬币的实验""" # 0表示反面，1表示正面 results = rng.binomial(n=n_tosses, p=0.5, size=n_experiments) proportions = results / n_tosses # 正面比例 return proportions # 不同样本量的实验 for n in [10, 50, 200, 1000]: props = coin_toss_experiment(n, 10000) print(f"投掷{n:>4d}次: 均值={props.mean():.4f} 标准差={props.std():.4f} 理论σ={np.sqrt(0.25/n):.4f}") # 大数定律演示：随着投掷次数增加，比例趋近于0.5 cumulative = np.cumsum(rng.binomial(1, 0.5, 10000)) / np.arange(1, 10001) print(f"投掷100次后正面比例: {cumulative[99]:.4f}") print(f"投掷1000次后正面比例: {cumulative[999]:.4f}") print(f"投掷10000次后正面比例: {cumulative[9999]:.4f}")

8.4 业务应用：库存管理模拟

# 模拟报童问题（News Vendor Problem） # 每天进货Q份报纸，成本c，售价p，未售出残值v # 需求量D服从泊松分布 rng = default_rng(42) c, p, v = 0.5, 2.0, 0.1 # 成本、售价、残值 demand_mean = 100 # 平均日需求 def simulate_profit(Q, n_days=10000): demand = rng.poisson(demand_mean, n_days) sold = np.minimum(demand, Q) leftover = Q - sold revenue = sold * p + leftover * v cost = Q * c profit = revenue - cost return profit.mean(), profit.std(), np.percentile(profit, 5) # 测试不同的进货量 for Q in [80, 90, 100, 110, 120]: mean_profit, std_profit, p5 = simulate_profit(Q) print(f"进货{Q:>3d}份: 日均利润={mean_profit:.2f} 标准差={std_profit:.2f} 95%VaR={p5:.2f}") # 寻找最优进货量 profits = [simulate_profit(Q)[0] for Q in range(70, 131)] opt_Q = 70 + np.argmax(profits) print(f"\n最优进货量: {opt_Q}份, 预期日利润: {max(profits):.2f}元")

九、核心要点总结

NumPy随机数模块要点速查

  • 两种API风格：旧风格 np.random.* (RandomState) 和新风格 rng = default_rng() (Generator)，推荐优先使用新风格API。
  • 种子设置： seed() 设置全局种子（旧风格），default_rng(seed) 创建带种子的Generator对象（新风格）。相同种子确保结果可复现。
  • 均匀分布： random() / uniform() / randint() / integers() 生成不同范围和类型的均匀分布随机数。
  • 正态分布： standard_normal() 生成标准正态分布，normal(loc,scale) 生成指定均值和标准差的正态分布。
  • 其他分布： 泊松、二项、Beta、Gamma、指数等数十种分布函数满足各类统计建模需求。
  • 随机抽样： choice() 实现有放回/无放回抽样，shuffle() 原地打乱数组，permutation() 返回打乱后的副本。
  • 蒙特卡洛应用： 随机数模块是实现蒙特卡洛模拟的核心工具，适用于π值估算、金融建模、风险管理等场景。
  • 性能优势： 新风格API比旧风格快2-5倍，推荐在性能敏感场景中优先使用。

迁移建议

如果你的项目正在使用旧风格API np.random.*，建议逐步迁移到新风格 Generator API。迁移步骤：1) 将全局的 np.random.seed() 替换为 default_rng(seed) 创建局部Generator对象。2) 将 np.random.rand()/randn() 替换为 rng.random()/standard_normal()。3) 将 np.random.randint() 替换为 rng.integers()。4) 其他函数大多保持同名直接替换。

十、进一步思考

掌握NumPy的随机数模块不仅是学习一个工具库，更是理解统计模拟和科学计算思维的重要一步。在实际应用中，以下几个方面值得深入探索：

深入方向

  • 并行随机数生成： 使用 Philox 或 PCG64.jumped() 方法在高性能计算中生成独立的随机数流，确保并行模拟的正确性。
  • 贝叶斯统计推断： 结合Beta分布等先验分布，使用蒙特卡洛方法进行后验采样，搭建完整的贝叶斯分析流程。
  • 随机过程模拟： 在金融工程和物理模拟领域，深入研究布朗运动、伊藤过程、马尔可夫链等随机过程的NumPy实现。
  • 与SciPy.stats的配合： NumPy负责随机数生成，SciPy.stats提供完整的概率分布统计函数（PDF、CDF、分位数等），两者配合使用可构建完整的统计分析流程。
  • 与深度学习的结合： 理解NumPy随机数机制有助于深入理解深度学习框架（如PyTorch、TensorFlow）中类似功能的设计思路和使用方法。

学习建议

建议通过动手实践来巩固所学知识：先创建独立的Generator对象熟悉新风格API，然后尝试用蒙特卡洛方法解决一个实际问题（如股票模拟、风险评估或物理模拟），并在实践中逐步掌握更多分布函数的使用方式。