Accès aux soins obstétricaux d'urgence sous blocage routier, Haïti

La question posée

Santé Commune Initiative est une ONG fictive spécialisée en santé reproductive. Deux ans plus tôt, elle avait identifié les communes haïtiennes où concentrer un programme de santé sexuelle et reproductive, à partir de la demande potentielle et de la densité des infrastructures.

Une femme en urgence obstétricale doit encore atteindre une structure capable de l'opérer, ce qui suppose des routes praticables et des hôpitaux ouverts. Le ciblage initial ne mesurait ni l'un ni l'autre.

Ce second travail mesure donc autre chose : quelle part de la demande obstétricale haïtienne se situe à moins d'une heure d'un hôpital, et ce que devient cette part quand des routes sont coupées ou que des structures ferment.

En un coup d'œil

Sous la définition la plus stricte de ce qu'est un hôpital, un peu moins de huit femmes sur dix vivent à moins d'une heure d'une structure. Le réseau tient cette couverture sur très peu de tronçons.

Les six sources

Toutes publiques. Les volumes ci-dessous sortent des journaux de la chaîne de traitement ; ils diffèrent parfois des chiffres annoncés par les producteurs.

Aucune ne se charge telle quelle. Le millésime COD-AB décide si la colonne des noms de commune s'appelle ADM2_FR ou adm2_name. Les exports OCHA alternent entre UTF-8 et CP1252 selon la date de génération, ce qui suffit à faire échouer une lecture. La population féminine de 15 à 49 ans n'existe pas comme colonne dans le COD-PS et doit être reconstituée en sommant sept tranches quinquennales. Le code détecte les noms de colonnes par motif plutôt que par égalité stricte, et tente trois encodages avant d'abandonner.

def lire_csv_tolerant(chemin: Path) -> pd.DataFrame:
    """Les exports OCHA alternent entre UTF-8 et CP1252 selon le millésime."""
    for encodage in ("utf-8", "cp1252", "latin-1"):
        try:
            return pd.read_csv(chemin, encoding=encodage)
        except UnicodeDecodeError:
            continue
    raise UnicodeDecodeError(f"encodage non reconnu pour {chemin}")

Nettoyer, qualifier, redresser

Trois opérations séparent les fichiers bruts des couches exploitables. Les choix faits à ce stade se retrouvent dans tous les chiffres qui suivent.

Qualifier les 2 073 points de santé

OpenStreetMap mélange dans un même fichier hôpitaux, cliniques, pharmacies, laboratoires et cabinets vétérinaires. Le tri se fait sur deux attributs qui ne concordent pas toujours, amenity et healthcare.

L'ordre des tests compte. L'exclusion passe en premier : une pharmacie reste hors périmètre même si elle porte par ailleurs une étiquette de soins, parce qu'aucune césarienne n'y sera pratiquée. Viennent ensuite trois niveaux, du plus sûr au plus incertain.

def classer(a: str, h: str) -> str:
    # Une pharmacie ou un laboratoire sort du périmètre même s'il porte par
    # ailleurs une étiquette de soins.
    if a in TAGS_HORS_PERIMETRE or h in TAGS_HORS_PERIMETRE:
        return "hors_perimetre"
    if h == TAG_SONUC:
        return "SONUC"        # les deux etiquettes concordent
    if a == TAG_SONUC:
        return "HOPITAL_NC"   # etablissement reel, nature incertaine
    if a in TAGS_SONUB or h in TAGS_SONUB:
        return "SONUB"        # premiere ligne, sans cesarienne
    return "hors_perimetre"

Ces 984 points de la catégorie intermédiaire sont l'héritage de la cartographie d'urgence de 2010 : dispensaires et centres de santé saisis comme hôpitaux par des contributeurs qui documentaient vite. L'étiquette healthcare, plus tardive et mieux tenue, sert donc de filtre de confirmation.

À 60 minutes, la couverture passe de 78,1 % avec les 198 structures confirmées à 94,7 % en y ajoutant les 984 incertaines ; à 30 minutes, l'écart atteint 33,6 points. Le choix d'une étiquette OpenStreetMap déplace le résultat plus que n'importe quel paramètre de modélisation. L'analyse retient la définition stricte et publie les deux autres.

