简介：PyVaR 是一套面向金融工程学习者与量化分析人员的风险价值计算代码包，围绕蒙特卡罗模拟与几何布朗运动（GBM）模型估算股票 VaR，并借助共享内存并行机制提升多核环境下的模拟效率。压缩包共 182 个文件、约 14.62MB，以 98 个 m 文件、29 个 txt 说明、11 个 f90 源码及若干 mex 二进制、mat 数据、py 脚本和 readme 文档为主，兼顾 MATLAB 调用、Fortran 数值计算与 Python 示例，目录层次便于按模块查阅。已有 858 人学习，适合需要理解 VaR 建模流程、并行加速思路与跨语言实现细节的读者。包内可看到蒙特卡罗路径生成、GBM 参数设置、示例脚本、使用说明与测试验证材料，能帮助读者把置信水平、持有期和价格分布等要素落到代码中，并观察并行任务拆分对计算速度的影响，为风险管理与投资决策场景提供可复用的实验基础。

1. 收盘前十分钟，风控要的那个数

下午两点五十，风控同事发来一句话：手上这只票，95% 置信度下，明天最多亏多少。要的不是“大概跌三个点”，而是一个带分布假设、能写进日报的数字——这就是 VaR（风险价值）。pyVaR 做的是这件事：用几何布朗运动（GBM）描述股价，用蒙特卡罗抽出几十万到上百万条未来路径，再从损益分布的尾部取分位数，顺带给出期望尾部损失 ES。整个计算走共享内存并行，进程池分块跑，路径结果落在同一块共享内存里，省掉子进程把大数组 pickle 回主进程的来回拷贝。它适合手里有一篮子标的、每天都要批量出一版风险测算，又不想为这点计算量去碰集群的人。模型口径和原版 MATLAB 脚本一致，我拿它对拍过数值，中间有几处对不上的坑，后面细讲。

2. 工程结构与运行入口：从一份 toml 配置到命令行跑通

2.1 技术选型：numpy 向量化打底，multiprocessing 做并行

这类计算天生适合向量化：生成标准正态随机数、乘波动率、累加、取终值，四步全是 C 层循环，numpy 不用写一行显式 for。真正卡住性能的从来不是算法，而是“一天几十次调用、每次几百万条路径”堆出来的墙钟时间。

选型上有三个明显的分岔口。第一，不要用 pandas 逐行模拟，DataFrame 的索引对齐和 NaN 语义在这个场景里纯属负担，价格矩阵一进 numpy 就该丢掉列名。第二，不要指望多线程，np.random.Generator.standard_normal和np.cumsum这两个热点都不释放 GIL，线程开 16 个也只是轮流排队，还额外多出上下文切换。第三，numba 或 Cython 确实能把单线程压到极限，但环境里多一套编译链，部署到别人的机器上就容易翻车，收益也不如直接多进程来得实在。

所以主干就是numpy出数值、multiprocessing.Pool出并行、shared_memory出结果回收。依赖上只需要numpy、pandas（读行情）、pyarrow（parquet）、scipy（可选，回测里的卡方分位数），没有重框架。

2.2 目录结构与模块职责

包名用全小写pyvar，避免和项目名大小写混用导致import踩坑。目录按“数据 → 模型 → 统计 → 调度”切，每个文件只干一件事：

pyvar/ ├── pyvar/ │ ├── __init__.py │ ├── data.py # 行情读取、清洗、GBM 参数估计 │ ├── gbm.py # GBM 路径/终值生成（单进程纯函数） │ ├── var.py # 分位数、VaR、ES、回测统计 │ ├── parallel.py # 共享内存 + 进程池调度 │ ├── config.py # toml 配置解析与参数校验 │ └── cli.py # 命令行入口 ├── configs/ │ └── default.toml ├── tests/ │ ├── test_gbm.py # 均值/方差收敛性测试 │ └── test_var.py # 与解析解、分位数方法对照 ├── requirements.txt └── README.md

