Python:使用 pandas 进行加权中位数算法

23

我有一个数据框长这个样子:

Out[14]:
    impwealth  indweight
16     180000     34.200
21     384000     37.800
26     342000     39.715
30    1154000     44.375
31     421300     44.375
32    1210000     45.295
33    1062500     45.295
34    1878000     46.653
35     876000     46.653
36     925000     53.476

我想要使用列impwealth中的频率权重来计算加权中位数,我的伪代码如下:

# Sort `impwealth` in ascending order 
df.sort('impwealth', 'inplace'=True)

# Find the 50th percentile weight, P
P = df['indweight'].sum() * (.5)

# Search for the first occurrence of `impweight` that is greater than P 
i = df.loc[df['indweight'] > P, 'indweight'].last_valid_index()

# The value of `impwealth` associated with this index will be the weighted median
w_median = df.ix[i, 'impwealth']

这种方法看起来有些笨重,而且我不确定它是否正确。在pandas参考资料中,我没有找到内置的方法来执行此操作。如何寻找加权中位数的最佳方法?


你确定你的伪代码是正确的吗?df['indweight'].sum() * (.5)将会得到一个约为219的值,而你的indweight值中没有一个超过该值。调用df['indweight'].median()将得到44.835,而使用mean()则会得到43.783。 - EdChum
我认为应该这样做...df['indweight'].sum() * (.5) 应该计算数据中落在第50个百分位以下的观测值数量,因为 indweight 是频率权重。因此,indweight 的均值和中位数超过其总和是有道理的。 - svenkatesh
@svenkatesh,你需要使用indweight.cumsum()而不是indweight本身。也许可以看一下我下面的答案。 - prooffreader
7个回答

27

如果你想要用纯pandas实现这个功能,这里有一种方法。它也不进行插值。(@svenkatesh,在你的伪代码中缺少累加和)

df.sort_values('impwealth', inplace=True)
cumsum = df.indweight.cumsum()
cutoff = df.indweight.sum() / 2.0
median = df.impwealth[cumsum >= cutoff].iloc[0]

这个数据的中位数为925000。


排序应该在indweight列上进行...所以第一行应该是df.sort_values('indweight', inplace=True) - Manuel F
@ManuelF,原始代码是正确的。需要对数据进行排序,而不是权重。 - xkudsraw

10

你尝试过wquantiles包吗?我以前从未使用过它,但它有一个加权中位数函数,似乎能给出至少合理的答案(你可能需要双重检查它是否使用了你期望的方法)。

In [12]: import weighted

In [13]: weighted.median(df['impwealth'], df['indweight'])
Out[13]: 914662.0859091772

2
就个人而言,我有点担心安装一个只需要几行代码就能完成的包,但如果你需要插值加权中位数,也许这是最好的方法。 - prooffreader

8
这个函数泛化了校对员的解决方案:
def weighted_median(df, val, weight):
    df_sorted = df.sort_values(val)
    cumsum = df_sorted[weight].cumsum()
    cutoff = df_sorted[weight].sum() / 2.
    return df_sorted[cumsum >= cutoff][val].iloc[0]

在这个例子中,它将是 weighted_median(df, 'impwealth', 'indweight')

5
您可以使用这个解决方案来使用numpy计算加权百分位数:Weighted percentile using numpy
def weighted_quantile(values, quantiles, sample_weight=None, 
                      values_sorted=False, old_style=False):
    """ Very close to numpy.percentile, but supports weights.
    NOTE: quantiles should be in [0, 1]!
    :param values: numpy.array with data
    :param quantiles: array-like with many quantiles needed
    :param sample_weight: array-like of the same length as `array`
    :param values_sorted: bool, if True, then will avoid sorting of
        initial array
    :param old_style: if True, will correct output to be consistent
        with numpy.percentile.
    :return: numpy.array with computed quantiles.
    """
    values = np.array(values)
    quantiles = np.array(quantiles)
    if sample_weight is None:
        sample_weight = np.ones(len(values))
    sample_weight = np.array(sample_weight)
    assert np.all(quantiles >= 0) and np.all(quantiles <= 1), \
        'quantiles should be in [0, 1]'

    if not values_sorted:
        sorter = np.argsort(values)
        values = values[sorter]
        sample_weight = sample_weight[sorter]

    weighted_quantiles = np.cumsum(sample_weight) - 0.5 * sample_weight
    if old_style:
        # To be convenient with numpy.percentile
        weighted_quantiles -= weighted_quantiles[0]
        weighted_quantiles /= weighted_quantiles[-1]
    else:
        weighted_quantiles /= np.sum(sample_weight)
    return np.interp(quantiles, weighted_quantiles, values)