Rattraper les cellules du littoral

Chaque cellule de population doit connaître sa commune, pour recevoir la bonne part de femmes en âge de procréer. La jointure spatiale échoue sur les cellules côtières : le trait de côte des limites administratives est généralisé, celui de la grille de population ne l'est pas, et des cellules habitées tombent quelques centaines de mètres en mer.

Les perdre reviendrait à effacer des populations littorales, précisément celles dont l'accès routier est souvent le plus contraint. Elles sont donc rattachées à la commune la plus proche par sjoin_nearest plutôt que jetées.

Redresser WorldPop sur le recensement

WorldPop date de 2020 et décrit la répartition de la population à l'intérieur de chaque commune ; le COD-PS fixe l'effectif communal de 2024 sans le localiser. Les cellules sont donc multipliées, commune par commune, par le facteur qui ramène leur somme sur le total 2024.

Ce facteur vaut 1,24 en médiane, soit quatre années de croissance démographique. Son étendue va de 0,10 à 2,68 : dans 31 communes, WorldPop sous-estime de plus de moitié ; dans 6, il surestime d'un tiers. La couche finale hérite là de l'incertitude du COD-PS. Le code journalise médiane et étendue à chaque exécution, et vérifie que le redressement boucle à 1 % près.

totaux_2020 = dedans.groupby("pcode")["pop_2020"].transform("sum")
cible = dedans["pcode"].map(communes.set_index("pcode")["pop_totale"])
facteur = np.where(totaux_2020 > 0, cible / totaux_2020, 0.0)
dedans["pop"] = dedans["pop_2020"] * facteur

total, attendu = dedans["pop"].sum(), communes["pop_totale"].sum()
exiger(
    abs(total - attendu) / attendu < 0.01,
    f"le redressement ne boucle pas : {total:,.0f} contre {attendu:,.0f}",
)

Après ces trois opérations : 8 259 cellules, 11,90 millions d'habitants, 3,25 millions de femmes de 15 à 49 ans.

D'une pile de polylignes à un graphe

Chaque polyligne OpenStreetMap connaît sa géométrie et ignore qu'elle en croise d'autres. Sans cette information, aucun plus court chemin ne se calcule.

Chaque paire de points consécutifs devient une arête, et les points de forme identiques deviennent un même sommet. Ce dédoublonnage soude le réseau : deux tronçons qui se rejoignent partagent alors un nœud, et un itinéraire peut passer de l'un à l'autre.

# Toutes les positions sauf le dernier point de chaque polyligne ouvrent une
# arête vers la position suivante.
debut_arete = np.ones(len(coords), dtype=bool)
debut_arete[offsets[1:] - 1] = False
idx_u = np.nonzero(debut_arete)[0]
idx_v = idx_u + 1

# Dédoublonnage des sommets au mètre.
cle = np.round(coords).astype(np.int64)
sommets, inverse = np.unique(cle, axis=0, return_inverse=True)
u = inverse.ravel()[idx_u]
v = inverse.ravel()[idx_v]

L'arrondi au mètre pèse moins que prévu

Le code arrondit les coordonnées avant de dédoublonner, au motif que la reprojection vers l'UTM 18 Nord passe par de l'arithmétique flottante et peut séparer de quelques millimètres deux extrémités identiques à la source.

Mesuré, l'effet est marginal. Sur 1 458 558 points de forme du réseau carrossable, le dédoublonnage exact fusionne 128 461 sommets, l'arrondi en ajoute 285. Soit 0,22 % des fusions.

La topologie vient donc d'OpenStreetMap, où les tronçons qui se croisent partagent déjà un nœud. L'arrondi rattrape ce que le flottant a décollé. Il reste utile, il coûte une ligne et évite des micro-coupures indétectables, mais le décrire comme ce qui rétablit la topologie serait faux.

Le graphe compte 1 329 812 sommets, 1 360 702 arêtes et 36 755 km de routes.