这种切法的好处是可测性：gbm.py和var.py都是纯函数，单元测试和串行对照都方便；parallel.py只负责“怎么切”，出问题九成是并发或内存问题，范围收得很窄。

2.3 环境安装与最小可运行命令

Python 3.10 以上，主要是shared_memory和np.random.Generator的行为在旧版本上有差异。虚拟环境建好之后，装依赖、装本包：

python -m venv .venv source .venv/bin/activate # Windows: .venv\Scripts\activate pip install numpy pandas pyarrow scipy pip install -e . # 以可编辑模式装进当前环境

配置文件用 toml，比 argparse 堆一屏参数清爽，也方便版本管理：

# configs/default.toml [market] prices = "data/prices.parquet" # 宽表：index=日期，columns=代码 lookback_days = 500 # 估计 mu/sigma 的回看窗口 confidence = 0.95 # VaR 置信度 horizon_days = 1 # 持有期 [model] steps_per_day = 1 # 每天拆几步，欧式一步即可 n_paths = 2000000 # 总路径数 [runtime] workers = 8 # 进程数，<= 物理核数 chunk_paths = 250000 # 每块路径数 seed = 20240101 # 固定种子，保证可复现

跑起来只要一行：

python -m pyvar.cli --config configs/default.toml --out reports/var.csv

输出是一张按标的排列的表，每行包含ticker、mu、sigma、var、es、n_paths、elapsed_s。跑通之后先别急着调大n_paths，先拿一两个标的和串行版本对一下数值，确认分位数口径没问题，再去压性能。

2.4 输入数据长什么样

行情必须是宽表，索引是交易日，列是标的代码，值是收盘价（复权后）。用长表也行，但每次 pivot 都是一笔额外开销，批量测算时没必要。读入后的第一件事是ffill：停牌日不填，收益率序列就会少一段，年化波动率直接偏小，VaR 自然偏乐观，这是最隐蔽的一类错误。

# pyvar/data.py import numpy as np import pandas as pd def load_prices(path: str) -> pd.DataFrame: """宽表读入：index=日期, columns=代码, values=复权收盘价""" df = pd.read_parquet(path) df = df.sort_index().ffill() # 停牌日沿用前值 df = df.dropna(how="any") # 上市晚的标的直接剔除，别补零 if (df <= 0).any().any(): raise ValueError("价格序列存在非正值，检查复权方式") return df def estimate_gbm_params(prices: pd.DataFrame, lookback: int = 500, ann: int = 252): """返回 {ticker: (mu, sigma)}，均为年化值""" px = prices.tail(lookback) log_ret = np.log(px / px.shift(1)).dropna() sigma = log_ret.std(ddof=1) * np.sqrt(ann) # 年化波动率 mu = log_ret.mean() * ann # 年化漂移 return {c: (float(mu[c]), float(sigma[c])) for c in prices.columns}

lookback决定参数稳定性：太短，sigma 跟着最近一波行情剧烈跳动；太长，又跟不上市场结构变化。500 个交易日是常见折中，约两年。ddof=1用样本标准差，和 MATLAB 的std默认行为一致，这点对后续对拍很重要。mu的估计噪声远大于sigma，实际做短期风险测算时，把mu直接设成无风险利率甚至 0、只让sigma从历史数据估，是更稳的做法——漂移项对一天的 VaR 影响几乎可以忽略，但估计误差会实打实传进结果里。

3. 核心算法落地：GBM 从参数到损益分布

3.1 两套写法：逐步模拟和一步到位，别用错场景

GBM 的解析解是S_T = S_0 * exp((mu - sigma^2/2)T + sigma*sqrt(T)*Z)，Z 是标准正态。很多入门实现会老老实实按天迭代，生成一个(n_paths, steps)的完整价格矩阵。这样做不是错，但内存占用是n_paths * steps * 8字节：100 万条路径、252 个交易日，就是 2GB，机器立马吃紧。

对纯欧式的 VaR 来说，只需要终值分布，中间路径一点用没有。一步到位只需要(n_paths, 1)，内存直接降两个数量级。所以工程里我一般写两条路：

