5

我正在尝试在 python 中复制一些 Matlab 代码。我找不到与 Matlab 函数完全等价的函数quantile。我发现最接近的是 python 的mquantiles.

Matlab 示例:

 quantile( [ 8.60789925e-05, 1.98989354e-05 , 1.68308882e-04, 1.69379370e-04],  0.8)

...给出:0.00016958

python中的相同示例:

scipy.stats.mstats.mquantiles( [8.60789925e-05, 1.98989354e-05, 1.68308882e-04, 1.69379370e-04], 0.8)

...给0.00016912

有谁知道如何精确复制 Matlab 的quantile功能?

4

3 回答 3

5

(在 More About => Algorithms 部分下)的文档给出了使用的确切算法。quantile这是一些 python 代码,它为平面数组的单个分位数执行此操作,使用瓶颈进行部分排序:

import numpy as np
import botteleneck as bn

def quantile(a, prob):
    """
    Estimates the prob'th quantile of the values in a data array.

    Uses the algorithm of matlab's quantile(), namely:
        - Remove any nan values
        - Take the sorted data as the (.5/n), (1.5/n), ..., (1-.5/n) quantiles.
        - Use linear interpolation for values between (.5/n) and (1 - .5/n).
        - Use the minimum or maximum for quantiles outside that range.

    See also: scipy.stats.mstats.mquantiles
    """
    a = np.asanyarray(a)
    a = a[np.logical_not(np.isnan(a))].ravel()
    n = a.size

    if prob >= 1 - .5/n:
        return a.max()
    elif prob <= .5 / n:
        return a.min()

    # find the two bounds we're interpreting between:
    # that is, find i such that (i+.5) / n <= prob <= (i+1.5)/n
    t = n * prob - .5
    i = np.floor(t)

    # partial sort so that the ith element is at position i, with bigger ones
    # to the right and smaller to the left
    a = bn.partsort(a, i)

    if i == t: # did we luck out and get an integer index?
        return a[i]
    else:
        # we'll linearly interpolate between this and the next index
        smaller = a[i]
        larger = a[i+1:].min()
        if np.isinf(smaller):
            return smaller # avoid inf - inf
        return smaller + (larger - smaller) * (t - i)

我只做了单分位数的一维案例,因为这就是我所需要的。如果您想要几个分位数,则可能值得进行完整排序;要按轴执行并且知道您没有任何 nan,您需要做的就是在排序中添加一个轴参数并将线性插值位矢量化。使用 nans 在每个轴上进行操作会有点棘手。

这段代码给出:

>>> quantile([ 8.60789925e-05, 1.98989354e-05 , 1.68308882e-04, 1.69379370e-04], 0.8)
0.00016905822360000001

matlab代码给出了0.00016905822359999999;不同的是3e-20。(小于机器精度)

于 2012-12-30T22:13:14.733 回答
4

您的输入向量只有 4 个值,这太少了,无法很好地近似底层分布的分位数。这种差异可能是 Matlab 和 SciPy 使用不同的启发式方法来计算欠采样分布的分位数的结果。

于 2012-12-16T17:47:41.363 回答
3

有点晚了,但是:

mquantiles 非常灵活。您只需要提供 alphap 和 betap 参数。在这里,由于 MATLAB 做了线性插值,所以需要将参数设置为 (0.5,0.5)。

In [9]: scipy.stats.mstats.mquantiles( [8.60789925e-05, 1.98989354e-05, 1.68308882e-04, 1.69379370e-04], 0.8, alphap=0.5, betap=0.5)

编辑:MATLAB 说它进行线性插值,但它似乎通过分段线性插值计算分位数,这相当于R中的 Type 5 分位数和 scipy 中的 (0.5, 0.5) 。

于 2015-04-20T12:45:29.617 回答