本文旨在解决 scipy 自定义连续随机变量中，昂贵常数（如 pdf 归一化常数和 cdf 积分常数）重复计算导致的性能问题。通过引入类级别的本地缓存机制，使用字典存储已计算的常数值，并以参数元组作为键，显著减少了重复计算，从而提升了自定义分布的评估效率。
在 SciPy 中定义自定义连续随机变量时，通常需要继承 scipy.stats.rv_continuous 类并实现 _pdf 和 _cdf 等核心方法。这些方法在计算概率密度函数和累积分布函数时，往往依赖于一些昂贵的、与分布参数相关的常数，例如 PDF 的归一化常数和 CDF 的积分常数。如果这些常数在每次评估 _pdf 或 _cdf 时都被重新计算，将会导致显著的性能瓶颈，尤其是在进行大量采样或统计分析时。
考虑一个自定义分布 Example_gen，其 _pdf 和 _cdf 方法依赖于两个昂贵的常数计算函数 _norm(a, b) 和 _C(a, b)：
from scipy.stats import rv_continuous
# 假设 N(a, b) 和 C(a, b) 是昂贵的常数计算函数
def N(a, b):
"""模拟昂贵的归一化常数计算"""
# 实际应用中可能涉及数值积分或其他复杂计算
import time
time.sleep(0.01) # 模拟耗时操作
return a + b + 1.0
def C(a, b):
"""模拟昂贵的积分常数计算"""
# 实际应用中可能涉及数值积分或其他复杂计算
import time
time.sleep(0.01) # 模拟耗时操作
return a - b + 0.5
# 假设 f(x, a, b) 是非归一化的PDF，F(x, a, b) 是其原函数
def f(x, a, b):
return x * a + b
def F(x, a, b):
return 0.5 * x**2 * a + b * x
class Example_gen(rv_continuous):
def _norm(self, a, b):
"""昂贵的归一化常数计算函数"""
return N(a, b)
def _C(self, a, b):
"""昂贵的积分常数计算函数"""
return C(a, b)
def _pdf(self, x, a, b):
return f(x, a, b) / self._norm(a, b)
def _cdf(self, x, a, b):
return (F(x, a, b) + self._C(a, b)) / self._norm(a, b)
Example = Example_gen()
# 示例：多次调用会重复计算 _norm 和 _C
# frozen_dist = Example(a=1, b=2)
# frozen_dist.pdf(0.5)
# frozen_dist.cdf(0.5)解决方案：本地缓存策略
为了避免重复计算这些昂贵的常数，我们可以采用本地缓存的策略。具体来说，可以在 Example_gen 类中定义类级别的字典来存储已经计算过的常数值。当需要某个常数时，首先检查缓存中是否存在对应参数的计算结果；如果存在，则直接返回缓存值；否则，执行昂贵的计算并将结果存入缓存，然后返回。
from scipy.stats import rv_continuous
import math
# 假设 N(a, b) 和 C(a, b) 保持不变，仍是昂贵的计算函数
# ... (N, C, f, F 函数定义同上) ...
class Example_gen(rv_continuous):
_n_cache = {} # 类级别的归一化常数缓存
_C_cache = {} # 类级别的积分常数缓存
def _norm(self, a, b):
"""昂贵的归一化常数计算函数，带有缓存机制"""
# 使用参数元组作为缓存键，对浮点数进行适当的四舍五入以避免精度问题
key = (round(a, 5), round(b, 5))
v = self._n_cache.get(key)
if v is None:
v = N(a, b) # 执行昂贵的计算
self._n_cache[key] = v
return v
def _C(self, a, b):
"""昂贵的积分常数计算函数，带有缓存机制"""
key = (round(a, 5), round(b, 5))
v = self._C_cache.get(key)
if v is None:
v = C(a, b) # 执行昂贵的计算
self._C_cache[key] = v
return v
def _pdf(self, x, a, b):
return f(x, a, b) / self._norm(a, b)
def _cdf(self, x, a, b):
return (F(x, a, b) + self._C(a, b)) / self._norm(a, b)
Example = Example_gen()
# 示例：使用缓存后的性能提升
# 第一次调用会计算并缓存常数，后续相同参数的调用将直接从缓存中获取
frozen_dist_1 = Example(a=1, b=2)
print("第一次调用 (a=1, b=2):")
import time
start_time = time.time()
frozen_dist_1.pdf(0.5)
frozen_dist_1.cdf(0.5)
print(f"耗时: {time.time() - start_time:.4f} 秒")
print("\n第二次调用 (a=1, b=2) - 应该更快:")
start_time = time.time()
frozen_dist_1.pdf(0.5)
frozen_dist_1.cdf(0.5)
print(f"耗时: {time.time() - start_time:.4f} 秒")
print("\n调用不同参数 (a=3, b=4) - 应该再次计算:")
frozen_dist_2 = Example(a=3, b=4)
start_time = time.time()
frozen_dist_2.pdf(0.5)
frozen_dist_2.cdf(0.5)
print(f"耗时: {time.time() - start_time:.4f} 秒")
print("\n再次调用 (a=3, b=4) - 应该更快:")
start_time = time.time()
frozen_dist_2.pdf(0.5)
frozen_dist_2.cdf(0.5)
print(f"耗时: {time.time() - start_time:.4f} 秒")注意事项与最佳实践
-
缓存键的生成：
- 浮点数精度： 由于浮点数运算的精度问题，直接使用浮点数元组作为字典键可能会导致相同逻辑值的参数被视为不同的键。例如 (1.0, 2.0) 和 (1.0000000000000001, 2.0) 可能被视为不同的键。因此，在生成缓存键时，对浮点数参数进行适当的四舍五入（如 round(a, 5)）是至关重要的，以确保具有相同有效数字的参数能够命中缓存。选择合适的舍入精度取决于实际应用中参数的精度要求。
- 参数顺序： 确保缓存键中参数的顺序始终一致，因为 (a, b) 和 (b, a) 是不同的键。
- 参数类型： 缓存键必须是不可变的（immutable），因此元组是理想的选择。
-
缓存的持久化：
在某些场景下，如果昂贵常数的计算结果需要在不同的程序运行会话之间保持，可以将缓存字典的内容序列化到文件（如 JSON 或 pickle）中。在程序启动时加载这些文件来初始化缓存，并在程序结束时将更新后的缓存写回文件。
-
例如，在类定义之外或类的 __init__ 方法中添加加载和保存逻辑：
import json # ... class Example_gen(rv_continuous): _n_cache = {} _C_cache = {} # 尝试从文件加载缓存 try: with open('n_cache.json', 'r') as f: _n_cache.update({eval(k): v for k, v in json.load(f).items()}) with open('C_cache.json', 'r') as f: _C_cache.update({eval(k): v for k, v in json.load(f).items()}) except FileNotFoundError: pass # 文件不存在，缓存为空 # ... (_norm, _C, _pdf, _cdf 方法) ... # 在程序退出前保存缓存 # import atexit # def save_caches(): # with open('n_cache.json', 'w') as f: # json.dump({str(k): v for k, v in Example_gen._n_cache.items()}, f) # with open('C_cache.json', 'w') as f: # json.dump({str(k): v for k, v in Example_gen._C_cache.items()}, f) # atexit.register(save_caches)请注意，使用 eval(k) 从字符串键转换回元组时需谨慎，确保键的来源是可信的。对于更复杂的数据结构，pickle 模块可能更合适。
-
缓存管理：
- 对于大多数自定义分布而言，常数计算的参数集合是有限且固定的，因此简单的字典缓存通常足够。
- 如果参数空间非常大，或者需要限制缓存的内存占用，可以考虑使用 functools.lru_cache 装饰器。然而，lru_cache 是基于函数调用的，并且默认是实例级别的（如果装饰的是实例方法），如果需要类级别的共享缓存，则需要将其应用于静态方法或类方法，并确保缓存键包含了所有相关参数。对于本例，类级别的字典更直接地实现了跨实例的常数共享。
总结
通过在 scipy.stats.rv_continuous 的子类中实现本地缓存机制，我们可以有效地预计算并存储那些昂贵的、依赖于分布参数的常数。这种方法显著减少了重复计算，从而大幅提升了自定义随机变量在进行 PDF、CDF 或其他统计函数评估时的性能。正确处理浮点数精度和缓存键的生成是确保缓存机制有效运行的关键。对于需要跨会话持久化缓存的场景，可以结合文件存储技术来进一步优化。