# pyvar/gbm.py import numpy as np def gbm_terminal_fast(S0, mu, sigma, T, n_paths, rng): """只求终值，欧式 VaR 用这个，内存 O(n_paths)""" z = rng.standard_normal(n_paths) return S0 * np.exp((mu - 0.5 * sigma ** 2) * T + sigma * np.sqrt(T) * z) def gbm_paths(S0, mu, sigma, T, steps, n_paths, rng): """完整路径矩阵，给障碍/路径依赖产品用，内存 O(n_paths * steps)""" dt = T / steps z = rng.standard_normal((n_paths, steps)) incr = (mu - 0.5 * sigma ** 2) * dt + sigma * np.sqrt(dt) * z log_paths = np.cumsum(incr, axis=1) # 沿时间维累加 return S0 * np.exp(log_paths)

两个函数在数学上等价，区别只在内存和用途。incr里(mu - 0.5*sigma**2)*dt是漂移修正项，少了这个减号，模拟出来的价格期望会系统性偏高，一年期上能差出几个百分点，属于那种“结果看着合理、但就是不对”的典型错误。np.cumsum(..., axis=1)的axis必须是 1，写成 0 就变成沿路径方向累加，出来是一片乱码数字。

T用年为单位，horizon_days=1时T = 1/252；steps只有需要中间路径时才大于 1，欧式场景固定 1。

3.2 参数估计：mu 和 sigma 从哪来，决定了结果可信度

估值口径必须统一。历史数据是日频收盘价，年化因子在 A 股常用 252，商品或加密资产会是 365，这个数字一改，sigma 差 20%，VaR 跟着差 20%。团队里最好在配置里写死，别让每个人自己拍。

mu的处理更微妙。用历史对数收益均值年化，500 天样本估出来的值动辄 ±30%，噪声比信号大得多。我一般的做法是：短期风险测算里把mu设成 0 或资金成本，只让sigma从数据来；只有在做长期压力测试、需要模拟价格中枢时才用历史漂移。这一条在 2.4 的代码里是显式暴露的，改配置就能切，不用动模型层。

还有一点容易忽略：GBM 假设对数收益服从正态，真实市场是尖峰厚尾，95% 分位数上还好，99% 上 GBM 会系统性低估尾部风险。所以项目里保留了分布可替换的位置，第 6 章会讲怎么换成 t 分布。

3.3 向量化路径生成：为什么不用循环

rng.standard_normal((n_paths, steps))一次调用就拿到整个随机数矩阵，内部是 C 的梅森旋转或 PCG64，比 Python 层循环快三到四个数量级。代价是内存峰值高，所以并行版本里要配合分块，一块一块生成，见第 4 章。

incr.sum(axis=1)和log_paths[:, -1]在数值上等价，前者少一次完整矩阵拷贝，只求终值时用sum。如果steps=1，两者一样。

一个细节：S0 * np.exp(...)里的乘法在数组外面，numpy 会广播，不用写成np.exp(...) * S0，没有性能差异，但前者语义更清楚——先算收益率，再折回价格。

3.4 从损益分布到 VaR 与 ES

拿到终值数组后，损益pnl = ST - S0，VaR 就是左尾分位数取正：

# pyvar/var.py import numpy as np def var_es(pnl: np.ndarray, confidence: float = 0.95): """pnl: 损益数组（正=赚，负=亏） 返回 (VaR, ES)，都取正数表示损失幅度""" q = np.percentile(pnl, (1 - confidence) * 100, method="linear") var = -q tail = pnl[pnl <= q] es = -tail.mean() if tail.size else var return float(var), float(es)

method="linear"是 numpy 1.22 之后的写法，之前叫interpolation，用旧版本会直接抛TypeError。这个参数是 VaR 对拍失败的头号嫌疑：MATLAB 的prctile默认插值方式、pandas 的quantile默认linear、Excel 的PERCENTILE又是另一种，同样一组数能差出千分之几，路径数少的时候这个差异肉眼可见。

