我有一个目标变量x
以及一些其他变量A
和B
。当要满足x
和A
的某些条件时,我想计算B
的平均值(以及其他统计信息)。一个现实的例子是,当太阳辐射(x
)和风速(A
)落入某些预定义的bin范围内时,通过一连串的测量值来计算平均气温(B
)。
我已经能够通过循环来完成此操作(请参见下面的示例),但是我了解到我应该避免在数据帧上循环。从我对该站点的研究中,我觉得可能存在使用pd.cut
或np.select
的更优雅/矢量化的解决方案,但坦率地说,我不知道该怎么做。
import pandas as pd
import numpy as np
n = 100
df = pd.DataFrame(
{
"x": np.random.randn(n),
"A": np.random.randn(n)+5,
"B": np.random.randn(n)+10
}
)
[df.head()
输出:
x A B
0 -0.585313 6.038620 9.909762
1 0.412323 3.991826 8.836848
2 0.211713 5.019520 9.667349
3 0.710699 5.353677 9.757903
4 0.681418 4.452754 10.647738
# define bin ranges
bins_A = np.arange(3, 8)
bins_B = np.arange(8, 13)
# prepare output lists
A_mins= []
A_maxs= []
B_mins= []
B_maxs= []
x_means= []
x_stds= []
x_counts= []
# loop over bins
for i_A in range(0, len(bins_A)-1):
A_min = bins_A[i_A]
A_max = bins_A[i_A+1]
for i_B in range(0, len(bins_B)-1):
B_min = bins_B[i_B]
B_max = bins_B[i_B+1]
# binning conditions for current step
conditions = np.logical_and.reduce(
[
df["A"] > A_min,
df["A"] < A_max,
df["B"] > B_min,
df["B"] < B_max,
]
)
# calculate statistics for x and store values in lists
x_values = df.loc[conditions, "x"]
x_means.append(x_values.mean())
x_stds.append(x_values.std())
x_counts.append(x_values.count())
A_mins.append(A_min)
A_maxs.append(A_max)
B_mins.append(B_min)
B_maxs.append(B_max)
binned = pd.DataFrame(
data={
"A_min": A_mins,
"A_max": A_maxs,
"B_min": B_mins,
"B_max": B_maxs,
"x_mean": x_means,
"x_std": x_stds,
"x_count": x_counts
}
)
[binned.head()
输出:
A_min A_max B_min B_max x_mean x_std x_count
0 3 4 8 9 0.971624 0.790972 2
1 3 4 9 10 0.302795 0.380102 3
2 3 4 10 11 0.447398 1.787659 5
3 3 4 11 12 0.462149 1.195844 2
4 4 5 8 9 0.379431 0.983965 4
我们将尝试将其保留在pandas / NumPy中,以便我们可以利用数据框方法或数组方法和ufunc,同时对其进行矢量化处理。当要解决复杂的问题或生成统计信息时,这使得扩展功能变得更加容易,这似乎就是这种情况。
现在,要解决该问题并使它靠近大熊猫,将是生成类似于分别在给定分箱A
和B
上对bins_A
和bins_B
的组合跟踪的中间ID或标签。为此,一种方法是分别对这两个数据使用searchsorted
-
tagsA = np.searchsorted(bins_A,df.A)
tagsB = np.searchsorted(bins_B,df.B)
现在,我们仅对边界内的情况感兴趣,因此需要屏蔽-
vm = (tagsB>0) & (tagsB<len(bins_B)) & (tagsA>0) & (tagsA<len(bins_A))
让我们将此掩码应用于原始数据帧-
dfm = df.iloc[vm]
为有效标签添加标签,它们代表A_mins
和B_min
等价物,因此将显示在最终输出中-
dfm['TA'] = bins_A[(tagsA-1)[vm]]
dfm['TB'] = bins_B[(tagsB-1)[vm]]
因此,我们标记的数据框已准备就绪,在将这两个标签分组后,可以将其用作describe-d
以获取通用统计信息-
df_out = dfm.groupby(['TA','TB'])['x'].describe()
示例运行使情况更清晰,同时与发布的解决方案进行了比较-
In [46]: np.random.seed(0)
...: n = 100
...: df = pd.DataFrame(
...: {
...: "x": np.random.randn(n),
...: "A": np.random.randn(n)+5,
...: "B": np.random.randn(n)+10
...: }
...: )
In [47]: binned
Out[47]:
A_min A_max B_min B_max x_mean x_std x_count
0 3 4 8 9 0.400199 0.719007 5
1 3 4 9 10 -0.268252 0.914784 6
2 3 4 10 11 0.458746 1.499419 5
3 3 4 11 12 0.939782 0.055092 2
4 4 5 8 9 0.238318 1.173704 5
5 4 5 9 10 -0.263020 0.815974 8
6 4 5 10 11 -0.449831 0.682148 12
7 4 5 11 12 -0.273111 1.385483 2
8 5 6 8 9 -0.438074 NaN 1
9 5 6 9 10 -0.009721 1.401260 16
10 5 6 10 11 0.467934 1.221720 11
11 5 6 11 12 0.729922 0.789260 3
12 6 7 8 9 -0.977278 NaN 1
13 6 7 9 10 0.211842 0.825401 7
14 6 7 10 11 -0.097307 0.427639 5
15 6 7 11 12 0.915971 0.195841 2
In [48]: df_out
Out[48]:
count mean std ... 50% 75% max
TA TB ...
3 8 5.0 0.400199 0.719007 ... 0.302472 0.976639 1.178780
9 6.0 -0.268252 0.914784 ... -0.001510 0.401796 0.653619
10 5.0 0.458746 1.499419 ... 0.462782 1.867558 1.895889
11 2.0 0.939782 0.055092 ... 0.939782 0.959260 0.978738
4 8 5.0 0.238318 1.173704 ... -0.212740 0.154947 2.269755
9 8.0 -0.263020 0.815974 ... -0.365103 0.449313 0.950088
10 12.0 -0.449831 0.682148 ... -0.436773 -0.009697 0.761038
11 2.0 -0.273111 1.385483 ... -0.273111 0.216731 0.706573
5 8 1.0 -0.438074 NaN ... -0.438074 -0.438074 -0.438074
9 16.0 -0.009721 1.401260 ... 0.345020 1.284173 1.950775
10 11.0 0.467934 1.221720 ... 0.156349 1.471263 2.240893
11 3.0 0.729922 0.789260 ... 1.139401 1.184846 1.230291
6 8 1.0 -0.977278 NaN ... -0.977278 -0.977278 -0.977278
9 7.0 0.211842 0.825401 ... 0.121675 0.398750 1.764052
10 5.0 -0.097307 0.427639 ... -0.103219 0.144044 0.401989
11 2.0 0.915971 0.195841 ... 0.915971 0.985211 1.054452
因此,如前所述,我们在A_min
和B_min
中有TA
和TB
,而相关的统计信息则在其他标题中捕获。请注意,这将是一个多索引数据框。如果我们需要捕获等效的数组数据,只需执行以下操作:df_out.loc[:,['count','mean','std']].values
用于统计,而np.vstack(df_out.loc[:,['count','mean','std']].index)
用于bin间隔开始。
或者,在不使用describe
的情况下捕获等效的统计数据,但是使用数据框方法,我们可以执行以下操作-
dfmg = dfm.groupby(['TA','TB'])['x']
dfmg.size().unstack().values
dfmg.std().unstack().values
dfmg.mean().unstack().values
替代#1:使用pd.cut
我们也可以按照问题中的建议使用pd.cut
来代替searchsorted
,因为它是自动处理越界越界的问题,使基本概念相同-
df['TA'] = pd.cut(df['A'],bins=bins_A, labels=range(len(bins_A)-1))
df['TB'] = pd.cut(df['B'],bins=bins_B, labels=range(len(bins_B)-1))
df_out = df.groupby(['TA','TB'])['x'].describe()
所以,这给了我们统计数据。对于A_min
和B_min
等效项,只需使用索引级别-
A_min = bins_A[df_out.index.get_level_values(0)]
B_min = bins_B[df_out.index.get_level_values(1)]
或使用某些网格方法-
mA,mB = np.meshgrid(bins_A[:-1],bins_B[:-1])
A_min,B_min = mA.ravel('F'),mB.ravel('F')
bincount
[我们可以再次以矢量化的方式利用np.bincount
获得所有这三个统计量度值,包括标准差-
lA,lB = len(bins_A),len(bins_B)
n = lA+1
x,A,B = df.x.values,df.A.values,df.B.values
tagsA = np.searchsorted(bins_A,A)
tagsB = np.searchsorted(bins_B,B)
t = tagsB*n + tagsA
L = n*lB
countT = np.bincount(t, minlength=L)
countT_x = np.bincount(t,x, minlength=L)
avg_all = countT_x/countT
count = countT.reshape(-1,n)[1:,1:-1].ravel('F')
avg = avg_all.reshape(-1,n)[1:,1:-1].ravel('F')
# Using numpy std definition for ddof case
ddof = 1.0 # default one for pandas std
grp_diffs = (x-avg_all[t])**2
std_all = np.sqrt(np.bincount(t,grp_diffs, minlength=L)/(countT-ddof))
stds = std_all.reshape(-1,n)[1:,1:-1].ravel('F')
sorting
来利用reduceat
方法-x,A,B = df.x.values,df.A.values,df.B.values
vm = (A>bins_A[0]) & (A<bins_A[-1]) & (B>bins_B[0]) & (B<bins_B[-1])
xm = x[vm]
tagsA = np.searchsorted(bins_A,A)
tagsB = np.searchsorted(bins_B,B)
tagsAB = tagsB*(tagsA.max()+1) + tagsA
tagsABm = tagsAB[vm]
sidx = tagsABm.argsort()
tagsAB_s = tagsABm[sidx]
xms = xm[sidx]
cut_idx = np.flatnonzero(np.r_[True,tagsAB_s[:-1]!=tagsAB_s[1:],True])
N = (len(bins_A)-1)*(len(bins_B)-1)
count = np.diff(cut_idx)
avg = np.add.reduceat(xms,cut_idx[:-1])/count
stds = np.empty(N)
for ii,(s0,s1) in enumerate(zip(cut_idx[:-1],cut_idx[1:])):
stds[ii] = np.std(xms[s0:s1], ddof=1)
要获得与pandas数据框样式化输出相同或相似的格式,我们需要调整形状。因此,它将是avg.reshape(-1,len(bins_A)-1).T
,依此类推。
如果您担心的是关于[[性能,则可以使用for循环,如果您使用numba,则可以进行较小的更改
这里有一个执行计算的函数。关键是calculate
使用numba,因此速度非常快。其余仅用于创建熊猫数据框:from numba import njit
def calc_numba(df, bins_A, bins_B):
""" wrapper for the timeit. It only creates a dataframe """
@njit
def calculate(A, B, x, bins_A, bins_B):
size = (len(bins_A) - 1)*(len(bins_B) - 1)
out = np.empty((size, 7))
index = 0
for i_A, A_min in enumerate(bins_A[:-1]):
A_max = bins_A[i_A + 1]
for i_B, B_min in enumerate(bins_B[:-1]):
B_max = bins_B[i_B + 1]
mfilter = (A_min < A)*(A < A_max)*(B_min < B)*(B < B_max)
x_values = x[mfilter]
out[index, :] = [
A_min,
A_max,
B_min,
B_max,
x_values.mean(),
x_values.std(),
len(x_values)
]
index += 1
return out
columns = ["A_min", "A_max", "B_min", "B_max", "mean", "std", "count"]
out = calculate(df["A"].values, df["B"].values, df["x"].values, bins_A, bins_B)
return pd.DataFrame(out, columns=columns)
性能测试使用
n = 1_000_000
以及相同的bins_A
和bins_B
,我们得到:
%timeit code_question(df, bins_A, bins_B) 15.7 s ± 428 ms per loop (mean ± std. dev. of 7 runs, 1 loop each) %timeit calc_numba(df, bins_A, bins_B) 507 ms ± 12.3 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
比问题代码快[<30
由于pandas
内置方法使用了类似的增强功能,因此很难击败numba的性能。
import pandas as pd
import numpy as np
n = 20
df = pd.DataFrame(
{
"x": np.random.randn(n),
"A": np.random.randn(n)+5,
"B": np.random.randn(n)+10
}
)
# define bin ranges
bins_A = np.arange(3, 8)
bins_B = np.arange(8, 13)
直到这里,我使用您的示例。然后我使用numpy]来介绍上下Bin边缘。
A_mins=bins_A[:-1] A_maxs=bins_A[1:] B_mins=bins_B[:-1] B_maxs=bins_B[1:]
以某种方式将它们放在一起,实际上您是在使用那些嵌套循环,我将自己限制为numpy,在这里我仍然可以精确地保留嵌套循环所要获得的结构。
A_mins_list=np.repeat(A_mins,len(B_mins)) A_maxs_list=np.repeat(A_maxs,len(B_mins)) B_mins_list=np.tile(B_mins,len(A_mins)) B_maxs_list=np.tile(B_maxs,len(A_mins))
新数据帧已使用bin信息初始化。
newdf=pd.DataFrame(np.array([A_mins_list,A_maxs_list,B_mins_list,B_maxs_list]).T,columns=['Amin','Amax','Bmin','Bmax'])
xvalues列在这里是最邪恶的,因为我必须使其成为一个numpy数组以适合数据框。然后,该子数组是一个numpy数组,还必须将其视为一个。请记住这一点,因为某些熊猫功能可能对此不起作用。在某些情况下,它必须是numpy函数。
newdf['xvalues']=newdf.apply(lambda row:np.array(df.x[(row.Amin<df.A) & (row.Amax>df.A) & (row.Bmin<df.B) & (row.Bmax>df.B)]),axis=1)
此外,您可以使用lambda函数执行任何所需的操作。就像我说的那样,这也许不是最有效的方法,但是代码有些清晰,只要您不需要数百万个条目的数据帧所需的最高性能,就可以通过[]
newdf['xmeans']=newdf.apply(lambda row: row.xvalues.mean(),axis=1) newdf['stds']=newdf.apply(lambda row: row.xvalues.std(),axis=1) newdf['xcounts']=newdf.apply(lambda row: row.xvalues.size,axis=1)
或您可能喜欢的任何东西。
使用cython可以避免使用lambda方式来显着提高性能,但是我不习惯cycyon,所以我宁愿将它留给专家...[另外请注意,如果您要尝试使用仅一个值的空数组或std,则可能会出现一些警告。如果需要,可以使用警告包将其抑制。