wake-up-neo.com

Finden Sie heraus, ob sich ein Punkt innerhalb einer konvexen Hülle für eine Menge von Punkten befindet, ohne die Hülle selbst zu berechnen

Wie lässt sich am einfachsten testen, ob sich ein Punkt P innerhalb einer konvexen Hülle befindet, die aus einer Menge von Punkten X besteht? 

Ich möchte einen Algorithmus, der in einem hochdimensionalen Raum arbeitet (z. B. bis zu 40 Dimensionen), der die konvexe Hülle nicht explizit berechnet. Irgendwelche Ideen?

42
dimi

Das Problem kann gelöst werden, indem ein realisierbarer Punkt eines linearen Programms ermittelt wird. Wenn Sie sich für die gesamten Details interessieren, anstatt nur eine LP in einen vorhandenen Solver einzubinden, würde ich empfehlen, Kapitel 11.4 in Boyd und Vandenberghes hervorragendes Buch über die konvexe Optimierung zu lesen.

Setzen Sie A = (X[1] X[2] ... X[n]), dh die erste Spalte ist v1, die zweite v2 usw.

Lösen Sie das folgende LP-Problem.

minimize (over x): 1
s.t.     Ax = P
         x^T * [1] = 1
         x[i] >= 0  \forall i

woher 

  1. x^T ist die Umsetzung von x
  2. [1] ist der all-1-Vektor.

Das Problem hat eine Lösung, wenn sich der Punkt in der konvexen Hülle befindet.

22
user1071136

Der Punkt liegt genau dann außerhalb der konvexen Hülle der anderen Punkte, wenn die Richtung aller Vektoren von diesen zu diesen anderen Punkten auf weniger als einer Hälfte eines Kreises/einer Kugel/Hypersphäre liegt.

Hier ist eine Skizze für die Situation von zwei Punkten, einem blauen innerhalb der konvexen Hülle (grün) und dem roten außerhalb:

 Sketch of the points

Für die rote gibt es Halbierungen des Kreises, so dass die Vektoren vom Punkt zu den Punkten auf der konvexen Hülle nur eine Hälfte des Kreises schneiden. Für den blauen Punkt ist es nicht möglich, einen solchen zu finden Halbierung.

20
Svante

Sie müssen die konvexe Hülle nicht selbst berechnen, da sie in mehrdimensionalen Räumen ziemlich störend erscheint. Es gibt eine bekannte Eigenschaft von konvexen Hüllen :

Jeder Vektor (Punkt) v in der konvexen Hülle von Punkten [v1, v2, .., vn] kann als sum(ki*vi) dargestellt werden, wobei 0 <= ki <= 1 und sum(ki) = 1. Dementsprechend wird kein Punkt außerhalb der konvexen Hülle eine solche Darstellung haben. 

Im m-dimensionalen Raum erhalten wir die Menge m linearer Gleichungen mit n unknowns.

edit
Ich bin nicht sicher über die Komplexität dieses neuen Problems im Allgemeinen, aber für m = 2 scheint es linear zu sein. Vielleicht korrigiert mich jemand mit mehr Erfahrung auf diesem Gebiet. 

9
Nikita Rybak

Obwohl der ursprüngliche Beitrag vor drei Jahren war, kann diese Antwort vielleicht immer noch hilfreich sein. Der Gilbert-Johnson-Keerthi (GJK) -Algorithmus ermittelt den kürzesten Abstand zwischen zwei konvexen Polytopen, von denen jedes als konvexe Hülle eines Generatorsatzes definiert ist - insbesondere muss die konvexe Hülle selbst nicht berechnet werden. In einem speziellen Fall, zu dem man gefragt wird, ist eines der Polytope nur ein Punkt. Warum nicht mit dem GJK-Algorithmus den Abstand zwischen P und der konvexen Hülle der Punkte X berechnen? Wenn dieser Abstand 0 ist, liegt P innerhalb von X (oder zumindest an seiner Grenze). Eine GJK-Implementierung in Octave/Matlab mit dem Namen ClosestPointInConvexPolytopeGJK.m sowie der zugehörige Code sind unter http://www.99main.com/~centore/MunsellAndKubelkaMunkToolbox/MunsellAndKubelkaMunkToolbox.html verfügbar. Eine einfache Beschreibung des GJK-Algorithmus finden Sie in Kap. 2 einer Arbeit unter http://www.99main.com/~centore/ColourSciencePapers/GJKinConstrainedLeastSquares.pdf . Ich habe den GJK-Algorithmus für einige sehr kleine Mengen X im 31-dimensionalen Raum verwendet und hatte gute Ergebnisse. Wie die Leistung von GJK mit den von anderen empfohlenen linearen Programmiermethoden verglichen wird, ist ungewiss (obgleich Vergleiche interessant wären). Das GJK-Verfahren vermeidet es, die konvexe Hülle zu berechnen oder die Hülle in Form linearer Ungleichungen auszudrücken, was beides zeitaufwändig sein kann. Hoffe, diese Antwort hilft.

2
Paul Centore