tail为空表示分位点正好落在最小值上，属于极小概率，兜底返回var避免nan。ES 是尾部均值，比 VaR 更能反映极端损失，两个一起输出，别只报 VaR。

注意：pnl的符号约定必须全项目统一。有的实现用收益率、有的用绝对金额，混着传进去会让 VaR 符号翻转，而且这种错在正收益占多数时看不出来，只有亏钱那天才暴露。

4. 共享内存并行：切分、回传与调优

4.1 为什么必须绕开 GIL，以及怎么切

numpy的部分 ufunc 会释放 GIL，但随机数生成和cumsum这类归约操作不在其列。实测下来，8 线程跑这段代码的加速比不到 1.3，还不如单线程。所以并行单元必须是进程。

切分策略就一个原则：每块要够大，摊平进程启动和内存映射的固定成本。我一般让chunk_paths >= 100000，workers取物理核数（超线程对纯浮点计算帮助有限，反而增加调度开销）。块数不用等于进程数，Pool.map会自动排队，块数多一点（比如 2 到 4 倍 workers）反而能减少长尾，防止某个进程分到的块特别慢把整体拖住。

但块数也不是越多越好。每多一块，就多一次共享内存 attach、一次随机数生成器初始化，单块 5 万条以下时，这部分开销能占到总时间的 40% 以上。

4.2 共享内存写回与可复现的种子管理

朴素的写法是让 worker 返回自己的terminal数组，主进程np.concatenate。问题是每次返回都要 pickle 一遍，100 万条 float64 是 8MB，8 个进程来回一趟就是几十兆的序列化开销，延迟还高。共享内存的思路是：主进程先划一块地，worker 直接把结果写进去，返回值只是一个状态码。

# pyvar/parallel.py import numpy as np from multiprocessing import Pool, shared_memory def _worker(job): name, n_total, start, stop, S0, mu, sigma, T, steps, seed = job n = stop - start rng = np.random.default_rng(seed) # 每个 worker 独立种子 dt = T / steps z = rng.standard_normal((n, steps)) incr = (mu - 0.5 * sigma ** 2) * dt + sigma * np.sqrt(dt) * z terminal = S0 * np.exp(incr.sum(axis=1)) # 只留终值，省内存 shm = shared_memory.SharedMemory(name=name) # attach 已有块 out = np.ndarray((n_total,), dtype=np.float64, buffer=shm.buf) out[start:stop] = terminal # 直接写进共享区 shm.close() # worker 只 close，不 unlink return 0 def simulate_terminal_parallel(S0, mu, sigma, T, steps, n_paths, workers=8, chunk=250_000, seed=0): shm = shared_memory.SharedMemory(create=True, size=n_paths * 8) out = np.ndarray((n_paths,), dtype=np.float64, buffer=shm.buf) starts = list(range(0, n_paths, chunk)) # 用 SeedSequence 派生子种子，保证可复现且互不重叠 seeds = [int(s.generate_state(1)[0]) for s in np.random.SeedSequence(seed).spawn(len(starts))] jobs = [(shm.name, n_paths, st, min(st + chunk, n_paths), S0, mu, sigma, T, steps, sd) for st, sd in zip(starts, seeds)] try: with Pool(workers) as pool: pool.map(_worker, jobs, chunksize=1) result = out.copy() # 拷出来再释放共享块 finally: del out # 先删 ndarray 引用 shm.close() shm.unlink() # 主进程负责销毁 return result

几个参数和顺序都不能乱。

size=n_paths * 8是按 float64 算的字节数，换 float32 要同步改成 4，两边不一致会直接抛ValueError: buffer is smaller than requested size。

SeedSequence(seed).spawn(n)是关键：如果用seed + i当每个 worker 的种子，相邻种子在 PCG64 下的前几个输出相关性强，等于人为制造了路径相关性，VaR 会偏小。spawn派生的是统计独立的子序列，可复现性也保住了——同一seed跑两遍，结果逐位相同。

写回时out[start:stop] = terminal用的是切片赋值，不产生中间数组；如果写成out = terminal，只是重新绑定名字，共享区里那片内存还是空的，主进程拿到的是一堆 0。这是我自己早期犯过的错，调了半天以为是并行没跑。

