我有一个CSV,它为每一行提供三个不同的x值的y值。当读到熊猫DataFrame时,它看起来是这样的:
5 10 20
0 -13.6 -10.7 -10.3
1 -14.1 -11.2 -10.8
2 -12.3 -9.4 -9.0也就是说,对于第0行,5的值为-13.6,在10为-10.7,在20为-10.3。这些值是形式中的算法的结果:
def calc(x, r, b, c, d):
if x < 10:
y = (x * r + b) / x
elif x >= 10 and x < 20:
y = ((x * r) + (b - c)) / x
else:
y = ((x * r) + (b - d)) / x
return y我希望找到每一行的r、b、c和d的值。我对每一种价值观都很了解。例如,对于每一行:r在np.arange(-.05,-.11,-.01)中,b在np.arange中(0,-20.05,-.05),c和d在np.arange中(0,85,5)。我还知道d是<= c。
目前,我正在用蛮力解决这个问题。对于每一行,我遍历r、b、c和d的每个组合,并测试三个x值处的值是否等于DataFrame中的已知值。这起作用,给了我几个组合,每一行基本上是相同的,除了四舍五入的差异。
问题是,当我需要对2,000+行运行它时,这种方法需要很长时间。我的问题是:有比迭代和测试每个组合更快的方法吗?我的理解是,这是一个约束满意问题,但在那之后,我不知道该缩小什么范围;有太多类型的约束满意度问题(似乎),我仍然迷失了方向(我甚至不确定这是一个这样的问题!)如果能帮助我找到正确的方向,我将不胜感激。
发布于 2016-07-19 21:00:00
我希望我正确地理解了这项任务。
如果你知道参数的分辨率/离散化,它看起来就像一个离散优化问题(一般是:hard),可以用CP-方法来解决。
但是,如果允许这些值是连续的(并重新制定公式),则如下所示:
这三个版本都可以有效地解决!
这里是一个非一般(可以很容易推广)的实现,(3)使用cvxpy来表示问题,用埃科斯来解决QP。这两个工具都是开源的。
代码
import numpy as np
import time
from cvxpy import *
from random import uniform
""" GENERATE TEST DATA """
def sample_params():
while True:
r = uniform(-0.11, -0.05)
b = uniform(-20.05, 0)
c = uniform(0, 85)
d = uniform(0, 85)
if d <= c:
return r, b, c, d
def calc(x, r, b, c, d):
if x < 10:
y = (x * r + b) / x
elif x >= 10 and x < 20:
y = ((x * r) + (b - c)) / x
else:
y = ((x * r) + (b - d)) / x
return y
N = 2000
sampled_params = [sample_params() for i in range(N)]
data_5 = np.array([calc(5, *sampled_params[i]) for i in range(N)])
data_10 = np.array([calc(10, *sampled_params[i]) for i in range(N)])
data_20 = np.array([calc(20, *sampled_params[i]) for i in range(N)])
data = np.empty((N, 3))
for i in range(N):
data[i, :] = [data_5[i], data_10[i], data_20[i]]
""" SOLVER """
def solve(row):
""" vars """
R = Variable(1)
B = Variable(1)
C = Variable(1)
D = Variable(1)
E = Variable(3)
""" constraints """
constraints = []
# bounds
constraints.append(R >= -.11)
constraints.append(R <= -.05)
constraints.append(B >= -20.05)
constraints.append(B <= 0.0)
constraints.append(C >= 0.0)
constraints.append(C <= 85.0)
constraints.append(D >= 0.0)
constraints.append(D <= 85.0)
constraints.append(D <= C)
# formula of model
constraints.append((1.0 / 5.0) * B + R == row[0] + E[0]) # alternate function form: b/x+r
constraints.append((1.0 / 10.0) * B - (1.0 / 10.0) * C == row[1] + E[1]) # alternate function form: b/x-c/x+r
constraints.append((1.0 / 20.0) * B - (1.0 / 20.0) * D == row[2] + E[2]) # alternate function form: b/x-d/x+r
""" Objective """
objective = Minimize(norm(E, 2))
""" Solve """
problem = Problem(objective, constraints)
problem.solve(solver=ECOS, verbose=False)
return R.value, B.value, C.value, D.value, E.value
start = time.time()
for i in range(N):
r, b, c, d, e = solve(data[i])
end = time.time()
print('seconds taken: ', end-start)
print('seconds per row: ', (end-start) / N)输出
('seconds taken: ', 20.620506048202515)
('seconds per row: ', 0.010310253024101258)https://stackoverflow.com/questions/38354320
复制相似问题