Les vitesses de parcours

Le poids d'une arête est un temps de parcours, calculé à partir d'une vitesse attribuée par classe de route. Après la définition de l'offre, ce paramètre est celui qui déplace le plus les résultats.

Les valeurs retenues décrivent des vitesses effectives de porte à porte, sensiblement plus basses que les vitesses réglementaires. Le réseau haïtien combine revêtement dégradé, franchissements de ravines non bitumés et traversées d'agglomération sans voie réservée. Les vitesses réglementaires produiraient des isochrones optimistes de 30 à 50 %, biais classique des analyses d'accessibilité fondées sur OpenStreetMap.

Dans l'aire métropolitaine de Port-au-Prince, ces vitesses tombent de 45 % supplémentaires. La zone cumule densité, marchés de rue occupant la chaussée et carrefours non régulés, et les vitesses observées en journée y descendent couramment sous 15 km/h sur les axes secondaires. La réduction s'applique si le milieu de l'arête tombe dans la zone, ce qui évite de découper les tronçons à la frontière.

# Congestion urbaine : appliquée si le milieu de l'arête tombe dans la ZMPP.
mx = (coords[idx_u, 0] + coords[idx_v, 0]) / 2
my = (coords[idx_u, 1] + coords[idx_v, 1]) / 2
dans_zmpp = shapely.contains_xy(zmpp, mx, my)
vitesses = np.where(dans_zmpp, vitesses * FACTEUR_CONGESTION_ZMPP, vitesses)

temps = longueurs / (vitesses * 1000.0 / 60.0)

Un jeu de vitesses ne se défend pas valeur par valeur mais par un test d'ensemble. Celui retenu : Port-au-Prince vers Saint-Marc, environ 95 km par la RN1, doit ressortir autour de deux heures, ce qui correspond aux temps rapportés hors période de blocage.

Les classes piétonnes sont exclues du graphe. Un transfert obstétrical n'emprunte ni sentier ni escalier.

La demande obstétricale

La couche de demande sort du redressement décrit plus haut, avec une pondération supplémentaire : chaque cellule est multipliée par la part communale des femmes de 15 à 49 ans, celle reconstituée en sommant les sept tranches quinquennales du COD-PS.

Deux arbitrages fixent la résolution. La grille WorldPop à 100 m est agrégée par blocs de dix, soit environ 1 km, parce que router depuis chaque cellule de 100 m serait coûteux et faussement précis au regard du réseau routier. Une cellule agrégée portant moins de dix habitants est écartée comme bruit de désagrégation.

Restent 8 259 cellules à router vers 198 hôpitaux.

Le calcul d'accessibilité

Résolu de front, le problème demande 8 259 recherches de plus court chemin sur un graphe de 1,3 million de sommets, à recommencer intégralement à chaque scénario simulé.

Le graphe étant non orienté, le temps d'une structure vers un nœud vaut aussi le temps de ce nœud vers la structure. On peut donc inverser le sens du calcul : partir simultanément des 198 hôpitaux et propager vers le réseau. Un seul parcours donne, en tout point, le temps vers la structure la plus proche.

matrice = graphe.matrice(aretes_coupees)
resultat = dijkstra(
    matrice,
    directed=False,
    indices=np.unique(sommets_offre),   # toutes les structures a la fois
    min_only=True,                      # on ne garde que le minimum
    return_predecessors=avec_chemins,
)
...
minutes = minutes_acces + distances[sommets_demande]

min_only=True fait la différence : SciPy ne conserve pour chaque nœud que la distance à la source la plus proche, au lieu de matérialiser une matrice de 198 lignes sur 1,3 million de colonnes.

Le dernier kilomètre

Une cellule de population tombe rarement sur un nœud routier. Le segment terminal, du centre de la cellule au nœud le plus proche, se parcourt à pied ou en moto-taxi sur piste. Il est compté à 4 km/h, et la distance à vol d'oiseau est multipliée par 1,3 puisqu'on ne marche pas en ligne droite.

