Adaptive spatial blocking for scalable clustering inference with applications to high-throughput spatial proteomics

Résumé technique : Blocage spatial adaptatif pour l'inférence de regroupement (clustering) scalable

Énoncé du problème
L'analyse des motifs ponctuels spatiaux, particulièrement dans la protéomique spatiale à haut débit, repose largement sur la fonction $K$ de Ripley pour évaluer le regroupement (clustering). Cependant, les méthodes existantes basées sur $K$, y compris le cadre $KAMP$ récemment proposé (Wrobel et Song, 2024), font face à de graves goulots d'étranglement computationnels lorsqu'elles sont appliquées à des jeux de données à grande échelle contenant des centaines de milliers à des millions de points. Les tests de permutation traditionnels nécessitent des calculs de distance en $O(n^2)$, et même l'approche analytique par moments de $KAMP$, bien que plus rapide, exige toujours le stockage des distances par paires pour tous les points, ce qui entraîne une utilisation prohibitive de la mémoire et une complexité computationnelle quadratique ($O(n^2)$) à mesure que les tailles d'échantillon augmentent. Les stratégies de sous-échantillonnage comme l'amincissement (thinning) indépendant offrent une solution partielle, mais entraînent une perte significative d'informations spatiales et de puissance statistique, ne résolvant pas fondamentalement la complexité quadratique.

Méthodologie
Pour répondre à ces défis de scalabilité, les auteurs proposent B-KAMP (Block-based KAMP), un cadre d'inférence de regroupement spatial basé sur un algorithme de blocage spatial adaptatif.

1. Algorithme de blocage spatial adaptatif :
   Le cœur de la méthode est un algorithme qui partitionne une grande fenêtre d'observation en une collection de blocs rectangulaires disjoints. L'algorithme opère en trois étapes :

   • Représentation par grille : Le domaine est discrétisé en une grille $R \times C$ où chaque cellule contient au plus $\sqrt{n}$ points.
   • Identification des blocs (Phases I & II) : En utilisant une approche basée sur l'histogramme inspirée du problème du « plus grand rectangle dans un histogramme », l'algorithme identifie de manière itérative des sous-régions rectangulaires satisfaisant des contraintes spécifiques :
     • Taille minimale : Les blocs doivent contenir au moins $\sqrt{m}$ points d'intérêt (POI) et un nombre suffisant de points de fond pour garantir une distribution de permutation non dégénérée.
     • Rapport d'aspect : Les blocs sont contraints par un rapport d'aspect maximal $\rho$ pour assurer la stabilité statistique (évitant les formes très allongées qui biaisent les estimateurs corrigés des effets de bord).
     • L'algorithme donne la priorité à la sélection des blocs ayant les plus petits nombres totaux de points afin de maximiser le nombre de blocs disjoints et de maintenir l'efficacité computationnelle.
   • Complément de bloc (Phase III) : Les cellules de grille résiduelles non assignées aux blocs sont fusionnées de manière itérative dans des blocs adjacents en étendant les limites des blocs vers le rectangle englobant minimal, à condition que le bloc résultant reste disjoint des autres et respecte les contraintes de rapport d'aspect. Cela minimise la perte de couverture spatiale.

2. Statistique de test basée sur les blocs :
   À l'intérieur de chaque bloc identifié $i$, la méthode calcule une statistique de la fonction $K$ standardisée $Z_i(r)$ en utilisant les moments analytiques dérivés du cadre $KAMP$ sous une distribution nulle de permutation de labels conditionnelle.

   • La statistique de test globale $Z(r)$ est définie comme une somme pondérée des statistiques par bloc : $Z(r) = \sum w_i Z_i(r)$.
   • Les poids $w_i$ sont choisis proportionnellement à $n_i/p_i$ (où $n_i$ est le nombre total de points et $p_i$ est la prévalence locale des POI dans le bloc) afin de donner plus de poids aux blocs présentant des signaux rares.
   • Sous la distribution nulle de permutation par bloc conditionnelle, les auteurs prouvent que lorsque le nombre de blocs $b \to \infty$, la statistique globale $Z(r)$ converge en distribution vers une loi normale standard $\mathcal{N}(0,1)$ via le théorème central limite de Lyapunov. Cela permet un calcul rapide des p-values sans réplicats de permutation empiriques.