finally里的顺序是del out→close→unlink。少了del，ndarray 的 buffer 还引用着共享区，unlink在某些平台上会抛BufferError。worker 里只能close，绝不能unlink——那会把还在被其他进程使用的块直接删掉。unlink只在主进程、且只做一次。

提示：Windows 上SharedMemory(name=...)的name不要带路径分隔符，用默认生成的随机名最稳；跨机器或跨容器的场景不要共享内存，改用np.memmap落到磁盘。

4.3 和 MATLAB parfor 版对拍，顺带把分块调优做掉

原版是 MATLAB 脚本，核心循环用parfor并行，口径和这里一致。对拍不能逐条路径比，因为rng(..., 'twister')和 numpy 的 PCG64 生成的随机序列完全不同步，要比的是统计量。

% 原版 MATLAB 核心片段（仅用于对拍口径） rng(20240101, 'twister'); Z = randn(nPaths, steps); incr = (mu - 0.5*sigma^2)*dt + sigma*sqrt(dt)*Z; ST = S0 * exp(sum(incr, 2)); % 注意沿第 2 维求和 varValue = -prctile(ST - S0, (1-confidence)*100);

对拍时固定S0、mu、sigma、T、nPaths，两边各跑 20 次取均值。我能接受的一致性标准是：var的相对差在 1% 以内，es在 2% 以内，且随nPaths增大差异单调收窄。如果差异不收敛，一定是口径问题，不是随机误差。

对不上的地方主要有三处，列出来备查：

调优部分我在 8 核测试机上跑了一组，固定 200 万条路径、steps=1：

8 进程之后开始回落，因为测试机只有 8 个物理核，超出的进程在抢资源。加速比到不了线性的 8 倍，剩下的损耗来自共享内存映射和进程启动，路径数越大这部分占比越小。要是把chunk_paths压到 5 万，8 进程的总耗时反而会涨到 3 秒以上。

5. 避坑与常见问题排查：五个我真实踩过的坑

5.1 进程开太多反而变慢

现象：把workers从 8 调到 32，总耗时从 2.1 秒涨到 3.8 秒，机器风扇狂转但 CPU 利用率上不去。

原因：物理核只有 8 个，多出来的 24 个进程在做时间片轮转，每个进程的共享内存 attach 和随机数生成器初始化都要完整走一遍，纯开销。另外每块路径数被摊薄到 6 万多条，固定成本占比飙升。

解决：workers上限锁死在multiprocessing.cpu_count()的物理核数，配置里加断言；chunk_paths保持 10 万以上。真要压榨，宁可开 8 个进程每块 50 万，也别开 16 个进程每块 12 万。

5.2 Windows 下启动就抛 RuntimeError

现象：Linux 上跑得好好的，换到 Windows 直接报An attempt has been made to start a new process before the current process has finished its bootstrapping phase。

原因：Linux 默认fork，子进程继承父进程内存；Windows 和 macOS 新版默认spawn，子进程要重新 import 主模块。如果Pool(...)写在模块顶层，import 时会再次执行到创建进程池那一行，无限递归，Python 用这个异常把它拦下来。

解决：所有进程池创建代码必须放在if __name__ == "__main__":下面，_worker函数定义在模块顶层而不是闭包里——闭包里的函数 spawn 时 pick 不到。实在绕不开，用multiprocessing.get_context("spawn")显式指定，行为在三个平台上就一致了。

5.3 共享内存泄漏与 BufferError

现象：跑几十次之后/dev/shm被占满，或者程序退出时报BufferError: cannot close exported pointers exist，严重时下一次运行直接FileExistsError。

原因：SharedMemory是操作系统级资源，进程正常退出不会自动清理，必须显式unlink。另一种情况是 ndarray 还引用着这块 buffer 就调了close，引用计数没归零，资源挂在那里。

解决：try/finally包住整个计算过程，finally里严格按del ndarray→shm.close()→shm.unlink()三步走。不要图省事只写close，也不要在 worker 里调unlink。调试阶段可以在启动时扫一遍/dev/shm，把残留的psm_*清掉。