Au-delà de 10 km à vol d'oiseau, la cellule est déclarée non rattachable. Lui attribuer une marche de plusieurs heures produirait un chiffre calculable qui ne décrit aucun trajet réel.

distance, sommet = graphe.arbre().query(np.column_stack([x, y]), k=1)
minutes = distance * FACTEUR_DETOUR / (VITESSE_ACCES_KMH * 1000.0 / 60.0)
minutes = np.where(distance > DISTANCE_MAX_RATTACHEMENT_M, np.inf, minutes)
return sommet.astype("int64"), minutes

Deux tronçons OpenStreetMap peuvent relier exactement les mêmes sommets. Si la matrice creuse additionne les doublons, le temps entre ces deux sommets devient la somme des deux tronçons au lieu du plus rapide, et rien dans les sorties ne le signale. Le code regroupe donc les arêtes parallèles et retient le minimum, une fois pour toutes, avant les simulations.

Charge et criticité

Protéger en priorité les tronçons les plus fréquentés produit un mauvais classement.

La démarche se fait en deux temps. D'abord une présélection par la charge : en remontant l'itinéraire de chaque cellule vers son hôpital, on cumule sur chaque tronçon la demande qui l'emprunte. Ensuite une mesure directe du dommage : pour chaque tronçon présélectionné, l'accessibilité du pays entier est recalculée sans lui.

for rang, idx in enumerate(candidats, start=1):
    masque = np.zeros(len(graphe.u), dtype=bool)
    masque[idx] = True
    minutes, _, _ = temps_vers_offre(
        graphe, sommets_sonuc, sommets_demande, minutes_acces, aretes_coupees=masque
    )
    decroche = couvert_avant & (minutes > SEUIL_REFERENCE)
    perdues = np.subtract(minutes, reference, out=np.zeros_like(reference), where=comparable)

Une contrainte d'espacement s'impose dès la présélection. Les arêtes sont les segments entre points de forme d'une polyligne, si bien que les cinquante tronçons les plus chargés du pays décrivent deux ou trois kilomètres du même axe. Sans espacement minimal, on évaluerait cent cinquante fois le même barrage. La règle impose 3 km entre deux candidats.

La configuration autorise 150 candidats. Après filtrage sur les communes exposées, les classes d'axes et la charge non nulle, puis espacement, 18 emplacements distincts subsistent et sont évalués un par un. Le plafond n'est jamais atteint : sur ce réseau, les points de coupure qui comptent se dénombrent en dizaines.

La corrélation de rang de Spearman entre les deux classements vaut 0,44. Le tronçon le plus chargé du pays porte 93 995 femmes et sa coupure coûte 52 033 minutes pondérées ; le plus coûteux en porte moitié moins, 47 243, pour 1 210 501 minutes, vingt-trois fois plus.

L'explication est géométrique. Une artère fréquentée d'un quartier maillé se contourne, sa coupure coûte quelques minutes. Un tronçon modeste qui constitue le seul franchissement d'une ravine coûte des heures. La charge compte les femmes dont l'itinéraire emprunte le tronçon, sans rien dire des détours possibles.

Le dommage est lui-même très concentré. Le tronçon le plus coûteux fait basculer 19 689 femmes au-delà d'une heure, soit 85 % de la perte totale des quinze coupures testées. À partir du sixième, chaque coupure supplémentaire n'ajoute presque rien.

Cinq scénarios

Trois scénarios portent sur le réseau routier, un sur l'offre de soins, un cumule les deux.

Séparer le choc de réseau du choc d'offre décide de l'affectation des moyens : réhabiliter des routes et maintenir des structures ouvertes ne se financent pas sur la même ligne budgétaire.

Les deux ne creusent pas les mêmes indicateurs. Fermer les 77 hôpitaux de l'aire métropolitaine fait chuter la couverture à 30 minutes de 56,9 % à 45,4 % et porte le temps médian de 19 à 32 minutes, tout en laissant la couverture à 60 minutes quasi intacte à 77,5 %. Les scénarios de réseau font l'inverse : sous l'encerclement, la part de la demande sans aucun accès routier passe de 4,8 % à 5,9 %. Pour ces femmes, il n'existe plus d'itinéraire du tout.

