7

您是否知道针对以下问题的快速/优雅的 Python/Scipy/Numpy 解决方案：您有一组 x、y 坐标以及相关值 w（所有 1D 数组）。现在将 x 和 y 放入 2D 网格（大小 BINSxBINS）并计算每个 bin 的 w 值的分位数（如中位数），这最终会产生一个具有所需分位数的 BINSxBINS 2D 数组。

使用一些嵌套循环很容易做到这一点，但我相信有一个更优雅的解决方案。

谢谢，马克

4

4 回答 4

5

这是我想出来的，我希望它有用。它不一定比使用循环更清洁或更好，但也许它会让你开始朝着更好的方向发展。

import numpy as np
bins_x, bins_y = 1., 1.
x = np.array([1,1,2,2,3,3,3])
y = np.array([1,1,2,2,3,3,3])
w = np.array([1,2,3,4,5,6,7], 'float')

# You can get a bin number for each point like this
x = (x // bins_x).astype('int')
y = (y // bins_y).astype('int')
shape = [x.max()+1, y.max()+1]
bin = np.ravel_multi_index([x, y], shape)

# You could get the mean by doing something like:
mean = np.bincount(bin, w) / np.bincount(bin)

# Median is a bit harder
order = bin.argsort()
bin = bin[order]
w = w[order]
edges = (bin[1:] != bin[:-1]).nonzero()[0] + 1
med_index = (np.r_[0, edges] + np.r_[edges, len(w)]) // 2
median = w[med_index]

# But that's not quite right, so maybe
median2 = [np.median(i) for i in np.split(w, edges)]

也看看 numpy.histogram2d

于 2012-04-24T22:39:16.720 回答
3

我只是想自己做这件事，听起来你想要命令“scipy.stats.binned_statistic_2d”，你可以找到给定箱的第三个参数的平均值、中值、标准偏差或任何定义的函数。

我意识到这个问题已经得到解答，但我相信这是一个很好的内置解决方案。

于 2014-11-19T12:24:03.200 回答
1

非常感谢您的代码。基于它，我找到了以下问题的解决方案（仅对您的代码进行了少量修改）：

import numpy as np
BINS=10
boxsize=10.0
bins_x, bins_y = boxsize/BINS, boxsize/BINS
x = np.array([0,0,0,1,1,1,2,2,2,3,3,3])
y = np.array([0,0,0,1,1,1,2,2,2,3,3,3])
w = np.array([0,1,2,0,1,2,0,1,2,0,1,2], 'float')

# You can get a bin number for each point like this
x = (x // bins_x).astype('int')
y = (y // bins_y).astype('int')
shape = [BINS, BINS]
bin = np.ravel_multi_index([x, y], shape)


# Median 
order = bin.argsort()
bin = bin[order]
w = w[order]
edges = (bin[1:] != bin[:-1]).nonzero()[0] + 1
median = [np.median(i) for i in np.split(w, edges)]

#construct BINSxBINS matrix with median values
binvals=np.unique(bin)
medvals=np.zeros([BINS*BINS])
medvals[binvals]=median
medvals=medvals.reshape([BINS,BINS])

print medvals
于 2012-04-25T21:28:57.317 回答
0

使用 numpy/scipy 它是这样的：

    import numpy as np
    import scipy.stats as stats

    x = np.random.uniform(0,200,100)
    y = np.random.uniform(0,200,100)
    w = np.random.uniform(1,10,100)

    h = np.histogram2d(x,y,bins=[10,10], weights=w,range=[[0,200],[0,200]])
    hist, bins_x, bins_y = h
    q = stats.mstats.mquantiles(hist,prob=[0.25, 0.5, 0.75])

    >>> q.round(2)
    array([ 512.8 ,  555.41,  592.73])

    q1 = np.where(hist<q[0],1,0)
    q2 = np.where(np.logical_and(q[0]<=hist,hist<q[1]),2,0)
    q3 = np.where(np.logical_and(q[1]<=hist,hist<=q[2]),3,0)
    q4 = np.where(q[2]<hist,4,0)

    >>>q1 + q2 + q3 + q4
    array([[4, 3, 4, 3, 1, 1, 4, 3, 1, 2],
   [1, 1, 4, 4, 2, 3, 1, 3, 3, 3],
   [2, 3, 3, 2, 2, 2, 3, 2, 4, 2],
   [2, 2, 3, 3, 3, 1, 2, 2, 1, 4],
   [1, 3, 1, 4, 2, 1, 3, 1, 1, 3],
   [4, 2, 2, 1, 2, 1, 3, 2, 1, 1],
   [4, 1, 1, 3, 1, 3, 4, 3, 2, 1],
   [4, 3, 1, 4, 4, 4, 1, 1, 2, 4],
   [2, 4, 4, 4, 3, 4, 2, 2, 2, 4],
   [2, 2, 4, 4, 3, 3, 1, 3, 4, 4]])

prob = [0.25, 0.5, 0.75] 是分位数设置的默认值，您可以更改或保留它。

于 2012-05-21T23:07:47.583 回答