3. Extension à l'analyse bivariée :
   Le cadre est étendu pour tester la colocalisation spatiale entre deux types de points (par exemple, les cellules plasmatiques et les macrophages) en utilisant un estimateur de fonction $K$ bivariée, avec une standardisation et une agrégation par blocs analogues.

Principales contributions

• Algorithme scalable : Le développement d'un algorithme de blocage spatial adaptatif en $O(RC)$ qui identifie efficacement des blocs rectangulaires valides pour des motifs ponctuels à grande échelle, évitant la complexité $O(n^2)$ des calculs par paires sur la fenêtre complète.
• Cadre statistique : La formulation d'un test de regroupement par blocs qui agrège les preuves locales via une somme pondérée, en exploitant l'asymptotique normale pour une inférence rapide.
• Implémentation pratique : La mise à disposition d'une stratégie de blocage polyvalente applicable à diverses méthodes de regroupement scalables, avec le code disponible via GitHub.

Résultats

• Études de simulation :

  • Efficacité computationnelle : B-KAMP passe à l'échelle linéairement avec le nombre de points et reste réalisable pour des jeux de données allant jusqu'à 100 000 points. En revanche, le $KAMP$ en fenêtre complète échoue en raison des contraintes de mémoire pour $n \ge 40 000$, et le $KAMP$ aminci ($KAMP_{0.25}$) reste coûteux en calcul. B-KAMP est systématiquement la méthode la plus rapide pour de grands $n$.
  • Puissance statistique : Bien que les méthodes en fenêtre complète atteignent généralement la puissance la plus élevée lorsque cela est possible, B-KAMP démontre une puissance compétitive, particulièrement aux rayons $r \approx 0.15$ (relatifs à la taille du bloc). Elle surpasse les méthodes fortement amincies ($KAMP_{0.25}$) et maintient un contrôle valide de l'erreur de Type I.
  • Robustesse : La méthode montre des performances robustes à travers diverses abondances de points et résolutions de grille.

• Application aux données réelles (Protéomique spatiale) :

  • Appliqué à un jeu de données de tissus intestinaux humains sains (environ 2,6 millions de cellules sur 65 images).
  • La méthode a détecté avec succès une forte agrégation spatiale des cellules plasmatiques et une colocalisation significative entre les cellules plasmatiques et les macrophages.
  • B-KAMP a fourni des preuves statistiques plus fortes (p-values plus faibles) que les bases de référence $KAMP$ aminci pour les grandes images et a parfaitement dimensionné là où les méthodes en fenêtre complète dépassaient la mémoire disponible.

Signification et affirmations
L'article affirme que B-KAMP offre un équilibre favorable entre la puissance statistique et l'efficacité computationnelle, permettant une inférence de regroupement spatial rigoureuse pour les données de protéomique spatiale à haut débit, ce qui était auparavant computationnellement prohibitif. Les auteurs positionnent la méthode non pas comme un remplacement des procédures en fenêtre complète pour les petits jeux de données, mais comme une alternative scalable nécessaire pour les images à grande échelle où les méthodes traditionnelles échouent.

Les auteurs maintiennent une position modeste concernant les limites :

• La méthode est un test par bloc conditionnel, ce qui signifie qu'elle n'utilise pas les relations spatiales à travers les frontières des blocs, ce qui peut atténuer les motifs de regroupement s'étendant sur plusieurs blocs et entraîner une certaine perte de puissance par rapport aux méthodes en fenêtre complète.
• La performance dépend de la résolution initiale de la grille et de la nature gloutonne (greedy) de l'étape de complément de bloc (Phase III).
• Des travaux futurs sont suggérés, incluant la sélection de grille pilotée par les données, l'optimisation globale des formes de blocs et des extensions aux données spatiales 3D.