采样自定义经验累积分布函数 (CDF) 的方法与实践

本文旨在详细阐述如何从用户定义的经验累积分布函数（cdf）中高效进行数据抽样。我们将探讨两种核心策略：一是基于分段常数插值的直接抽样，二是利用平滑插值（如线性或样条插值）进行更为精细的抽样。文章将结合python中的numpy和scipy库，提供清晰的示例代码，帮助读者掌握逆变换抽样原理在实际应用中的实现，从而有效生成符合特定cdf分布的随机样本。

在统计分析和模拟中，我们经常需要从一个已知的概率分布中生成随机样本。当这个分布不是标准参数分布，而是通过一系列数据点定义的经验累积分布函数（Empirical CDF）时，传统的抽样方法可能不再适用。此时，逆变换抽样（Inverse Transform Sampling）原理成为一种强大的工具，它允许我们通过均匀分布的随机数来生成符合任意给定CDF的样本。

理解逆变换抽样原理

逆变换抽样基于以下核心思想：如果 $U$ 是一个服从 $[0, 1]$ 上均匀分布的随机变量，而 $F(x)$ 是某个随机变量 $X$ 的累积分布函数，那么 $X = F^{-1}(U)$ 将服从 $F(x)$ 所定义的分布。在处理经验CDF时，我们实际上是利用插值方法来近似 $F^{-1}(U)$。

准备自定义CDF数据

首先，我们需要定义一个自定义的经验CDF。这通常由一系列 (x, cdf) 对构成，其中 x 是某个值，cdf 是该值对应的累积概率。

import pandas as pd
import numpy as np
from scipy.interpolate import interp1d

# 定义自定义CDF数据
# 'x' 代表随机变量的值
# 'cdf' 代表累积概率 F(x)
cdf_data = pd.DataFrame.from_dict(
    {'x': [10e6, 20e6, 50e6, 100e6, 250e6],
     'cdf': [0.4, 0.6, 0.7, 0.8, 1]}
)

print("自定义CDF数据：")
print(cdf_data)

方法一：直接分段常数插值抽样（无平滑）

这种方法假设CDF在给定数据点之间是分段常数或分段线性的，但在逆变换时，numpy.interp 默认采用线性插值，这在处理CDF逆函数时通常足够且直观。它不进行额外的平滑处理，直接根据CDF数据点的分布进行抽样。

原理：

1. 生成一组在 $[0, 1]$ 区间内均匀分布的随机数。这些随机数可以被视为累积概率 $U$。
2. 使用 numpy.interp 将这些均匀分布的随机数映射到CDF的 x 值。numpy.interp(x_new, xp, fp) 函数的作用是，给定 xp 和 fp（即已知点的x和y坐标），计算 x_new 对应的新y值。在这里，xp 是CDF的累积概率（cdf_data['cdf']），fp 是对应的 x 值（cdf_data['x']），x_new 则是我们生成的均匀随机数。

代码示例：

# 生成10000个在0到1之间均匀分布的随机数
num_samples = 10000
uniform_samples = np.random.uniform(0, 1, num_samples)

# 使用numpy.interp进行抽样
# xp (x-coordinates of the data points) 是 cdf 值
# fp (y-coordinates of the data points) 是 x 值
samples_method1 = np.interp(uniform_samples, cdf_data['cdf'], cdf_data['x'])

print("\n使用numpy.interp进行直接抽样（前10个样本）：")
print(samples_method1[:10])

输出示例：

[1.12425511e+08 7.90655259e+07 8.86092443e+07 ... 3.05881874e+07
 2.35671711e+08 1.00000000e+07]

这种方法简单高效，适用于大多数不需要精细平滑的场景。

法语写作助手

法语助手旗下的AI智能写作平台，支持语法、拼写自动纠错，一键改写、润色你的法语作文。

31

查看详情

方法二：平滑插值抽样（例如样条插值）

当需要CDF在数据点之间表现出更平滑的过渡时，例如，如果原始数据是连续的，或者我们希望模拟一个连续的潜在分布，可以使用 scipy.interpolate.interp1d。这个函数允许我们指定不同的插值类型（kind参数），如线性（'linear'）、二次（'quadratic'）、三次样条（'cubic'）等。

原理：

