ohtli 0.9.0

Cálculo y visualización de β-esqueletos para el modelado de redes de conexión entre sitios arqueológicos.

Ohtli

Librería para el cálculo y visualización de β-esqueletos para el modelado de redes de conexión entre sitios arqueológicos.

Instalación

pip install ohtli

O en modo desarrollo:

git clone https://github.com/tu-usuario/ohtli.git
cd ohtli
pip install -e ".[dev]"

Uso básico

import numpy as np
from ohtli import adjacency_matrix, build_graph, plot_skeleton
from ohtli import relative_asymmetry, control_value, plot_metric_bars
from ohtli.datasets import load_metro_cdmx

coords = np.array([[0,0],[1,0],[2,0],[1,1]], dtype=float)
names  = ["A", "B", "C", "D"]

# Matriz de adyacencia
adj = adjacency_matrix(coords, beta=1.0)

# Grafo de NetworkX
G = build_graph(coords, beta=2.0, labels=names)

# Visualización
ax, G = plot_skeleton(coords, beta=1.0, labels=names)

# Métricas
ra = relative_asymmetry(G)
cv = control_value(G)
plot_metric_bars(ra, metric_name="Asimetría relativa")

metro  = load_metro_cdmx()
coords = metro["coords"]
names  = metro["names"]
adj    = metro["adj"]
lines  = metro["lines"]

Estructura del paquete

ohtli/
├── geometry.py   # cálculo de lunas y vecindades
├── brute_force.py   # algoritmo fuerza bruta
├── delaunay.py # algoritmos deulaunay
├── datasets.py # datos del metro de la Ciudad de México
├── metrics.py    # asimetría relativa y valor de control
└── viz.py        # visualizaciones

El contraejemplo de Toussaint al algoritmo de Urquhart

Contexto

El algoritmo de Urquhart (1980) propone calcular el grafo de vecindad relativa (RNG, β = 2) en O(n log n) evaluando, para cada arista de la triangulación de Delaunay, únicamente los dos vértices que completan los triángulos que comparten esa arista.

Toussaint (1980) demostró que este algoritmo no es correcto para el RNG mediante el contraejemplo descrito a continuación.

---

El contraejemplo (Toussaint, Electronics Letters 1980)

Considera cinco puntos a, b, c, d, e con la siguiente disposición:

         d
        /
       / 
      a ─────────── b
                   /
                  c
    
    (e está dentro de la luna de (a,b),
     pero fuera del campo visual de Urquhart)

Condiciones geométricas precisas:

• a y b tienen la misma coordenada y (están en la misma altura).
• a y d tienen la misma coordenada x; b y c tienen la misma coordenada x.
• d tiene coordenada y ligeramente mayor que c, lo que hace que a y c sean vecinos de Voronoi y, por tanto, vecinos de Delaunay.
• e se encuentra dentro de la luna LUNE(a, b) con β = 2.
• c y d se encuentran fuera de la luna LUNE(a, b).

La triangulación de Delaunay de estos cinco puntos incluye la arista (a, b). Los únicos vecinos de Delaunay comunes a a y b son c y d.

El algoritmo de Urquhart, al evaluar la arista (a, b), solo prueba c y d. Como ambos están fuera de la luna, concluye que (a, b) ∈ RNG. Pero e — que sí está dentro de la luna — nunca es evaluado, y por tanto la arista se incluye incorrectamente en el grafo resultante.

Resultado más fuerte

Toussaint demuestra que el problema no se resuelve ampliando el conjunto de vecinos probados: cualquier algoritmo que decida la pertenencia de una arista al RNG basándose únicamente en subconjuntos de los vecinos de Delaunay de sus dos vértices está condenado a fallar. Siempre es posible construir una configuración donde el punto invalidador no sea vecino de Delaunay de ninguno de los dos vértices.

---

¿Qué calcula entonces el algoritmo de Urquhart?

El propio Urquhart reconoció el error en su respuesta (Electronics Letters, 1980) y precisó que su algoritmo calcula un grafo intermedio U tal que:

MST ⊆ RNG ⊆ U ⊆ GG ⊆ DT

donde GG es el Grafo de Gabriel. El algoritmo de Urquhart es exacto para el Grafo de Gabriel (β = 1) pero solo aproximado para el RNG (β = 2) y cualquier otro β ≠ 1.

---

Verificación estadística con ohtli

El siguiente código muestra que urquhart_approximate produce aristas extra (falsas) en comparación con beta_skeleton_delaunay (exacto) para β ≠ 1:

import numpy as np
from ohtli.delaunay import beta_skeleton_delaunay, urquhart_approximate

np.random.seed(42)
extra_edges = 0
trials = 100

for _ in range(trials):
    coords = np.random.rand(20, 2) * 10
    G_exact  = beta_skeleton_delaunay(coords, beta=2.0)
    G_approx = urquhart_approximate(coords, beta=2.0, warn=False)
    extra_edges += len(set(G_approx.edges()) - set(G_exact.edges()))

print(f"Aristas falsas en {trials} conjuntos: {extra_edges} total")
print(f"Promedio por conjunto: {extra_edges / trials:.2f}")
# Resultado típico: ~0.25 aristas falsas por conjunto de 20 puntos

Para β = 1 (Grafo de Gabriel) el resultado es exacto:

for _ in range(trials):
    coords = np.random.rand(20, 2) * 10
    G_exact  = beta_skeleton_delaunay(coords, beta=1.0)
    G_approx = urquhart_approximate(coords, beta=1.0, warn=False)
    assert set(G_exact.edges()) == set(G_approx.edges()), "Diferencia inesperada"

print("Gabriel: exacto en todos los casos ✓")

---

Cuándo usar cada función

Función	Complejidad	β soportado	Corrección
beta_skeleton_delaunay	O(n²)	β ≥ 1	Exacta
urquhart_approximate	O(n log n)	β ≥ 1	Exacta solo para β = 1 (GG)
build_graph (fuerza bruta)	O(n³)	cualquier β	Exacta

Referencia

• Kirkpatrick, D. G. y Radke, J. D. (1985). A Framework for Computational Morphology. Computational Geometry, pp. 217-248.
• Urquhart, R. B. (1980). Algorithms for Computation of Relative Neighbourhood Graph. Electronics Letters, 16(14), 556–557.
• Toussaint, G. T. (1980). Comment: Algorithms for Computing Relative Neighbourhood Graph. Electronics Letters, 16(22), 860.
• Urquhart, R. B. (1980). Reply. Electronics Letters, 16(22), 860–861.
• Toussaint, G. T. (1980). The Relative Neighbourhood Graph of a Finite Planar Set. Pattern Recognition, 12(4), 261–268.