2

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

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

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

我的问题是如何mquantiles()产生同样的打破惯例。我注意到这个函数有两个参数，alphap并且betap（文档调用它们alpha，beta但是你需要额外的'p'才能让它工作，至少我这样做......如果我只使用'alpha'和' beta' 与 Python 2.7.1 和 SciPy 0.10.0）。但即使在 SciPy 文档中，我也看不出这些参数的组合是否会在平坦的 CDF 范围内产生平均值。

我看到了看起来像计算为这个范围的中位数或众数的选项，但不是平均值（也不清楚这些具有 alpha 和 beta 的 SciPy 中位数/众数选项是否被计算为观测值的中位数/众数或产生平坦 CDF 值的范围。）

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

4

1 回答 1

7

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 匹配上分位数的界限大于或等于该分位数中的所有值。

于 2012-07-07T05:17:27.453 回答