python - 获取 SciPy 分位数以匹配 Stata xtile 函数

Question

我继承了一些旧的 Stata 代码 (Stata11)，它使用 xtile 函数按分位数对向量中的观测值进行分类(在这种情况下，只是标准的 5 个五分位数、20%、40 %、60%、80%、100%)。

我正在尝试用 Python 复制一段代码，我正在使用 SciPy.stats.mstats 函数 mquantiles() 进行计算。

据我从 Stata 文档和在线搜索中得知，Stata xtile 方法试图反转数据的经验 CDF，并使用所有观测值的等加权平均值CDF 是平坦的，可以制作切割点。这似乎是一种非常糟糕的分位数分类方法，但事实就是如此，我相信在某些情况下这是正确的做法。

我的问题是如何使 mquantiles() 产生相同类型的破坏约定。我注意到这个函数有两个参数，alphap 和 betap(文档称它们为 alpha 和 beta 但你需要额外的'p'才能让它工作，至少我这样做......如果我只使用'alpha'和'beta'与Python 2.7.1和SciPy 0.10.0，我会收到错误消息)。但即使在 SciPy 文档中，我也看不到这些参数的组合是否会产生平坦 CDF 范围内的平均值。

我看到了计算为该范围的中值或模式的选项，但不是意思(也不清楚这些带有 alpha 和 beta 的 SciPy 中值/模式选项是否被计算为 < em>观察值 或产生平坦 CDF 值的范围。)

任何帮助消除这些不同选项的歧义并找到一些帮助我在 Python 中重新创建 Stata 约定的文档都会很棒。请避免只说“编写您自己的分位数函数”的答案。首先，这无助于我理解 Stata 或 SciPy 的约定，其次，鉴于这些数值库，编写我自己的分位数函数应该是最后的选择。我当然可以做到，但如果我需要的话，这对周围的一切都会很糟糕。

Answer 1

scipy.stats.mquantiles 文档很差，而且在某些地方错误，现在已修复，这样可能会有帮助... http://docs.scipy.org/scipy/docs/scipy.stats.mstats_basic.mquantiles/ .当您指出 alpha/beta、alphap/betap 差异时，该过程就开始了。谢谢。

mquantiles的实现遵循R。

最大的区别在于 R 有 9 种离散类型，因为 scipy.stats.mquantiles 从 'alphap' 和 'betap' 计算 'm'，scipy 有连续范围的“类型”(因为缺少更好的词)。

我承认我不了解所涉及统计数据的所有来龙去脉，因此我决定进行强力评估。我在 http://www.biostat.sdu.dk/~biostat/StataReferenceManual/StataRef.pdf 找到了一个 xtile 示例并能够将结果与 alphap=0.5 和 betap=0.5(分段线性)相匹配。不是确定的也不是详尽无遗的，但我现在拥有的一切。

In [1]: import scipy.stats as st

In [9]: st.mstats.mquantiles([23,56,67,123,99,17],prob=[0.5],alphap=0.5,betap=.5)
Out[9]: array([ 61.5])

In [10]: st.mstats.mquantiles([23,56,67,123,99,17],prob=[0.33,0.66],alphap=0.5,betap=.5)
Out[10]: array([ 38.84,  81.72])

In [11]: st.mstats.mquantiles([23,56,67,123,99,17],prob=[0.25,0.5,0.75],alphap=0.5,betap=.5)
Out[11]: array([ 23. ,  61.5,  99. ])

最后一个有点问题，因为两个划分点恰好位于数据集中的值上。 Stata/xtile(至少在我发现的例子中)没有给出分位数的分割点，而是给出了分位数本身。给定排序后的数据集 [17,23,56,67,99,123]，Stata/xtile 给出的分类为 [1,1,2,3,3,4] 这意味着对于 scipy.stat.mquantiles 匹配上分位数的边界大于或等于该分位数中的所有值。