Scipy 使用分组依据 pandas 数据框最小化
Scipy minimize with pandas dataframe with group by
我有一个数据框(下面的示例 df)并试图最小化它的成本函数。
GrpId = ['A','A','A','A','A','A','B','B','B','B','B','B','B']
col1 = [69.1,70.5,71.4,72.8,73.2,74.2,208.0,209.2,210.2,211.0,211.2,211.7,212.5]
col2 = [2,3.1,1.1,2.1,6.0,1.1,1.2,1.3,3.1,2.9,5.0,6.1,3.2]
d = {'GrpId':GrpId,'col1':col1,'col2':col2}
df1 = pd.DataFrame(d)
下面是最小化和成本函数。
col1_const=[0,0,0,0,60.0,0,0,0]
col2_const=[0,0,0,0,0,100.0,0,0]
def main(type1,type2,type3,df):
vall0=[type1,type2,type3]
res=minimize(cost_fun, vall0, args=(df), method = 'SLSQP', tol=0.01)
[type1,type2,type3]=res.x
return type1,type2,type3
def cost_fun(v, df):
df['col1_res'][i] = np.where((df['col1'][i]!=np.nan), ((1/0.095)*(np.sqrt(df['col1'][i])-np.sqrt(col1_const[4]*(0.1*v[1]+v[2])**2)))**2 ,0)
df['col2_res'][i] = np.where((df['col2'][i]!=np.nan), ((1/0.12)*(np.sqrt(df['col2'][i])-np.sqrt(col2_const[5]*(0.1*v[0]+v[2])**2)))**2 ,0)
res=0.5*np.sqrt(df['col1_res'][i]+df['col2_res'][i])
return res
然后我在循环中迭代这个函数,如下所示,它工作但需要大量时间和内存,
df1['type1']=np.nan
df1['type2']=np.nan
df1['type3']=np.nan
df1['col3']=np.nan
df1['col1_res']=np.nan
df1['col2_res']=np.nan
for i in range(len(df1.GrpId)):
if i==0:
df1['type1'][i], df1['type2'][i], df1['type3'][i]= main(0.125, 0.125, 0.125,df1)
else:
df1['type1'][i], df1['type2'][i], df1['type3'][i]= main(df1['type1'][i-1], df1['type2'][i-1], df1['type3'][i-1],df1)
df1['col3'][i]=df1['type1'][i]+df1['type2'][i]
请注意,我有更大的数据框,有更多的行和列,对于这个问题,我刚刚创建了一个示例 code/case。
我的问题是,
- 如何在不迭代的情况下做同样的事情
col1_const[4]
值将根据组(按 GrpId 分组)更改 - 我有另一个函数来计算每个组的 col1_const[4] 值。在这种情况下,如何按组将此值传递给 cost_fun。
首先,我认为没有必要在 objective 函数中检查 != np.nan
。相反,您可以清理数据框并将所有 np.nan
替换为零。 objective 函数在优化例程中被多次调用,因此应尽可能高效快速地编写。因此,我们删除了 np.where
的调用。另请注意,依赖索引变量 i
在外部范围内已知的事实是不好的做法,并且会使代码难以阅读。我会推荐这样的东西:
col1 = df1.col1.values[df1.col1.values != np.nan]
col2 = df1.col2.values[df1.col2.values != np.nan]
col1_const = np.array([0,0,0,0,60.0,0,0,0])
col2_const = np.array([0,0,0,0,0,100.0,0,0])
def cost_fun1(v, *args):
i, col1, col2, col1_const, col2_const = args
col1_res = ((1/0.095)*(np.sqrt(col1[i]) - np.sqrt(col1_const[4]*(0.1*v[1]+v[2])**2)))**2
col2_res = ((1/0.12)*(np.sqrt(col2[i]) - np.sqrt(col2_const[5]*(0.1*v[0]+v[2])**2)))**2
return 0.5*np.sqrt(col1_res + col2_res)
接下来,更重要的是,您要解决多个优化问题,而不是解决一个 large-scale 优化问题。从数学上讲,因为您的 objective 函数保证为正,所以您可以按照与 相同的方式重新表述您的问题。然后,cost_fun2
基本上 returns 所有索引 i 的所有 cost_fun1
的总和。使用一点重塑魔法,函数看起来几乎相同:
def cost_fun2(vv, *args):
col1, col2, col1_const, col2_const = args
v = vv.reshape(3, col1.size)
col1_res = ((1/0.095)*(np.sqrt(col1) - np.sqrt(col1_const[4]*(0.1*v[1]+v[2])**2)))**2
col2_res = ((1/0.12)*(np.sqrt(col2) - np.sqrt(col2_const[5]*(0.1*v[0]+v[2])**2)))**2
return np.sum(0.5*np.sqrt(col1_res + col2_res))
然后,我们简单地解决问题,然后将解决方案的值写入数据框:
from scipy.optimize import minimize
# initial guess
x0 = np.ones(3*col1.size)
# solve the problem
res = minimize(lambda vv: cost_fun2(vv, col1, col2, col1_const, col2_const), x0=x0, method="trust-constr")
# write to dataframe
type1_vals, type2_vals, type3_vals = np.split(res.x, 3)
df1['type1'] = type1_vals
df1['type2'] = type2_vals
df1['type3'] = type3_vals
如果数据框中需要 col1_res
和 col2_res
,可以直接修改 objective 函数。
最后但同样重要的是,根据数据帧的大小,强烈建议将精确的 objective 梯度传递给 scipy.optimize.minimize
以获得良好的收敛性能。目前,梯度是通过有限差分来近似的,这很慢并且容易出现舍入误差。
我有一个数据框(下面的示例 df)并试图最小化它的成本函数。
GrpId = ['A','A','A','A','A','A','B','B','B','B','B','B','B']
col1 = [69.1,70.5,71.4,72.8,73.2,74.2,208.0,209.2,210.2,211.0,211.2,211.7,212.5]
col2 = [2,3.1,1.1,2.1,6.0,1.1,1.2,1.3,3.1,2.9,5.0,6.1,3.2]
d = {'GrpId':GrpId,'col1':col1,'col2':col2}
df1 = pd.DataFrame(d)
下面是最小化和成本函数。
col1_const=[0,0,0,0,60.0,0,0,0]
col2_const=[0,0,0,0,0,100.0,0,0]
def main(type1,type2,type3,df):
vall0=[type1,type2,type3]
res=minimize(cost_fun, vall0, args=(df), method = 'SLSQP', tol=0.01)
[type1,type2,type3]=res.x
return type1,type2,type3
def cost_fun(v, df):
df['col1_res'][i] = np.where((df['col1'][i]!=np.nan), ((1/0.095)*(np.sqrt(df['col1'][i])-np.sqrt(col1_const[4]*(0.1*v[1]+v[2])**2)))**2 ,0)
df['col2_res'][i] = np.where((df['col2'][i]!=np.nan), ((1/0.12)*(np.sqrt(df['col2'][i])-np.sqrt(col2_const[5]*(0.1*v[0]+v[2])**2)))**2 ,0)
res=0.5*np.sqrt(df['col1_res'][i]+df['col2_res'][i])
return res
然后我在循环中迭代这个函数,如下所示,它工作但需要大量时间和内存,
df1['type1']=np.nan
df1['type2']=np.nan
df1['type3']=np.nan
df1['col3']=np.nan
df1['col1_res']=np.nan
df1['col2_res']=np.nan
for i in range(len(df1.GrpId)):
if i==0:
df1['type1'][i], df1['type2'][i], df1['type3'][i]= main(0.125, 0.125, 0.125,df1)
else:
df1['type1'][i], df1['type2'][i], df1['type3'][i]= main(df1['type1'][i-1], df1['type2'][i-1], df1['type3'][i-1],df1)
df1['col3'][i]=df1['type1'][i]+df1['type2'][i]
请注意,我有更大的数据框,有更多的行和列,对于这个问题,我刚刚创建了一个示例 code/case。
我的问题是,
- 如何在不迭代的情况下做同样的事情
col1_const[4]
值将根据组(按 GrpId 分组)更改 - 我有另一个函数来计算每个组的 col1_const[4] 值。在这种情况下,如何按组将此值传递给 cost_fun。
首先,我认为没有必要在 objective 函数中检查 != np.nan
。相反,您可以清理数据框并将所有 np.nan
替换为零。 objective 函数在优化例程中被多次调用,因此应尽可能高效快速地编写。因此,我们删除了 np.where
的调用。另请注意,依赖索引变量 i
在外部范围内已知的事实是不好的做法,并且会使代码难以阅读。我会推荐这样的东西:
col1 = df1.col1.values[df1.col1.values != np.nan]
col2 = df1.col2.values[df1.col2.values != np.nan]
col1_const = np.array([0,0,0,0,60.0,0,0,0])
col2_const = np.array([0,0,0,0,0,100.0,0,0])
def cost_fun1(v, *args):
i, col1, col2, col1_const, col2_const = args
col1_res = ((1/0.095)*(np.sqrt(col1[i]) - np.sqrt(col1_const[4]*(0.1*v[1]+v[2])**2)))**2
col2_res = ((1/0.12)*(np.sqrt(col2[i]) - np.sqrt(col2_const[5]*(0.1*v[0]+v[2])**2)))**2
return 0.5*np.sqrt(col1_res + col2_res)
接下来,更重要的是,您要解决多个优化问题,而不是解决一个 large-scale 优化问题。从数学上讲,因为您的 objective 函数保证为正,所以您可以按照与 cost_fun2
基本上 returns 所有索引 i 的所有 cost_fun1
的总和。使用一点重塑魔法,函数看起来几乎相同:
def cost_fun2(vv, *args):
col1, col2, col1_const, col2_const = args
v = vv.reshape(3, col1.size)
col1_res = ((1/0.095)*(np.sqrt(col1) - np.sqrt(col1_const[4]*(0.1*v[1]+v[2])**2)))**2
col2_res = ((1/0.12)*(np.sqrt(col2) - np.sqrt(col2_const[5]*(0.1*v[0]+v[2])**2)))**2
return np.sum(0.5*np.sqrt(col1_res + col2_res))
然后,我们简单地解决问题,然后将解决方案的值写入数据框:
from scipy.optimize import minimize
# initial guess
x0 = np.ones(3*col1.size)
# solve the problem
res = minimize(lambda vv: cost_fun2(vv, col1, col2, col1_const, col2_const), x0=x0, method="trust-constr")
# write to dataframe
type1_vals, type2_vals, type3_vals = np.split(res.x, 3)
df1['type1'] = type1_vals
df1['type2'] = type2_vals
df1['type3'] = type3_vals
如果数据框中需要 col1_res
和 col2_res
,可以直接修改 objective 函数。
最后但同样重要的是,根据数据帧的大小,强烈建议将精确的 objective 梯度传递给 scipy.optimize.minimize
以获得良好的收敛性能。目前,梯度是通过有限差分来近似的,这很慢并且容易出现舍入误差。