Historiquement la géostatistique s'est développée dans des espaces Euclidiens qui permettent de définir facilement les notions de stationnarité et d'isotropie. Cependant les phénomènes naturels montrent les limites de ce formalisme : la projection sur une carte à deux dimensions de phénomènes définis sur le globe terrestre élimine leurs périodicités ; un paysage résultant d'une histoire géologique complexe présente des directions de continuités variables dans l'espace. Ces géométries complexes peuvent être décrites grâce à des variétés Riemanniennes, qui localement ressemblent à un espace Euclidien mais avec une métrique variant à l'intérieur du domaine considéré. L'approche SPDE (Lindgren et al 2011, Carrizo et al 2018) a été étendue à des fonctions aléatoires généralisées sur ces espaces (Pereira et Desassis 2018 et 2019) : un opérateur pseudo-différentiel, spécifié par une fonction positive de carré sommable, est appliqué à un bruit blanc défini sur la variété. Les techniques numériques de résolution par éléments finis des EDP sont utilisées pour résoudre les problèmes opérationnels d'estimation par krigeage et de simulations conditionnelles. Elles permettent de traiter de gros jeux de données sur des grilles de grande taille ainsi que de tenir compte d'anisotropies variables. Ces techniques sont illustrées sur des cas synthétiques pour différentes géométries. L'inférence des paramètres de ces modèles est l'objet de travaux en cours. Ces modèles définis sur ces géométries non Euclidienne permettent d'intégrer des connaissances a priori sur la géométrie des phénomènes étudiés.