调用 weighted_quantile(df.impwealth, quantiles=0.5, df.indweight) 函数。


你介意解释一下以下代码行的作用吗:weighted_quantiles = np.cumsum(sample_weight) - 0.5 * sample_weight - Jake
假设我们有值 [3, 10, 12] 和相关权重 [0.2, 0.5, 0.3](按先前的行排序)。 np.cumsum 将产生 [0.2, 0.7, 1.0],但实际上这些是相关分位数的右边缘。为了使它们居中,我们从每个桶中减去一半的权重,得到 [0.1, 0.45, 0.85]。然后我们从中插值得到加权分位数。 - Max Ghenis
非常感谢。还有一个问题(如果这很愚蠢,对不起),为什么我们要将分位数居中? - Jake
假设你只有两个值(3, 4),它们各自的分位数应该是(0, 0.5)(0.5, 1)(0.25, 0.75)还是(0, 1)?前两个存在问题,因为它们是不对称的。第三个是这个函数的默认设置,第四个是numpy.percentile所做的,可以通过在此处激活old_style=True参数来实现。默认设置的优点是,如果从分位数中进行抽样,你有一个非零的机会得到观察到的值,例如,在分位数0-0.25中,它将是3。然而,那种梯形分布可能比old_style平坦分布更难以理解。 - Max Ghenis

4
你可以使用我编写的这个函数来实现相同的目的。
注意:weighted 在最后使用插值法选择0.5分位数(您可以自行查看代码)。
我的编写的函数只返回一个边界0.5权重。
import numpy as np

def weighted_median(values, weights):
    ''' compute the weighted median of values list. The 
weighted median is computed as follows:
    1- sort both lists (values and weights) based on values.
    2- select the 0.5 point from the weights and return the corresponding values as results
    e.g. values = [1, 3, 0] and weights=[0.1, 0.3, 0.6] assuming weights are probabilities.
    sorted values = [0, 1, 3] and corresponding sorted weights = [0.6,     0.1, 0.3] the 0.5 point on
    weight corresponds to the first item which is 0. so the weighted     median is 0.'''

    #convert the weights into probabilities
    sum_weights = sum(weights)
    weights = np.array([(w*1.0)/sum_weights for w in weights])
    #sort values and weights based on values
    values = np.array(values)
    sorted_indices = np.argsort(values)
    values_sorted  = values[sorted_indices]
    weights_sorted = weights[sorted_indices]
    #select the median point
    it = np.nditer(weights_sorted, flags=['f_index'])
    accumulative_probability = 0
    median_index = -1
    while not it.finished:
        accumulative_probability += it[0]
        if accumulative_probability > 0.5:
            median_index = it.index
            return values_sorted[median_index]
        elif accumulative_probability == 0.5:
            median_index = it.index
            it.iternext()
            next_median_index = it.index
            return np.mean(values_sorted[[median_index, next_median_index]])
        it.iternext()

    return values_sorted[median_index]
#compare weighted_median function and np.median
print weighted_median([1, 3, 0, 7], [2,3,3,9])
print np.median([1,1,0,0,0,3,3,3,7,7,7,7,7,7,7,7,7])

加权中位数函数与已接受的答案在代码上看起来非常相似,但在结尾处不进行插值。 - Ash

2
你可以使用 robustats 库来计算加权中位数:
import numpy as np
import robustats # pip install robustats


# Weighted Median
x = np.array([1.1, 5.3, 3.7, 2.1, 7.0, 9.9])
weights = np.array([1.1, 0.4, 2.1, 3.5, 1.2, 0.8])

weighted_median = robustats.weighted_median(x, weights)

print("The weighted median is {}".format(weighted_median))

0

有一个weightedstats包,可以通过condapip两种方式获得,它可以执行weighted_median操作。

假设您正在使用conda,从终端(Mac/Linux)或Anaconda提示符(Win)中输入:

conda activate YOURENVIRONMENT
conda install -c conda-forge -y weightedstats

-y 的意思是“不要询问我是否确认更改,直接执行”)

然后在您的Python代码中:

import pandas as pd
import weightedstats as ws

df = pd.read_csv('/your/data/file.csv')
ws.weighted_median(df['values_col'], df['weights_col'])

我不确定它是否适用于所有情况,但我刚刚对一些简单数据进行了比较,与R包matrixStats中的weightedMedian()函数进行了比较,并且两者得到了相同的结果。


顺便提一下,使用weightedstats,您也可以计算weighted_mean(),虽然这也可以通过NumPy实现:

np.average(df['values_col'], weights=df['weights_col'])

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接