Algoritmo para encontrar la masa agregada de "barra de granola" -Estructuras similares?

Soy un investigador de ciencias planetarias y un proyecto en el que estoy trabajando es la simulación de cuerpos N de los anillos de Saturno. El objetivo de este estudio en particular es observar cómo las partículas se agrupan bajo su propia gravedad y medir la masa agregada de los grupos frente a la velocidad media de todas las partículas en la célula. Estamos tratando de averiguar si esto puede explicar algunas observaciones hechas por la nave espacial Cassini durante el solsticio de verano de Saturno cuando se vieron grandes estructuras proyectando sombras en los anillos casi de borde. A continuación se muestra una captura de pantalla de cómo se ve cualquier paso de tiempo dado. (Cada partícula tiene 2 m de diámetro y la celda de simulación en sí tiene unos 700 m de diámetro).

El código que estoy usando ya escupe la velocidad media en cada paso de tiempo. Lo que necesito hacer es encontrar una manera de determinar la masa de partículas en los grupos y NO las partículas perdidas entre ellos. Conozco la posición, masa, tamaño, etc. de cada partícula, pero no sé fácilmente que, digamos, las partículas 30,000-40,000 junto con 102,000-105,000 forman una hebra que para el ojo humano es obvia.

Por lo tanto, el algoritmo que necesito escribir debería ser un código con la menor cantidad posible de parámetros ingresados ​​por el usuario (para replicabilidad y objetividad) que atraviese todas las posiciones de las partículas, descubra qué partículas pertenecen a los grupos y luego calcule masa. Sería genial si pudiera hacerlo para "cada" grupo / filamento en lugar de todo sobre la celda, pero no creo que realmente lo necesite para separarlos.

Lo único en lo que estaba pensando era en hacer algún tipo de cálculo de distancia N ² donde calcularía la distancia entre cada partícula y si, por ejemplo, las 100 partículas más cercanas estuvieran dentro de una cierta distancia, esa partícula se consideraría parte de un racimo. Pero eso parece bastante descuidado y esperaba que ustedes, programadores y CS, conocieran una solución más elegante.

Editado con mi solución: Lo que hice fue tomar una especie de enfoque de vecino / clúster más cercano y hacer primero la implementación rápida y sucia de N ² . Entonces, tome cada partícula, calcule la distancia a todas las demás partículas, y el umbral para un grupo o no fue si había N partículas dentro de la distancia d (dos parámetros que deben establecerse a priori , desafortunadamente, pero como dijeron algunos respuestas / comentarios, no iba a salirse con la suya al no tener algunos de esos).

Luego lo aceleré no ordenando distancias, sino simplemente haciendo una búsqueda de orden N e incrementando un contador para las partículas dentro de d , y eso aceleró las cosas por un factor de 6. Luego agregué un "árbol de programadores estúpido" (porque sé casi nada sobre códigos de árbol). Divido la celda de simulación en un número establecido de cuadrículas (mejores resultados cuando el tamaño de cuadrícula ≈7 d ) donde la cuadrícula principal se alinea con la celda, una cuadrícula se compensa por la mitad en x e y , y las otras dos se compensan por 1/4 pulg. ± x y ± y . El código luego divide las partículas en las cuadrículas, luego cada partícula N solo debe tener distancias calculadas a las otras partículas en esa celda.

Teóricamente, si este fuera un árbol real, debería obtener el orden N * log ( N ) en lugar de N ² velocidades. Llegué a algún lugar entre los dos, donde para un subconjunto de 50,000 partículas obtuve un aumento de 17x en la velocidad, y para una celda de 150,000 partículas, obtuve un aumento de 38x en la velocidad. 12 segundos para el primero, 53 segundos para el segundo, 460 segundos para una celda de 500,000 partículas. Esas son velocidades comparables a cuánto tiempo tarda el código en ejecutar la simulación 1 paso adelante, por lo que es razonable en este punto. Ah, y está completamente enhebrado, por lo que necesitará tantos procesadores como pueda.

Mi primera sugerencia es dividir su problema en dos problemas: primero, descubra lo que quiere y luego descubra cómo obtener eficientemente lo que desea. No puede obtener eficientemente algo que aún no ha definido. Pondré algunas ideas en esta respuesta que pueden ayudarlo a encontrar esta definición. Le sugiero que haga una implementación ineficiente de las ideas que le gustan primero, aplíquelas a algunos conjuntos de datos no demasiado grandes, evalúe los resultados a mano, adapte su definición y repita (posiblemente haciendo otra pregunta aquí), hasta que esté satisfecho con tu definición Después de eso, le sugiero que haga otra pregunta sobre cómo calcular eficientemente el resultado de su definición (si aún necesita ayuda).

Entonces, veamos qué correspondería con nuestra idea intuitiva de un 'capítulo'. Sus hebras parecen consistir en puntos distribuidos de manera más o menos uniforme, aunque debe verificar esto haciendo una imagen ampliada (del conjunto de datos original): la resolución de su imagen es demasiado baja para decir con certeza que los puntos realmente están distribuidos de manera aproximadamente uniforme . Asumiré que son para esta respuesta.

Una idea inicial puede ser mirar al vecino más cercano de cada punto. Seleccionemos un punto X, llamemos a su vecino más cercano Y y establezcamos D como la distancia entre X e Y. Luego miramos el círculo C alrededor de X con radio D * A, donde A es un parámetro de ajuste, digamos A = 3. Si X es parte de una cadena, esperamos que para cada punto Z en C, la distancia de Z a su vecino W más cercano sea aproximadamente la misma que D. Si es significativamente más corto, diga más que A (o tal vez algún otro parámetro B) entonces X aparentemente está cerca de puntos que están mucho más cerca uno del otro que X, por lo que X probablemente no sea parte de una cadena.

Sin embargo, este criterio no está completo. Solo da un criterio para detectar un 'borde' entre áreas densas con puntos y áreas menos densas con puntos. Todavía tenemos que agrupar los puntos en hebras.

Hay una característica en su imagen que muestra que esto no es simple. En la esquina inferior derecha de la imagen, hay un área relativamente grande con muchos puntos extraviados. Estos puntos perdidos están distribuidos de manera más o menos uniforme, por lo que si elimináramos todos los puntos en el filamento que lo rodea (y todos los demás puntos), ¡esperaríamos que cualquier algoritmo de detección de filamentos marque este conjunto de puntos perdidos como un filamento! Por lo tanto, debemos tener cuidado al hacer nuestros grupos.

Una idea puede ser hacer lo siguiente. Vamos a hacer un gráfico en estos puntos, donde los vértices son los puntos y los bordes significan que dos puntos tienen una densidad similar. Para cada punto, verificamos el criterio anterior. Si sale, conectamos X con un borde a todos los puntos en C. Si no sale, no agregamos ningún borde, y marcamos X como 'perdido'. Después de hacer esto para cada punto, consideramos el conjunto de componentes conectados. Deben consistir en un solo componente conectado (en el caso de su imagen, pero otros conjuntos de datos pueden tener múltiples) que consta de todos los puntos en los filamentos, más (potencialmente muchos) más componentes que consisten en puntos dispersos individuales y estos 'filamentos perdidos'. Sin embargo, estos filamentos extraviados tienen puntos que han sido marcados como 'extraviados', por lo que simplemente puede ignorar cualquier componente que contenga un punto que haya sido marcado como 'extraviados'.

Un peligro de esta idea es que puede tener una característica en la que la densidad de un filamento disminuye progresivamente a medida que se mueve a lo largo del filamento, hasta que la densidad es tan baja que es solo un conjunto de puntos perdidos. Como nuestro criterio es 'local', es posible que no detecte esto y marque estos puntos perdidos como parte de la cadena. No estoy seguro de si esto será un problema: supongo que la mayoría de los puntos perdidos deberían ser captados por el criterio, ya que los cambios en la densidad parecen bastante abruptos en su imagen.

