Spatial Proportional Hazards Model with Differential Regularization

1. Problématique et Contexte

L'article s'intéresse à l'analyse de survie dans un contexte spatial. Le modèle de risques proportionnels (PH) de Cox est la référence standard, supposant une relation log-linéaire entre les covariables et la fonction de risque. Cependant, cette hypothèse est souvent insuffisante pour capturer les effets spatiaux complexes lorsque les données sont associées à des localisations géographiques précises dans un domaine bidimensionnel.

Les défis principaux identifiés sont :

• Limites des approches existantes : Les méthodes non paramétriques classiques (comme les lissages par noyaux ou les splines unidimensionnelles) peinent à s'adapter aux domaines spatiaux bidimensionnels, en particulier ceux aux frontières irrégulières, aux formes non convexes ou aux trous internes. Elles souffrent souvent de biais aux frontières, d'instabilité ou de coûts computationnels prohibitifs.
• Géométrie du domaine : La plupart des modèles spatiaux existants supposent que la dépendance spatiale est uniquement fonction de la distance euclidienne, négligeant la géométrie réelle du domaine (ex: côtes, îles, zones urbaines complexes).
• Traitement du risque de base : De nombreuses approches bayésiennes ou à vraisemblance composite imposent une spécification paramétrique du risque de base, perdant ainsi l'avantage de l'approche de Cox qui traite ce risque comme un paramètre de nuisance non spécifié.

L'objectif est donc de développer un modèle PH non paramétrique capable de modéliser un effet spatial lisse sur des domaines irréguliers tout en restant dans le cadre de la vraisemblance partielle de Cox.

2. Méthodologie Proposée

Les auteurs proposent un modèle où la fonction de risque conditionnelle $\lambda(t | x, p)$ est définie par :
$$\lambda(t | x, p) = \lambda_0(t) \exp(x^\top \beta_0 + h_0(p))$$
où $h_0(p)$ est une fonction inconnue et lisse représentant l'effet spatial en la position $p$.

A. Estimation par Vraisemblance Partielle Pénalisée

L'estimation repose sur la maximisation d'une fonctionnelle de vraisemblance partielle pénalisée :
$$Q_n(\beta, h) = \frac{1}{n} \sum_{i=1}^n \delta_i \left( x_i^\top \beta + h(p_i) - \log \sum_{j=1}^n I(Y_j \ge Y_i) e^{x_j^\top \beta + h(p_j)} \right) - \frac{\lambda}{2} \int_\Omega (\Delta h(p))^2 dp$$

• Le terme de pénalité utilise l'opérateur Laplacien $\Delta$ (courbure locale) pour imposer la régularité de la fonction spatiale $h$.
• La contrainte de centrage $\int_\Omega h(p) dp = 0$ assure l'identifiabilité du modèle.

B. Approximation par Éléments Finis (FEM) et Méthode des Tamis (Sieve)

Pour rendre le problème d'optimisation infiniment dimensionnel traitable numériquement, les auteurs utilisent la méthode des éléments finis (FEM) :

• Maillage : Le domaine $\Omega$ est discrétisé via une triangulation.
• Espace fonctionnel : Contrairement aux éléments linéaires standards ($C^0$), l'approche utilise des éléments finis $C^1$ (par exemple, l'élément d'Argyris de degré 5). Cela est nécessaire car le terme de pénalité implique le carré du Laplacien, qui nécessite que la fonction soit deux fois dérivable ($h \in H^2(\Omega)$).
• Méthode des Tamis (Sieve Method) : L'espace de recherche est restreint à une suite de sous-espaces de dimension finie $H_n$ qui deviennent denses dans l'espace infini lorsque la taille du maillage $\eta(n)$ tend vers zéro. Cela permet d'établir des propriétés asymptotiques rigoureuses.

C. Optimisation Numérique

Le problème résultant est une maximisation concave d'une fonction différentiable. Les auteurs proposent d'utiliser l'algorithme BFGS (quasi-Newton) avec des gradients analytiques, ce qui assure une convergence rapide et efficace, même pour des dimensions élevées.

3. Contributions Clés et Résultats Théoriques

• Identifiabilité : Le théorème 1 prouve que le modèle est unique sous des hypothèses standard (censure indépendante, support non dégénéré, etc.).
• Propriétés Asymptotiques :
  • Consistance : Le théorème 2 établit que l'estimateur de tamis converge en probabilité vers le vrai paramètre à un taux contrôlé dépendant de la taille de l'échantillon $n$, de la dimension du tamis $K(n)$ et du paramètre de lissage $\lambda$.
  • Normalité Asymptotique : Le théorème 3 démontre que la composante paramétrique $\hat{\beta}_n$ suit une loi normale asymptotique, permettant des inférences statistiques standard (intervalles de confiance, tests) sur les coefficients de régression, malgré la présence de l'effet spatial non paramétrique.
• Gestion des Domaines Irréguliers : L'utilisation de la FEM permet de traiter nativement des géométries complexes (frontières courbes, formes non convexes) sans biais de frontière significatif, contrairement aux méthodes de lissage par noyaux.

4. Résultats Empiriques

A. Simulations

Des simulations sur un domaine en forme de fer à cheval montrent que :

• Le modèle proposé surpasse nettement le modèle PH standard (qui ignore l'espace) et le modèle PH avec splines (GAM) en termes de biais et d'erreur quadratique moyenne (MSE) pour la récupération de l'effet spatial.
• Les estimateurs des coefficients $\beta$ sont non biaisés et les probabilités de couverture des intervalles de confiance sont proches des valeurs théoriques (95%).
• La performance s'améliore avec la taille de l'échantillon et diminue avec le taux de censure.

B. Applications Réelles

1. Temps de réponse des ambulances à San Francisco :
   • Le modèle capture des disparités géographiques dans les temps d'arrivée, même après ajustement pour les covariables (météo, heure, type d'appel).
   • Il révèle des zones à risque élevé (réponse rapide) et faible, corrélées à la configuration du réseau routier et à la topographie, y compris des gradients spécifiques sur des îles comme Treasure Island.
2. Données sismiques crowdsourcées (Campi Flegrei, Italie) :
   • Analyse des temps de déclenchement des smartphones lors d'un séisme.
   • Le modèle identifie un gradient Est-Ouest dans l'amplification des ondes sismiques, cohérent avec la géologie locale, ce qui n'aurait pas été détecté par un modèle purement basé sur la distance à l'épicentre.

5. Signification et Perspectives

Cet article représente une avancée significative en statistique spatiale et en analyse de survie :

• Flexibilité Géométrique : Il offre une solution robuste pour modéliser des effets spatiaux sur des domaines réels et complexes, là où les méthodes traditionnelles échouent.
• Cadre Théorique Solide : L'intégration de la méthode des tamis avec la pénalisation différentielle permet de prouver des propriétés asymptotiques rigoureuses, comblant un vide théorique pour les modèles PH spatiaux non paramétriques.
• Efficacité Computationnelle : L'approche par éléments finis $C^1$ combinée à l'optimisation BFGS rend le modèle applicable à des jeux de données de grande taille.

Les auteurs suggèrent des pistes futures pour l'extension de la méthode aux données agrégées (aréales), l'amélioration des diagnostics de validité de l'hypothèse de risques proportionnels en présence d'effets spatiaux, et l'incorporation d'effets spatio-temporels dynamiques.

En résumé, cette méthode fournit un outil puissant pour extraire des structures spatiales latentes complexes dans des données de survie, offrant une meilleure interprétabilité et une précision supérieure par rapport aux approches existantes.