Dans le scénario combiné, 118 721 femmes de 15 à 49 ans passent au-delà d'une heure. Les pertes suivent le corridor de la RN1 et touchent surtout Gressier, avec 72 % de la demande communale, Croix-des-Bouquets à 60 %, Cabaret à 56 % et Ganthier à 55 %.

Le choix du seuil

Le seuil de 120 minutes est le standard de la littérature sur l'accès aux soins obstétricaux d'urgence complets. Il a été écarté comme indicateur principal, pour une raison mesurée sur ces données.

À 120 minutes, la couverture de départ atteint 94,1 %. Fermer 77 hôpitaux sur 198 la fait passer à 94,0 %. Un indicateur qui perd un dixième de point quand on ferme 39 % du parc hospitalier ne mesure plus rien. Le suivre reviendrait à conclure que la fermeture n'a pas eu lieu.

Le seuil de référence est donc de 60 minutes, également mieux fondé cliniquement : l'hémorragie du post-partum, première cause de mortalité maternelle, tue en une à deux heures sans prise en charge.

Les seuils de 30 et 120 minutes restent publiés dans tous les tableaux. Ils permettent la comparaison avec les études existantes, et le seuil court est le seul à capter l'allongement des trajets provoqué par les fermetures, que 60 minutes manque entièrement.

Limites

La qualification SONU reste un proxy

OpenStreetMap localise les structures sans rien dire de leur capacité obstétricale : ni bloc opératoire fonctionnel, ni réserve de sang, ni équipe de nuit. Les couvertures mesurent donc la distance-temps à des structures étiquetées hôpitaux, sans garantie qu'une prise en charge chirurgicale y soit possible.

Les points de contrôle sont des scénarios

La règle de placement des points de contrôle est un choix de modélisation parmi d'autres possibles. Une version d'ACLED avec coordonnées s'intégrerait à la même chaîne en remplaçant le seul fichier de points.

Le réseau est traité comme non orienté

Les sens uniques ne sont pas modélisés. À l'échelle des déplacements intercommunaux, cette simplification pèse moins que l'incertitude sur les vitesses, mais elle reste une limite, d'autant que l'attribut oneway est inégalement renseigné dans l'OpenStreetMap haïtien.

La capacité des structures n'est pas modélisée

L'analyse mesure un temps d'accès. Elle ignore la capacité d'une structure à recevoir une patiente de plus : un établissement à 20 minutes compte comme accessible quels que soient ses lits libres, son personnel ou la demande simultanée.

Reproduire

Six étapes, chacune écrivant ses sorties sur disque et lisant celles de la précédente : la chaîne reprend à mi-parcours sans refaire les étapes déjà passées.

python -m venv .venv
python -m pip install -r requirements.txt

python run_all.py                 # chaine complete
python run_all.py --depuis 3      # reprend apres la construction du graphe
python run_all.py --etape 4       # relance les seuls scenarios

Sources déjà téléchargées, la chaîne complète prend une vingtaine de minutes, dont une dizaine pour la seule étape 4 au premier lancement, puisqu'elle recalcule l'accessibilité nationale pour chaque tronçon candidat. Elle descend à trois minutes quand les couches sont en cache disque.

L'étape 5 écrit sept cartes, huit tableaux et une carte interactive. Cette dernière rassemble en un seul fichier les couches que les cartes statiques présentent séparément : déficit d'accès par commune, tronçons les plus chargés, cellules décrochées dans le scénario combiné, hôpitaux et points de contrôle simulés, chacune activable indépendamment.

Tous les paramètres de modélisation vivent dans src/config.py, jamais en dur dans les scripts, et chacun y porte sa justification. Un paramètre non justifié est un paramètre non défendable devant un bailleur.

Les quatre figures propres à cette page sortent d'un carnet séparé, notebooks/documentation_figures.ipynb, qui relit les sorties de la chaîne sans rien recalculer.