Contoh titik yang diberikan secara stokastik dalam ruang 3D dengan jarak tetangga terdekat minimum dan kepadatan maksimum
Saya memiliki n
poin dalam ruang 3D. Saya ingin mengambil sampel secara stokastik subset titik dengan semua jarak tetangga terdekat lebih besar dari r
. Ukuran subset m
tidak diketahui, tetapi saya ingin titik sampelnya sepadat mungkin.
Ada pertanyaan serupa, tetapi semuanya tentang menghasilkan poin, bukan pengambilan sampel dari poin yang diberikan.
Hasilkan titik acak dalam ruang 3D dengan jarak tetangga terdekat minimum
Hasilkan titik acak 3-d dengan jarak minimum di antara masing-masing titik?
Katakanlah saya memiliki 300 titik 3D acak,
import numpy as np
n = 300
points = np.random.uniform(0, 10, size=(n, 3))
Apa cara tercepat untuk mendapatkan subset m
poin dengan jarak tetangga terdekat minimum r = 1
sambil memaksimalkan m
?
Jawaban
Mungkin ada skema perkiraan bikriteria yang efisien, tetapi mengapa repot-repot ketika pemrograman bilangan bulat rata-rata begitu cepat?
import numpy as np
n = 300
points = np.random.uniform(0, 10, size=(n, 3))
from ortools.linear_solver import pywraplp
solver = pywraplp.Solver.CreateSolver("SCIP")
variables = [solver.BoolVar("x[{}]".format(i)) for i in range(n)]
solver.Maximize(sum(variables))
for j, q in enumerate(points):
for i, p in enumerate(points[:j]):
if np.linalg.norm(p - q) <= 1:
solver.Add(variables[i] + variables[j] <= 1)
solver.EnableOutput()
solver.Solve()
print(len([i for (i, variable) in enumerate(variables) if variable.SolutionValue()]))
Ini bukan bagian yang besar secara optimal, tetapi harus dekat sementara tidak membutuhkan waktu yang lama, gunakan KDTree
untuk mengoptimalkan penghitungan jarak:
from scipy.spatial import KDTree
import numpy as np
def space_sample(n = 300, low = 0, high = 10, dist = 1):
points = np.random.uniform(low, high, size=(n, 3))
k = KDTree(points)
pairs = np.array(list(k.query_pairs(dist)))
def reduce_pairs(pairs, remove = []): #iteratively remove the most connected node
p = pairs[~np.isin(pairs, remove).any(1)]
if p.size >0:
count = np.bincount(p.flatten(), minlength = n)
r = remove + [count.argmax()]
return reduce_pairs(p, r)
else:
return remove
return np.array([p for i, p in enumerate(points) if not(i in reduce_pairs(pairs))])
subset = space_sample()
Menghapus node yang paling terhubung secara berulang tidak optimal (menyimpan sekitar 200 dari 300 poin), tetapi kemungkinan merupakan algoritme tercepat yang mendekati optimal (satu-satunya hal yang lebih cepat adalah menghapus secara acak). Anda mungkin bisa @njit
reduce_pairs
membuat bagian itu lebih cepat (akan coba jika saya punya waktu nanti).
Menguji jawaban @David Eisenstat dengan 30 poin yang diberikan:
from ortools.linear_solver import pywraplp
import numpy as np
def subset_David_Eisenstat(points, r):
solver = pywraplp.Solver.CreateSolver("SCIP")
variables = [solver.BoolVar("x[{}]".format(i)) for i in range(len(points))]
solver.Maximize(sum(variables))
for j, q in enumerate(points):
for i, p in enumerate(points[:j]):
if np.linalg.norm(p - q) <= r:
solver.Add(variables[i] + variables[j] <= 1)
solver.EnableOutput()
solver.Solve()
indices = [i for (i, variable) in enumerate(variables) if variable.SolutionValue()]
return points[indices]
points = np.array(
[[ 7.32837882, 12.12765786, 15.01412241],
[ 8.25164031, 11.14830379, 15.01412241],
[ 8.21790113, 13.05647987, 13.05647987],
[ 7.30031002, 13.08276009, 14.05452502],
[ 9.18056467, 12.0800735 , 13.05183844],
[ 9.17929647, 11.11270337, 14.03027534],
[ 7.64737905, 11.48906945, 15.34274827],
[ 7.01315886, 12.77870699, 14.70301668],
[ 8.88132414, 10.81243313, 14.68685022],
[ 7.60617372, 13.39792166, 13.39792166],
[ 8.85967682, 12.72946394, 12.72946394],
[ 9.50060705, 11.43361294, 13.37866092],
[ 8.21790113, 12.08471494, 14.02824481],
[ 7.32837882, 12.12765786, 16.98587759],
[ 8.25164031, 11.14830379, 16.98587759],
[ 7.30031002, 13.08276009, 17.94547498],
[ 8.21790113, 13.05647987, 18.94352013],
[ 9.17929647, 11.11270337, 17.96972466],
[ 9.18056467, 12.0800735 , 18.94816156],
[ 7.64737905, 11.48906945, 16.65725173],
[ 7.01315886, 12.77870699, 17.29698332],
[ 8.88132414, 10.81243313, 17.31314978],
[ 7.60617372, 13.39792166, 18.60207834],
[ 8.85967682, 12.72946394, 19.27053606],
[ 9.50060705, 11.43361294, 18.62133908],
[ 8.21790113, 12.08471494, 17.97175519],
[ 7.32837882, 15.01412241, 12.12765786],
[ 8.25164031, 15.01412241, 11.14830379],
[ 7.30031002, 14.05452502, 13.08276009],
[ 9.18056467, 13.05183844, 12.0800735 ],])
Jika jarak minimum yang diharapkan adalah 1:
subset1 = subset_David_Eisenstat(points, r=1.)
print(len(subset1))
# Output: 18
Periksa jarak minimum:
from scipy.spatial.distance import cdist
dist = cdist(subset1, subset1, metric='euclidean')
# Delete diagonal
res = dist[~np.eye(dist.shape[0],dtype=bool)].reshape(dist.shape[0],-1)
print(np.min(res))
# Output: 1.3285513450926985
Ubah jarak minimum yang diharapkan menjadi 2:
subset2 = subset_David_Eisenstat(points, r=2.)
print(len(subset2))
# Output: 10
Periksa jarak minimum:
from scipy.spatial.distance import cdist
dist = cdist(subset2, subset2, metric='euclidean')
# Delete diagonal
res = dist[~np.eye(dist.shape[0],dtype=bool)].reshape(dist.shape[0],-1)
print(np.min(res))
# Output: 2.0612041004376223