评估凸组合的权重
Evaluate the weights of a convex combination
我正在使用 scipy.spatial.ConvexHull API 来评估一组点的凸包并且效果很好。给定以下代码:
p1 = points[11]
hull = ConvexHull(points)
vertices = [points[i] for i in hull.vertices]
如何计算 vertices
的凸组合的系数(权重)等于 p1
?
非常感谢,
莫舍
注意:如果顶点数大于d+1
,其中d
是维度,那么组合将不是 独一无二。在答案的其余部分,为了简单起见,我假设 d=2
。
输入:
- 顶点
v0 = (x0, y0), v1 = (x1, y1), ..., vn = (xn, yn)
;
- 一个点
p1 = (x,y)
;
输出:
- 组合
a0, a1, ..., an
;
这样:
x = a0 * x0 + a1 * x1 + ... + an * xn
;
y = a0 * y0 + a1 * y1 + ... + an * yn
;
1 = a0 + a1 + ... + an
;
- 对于所有 i,
ai >= 0
。
这是一个包含三个未知线性方程组 (a0, a1, ..., an)
以及 n+1 个线性不等式的系统。我们可以使用 scipy's linear programming module scipy.optimize.linprog
求解这个系统
示例解决方案:
import random # generate points
import scipy.spatial as scispa # ConvexHull
import scipy.optimize as sciopt # linprog
from scipy.sparse import identity # identity matrix
points = [(random.randrange(0,100), random.randrange(0,100)) for _ in range(20)]
p1 = points[11]
hull = scispa.ConvexHull(points)
vertices = [points[i] for i in hull.vertices]
c = [1 for _ in vertices]
A = [[x for x,y in vertices], [y for x,y in vertices], [1 for _ in vertices]]
b = [p1[0], p1[1], 1]
s = sciopt.linprog(c, A_eq = A, b_eq = b)
>>> s.x
array([0.13393774, 0.06470577, 0.07367599, 0.09523271, 0.18924727,
0.26909487, 0.17410566])
>>> p1
(36, 57)
>>> [sum(a*xi for a,(xi,_) in zip(s.x, vertices)), sum(a*yi for a,(_,yi) in zip(s.x, vertices))]
[36.00000000719907, 57.00000000671608]
重要说明:我开始这个回答时警告说,如果平面中有超过 3 个顶点,或者维度 d 中有超过 d+1 个顶点,那么会有是我们方程组的多于 1 个解。上面的代码只有 returns 一种解决方案,意味着做出了任意选择。您可以控制这个任意选择:scipy.optimize.linprog
是一个优化工具;这里它最小化我们作为参数给出的向量 c
和解向量 s.x
的点积。这里我把c
的所有系数都设置为1,意思是linprog
会找到一个解中系数之和最小的解;尝试提供不同的向量 c
,你会得到不同的解决方案。
我正在使用 scipy.spatial.ConvexHull API 来评估一组点的凸包并且效果很好。给定以下代码:
p1 = points[11]
hull = ConvexHull(points)
vertices = [points[i] for i in hull.vertices]
如何计算 vertices
的凸组合的系数(权重)等于 p1
?
非常感谢, 莫舍
注意:如果顶点数大于d+1
,其中d
是维度,那么组合将不是 独一无二。在答案的其余部分,为了简单起见,我假设 d=2
。
输入:
- 顶点
v0 = (x0, y0), v1 = (x1, y1), ..., vn = (xn, yn)
; - 一个点
p1 = (x,y)
;
输出:
- 组合
a0, a1, ..., an
;
这样:
x = a0 * x0 + a1 * x1 + ... + an * xn
;y = a0 * y0 + a1 * y1 + ... + an * yn
;1 = a0 + a1 + ... + an
;- 对于所有 i,
ai >= 0
。
这是一个包含三个未知线性方程组 (a0, a1, ..., an)
以及 n+1 个线性不等式的系统。我们可以使用 scipy's linear programming module scipy.optimize.linprog
示例解决方案:
import random # generate points
import scipy.spatial as scispa # ConvexHull
import scipy.optimize as sciopt # linprog
from scipy.sparse import identity # identity matrix
points = [(random.randrange(0,100), random.randrange(0,100)) for _ in range(20)]
p1 = points[11]
hull = scispa.ConvexHull(points)
vertices = [points[i] for i in hull.vertices]
c = [1 for _ in vertices]
A = [[x for x,y in vertices], [y for x,y in vertices], [1 for _ in vertices]]
b = [p1[0], p1[1], 1]
s = sciopt.linprog(c, A_eq = A, b_eq = b)
>>> s.x
array([0.13393774, 0.06470577, 0.07367599, 0.09523271, 0.18924727,
0.26909487, 0.17410566])
>>> p1
(36, 57)
>>> [sum(a*xi for a,(xi,_) in zip(s.x, vertices)), sum(a*yi for a,(_,yi) in zip(s.x, vertices))]
[36.00000000719907, 57.00000000671608]
重要说明:我开始这个回答时警告说,如果平面中有超过 3 个顶点,或者维度 d 中有超过 d+1 个顶点,那么会有是我们方程组的多于 1 个解。上面的代码只有 returns 一种解决方案,意味着做出了任意选择。您可以控制这个任意选择:scipy.optimize.linprog
是一个优化工具;这里它最小化我们作为参数给出的向量 c
和解向量 s.x
的点积。这里我把c
的所有系数都设置为1,意思是linprog
会找到一个解中系数之和最小的解;尝试提供不同的向量 c
,你会得到不同的解决方案。