使用NumPy基于已知分布判断新观测所属状态及反查分位数的方法

你不需要自己手写暴力遍历或者二分查找逻辑，NumPy提供了完全适配这个需求的底层高性能实现：numpy.searchsorted。

核心逻辑 #

经验分位数的本质是：对于给定的样本分布，某一观测值对应的分位值，等于样本中小于/不大于该观测值的样本数量占总样本数的比例。而numpy.searchsorted的原生功能，就是在已经排序的一维数组中，通过C语言实现的二分查找，返回目标值插入数组时的合法位置索引——这个索引值和总样本数的比值，就是你要计算的经验分位值。
这个实现的时间复杂度是严格的O(logn)，相比Python层面手写的二分查找、暴力遍历，15000样本规模下单次计算的耗时可以从毫秒级降到微秒级，也是NumPy社区处理这类问题的标准写法。

具体使用步骤 #

• 一次性完成预排序：因为searchsorted要求输入的样本数组是已排序的，所以在第一次计算前，把S1、S2的样本数组分别排序即可，后续所有新观测值的计算都不需要重复排序：

  # 假设s1_samples、s2_samples是你存储的两个NumPy数组
  s1_sorted = np.sort(s1_samples)
  s2_sorted = np.sort(s2_samples)
  

• 单观测值分位计算：拿到新观测值x_new后，直接调用函数计算两个分布下的对应分位：

  # 通过side参数匹配你之前使用numpy.quantile的分位对齐规则
  # side='left'表示：将目标值插入到所有等于它的元素的最左侧位置
  rank_s1 = np.searchsorted(s1_sorted, x_new, side='left')
  quantile_s1 = rank_s1 / len(s1_sorted)
  
  rank_s2 = np.searchsorted(s2_sorted, x_new, side='left')
  quantile_s2 = rank_s2 / len(s2_sorted)
  

• 批量观测值计算：如果后续有大量新观测值需要判断归属，不需要写循环，直接把存储新观测值的数组作为输入传给searchsorted，就可以一次性返回所有观测值对应的分位结果，完全兼容NumPy的向量化操作习惯。

对齐numpy.quantile的规则说明 #

如果你之前已经习惯了numpy.quantile的分位计算逻辑，只需要根据需要调整searchsorted的side参数即可：

• side='left'的计算结果，和分位定义中“观测值严格小于当前值的样本占比”的规则完全一致
• side='right'的计算结果，和分位定义中“观测值小于等于当前值的样本占比”的规则完全一致

调整参数后，计算得到的分位值和原生numpy.quantile的计算逻辑是完全自洽的，不会出现规则偏差。

归属判断的注意点 #

因为两个分布存在重叠，拿到两个分位值后，不需要套用经典分布的先验假设，只需要根据对两个状态的业务认知设定判断阈值即可：比如当某一观测值在S1分布中处于90%分位、在S2分布中处于20%分位时，就可以将其判定为更有可能归属S1状态。

注意：如果两个分布的样本量有差异，不要直接对比分位的绝对值，可以结合你之前通过numpy.quantile计算得到的两个分布的分位区间边界来设定判断规则，避免样本量偏差带来的判断误差。

内容的提问来源于stack exchange，提问作者johnnyp