您是否知道以下問題的快速/優雅的Python/Scipy/Numpy解決方案: 您有一組帶有關聯值w(所有1D數組)的x,y座標。現在bin x和y到一個2D網格(大小爲BINSxBINS)並計算每個bin的w值的分位數(如中位數),最終得到具有所需分位數的BINSxBINS 2D數組。Python中的分位數/中位數/二進制分組
這很容易做一些嵌套循環,但我相信有一個更優雅的解決方案。
謝謝, 馬克
您是否知道以下問題的快速/優雅的Python/Scipy/Numpy解決方案: 您有一組帶有關聯值w(所有1D數組)的x,y座標。現在bin x和y到一個2D網格(大小爲BINSxBINS)並計算每個bin的w值的分位數(如中位數),最終得到具有所需分位數的BINSxBINS 2D數組。Python中的分位數/中位數/二進制分組
這很容易做一些嵌套循環,但我相信有一個更優雅的解決方案。
謝謝, 馬克
這是我想出了,我希望它是非常有用的。它不一定比使用循環更清潔或更好,但也許它會讓你開始朝更好的方向發展。
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
非常感謝您的代碼。在此基礎上,我發現我的問題的下述溶液(只有代碼的一小的修改):
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
隨着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]])
概率= [0.25,0.5 ,0.75]是分位數設置的默認值,您可以更改它或將其保留。
不幸的是,這是行不通的。 'mquantiles'只適用於沒有裝箱的數據。 – imsc 2012-08-28 13:21:40
我只是試圖自己做這個,聽起來像你想從你的命令「scipy.stats.binned_statistic_2d」可以找到平均值,中位數,標準偏差或任何定義的函數的第三個參數給定箱。
我意識到這個問題已經得到解答,但我相信這是一個很好的解決方案。
謝謝,這看起來已經不錯了。但似乎中間部分並不完全正確。例如訂單在錯誤的行中定義。 – Mark 2012-04-24 22:59:31
是的,你是對的。 – 2012-04-24 23:16:43
使用'np.histogram2d'你可以[非常容易地](https://stackoverflow.com/a/12588656/60982)按平均值進行分箱,但不是按中值計算。 – letmaik 2014-05-23 21:13:09