
CuPy 统计函数参考指南顺序统计量、均值方差、相关性与直方图的 GPU 实现与用法【免费下载链接】cupyNumPy SciPy for GPU项目地址: https://gitcode.com/GitHub_Trending/cu/cupyCuPy 是一套面向 GPU 的 NumPy/SciPy 兼容数值计算库本文以其 API 参考文档 statistics.rst 为骨架系统讲解 CuPy 统计子模块cupy._statistics提供的四大类函数顺序统计量、均值与方差、相关性分析、直方图统计。读完本文你将掌握每个函数的参数语义、与 NumPy 的兼容性边界、设备同步注意事项以及它们背后的 CUDA Kernel 与 CUB 加速实现原理能够直接在自己的 GPU 计算任务中正确选用并写出高性能统计代码。一、统计模块在 CuPy 中的定位CuPy 的公共 API 设计与 NumPy 保持一致。在 statistics.rst 中全部统计函数被划分为四个分组与 NumPy 官方文档的routines.statistics章节一一对应Order statistics顺序统计量ptp、percentile、quantileAverages and variances均值与方差median、average、mean、std、var、nanmedian、nanmean、nanstd、nanvarCorrelating相关性corrcoef、correlate、covHistograms直方图histogram、histogram2d、histogramdd、bincount、digitize从源码结构看这些函数分散在 cupy/_statistics 目录的四个文件中order.py顺序统计量、meanvar.py均值方差、correlation.py相关性、histogram.py直方图并在 cupy/init.py 中被统一导出到cupy命名空间顶层因此用法与 NumPy 完全一致cupy.mean(x)、cupy.histogram(x)等。此外amin/amax/nanmin/nanmax也在该模块中实现并作为cupy.min/cupy.max导出。需要特别留意的是statistics.rst 中以#注释掉的条目nanpercentile、nanquantile、histogram_bin_edges表示当前版本尚未提供这些 API使用前应通过hasattr(cupy, nanpercentile)等方式确认避免导入错误。二、顺序统计量ptp、percentile 与 quantile2.1 ptp峰值到峰值的取值范围ptppeak-to-peak返回数组沿指定轴的最大值与最小值之差import cupy as cp a cp.array([[4, 9, 2], [10, 5, 13], [1, 14, 7]]) cp.ptp(a) # 13整个数组的 max - min cp.ptp(a, axis0) # array([ 9, 9, 11]) cp.ptp(a, axis1) # array([7, 8, 13])参数axis缺省时作用于展平后的数组out可指定输出数组keepdimsTrue时被归约的轴保留为长度 1 的轴。从 order.py 的源码看ptp直接委托给ndarray.ptp方法。注意当某个归约切片内含有 NaN 时对应的ptp结果也是 NaN文档同时注明若启用 cuTENSOR 加速器含 NaN 的归约轴输出值可能被折叠collapsed这是加速路径与逐元素路径的已知差异。2.2 percentile 与 quantile分位数计算的完整参数percentile计算 q-th 百分位数q 取值范围 0100quantile计算 q-th 分位数q 取值范围 01。两者共享同一个底层实现_quantile_unchecked区别仅在于百分位数先把q除以 100再校验取值区间cp.percentile(a, 50) # 中位数 cp.percentile(a, [25, 50, 75]) # 四分位数q 支持 tuple/list/ndarray cp.quantile(a, 0.5) cp.quantile(a, 0.5, axis0) cp.quantile(a, 0.5, methodhigher, keepdimsTrue)核心参数如下参数说明q分位点集合percentile要求在 [0, 100]quantile要求在 [0, 1]越界抛出ValueErroraxis沿哪条轴可传 int 或 tuple计算默认展平out输出数组overwrite_input为True时允许中间计算原地修改输入a以节省内存函数返回后输入内容不可预期method分位点落在两个数据点之间时的插值方法默认linearkeepdims为True时保留被归约轴为长度 1interpolationmethod的已废弃旧名传入时会发出DeprecationWarning并自动映射method 插值方法的完整集合从 order.py 的分发逻辑看支持linear默认、lower、higher、midpoint、nearest、inverted_cdf以及基于 Hyndman Fan (1996) 连续插值法的hazenHF type 5、weibulltype 6、median_unbiasedtype 8、normal_unbiasedtype 9。后四者通过源码顶部的_QUANTILE_PARAMS字典映射为对应的(alpha, beta)参数对_QUANTILE_PARAMS { hazen: (0.5, 0.5), weibull: (0, 0), median_unbiased: (1/3, 1/3), normal_unbiased: (3/8, 3/8), }尚未实现的方法averaged_inverted_cdf、closest_observation、interpolated_inverted_cdfNumPy 1.22 新增当前会直接抛出ValueError测试文件 test_order.py 中也将这三个方法注释为TODO(takagi) Not implemented。传入未知方法名同样抛ValueError。GPU 上的实现方式_quantile_unchecked先把归约轴上的数据排序ap.sort(axisaxis)将归约轴搬到最后一维并展平然后按插值方法把q换算成虚拟下标indices当需要线性插值时会现场编译一个名为cupy_percentile_weightnening的ElementwiseKernel见 order.py对排序后的数组按idx_below floor(idx)、weight_above idx - idx_below做加权插值。因此分位数计算在 GPU 上等价于排序 插值核对大型数组有良好吞吐。三、均值与方差mean / std / var / average / median 及 NaN 安全变体3.1 mean / std / var直接的 ndarray 方法委托mean、std、var在 meanvar.py 中直接委托给ndarray的同名方法共享以下参数语义axisint、int 序列或None默认对展平数组计算dtype指定计算与输出的数据类型如用float32累积out输出数组keepdims保留被归约轴std/var额外支持ddofdelta degrees of freedom默认为 0除以 N设为 1 时除以N-1得到无偏样本方差。x cp.arange(6, dtypecp.float32).reshape(2, 3) cp.mean(x) # 2.5 cp.mean(x, axis1) # array([1., 4.]) cp.var(x, ddof1) cp.std(x, axis0, keepdimsTrue)3.2 average支持权重的加权平均average在 meanvar.py 中实现返回沿轴的加权平均cp.average(cp.array([1, 2, 3, 4]), weightscp.array([4, 3, 2, 1])) # 2.0 data cp.arange(6).reshape((3, 2)) cp.average(data, axis1, weightscp.array([1./4, 3./4]))关键行为weightsNone时退化为mean权重数组与a形状不同时要求axis必须显式指定且权重为一维数组其长度须等于该轴长度所有权重之和为 0 时抛出ZeroDivisionError整数或布尔输入会通过numpy.promote_types提升到至少float64的result_dtypereturnedTrue时返回(avg, scl)元组其中scl是权重和注意提供weights时该函数可能触发设备同步源码中以cupy.any(scl 0.0)注释# synchronize!标出混合 CPU/GPU 流水线中需留意。3.3 median 与 nanmedianmedian委托给 Cython 层的_statistics._median见 meanvar.py支持axis含多轴、out、overwrite_input、keepdims。nanmedian有一个值得注意的分流当输入 dtype 属于浮点或复数efdFD时走专门的_nanmedian实现否则直接复用median整数数组不含 NaN无需特殊处理。cp.median(cp.array([[10, 7, 4], [3, 2, 1]])) # 3.5 cp.median(cp.array([[10, 7, 4], [3, 2, 1]]), axis0) # array([6.5, 4.5, 2.5]) cp.nanmedian(cp.array([1.0, cp.nan, 3.0])) # 1.03.4 NaN 安全变体nanmean / nanvar / nanstd / nanmin / nanmax这一组函数忽略 NaN 值参与计算。NaN 安全归约的核心实现位于 Cython 文件 cupy/_core/_routines_statistics.pyx例如_nanmean_func通过create_reduction_func生成归约核_nanvar_core使用ReductionKernel其归约逻辑在模板函数nanvar_impl中实现——对 NaN 元素贡献 0对有效元素贡献(x - mean) * (x - mean)最后除以max(_count - ddof, 0LL)见 _routines_statistics.pyx 附近的实现并提供了复数版本的专用内核_nanvar_core_complex64/complex128。Python 层的分流规则meanvar.pynanmean/nanvar/nanstd对整数与布尔 dtypedtype.kind in biu直接退化为普通mean/var/std避免无谓的 NaN 检查nanmin/nanmaxorder.py先调用 Cython 层归约再通过content.isnan(res).any()检测是否存在全 NaN 切片若有则发出RuntimeWarning(All-NaN slice encountered)并返回 NaN。cp.nanmean(cp.array([[1, cp.nan], [3, 4]])) # 2.666... cp.nanstd(cp.array([1.0, cp.nan, 2.0, 4.0]), ddof1) cp.nanmin(cp.array([1.0, cp.nan, 3.0])) # 1.0设备同步警告nanmin/nanmax以及带权重的average、histogram系列文档明确标注 This function may synchronize the device因为全 NaN 切片检测需要把结果取回主机判断。在追求极致吞吐的循环中应评估其影响。四、相关性分析corrcoef / cov / correlate4.1 cov协方差矩阵cov在 correlation.py 中实现计算协方差矩阵x cp.array([[0, 2], [1, 1], [2, 0]]).T cp.cov(x) # 2x2 协方差矩阵 cp.cov(x, rowvarFalse)参数说明y额外的变量观测集合与a按行拼接rowvarTrue默认表示每行是一个变量、每列是一次观测False则转置解释biasFalse时按(N-1)归一化无偏估计biasTrue时按N归一化ddof非None时覆盖bias隐含的默认值ddof1无偏、ddof0简单平均fweights整数频率权重与aweights观测向量权重均要求为cupy.ndarray且长度须等于观测数两者同时给出时按w * aweights相乘dtype缺省时结果至少为float64精度通过numpy.promote_types与float64提升。实现上cov先做均值中心化X - X.mean(axis1)[:, None]再通过矩阵乘法X.dot(X_T.conj()) / fact得到协方差矩阵其中fact由ddof/bias与权重共同决定当自由度fact 0时发出RuntimeWarning并置 0。注意corrcoef中的bias与ddof参数已废弃传入会触发DeprecationWarning且不生效。4.2 corrcoef皮尔逊相关系数corrcoef在 correlation.py 中基于cov实现先求协方差矩阵取其对角线的实部开方得到标准差再对矩阵逐项归一化out / stddev[:, None]; out / stddev[None, :]最后把实部以及复数的虚部clip 到[-1, 1]区间cp.corrcoef(x) cp.corrcoef(x, y) # y 为额外变量集4.3 correlate一维互相关correlate计算两个一维序列的离散互相关mode支持valid默认、same、full。实现上先通过_choose_conv_method在直接卷积与 FFT 卷积之间做选择然后分别调用_dot_convolve或_fft_convolve核心等价于convolve(a, v[::-1])correlation.py。输入为空或非一维数组时会抛出ValueError。cp.correlate(cp.array([1, 2, 3]), cp.array([0, 1, 0.5]), modefull)五、直方图histogram 系列、bincount 与 digitize5.1 histogram一维直方图histogram在 histogram.py 中实现返回(hist, bin_edges)元组cp.histogram(cp.arange(10), bins5) cp.histogram(x, binscp.array([0., 2., 4., 6., 8., 10.])) # 显式 bin 边界 cp.histogram(x, range(0, 10), densityTrue) cp.histogram(x, bins10, weightsw) # 加权直方图参数语义与 NumPy 对齐bins整数等宽分箱个数须为正数或一维数组单调递增的 bin 边界字符串形式的 bin 算法如auto、fd当前抛出NotImplementedError源码注释only integer and array bins are implementedrange(min, max)二元组缺省为(x.min(), x.max())边界外取值被忽略且range同时影响自动分箱宽度densityTrue返回概率密度bin_count / sample_count / bin_volumeweights与x同形状的权重数组支持可安全转换到 float/complex 的 dtype如Decimal等对象 dtype 不支持复数输入抛出NotImplementedError布尔输入会先告警并转为uint8。GPU 实现亮点直方图的核心是二分查找分箱 原子累加。源码中预编译了两个ElementwiseKernel_histogram_kernel与_weighted_histogram_kernelhistogram.py每个元素在 GPU 上通过 while 循环二分定位 bin再用atomicAdd累加计数。当启用 CUB 加速器且数据量不超过0x7fffffff时优先调用cub.cub_histogramcupy/cuda/cub.pyx 提供绑定并针对 NumPy 的最后一 bin 右闭语义做修正整数 bin 上界 1浮点 bin 上界用cupy.nextafter外推一位HIP 平台上计数临时用uint64。CUB 路径失败时自动回退到 CuPy 内核。该函数与histogramdd/histogram2d均在文档中标注 may synchronize the device。5.2 histogram2d 与 histogramdd多维直方图histogram2d(x, y, bins10, rangeNone, weightsNone, densityNone)二维直方图返回(H, xedges, yedges)。bins可为单个整数两个维度共用、二元组各维分箱数或一维数组序列各维边界bins为cupy.ndarray时表示两维共用同一组边界。histogramdd(sample, bins10, rangeNone, weightsNone, densityFalse)D 维直方图返回(H, edges)。sample为(N, D)数组每行一个 D 维坐标点或(X, Y, Z, ...)形式的坐标序列内部用cupy.stack组装bins可为标量、长度为 D 的序列或各维边界数组序列。histogramdd的实现流程histogram.py逐维用linspace生成等宽边界或用给定数组作边界 → 用cupy.searchsorted(edges[i], sample[:, i], sideright)计算每个样本落入的 bin 编号源码注明刻意避开cupy.digitize以规避 NumPy gh-11022 问题→ 对恰好落在右边界上的样本编号 -1 修正 →cupy.ravel_multi_index展平为扁平索引 →cupy.bincount(xy, weights, minlength...)统计 → reshape 后裁掉每个维度的离群 bin首尾各一densityTrue时再除以各维 bin 宽度与总样本数。5.3 bincount非负整数计数bincount(x, weightsNone, minlengthNone)统计非负整数数组每个取值的出现次数输出长度等于max(cupy.max(x) 1, minlength)cp.bincount(cp.array([0, 1, 1, 3, 2, 1, 7])) # array([1, 3, 1, 1, 0, 0, 0, 1]) cp.bincount(cp.array([0, 1, 1, 3]), weightscp.array([0.5, 0.25, 0.75, 1.0]))约束x必须是一维非负整数数组浮点 dtype 抛TypeError出现负数抛ValueErrorweights形状须与x相同minlength须非负。计数内核_bincount_kernel同样使用atomicAdd无权重时优先走 CUB 路径HIP 平台除外有权重时使用_bincount_with_weight_kernel并以float64累积。空输入返回numpy.intp类型的零数组。该函数会触发设备同步需要取回x.max()确定输出长度。5.4 digitize样本所属 bin 的下标digitize(x, bins, rightFalse)返回x中每个值所属区间的下标结果形状与x相同bins须为一维单调数组x cp.array([0.2, 6.4, 3.0, 1.6]) bins cp.array([0.0, 1.0, 2.5, 4.0, 10.0]) cp.digitize(x, bins) # array([1, 4, 3, 2])rightTrue时区间左开右闭即bins[i-1] x bins[i]否则左闭右开bins[i-1] x bins[i]。实现上直接转发给排序模块的_searchsortedhistogram.py。与 NumPy 的一个刻意差异在文档注释中注明为避免设备同步digitize不会对bins的单调性做校验非单调数组的结果是未定义的。六、性能相关加速器CUB与设备同步清单统计模块的性能与行为受两个机制影响1. 例程加速器routine accelerators。histogram与bincount会遍历_accelerator.get_routine_accelerators()由CUPY_ACCELERATORS环境变量或cupy._core.set_routine_accelerators配置优先调用 CUB 的DeviceHistogram绑定 cub_histogram失败或未启用时回退到 CuPy 自研atomicAdd内核。CUB 路径对输入规模有 0x7fffffff的限制源码TODO(leofang): support 2^31 elements in x?。2. 设备同步点。以下函数在特定条件下会把设备数据取回主机或执行同步密集调用时需权衡函数同步触发条件nanmin/nanmax需要检测是否存在全 NaN 切片isnan(...).any()average提供weights且需要校验权重和是否为 0histogram/histogram2d/histogramdd/bincount需要从设备取回范围、最大值或 bin 边界等标量如x.min()、x.max()、int(bins)percentile/quantileq为cupy.ndarray时校验取值区间需要同步七、如何查阅与验证文档、源码与测试的对应关系如果你需要深入某个函数的边界行为可以按参考文档 → Python 封装 → Cython 内核 → 测试的链路查阅API 参考docs/source/reference/statistics.rst 列出全部函数及其分组并链接到各函数生成的独立页面Python 封装层cupy/_statistics 下的order.py、meanvar.py、correlation.py、histogram.py包含完整的 docstring、参数校验与 NumPy 兼容逻辑Cython 归约内核cupy/_core/_routines_statistics.pyx 定义了cupy_min/cupy_max归约函数以及nanmean/nanvar/nanstd等 NaN 安全归约内核测试套件tests/cupy_tests/statistics_tests 下的test_order.py、test_meanvar.py、test_correlation.py、test_histogram.py覆盖了各函数的 dtype 组合、axis/keepdims/method组合与异常路径是验证预期行为的第一手资料。例如 test_order.py 中列出了percentile/quantile全部受测的method取值可作为哪些插值方法可用的权威清单。八、使用建议与已知限制与 NumPy 的双向互操作以上函数均接受并返回cupy.ndarray可与 NumPy 数组通过cupy.asarray/numpy.asarray互转多数函数保持与 NumPy 相同的 dtype 提升规则如cov/corrcoef至少float64精度。已知未实现项nanpercentile、nanquantile、histogram_bin_edges文档中已注释histogram的字符串 bin 算法percentile/quantile的averaged_inverted_cdf、closest_observation、interpolated_inverted_cdf三种 methodcov的对象 dtype 权重与负数权重不报错的性能优化行为也与 NumPy 略有差异。同步代价意识在 GPU 流水线中尽量一次性提交大批量数据避免在小数组上频繁调用会同步的统计函数以免多次 kernel 启动与设备同步吃掉并行收益。总而言之CuPy 的统计模块在 API 层面完整对齐 NumPy在实现层面则针对 GPU 特性做了排序 插值内核、原子累加直方图、CUB 加速回退与复数/NaN 专用归约内核等深度优化。将 statistics.rst 与上述源码、测试配合使用你就能在 GPU 上写出既符合 NumPy 习惯又充分释放硬件性能的统计计算代码。【免费下载链接】cupyNumPy SciPy for GPU项目地址: https://gitcode.com/GitHub_Trending/cu/cupy创作声明:本文部分内容由AI辅助生成(AIGC),仅供参考