Cálculos Raster son muy eficaces para este tipo de cálculos, especialmente cuando usted tiene el tiempo de viaje de los datos.
Las entradas son los rásteres de la distancia a las ciudades, una trama por la ciudad. Por ejemplo, si usted está utilizando la distancia Euclídea, en primer lugar, calcular la distancia Euclidiana de la cuadrícula para cada ciudad. De lo contrario, calcular la "distancia" de cualquier manera que pueda, tal vez usando CostDistance o tal vez rasterizar una red de viajes de cálculo de tiempo.
El área asociada a una ciudad consiste de todos los puntos en los que la población de la ciudad, dividida por el cuadrado de la distancia, es más grande que la comparables de la población/pie^2 valores para todas las demás ciudades. Para prepararse para esto, utilice su distancia rejillas para calcular estos habitantes/pie^2 valores: que es un simple "mapa de álgebra de la operación". Ahora estás casi un hecho: la más alta Posición de comando crea una versión codificada de el resultado que usted desea. Su argumento es una lista de su población/pie^2 rejillas de y para los valores codificados, 1 se utiliza para la primera red en la lista, 2 para el segundo, y así sucesivamente.
Ejemplo
Los comentarios indican que esto puede ser un gran problema: 2300 o así que los centros, tal vez cubriendo un enorme (en todo el país?) de la zona, por lo que requiere una gran cuadrícula. No intente esto con Spatial Analyst: como un propósito general SIG raster, es demasiado lento.
Aquí hay un ejemplo R que implican una décima parte de los centros en un 1000 por 1667 cuadrícula. Tardó diez segundos para calcular. (Que sigue siendo un poco lento, pero es al menos un orden de magnitud más rápido que se podía hacer en ArcGIS 10.x.)
![Map]()
Los números en este mapa de identificar los centros, con valores de 1 a 230. Los puntos negros muestran las ubicaciones de los centros (con zonas en proporción a su población). Debido a que este está basado en la distancia Euclídea, los límites entre "territorios" son todas las partes de arcos circulares.
El cálculo se efectúa mediante el establecimiento de dos redes, una para el valor máximo de "gravedad" visto hasta ahora en cada celda y otro para el actual "dueño" de la célula, y el cálculo de una nueva cuadrícula para cada centro. Identifica todas las células donde el nuevo centro de la "gravedad" excede el máximo actual y las actualizaciones de la [gravedad] cuadrícula con el nuevo centro del valor y la [titular] cuadrícula con el identificador del nuevo propietario.
El tiempo de cálculo es directamente proporcional con el número de centros, por lo 2300 centros requeriría de dos minutos de cálculo. También es directamente proporcional con el número de filas y el número de columnas. Una descripción más detallada de la cuadrícula de, digamos, 10.000 filas por 16,667 columnas sería, por tanto, requieren alrededor de 160 minutos, o tres horas.
Para lograr estos moderadamente rápido velocidades, el algoritmo necesita memoria RAM para almacenar tres rejillas a la vez: [gravedad], [propietario], y los valores de cada centro a su vez (que son re-utilizados para cada centro). Suponiendo 16 bytes por valor de (doble precisión con algunos R sobrecarga), el cálculo de 10.000 por 16,667 rejilla requieren de 8 GB, que me confirmó con las pruebas. (Que cubra todo el territorio continental de estados unidos en una de 300 metros de resolución.) En esta resolución, el cálculo se toma 5 segundos por el centro. Extrapolando a 2300 centros de da tres horas.
Estar preparado para el contra-intuitiva de los resultados. (No estoy diciendo que esté mal, pero podría ser sorprendente.) Por ejemplo, cuando sólo dos ciudades están involucrados, el área asociada con el más pequeño de la ciudad será un perfecto disco: pero no va a estar centrada en la propia ciudad. Estará rodeado de un (infinito) área asociada con el más grande de la ciudad. Que significa que una pluralidad de personas que viven lo suficientemente lejos en el otro lado de la pequeña ciudad se prevé viajar a la gran ciudad, incluso cuando que consiste en pasar de la derecha a través de la pequeña ciudad.
Aquí es el R código utilizado para producir el ejemplo. Para utilizarlo en la práctica usted desea leer el centro de coordenadas (tal vez de un shapefile) y, después, escriba la asignación de cuadrícula en el disco para su posterior tratamiento y mejorar los resultados. Los sencillos procedimientos se ilustran (con código de trabajo) en otras R-basado en hilos aquí.
#
# Create population centers and their populations.
#
set.seed(17)
x.min <- 0; x.max <- 5000000 # Extent of easting coordinates
y.min <- 0; y.max <- 3000000 # Extent of northing coordinates
n <- 230 # Number of centers
center <- matrix(runif(n*2, min=c(x.min,y.min), max=c(x.max, y.max)),
ncol=2, nrow=n, byrow=TRUE)
pop <- rgamma(n, 5, 10^-2 / (4*n)) # Total population around 110M people
#
# Create row and column coordinate arrays.
#
n.rows <- 10^3
cellsize <- (y.max - y.min) / n.rows
n.cols <- ceiling((x.max - x.min) / cellsize)
x.max <- x.min + n.cols * cellsize # Assures square cells are used
y.0 <- seq(y.max-cellsize/2, y.min+cellsize/2, length.out=n.rows) #
x.0 <- seq(x.min+cellsize/2, x.max-cellsize/2, length.out=n.cols)
#
# For each center, compute a "gravity" grid (based on population times inverse
# squared Euclidean distance) and accumulate the maxima.
#
system.time(
{
#
# It is slightly more efficient to process the larger
# centers first.
#
i <- order(pop, decreasing=TRUE)
pop <- pop[i]
center <- center[i, , drop=FALSE]
#
# Initialize the output.
#
owner <- matrix(0, n.rows, n.cols)
gravity.max <- matrix(0, n.rows, n.cols)
#
# Loop over the centers, updating the output.
# (This is a "highest position" calculation.)
#
for (i in 1:n) {
r <- pop[i] / outer((y.0 - center[i,2])^2, (x.0 - center[i,1])^2, "+")
update <- which(r >= gravity.max)
gravity.max[update] <- r[update]
owner[update] <- i
}
})
#
# Plot the results.
#
library(raster)
par(mfrow=c(1,1))
plot(raster(log(gravity.max), xmn=x.min, xmx=x.max, ymn=y.min, ymx=y.max),
main="Log Gravity")
plot(raster(owner, xmn=x.min, xmx=x.max, ymn=y.min, ymx=y.max),
main="Gravity-based allocation")
points(center, pch=19, cex = 20*sqrt(pop/n/max(pop)))
Necesito para delimitar los límites de la ciudad basándose en los centros de'magnitude. CenterA de empleo/CenterB del empleo = (Distancia de centerA a punto de ruptura E/Distancia de centerB para romper el punto E)^2. AB=AE+BE. 