Ich hatte das gleiche Problem mit 16 Dimensionen. Da auch qhull nicht richtig funktionierte, da zu viele Gesichter erzeugt werden mussten, entwickelte ich einen eigenen Ansatz, indem ich testete, ob zwischen dem neuen Punkt und den Referenzdaten eine trennende Hyperebene gefunden werden kann (ich nenne dies "HyperHull";). . 

Das Problem, eine trennende Hyperebene zu finden, kann in ein konvexes quadratisches Programmierproblem umgewandelt werden (siehe: SVM ). Ich habe dies in Python mit cvxopt mit weniger als 170 Codezeilen (einschließlich E/A) durchgeführt. Der Algorithmus arbeitet ohne Änderung in einer beliebigen Dimension, auch wenn das Problem besteht, dass je höher die Dimension ist, desto höher die Anzahl der Punkte auf der Hülle (siehe: Auf der konvexen Hülle von zufälligen Punkten in einem Polytop ). Da der Rumpf nicht explizit konstruiert, sondern nur geprüft wird, ob sich ein Punkt im Inneren befindet oder nicht, hat der Algorithmus sehr große Vorteile in höheren Dimensionen im Vergleich zu z. schneller Rumpf.

Dieser Algorithmus kann "natürlich" parallelisiert werden, und die Geschwindigkeit sollte der Anzahl der Prozessoren entsprechen.

2
Michael Hecht

Sind Sie bereit, eine heuristische Antwort zu akzeptieren, die normalerweise funktionieren sollte, aber nicht garantiert ist? Wenn Sie dann sind, könnten Sie diese zufällige Idee versuchen.

Sei f(x) der Würfel der Entfernung zu P mal der Anzahl der Dinge in X, minus der Summe der Würfel der Entfernung zu allen Punkten in X. Beginne irgendwo zufällig und benutze einen Hügel Steigern Sie den Algorithmus, um f(x) für x in einer Kugel zu maximieren, die sehr weit von P entfernt ist. Abgesehen von degenerierten Fällen, wenn P nicht in der konvexen Hülle ist, sollte dies eine sehr gute Wahrscheinlichkeit haben, die Normalen zu finden Hyperebene, von der sich P auf einer Seite befindet, und alles in X ist auf der anderen Seite von.

1
btilly

Eine Überprüfung, bei der geprüft wird, ob sich ein Punkt in einem Rumpfraum befindet, wobei scipy.optimize.minimize verwendet wird.

Basierend auf der Antwort von Benutzer1071136.

Wenn Sie die konvexe Hülle berechnen, geht es viel schneller. Ich habe also ein paar Zeilen für die Leute hinzugefügt, die dies tun möchten. Ich wechselte von Graham Scan (nur 2D) zum Scipy-Qhull-Algorithmus. 

scipy.optimize.minimiere Dokumentation:
https://docs.scipy.org/doc/scipy/reference/optimize.nonlin.html

import numpy as np
import scipy.optimize
import matplotlib.pyplot as plt
from scipy.spatial import ConvexHull


def hull_test(P, X, use_hull=True, verbose=True, hull_tolerance=1e-5, return_hull=True):
    if use_hull:
        hull = ConvexHull(X)
        X = X[hull.vertices]

    n_points = len(X)

    def F(x, X, P):
        return np.linalg.norm( np.dot( x.T, X ) - P )

    bnds = [[0, None]]*n_points # coefficients for each point must be > 0
    cons = ( {'type': 'eq', 'fun': lambda x: np.sum(x)-1} ) # Sum of coefficients must equal 1
    x0 = np.ones((n_points,1))/n_points # starting coefficients

    result = scipy.optimize.minimize(F, x0, args=(X, P), bounds=bnds, constraints=cons)

    if result.fun < hull_tolerance:
        hull_result = True
    else:
        hull_result = False

    if verbose:
        print( '# boundary points:', n_points)
        print( 'x.T * X - P:', F(result.x,X,P) )
        if hull_result: 
            print( 'Point P is in the hull space of X')
        else: 
            print( 'Point P is NOT in the hull space of X')

    if return_hull:
        return hull_result, X
    else:
        return hull_result

Testen Sie einige Beispieldaten:

n_dim = 3
n_points = 20
np.random.seed(0)

P = np.random.random(size=(1,n_dim))
X = np.random.random(size=(n_points,n_dim))

_, X_hull = hull_test(P, X, use_hull=True, hull_tolerance=1e-5, return_hull=True)

Ausgabe:

# boundary points: 14
x.T * X - P: 2.13984259782e-06
Point P is in the hull space of X

Visualisieren Sie es:

rows = max(1,n_dim-1)
cols = rows
plt.figure(figsize=(rows*3,cols*3))
for row in range(rows):
    for col in range(row, cols):
        col += 1
        plt.subplot(cols,rows,row*rows+col)
        plt.scatter(P[:,row],P[:,col],label='P',s=300)
        plt.scatter(X[:,row],X[:,col],label='X',alpha=0.5)
        plt.scatter(X_hull[:,row],X_hull[:,col],label='X_hull')
        plt.xlabel('x{}'.format(row))
        plt.ylabel('x{}'.format(col))
plt.tight_layout()

 Visualization of hull test

1
cnash