我试图用cvxpy
找到凸似然函数的最优值。变量表示分段线性函数在预先指定的x值处的y值,x值表示一天内的时间。我想知道是否有经验的人对如何正确实现下面的似然函数有意见
The log-likelihood function is defined like this
The function ^{
Here is a link to the full research paper that I am trying to implement,但我怀疑你可能需要某种大学/研究人员的权限来阅读它
这是我写的代码。下面是我目前的错误和堆栈跟踪,但我怀疑还有更大的概念错误-但我现在很迷茫,所以请问一些不清楚的问题,我将尝试进一步说明
def ll(values, times, y):
N = pd.DataFrame({
'count_at_time': values,
'timestamp': pd.to_datetime(times.flatten()),
})
N = N.set_index('timestamp')
dates = N.index.normalize().drop_duplicates().strftime('%Y-%m-%d').to_list()
f = np.vectorize(lambda date: sum_term(N, y, date))
terms = f(dates)
n = len(dates)
ll = (1.0/n) * cp.sum(terms)
return ll
# N[i](s) is the cumulative arrivals on day i at time s
def sum_term(N, y, date):
t = Intervals.get(date)
# hack: add 1 second to simulate open interval
t0 = t[0].start_time.to_numpy().astype(int) + 1000000000
tp = t[-1].end_time.to_numpy().astype(int)
def t1(s):
v1 = rate(y, t, s)
v2 = at(N, date, s)
return cp.log(v1 * v2)
xvals = np.linspace(t0, tp, num=50, endpoint=False)
yvals = list(map(lambda s: t1(s), xvals))
v1 = np.trapz(yvals, x=xvals)
yvals2 = list(map(lambda s: rate(y, t, s), xvals))
v2 = np.trapz(yvals2, x=xvals)
return v1 - v2
def rate(y, t, s):
thing = t.start_time.astype(np.int64).values
# end of day hack
thing = np.append(thing, t[23].end_time.to_numpy().astype(int))
i = np.searchsorted(thing, s, side='left')
m = y[i-1]
k = (y[i]- y[i-1]) / (thing[i] - thing[i-1])
x = (s - thing[i-1]) # difference between time s and time at interval point t[i-1]
return m + k * x
# data is:
# - arrival time grouped by day
# - count of arrivals at that time that day
data = ExperimentsService.run().data
N = data.count_at_time.values
X = np.array(data.index.astype(np.int64))
print("X", X, type(X), np.shape(X))
y = cp.Variable(25)
print(y.shape)
objective = cp.Minimize(-ll(N, X, y))
constraints = [ y >= 0, y[0] == y[-1] ]
prob = cp.Problem(objective, constraints)
val = prob.solve()
solution = y.value
print("Optimal value", val)
print("Optimal variables", solution)
Traceback (most recent call last):
File "convex.py", line 132, in <module>
objective = cp.Minimize(-ll(N, X, y))
File "convex.py", line 34, in ll
ll = (1.0/n) * cp.sum(terms)
File "/Users/aron/.local/share/virtualenvs/project-nFRdcLaN/lib/python3.7/site-packages/cvxpy/atoms/affine/sum.py", line 109, in sum
return Sum(expr, axis, keepdims)
File "/Users/aron/.local/share/virtualenvs/project-nFRdcLaN/lib/python3.7/site-packages/cvxpy/atoms/affine/sum.py", line 39, in __init__
super(Sum, self).__init__(expr, axis=axis, keepdims=keepdims)
File "/Users/aron/.local/share/virtualenvs/project-nFRdcLaN/lib/python3.7/site-packages/cvxpy/atoms/axis_atom.py", line 31, in __init__
super(AxisAtom, self).__init__(expr)
File "/Users/aron/.local/share/virtualenvs/project-nFRdcLaN/lib/python3.7/site-packages/cvxpy/atoms/atom.py", line 43, in __init__
self.args = [Atom.cast_to_const(arg) for arg in args]
File "/Users/aron/.local/share/virtualenvs/project-nFRdcLaN/lib/python3.7/site-packages/cvxpy/atoms/atom.py", line 43, in <listcomp>
self.args = [Atom.cast_to_const(arg) for arg in args]
File "/Users/aron/.local/share/virtualenvs/project-nFRdcLaN/lib/python3.7/site-packages/cvxpy/expressions/expression.py", line 446, in cast_to_const
return expr if isinstance(expr, Expression) else cvxtypes.constant()(expr)
File "/Users/aron/.local/share/virtualenvs/project-nFRdcLaN/lib/python3.7/site-packages/cvxpy/expressions/constants/constant.py", line 44, in __init__
self._value = intf.DEFAULT_INTF.const_to_matrix(value)
File "/Users/aron/.local/share/virtualenvs/project-nFRdcLaN/lib/python3.7/site-packages/cvxpy/interface/numpy_interface/ndarray_interface.py", line 50, in const_to_matrix
return result.astype(numpy.float64)
ValueError: setting an array element with a sequence.
目前没有回答
相关问题 更多 >
编程相关推荐