Wygeneruj losowe punkty rozproszone jak miasta? - algorytm losowy, analiza skupień

Jak można wygenerować powiedzmy 1000 losowych punktów z rozkładem takim jak ten miasta w np. Ohio?
„Obawiam się, że nie mogę precyzyjnie zdefiniować„ rozproszonych miast ”; równomiernie rozmieszczone centra + małe chmury Gaussa są łatwe, ale ad hoc.
Dodano: Musi istnieć rodzina dystrybucji 2d z parametrem grupowania, który można zmieniać, aby dopasować dany zestaw punktów?

Odpowiedzi:

Może spojrzysz na Waltera Christallera Teoria miejsc centralnych. Chyba gdzieś musi być jakiś generator, albo sam możesz ugotować.

Zacznij od modelu elementów wodnych w twoimobszar docelowy (lub stwórz jeden, jeśli jest to miejsce wyimaginowane), a następnie połącz miasta w pobliżu skrzyżowań rzek, wzdłuż brzegów jezior, skrzyżowań jezioro-rzeka, a następnie wymyśl drogi łączące te duże miasta. Autostrady w rozsądnych odstępach, wolą być blisko skrzyżowań na autostradach, a teraz posypać kilka małych miast pustymi przestrzeniami.

Klastry Gaussa z rozmiarami klastrów Poissona działają dość dobrze.

Problem: generuj losowe punkty, które skupiają się mniej więcej jak dane miasta, np. W USA.

Podproblemy:
a) opisz klastry z rzędami liczb, tak aby „klaster A był jak klaster B” upraszcza do „klastrów (A) jest jak„ klaster zbiorczy (B) ”.
Uruchamianie N = 100, a następnie 1000 punktów przez fcluster poniżej, z ncluster = 25, daje

N 100 ncluster 25: 22 + 3  r 117
sizes: av 4     10   9   8   7   6   6   5   5   4   4   4 ...
radii: av 117  202 198 140 134  64  62  28 197 144 148 132 ...

N 1000 cluster 25: 22 + 3  r 197
sizes: av 45  144 139 130  85  84  69  63  43  38  33  30  ...
radii: av 197  213 279 118 146 282 154 245 212 243 226 235 ...

b) znajdź kombinację losowych generatorów z 2 lub 3 parametrami które mogą być zmieniane w celu generowania różnych klastrów.
Klastry Gaussa z rozmiarami klastrów Poissona mogą dość dobrze pasować do grupowania miast:

def randomclusters( N, ncluster=25,  radius=1, box=box ):
""" -> N 2d points: Gaussian clusters, Poisson cluster sizes """
pts = []
lam = eval( str( N // ncluster ))
clustersize = lambda: np.random.poisson(lam - 1) + 1
# poisson 2:  14  27  27  18   9   4  %
# poisson 3:   5  15  22  22  17  10  %
while len(pts) < N:
u = uniformrandom2(box)
csize = clustersize()
if csize == 1:
pts.append( u )
else:
pts.extend( inbox( gauss2( u, radius, csize )))
return pts[:N]


# Utility functions --

import scipy.cluster.hierarchy as hier

def fcluster( pts, ncluster, method="average", criterion="maxclust" ):
""" -> (pts, Y pdist, Z linkage, T fcluster, clusterlists)
ncluster = n1 + n2 + ... (including n1 singletons)
av cluster size = len(pts) / ncluster
"""
# Clustering is pretty fast:
# sort pdist, then like Kruskal"s MST, O( N^2 ln N )
# Many metrics and parameters are possible; these satisfice.
pts = np.asarray(pts)
Y = scipy.spatial.distance.pdist( pts )  # N*(N-1)/2
Z = hier.linkage( Y, method )  # N-1, like mst
T = hier.fcluster( Z, ncluster, criterion=criterion )
clusters = clusterlists(T)
return (pts, Y, Z, T, clusters)

def clusterlists(T):
""" T = hier.fcluster( Z, t ) e.g. [a b a b c a]
-> [ [0 2 5] [1 3] ] sorted by len, no singletons [4]
"""
clists = [ [] for j in range( max(T) + 1 )]
for j, c in enumerate(T):
clists[c].append( j )
clists.sort( key=len, reverse=True )
n1 = np.searchsorted(  map( len, clists )[::-1], 2 )
return clists[:-n1]

def radius( x ):
""" rms |x - xmid| """
return np.sqrt( np.mean( np.var( x, axis=0 )))
# * 100  # 1 degree lat/long ~ 70 .. 111 km

W Javie jest to zapewnione przez new Random().nextGaussian(). Ponieważ źródło java jest dostępne, możesz na nie spojrzeć:

synchronized public double nextGaussian() {
// See Knuth, ACP, Section 3.4.1 Algorithm C.
if (haveNextNextGaussian) {
haveNextNextGaussian = false;
return nextNextGaussian;
} else {
double v1, v2, s;
do {
v1 = 2 * nextDouble() - 1; // between -1 and 1
v2 = 2 * nextDouble() - 1; // between -1 and 1
s = v1 * v1 + v2 * v2;
} while (s >= 1 || s == 0);
double multiplier = StrictMath.sqrt(-2 * StrictMath.log(s)/s);
nextNextGaussian = v2 * multiplier;
haveNextNextGaussian = true;
return v1 * multiplier;
}
}

Tworzenie 30000 domów za pomocą

x = r.nextGaussian() * rad/4 + rad;
y = r.nextGaussian() * rad/4 + rad;

daje to piękne miasto: