Multicellular simulations with shape and volume constraints using optimal transport

Resumen Técnico: Simulaciones Multicelulares con Restricciones de Forma y Volumen Mediante Transporte Óptimo

Enunciado del Problema
Los sistemas vivos, desde componentes intracelulares hasta tejidos macroscópicos y colonias bacterianas, exhiben dinámicas complejas impulsadas por la exclusión de volumen, las interacciones de forma y la deformabilidad. Comprender cómo estas restricciones conducen a la autoorganización macroscópica es un desafío fundamental en biología del desarrollo y biofísica. Los modelos computacionales existentes enfrentan una compensación entre precisión física, costo numérico y manejabilidad matemática. Los enfoques comunes incluyen:

• Sistemas de partículas puntuales: Matemáticamente simples pero a menudo carecen de una representación realista de la forma.
• Modelos de campo de fase y nivel de conjunto: Capaces de representar formas arbitrarias, pero pueden ser computacionalmente costosos y difíciles para interacciones de contacto.
• Autómatas celulares y modelos de red: Discretizan el espacio en vóxeles, limitando potencialmente la flexibilidad geométrica.
• Modelos de Voronoi y de vértices: Eficientes para estructuras similares a tejidos, pero a menudo dependen de restricciones suaves para el volumen o carecen de control sobre deformaciones de forma específicas.

Existe una brecha significativa en la modelización de sistemas con restricciones de volumen estrictas, formas dinámicas arbitrarias y altos conteos de partículas (de decenas a decenas de miles) dentro de un marco unificado que apoye tanto la deformabilidad individual como las interacciones colectivas.

Metodología
Los autores introducen un marco de simulación a mesoescala fundamentado en la teoría del transporte óptimo, utilizando específicamente teselaciones de Laguerre.

1. Modelo Estático (Teselación de Laguerre):
   El espacio ocupado por cada partícula $i$ se modela como una celda de Laguerre $L_i$, definida como el conjunto de puntos $x$ que satisfacen:
   $$L_i = \{x \in \Omega, c(x, x_i) - w_i \leq c(x, x_j) - w_j \text{ para todo } j\}$$
   Aquí, $x_i$ es la posición de la partícula, $v_i$ es su volumen objetivo, y $c(x, x_i)$ es una función de costo que codifica las propiedades de forma y deformabilidad de la partícula. Los pesos $w_i$ (potenciales de Kantorovich) se determinan de manera única para satisfacer la restricción estricta de volumen $|L_i| = v_i$.

   • Esta formulación generaliza los diagramas de Voronoi (donde $w_i=0$ y $c$ es la distancia euclidiana) para permitir límites curvos y formas no esféricas mediante la elección de funciones de costo arbitrarias (por ejemplo, costos de potencia, costos anisotrópicos).
   • El espacio vacío se modela como una partícula adicional con costo cero.

2. Modelo Dinámico:
   La evolución de las posiciones de las partículas $x_i$ sigue una ecuación de descenso de gradiente de primer orden:
   $$\dot{x}_i = -\tau_i \nabla_{x_i} T_c + F_i^{\text{point}} + F_i^{\text{surf}} + \dot{\xi}_i$$

   • Fuerza de Incompresibilidad: El término $-\nabla_{x_i} T_c$ actúa como una fuerza repulsiva derivada del costo total de transporte $T_c$. Asegura que las partículas no se superpongan y hace cumplir naturalmente las restricciones de volumen. A diferencia de las fuerzas por pares, esta fuerza depende de la geometría global de todas las superficies de contacto.
   • Fuerzas Externas/Interacción ($F^{\text{point}}$): Incluye autopropulsión, gravedad o atracción, similar a los modelos de partículas puntuales.
   • Fuerzas Superficiales ($F^{\text{surf}}$): Modela efectos de presión y tensión superficial a lo largo de las interfaces, permitiendo fenómenos como la clasificación celular y la adhesión.
   • Dinámica de Forma: Las formas pueden evolucionar mediante costos parametrizados (por ejemplo, elipsoides definidos por matrices de covarianza) o funciones de potencial de nivel de conjunto que responden a campos externos (por ejemplo, quimiotaxis).

3. Implementación Numérica:

   • El problema se resuelve como un problema de transporte óptimo semi-discreto en una red discreta (vóxeles).
   • Los autores emplean una estrategia multiescala de dividir y conquistar y un solucionador cuasi-Newton (L-BFGS-B) para calcular los potenciales de Kantorovich $w_i$ de manera eficiente.
   • La implementación está masivamente paralelizada en GPUs utilizando las bibliotecas KEOPS y GEOMLOSS, permitiendo simulaciones de hasta 50.000 partículas en 3D con un rendimiento de vanguardia.

Contribuciones Clave

• Marco Unificado: El método unifica conceptos de teselaciones de Voronoi, métodos de nivel de conjunto y modelos de vértices bajo el paraguas del transporte óptimo. Apoya formas arbitrarias, restricciones de volumen estrictas y mecanismos de interacción complejos dentro de un solo marco.
• Ejecución Estricta del Volumen: A diferencia de muchos modelos existentes que utilizan restricciones suaves o fuerzas de relajación, este enfoque garantiza matemáticamente las restricciones de volumen en todo momento mediante la optimización de los potenciales de Kantorovich.
• Escalabilidad: El solucionador acelerado por GPU permite la simulación de sistemas a gran escala (hasta 50.000 partículas) en 3D, lo cual es computacionalmente prohibitivo para muchos métodos tradicionales de elementos finitos o de campo de fase.
• Versatilidad en la Modelización: El marco modela con éxito diversos fenómenos, incluyendo:
  • Transiciones de fase fluido-sólido en partículas brownianas activas impulsadas por la deformabilidad.
  • Orden orientacional emergente en partículas enjambres con forma de varilla únicamente debido a la exclusión de volumen.
  • Migración quimiotáctica impulsada por deformaciones de forma.
  • Clasificación celular y fenómenos de englobamiento en agregados heterogéneos 3D, reproduciendo resultados consistentes con la Hipótesis de la Tensión Interfacial Diferencial (DITH).

Resultados
Los autores validan el marco a través de varias pruebas de referencia:

• Crecimiento y División: Un agregado 3D crece desde una sola célula hasta 50.000 células mediante división, llenando un dominio en aproximadamente un día de cómputo.
• Transiciones de Fase: En un sistema de partículas brownianas activas, el parámetro de deformabilidad $\alpha$ impulsa una transición de una fase fluida (comportamiento difusivo) a una fase sólida (atascada, empaquetamiento hexagonal). El modelo captura la difusión efectiva que desaparece y la saturación del desplazamiento cuadrático medio en la fase sólida.
• Alineación de Forma de Varilla: Las simulaciones de partículas autopropulsadas con forma de varilla demuestran la aparición espontánea de orden orientacional de largo alcance debido a deformaciones inducidas por colisiones, consistente con las teorías de física de polímeros.
• Clasificación Celular 3D: El modelo reproduce diversos patrones de clasificación (englobamiento, separación, tablero de ajedrez) en mezclas 3D de dos tipos de células con diferentes tensiones superficiales y suavidad, coincidiendo con predicciones teóricas y observaciones experimentales.
• Dinámica de Multitudes: El modelo simula la evacuación de multitudes, mostrando que las partículas deformables (tipo fluido) evacuan de manera segura, mientras que las partículas de esfera dura (tipo sólido) forman arcos estables y se atascan en las salidas.

Importancia y Afirmaciones
El artículo afirma introducir un enfoque novedoso para modelar sistemas multicelulares que cierra la brecha entre descripciones físicas detalladas y eficiencia computacional.

• Validación: La contribución principal es la validación de este nuevo marco frente a modelos computacionales clásicos, demostrando que puede reproducir resultados establecidos (por ejemplo, clasificación celular, transiciones de fase) mientras ofrece un manejo superior de las restricciones de volumen y la complejidad 3D.
• Relevancia Biológica: Los autores destacan el potencial del método para la modelización in silico de procesos de desarrollo complejos, como el desarrollo del blastocisto, donde tanto las deformaciones celulares individuales como los comportamientos de agregados densos son críticos.
• Extensión Teórica: El trabajo sugiere una nueva perspectiva de transporte óptimo para métodos existentes (partículas puntuales, nivel de conjunto, vértices) y proporciona una vía hacia la derivación de modelos continuos de grano grueso (aproximaciones de campo medio) para sistemas con interacciones de volumen, una pregunta ampliamente abierta en biología matemática.
• Accesibilidad: Los autores proporcionan una implementación en Python disponible gratuitamente, haciendo el método accesible para futuras investigaciones en física computacional y biología.

Los autores mantienen una postura modesta, señalando que, aunque el marco es nuevo, el enfoque inmediato fue la validación frente a modelos existentes. Identifican direcciones futuras, como extender la teoría de vértices activos a las teselaciones de Laguerre, refinar la derivación de límites continuos y mejorar el manejo de costos dinámicos complejos e interacciones superficiales para simulación en tiempo real.