Si se produce este problema, puede probar una alternativa a simplemente tomar los componentes conectados. Para cada punto X, calculamos la distancia a su vecino más cercano D (X). Comenzamos en el punto con un mínimo de D (X) y realizamos un BFS (o DFS , el orden no importa). Agregamos cualquier punto Y cuya D (Y) no sea mucho mayor que la D (X) (por un factor ajustable) con la que comenzamos. Si encontramos un punto Y que tiene D (Y) demasiado grande, eliminamos el borde (X, Y), marcamos Y como 'extraviados' y actuamos como si nunca visitáramos Y en nuestro BFS. Si se ajusta correctamente, esto debería evitar el problema que describo anteriormente.

Una idea alternativa para solucionar este problema es un poco más local: podría hacer un BFS y realizar un seguimiento de la D (X) más baja (uso D (X) como una medida de la densidad alrededor de un punto) encontrada como máximo 10 BFS-pasos antes, y si encontramos una Y que tiene D (Y) que es mucho más grande que esta D (X), hacemos lo mismo que la otra solución (potencial) que ofrecí.

Como descargo de responsabilidad: todas las ideas anteriores que pensé en este momento, realmente no sé si este problema en particular se ha estudiado antes, por lo que es posible que esté surgiendo tonterías. Simplemente pruebe las ideas (ya sean mías o propias) que le parezcan razonables y descubra si realmente funcionan, y solo entonces concéntrese en implementarlas de manera eficiente.

Con la descomposición modular , puede crear un árbol que contendrá todas las partículas como hojas y los nodos superiores los agruparán. En función de ese árbol, puede definir medidas que se apliquen a cada nodo desde la raíz hasta las hojas hacia abajo. Usted detiene este recorrido hacia abajo cuando las mediciones alcanzan los umbrales definidos por el usuario. Una de esas medidas puede ser la densidad del casco convexo de todas las partículas en un grupo.

Creo que buscas un algoritmo de agrupación de aprendizaje automático.

Esta página del kit de herramientas Python SciKit Learn tiene imágenes que sugieren que el algoritmo DBSCAN (Wikipedia) podría ser lo que estás buscando. Parece ideal ya que su parámetro de entrada es el tamaño del vecindario, mientras que la mayoría de los otros algoritmos de agrupación desean la cantidad de agrupaciones, lo que no sabría de antemano.

"Un algoritmo basado en la densidad para descubrir clústeres en grandes bases de datos espaciales con ruido" Ester, M., HP Kriegel, J. Sander y X. Xu, en las actas de la 2da Conferencia internacional sobre descubrimiento de conocimiento y minería de datos, Portland, OR , AAAI Press, págs. 226–231. 1996

He estado pensando en este problema. No soy un experto en física, así que tengan paciencia conmigo.

Parece que no es la distancia entre las partículas lo que cuenta para determinar los grupos. Es si los campos de gravedad se superponen o no.

Tome una partícula P y determine qué otras partículas tienen campos de gravedad superpuestos.

Luego toma uno de esos y haz lo mismo. Su objetivo no es encontrar todas las partículas en el grupo, sino encontrar sus límites.

Repita esto hasta que se encuentren todos los grupos.

Ahora regrese y determine la masa de los grupos. Habrás eliminado las partículas perdidas, y puedes usar los límites del grupo para encontrar la masa.

No estoy seguro de si esto ayuda, pero es todo lo que se me ocurre.

Al final de cada paso de tiempo, podría convertir los datos en un gráfico, calcular el árbol de expansión mínimo y luego comenzar a eliminar los bordes que exceden un cierto umbral. Eso debería darle grupos y una forma fácil de enumerar a través de las partículas en cada grupo.