1. 首先，创建一个插值函数，该函数能够根据给定的CDF累积概率 p 返回对应的 x 值。这本质上是构建CDF的逆函数 $F^{-1}(p)$。
2. 然后，生成一组在 $[0, 1]$ 区间内均匀分布的随机数作为输入 p。
3. 将这些随机数输入到构建的插值函数中，即可得到平滑抽样后的 x 值。

代码示例：

# 构建CDF的逆函数插值器
# x轴是cdf值，y轴是x值
# kind参数可以选择插值类型，例如'linear', 'quadratic', 'cubic'
# fill_value='extrapolate' 允许在定义域外进行外推，但通常建议限制在[0,1]范围内
# 或者使用 fill_value=(cdf_data['x'].iloc[0], cdf_data['x'].iloc[-1]) 来限制边界
cdf_inverse_interpolator = interp1d(
    cdf_data['cdf'],
    cdf_data['x'],
    kind='cubic', # 可以尝试 'linear', 'quadratic', 'cubic'
    bounds_error=False,
    fill_value=(cdf_data['x'].iloc[0], cdf_data['x'].iloc[-1]) # 限制边界值
)

# 生成10000个在0到1之间均匀分布的随机数
uniform_samples_smooth = np.random.uniform(0, 1, num_samples)

# 使用插值器进行抽样
samples_method2 = cdf_inverse_interpolator(uniform_samples_smooth)

print("\n使用scipy.interpolate.interp1d进行平滑抽样（三次样条，前10个样本）：")
print(samples_method2[:10])

# 我们可以比较不同kind参数的效果
# 例如，使用线性插值
cdf_inverse_linear_interpolator = interp1d(
    cdf_data['cdf'],
    cdf_data['x'],
    kind='linear',
    bounds_error=False,
    fill_value=(cdf_data['x'].iloc[0], cdf_data['x'].iloc[-1])
)
samples_linear = cdf_inverse_linear_interpolator(uniform_samples_smooth)
print("\n使用scipy.interpolate.interp1d进行线性抽样（前10个样本）：")
print(samples_linear[:10])

不同 kind 参数的影响：

• 'linear'：在数据点之间进行线性插值，结果会比较“尖锐”，适合分段线性的近似。
• 'quadratic'：使用二次多项式插值，曲线会更平滑。
• 'cubic'：使用三次样条插值，通常能提供非常平滑的曲线，更适合模拟连续且变化平缓的分布。但需要至少4个数据点。

选择哪种 kind 取决于你对CDF平滑度的假设以及数据的性质。

注意事项

1. CDF的单调性：累积分布函数必须是非降的。在准备数据时，请确保 cdf 列的值是递增的。如果数据不满足这个条件，插值结果将不准确。
2. 边界处理：numpy.interp 和 scipy.interpolate.interp1d 在处理输入值超出 xp 或插值器定义域时有不同的行为。
   • numpy.interp 默认会使用边界值进行外推。
   • scipy.interpolate.interp1d 在 bounds_error=True 时会抛出错误，在 bounds_error=False 时会使用 fill_value。在CDF抽样中，均匀随机数总是在 [0, 1] 之间，因此如果CDF的 cdf 范围是 [0, 1]，通常不会有超出范围的问题。但设置 fill_value 可以提供更稳健的行为。
3. 插值方法的选择：
   • numpy.interp 适用于快速、简单的线性插值场景，尤其是在数据点较多且不需要高度平滑时。
   • scipy.interpolate.interp1d 提供了更灵活的插值类型选择，适用于需要更精细控制平滑度的场景。
4. 抽样数量：生成的样本数量越多，样本分布越能逼近原始CDF。在实际应用中，通常需要生成数千甚至数万个样本才能获得具有统计意义的结果。

总结

从自定义经验CDF中抽样是数据科学和模拟中的一项基本任务。通过逆变换抽样原理，结合Python中的 numpy.interp 和 scipy.interpolate.interp1d，我们可以灵活地实现这一目标。numpy.interp 提供了一种快速直接的线性插值方法，适用于大多数场景；而 scipy.interpolate.interp1d 则通过其 kind 参数，允许用户选择不同程度的平滑插值，如线性、二次或三次样条，以更好地拟合潜在的连续分布。根据具体的应用需求和对分布平滑度的假设，选择合适的插值策略至关重要。

以上就是采样自定义经验累积分布函数 (CDF) 的方法与实践的详细内容，更多请关注php中文网其它相关文章！