5.4 VaR 数值和 MATLAB 版对不上

现象：同一组参数，MATLAB 出 0.0231，Python 出 0.0227，差 1.7%，而且加大路径数差异不收敛。

原因：查下来是三处叠加。一是np.percentile的method没指定，默认linear和 MATLABprctile的插值方式不同；二是年化因子，一边用 252 一边用 250；三是最隐蔽的——MATLAB 那边用的是简单收益，Python 这边估参数用的是对数收益，两条链路的前提假设不一样。

解决：把插值方法、年化因子、收益类型三件事写进配置并打印到日志里，对拍前先核这三项。年化因子对 sigma 的影响是线性的，收益类型的影响在波动率大的标的上能到 5% 以上，不是小数。

5.5 路径矩阵把内存吃满

现象：n_paths调到 500 万，机器卡死，dmesg里能看到 OOM killer。

原因：gbm_paths一次性生成(5000000, 252)的 float64 矩阵，将近 10GB，再加上incr和log_paths两个中间数组，峰值能到 30GB。

解决：欧式场景直接换gbm_terminal_fast，内存降到 40MB。确实要中间路径的，按chunk_paths分块生成、算完即弃，或者用float32把内存砍半——代价是尾部精度损失，99% 分位数上要谨慎。别指望 swap 能救，这类大数组一换出就是灾难。

6. 二次开发：加对偶变量、换厚尾分布、做回测验证

拿到项目后最值得动的一刀是方差缩减。蒙特卡罗的误差按1/sqrt(n)收敛，路径数翻 4 倍才能把标准误减半，代价是 4 倍时间。对偶变量法不增加路径数，只让一半路径用Z、另一半用-Z，同样的计算量下标准误能降 30% 到 50%。

# 在 gbm_terminal_fast 基础上改造 def gbm_terminal_antithetic(S0, mu, sigma, T, n_paths, rng): half = n_paths // 2 z_half = rng.standard_normal(half) z = np.concatenate([z_half, -z_half]) # 对偶配对 return S0 * np.exp((mu - 0.5 * sigma ** 2) * T + sigma * np.sqrt(T) * z)

n_paths必须是偶数，奇数时最后一条会被丢掉；half用整除，在配置校验里加一条if n_paths % 2: raise。

换厚尾分布也简单，把标准正态换成自由度 5 到 7 的 t 分布再标准化，尾部就会明显变厚，99% VaR 会往上走一截：

from scipy import stats z = stats.t.rvs(df=5, size=n_paths, random_state=rng) z = z / np.sqrt(5 / 3) # t(5) 方差是 df/(df-2)=5/3，标准化到方差 1

不标准化的话，波动率和 t 分布的方差叠在一起，sigma的物理含义就变了，VaR 会莫名偏大，又是一个隐蔽的坑。

验证方法别只看单次结果，要跑重复实验：固定参数跑 200 次，统计var的均值和标准误，比较改造前后的收敛速度；再拿真实历史数据做 Kupiec POF 回测，看例外率是否落在置信区间里：

def kupiec_pof(exceptions: int, n: int, confidence: float): """返回 LR 统计量，与卡方分布 3.84(df=1, 95%) 比较""" p = 1 - confidence ph = exceptions / n return -2 * (np.log((1 - p) ** (n - exceptions) * p ** exceptions) - np.log((1 - ph) ** (n - exceptions) * ph ** exceptions))

exceptions是实际损失超过 VaR 的天数。LR 超过 3.84 说明模型在统计上不可接受，例外率偏高就是低估风险，偏低则是模型过于保守——后者不一定是好事，资本占用会虚高。这套实验我一般跑 250 个交易日，样本太短统计量没意义。

从那以后我每次改模型，都会先跑一遍对拍脚本、再跑一遍回测，两关都过才敢把新版本挂到日报上，光看单次数字顺眼就上线，迟早翻车。希望帮到你。

本文还有配套的精品资源，点击获取