Half-tapering strategy for conditional simulation with large datasets

Gaussian conditional realizations are routinely used for risk assessment and planning in a variety of Earth sciences applications. Conditional realizations can be obtained by first creating unconditional realizations that are then post-conditioned by kriging. Many efficient algorithms are available for the first step, so the bottleneck resides in the second step. Instead of doing the conditional simulations with the desired covariance (F approach) or with a tapered covariance (T approach), we propose to use the taper covariance only in the conditioning step (Half-Taper or HT approach). This enables to speed up the computations and to reduce memory requirements for the conditioning step but also to keep the right short scale variations in the realizations. A criterion based on mean square error of the simulation is derived to help anticipate the similarity of HT to F. Moreover, an index is used to predict the sparsity of the kriging matrix for the conditioning step. Some guides for the choice of the taper function are discussed. The distributions of a series of 1D, 2D and 3D scalar response functions are compared for F, T and HT approaches. The distributions obtained indicate a much better similarity to F with HT than with T.Comment: 39 pages, 2 Tables and 11 Figure

New Developments in Covariance Modeling and Coregionalization for the Study and Simulation of Natural Phenomena

RÉSUMÉ La géostatistique s’intéresse à la modélisation des phénomènes naturels par des champs aléatoires univariables ou multivariables. La plupart des applications utilisent un modèle stationnaire pour représenter le phénomène étudié. Il est maintenant reconnu que ce modèle n’est pas assez flexible pour représenter adéquatement un phénomène naturel montrant des comportements qui varient considérablement dans l’espace (un exemple simple de cette hétérogénéité est le problème de l’estimation de l’épaisseur du mort-terrain en présence d’affleurements). Pour le cas univariable, quelques modèles non-stationnaires ont été développés récemment. Toutefois, ces modèles n’ont pas un support compact, ce qui limite leur domaine d’application. Il y a un réel besoin d’enrichir la classe des modèles non-stationnaires univariable, le premier objectif poursuivi par cette thèse.----------ABSTRACT Geostatistics focus on modeling natural phenomena by univariate or multivariate spatial random fields. Most applications rely on the choice of a stationary model to represent the studied phenomenon. It is now acknowledged that this model is not flexible enough to adequately represent a natural phenomenon showing behaviors that vary substantially in space (a simple example of such heterogeneity is the problem of estimating overburden thickness in the presence of outcrops). For the univariate case, a few non-stationary models were developed recently. However, these models do not have compact support, which limits in practice their range of application. There is a definite need to enlarge the class of univariate non-stationary models, a first goal pursued by this thesis

Multivariate Mat\'ern Models -- A Spectral Approach

The classical Mat\'ern model has been a staple in spatial statistics. Novel data-rich applications in environmental and physical sciences, however, call for new, flexible vector-valued spatial and space-time models. Therefore, the extension of the classical Mat\'ern model has been a problem of active theoretical and methodological interest. In this paper, we offer a new perspective to extending the Mat\'ern covariance model to the vector-valued setting. We adopt a spectral, stochastic integral approach, which allows us to address challenging issues on the validity of the covariance structure and at the same time to obtain new, flexible, and interpretable models. In particular, our multivariate extensions of the Mat\'ern model allow for time-irreversible or, more generally, asymmetric covariance structures. Moreover, the spectral approach provides an essentially complete flexibility in modeling the local structure of the process. We establish closed-form representations of the cross-covariances when available, compare them with existing models, simulate Gaussian instances of these new processes, and demonstrate estimation of the model's parameters through maximum likelihood. An application of the new class of multivariate Mat\'ern models to environmental data indicate their success in capturing inherent covariance-asymmetry phenomena

Spatial interpolation

The theory and practical application of techniques of statistical interpolation are studied in this thesis, and new developments in multivariate spatial interpolation and the design of sampling plans are discussed. Several applications to studies in soil science are presented.Sampling strategies for collecting spatial data (both the number of observations and their location in the region to be studied) are discussed. It is shown that nested sampling is unsuitable if data are to be collected in an area in order to determine the spatial semivariogram, because semivariogram values for only a few distances are obtained. Furthermore, grid sampling is preferable to nested sampling if spatial interpolation is intended. Sequential sampling is advantegeous if the mean of a variable within an area is to be estimated and there is no spatial correlation. Sequential sampling requires only about 30% of the number of observations required by standard sampling schemes.In this thesis universal kriging and cokriging (that is kriging and cokriging in the presence of a trend) are formulated in terms of regression procedures. Universal kriging is a special case of universal cokriging. Multivariate increments are extensions of univariate increments and of multivariate stationary variables. Conditions are formulated which permissible polynomial pseudo-covariance and pseudo- crosscovariance functions describing the spatial structure of the variables (or their increments) and their interaction, respectively, have to obey. The coefficients of these functions have been estimated by using the restricted maximum likelihood (REML) method. A practical application of universal cokriging is described.The application of spatial interpolation in soil science is examined. One of the studies described investigates the problems of scale and the use of soil survey information on moisture deficits caused by groundwater extraction in the Mander area in the Netherlands. In another study the available water and the infiltration rates on terraces of the Allier river in the Limagne area in France are investigated. Cokriging becomes more precise as compared to kriging (and there is a concomitant reduction in costs) if the predictand is strongly correlated with the covariable. This is particularly true if the sampling of the covariable is denser than that of the predictand and the costs of sampling of the covariable are much less than those of sampling the predictand. Stratification of the survey area, e.g. by means of soil map delineations increases the precision of predictions when applying cokriging. An obvious gain in precision is achieved for homogeneous soil units, where the measured values are relatively small, and there is no spatial structure. Also, the use of cokriging permits fewer observations as compared to kriging, if a certain predescribed precision of predictions is defined. When simulation calculation models are used, e.g. to obtain values for moisture deficits, one should first calculate model results for every observation point, and then interpolate, rather than interpolate the input variables, and then calculate model results.</TT

Caractérisation stochastique des hétérogénéités hydrostratigraphiques appliquée à la modélisation hydrogéologique régionale

Résumé L’une des principales sources d’incertitude des modèles hydrogéologiques porte sur les hétérogénéités hydrostratigraphiques. Ces dernières jouent un rôle majeur sur l’écoulement de l’eau souterraine et le transport des contaminants. Pourtant la caractérisation hydrogéologique inclut très rarement l’évaluation de l’effet de l’incertitude des hétérogénéités. Cette situation s’explique en partie par la difficulté à reproduire, de manière satisfaisante, ces caractéristiques physiques en modèles équivalents. Peu ou pas de méthodologie complète et efficace permet d’effectuer la modélisation hydrogéologique régionale dans un cadre stochastique tout en quantifiant l’incertitude des hétérogénéités inter et intra unités. La difficulté de simuler efficacement les hétérogénéités inter unités de ce type de contexte est un frein à la modélisation hydrogéologique stochastique régionale. Également, la difficulté associée à la paramétrisation stochastique régionale intra unités (conductivité hydraulique ou K ) de ces modèles représente aussi un obstacle. Ces difficultés constituent une problématique majeure quand vient le besoin de quantifier adéquatement l’incertitude des modèles. Les méthodes de simulation géostatistique moderne, telles que la méthode plurigaussienne et multipoint, permettent de produire des modèles équivalents des hétérogénéités inter unités, appelés réalisations. Malgré leurs récents développements, ces méthodes conservent différentes lacunes telles que le manque de réalisme géologique dans certains contextes. Par exemple, la méthode plurigaussienne ne peut facilement incorporer des transitions asymétriques entre unités, de même que la méthode multipoint qui éprouve aussi des difficultés avec la non- stationnarité ainsi que les forts contrastes d’épaisseurs. Par contre, la méthode MCP peut s’affranchir de ces difficultés notamment grâce au transiogramme et sa propriété de forçage 0/1 des probabilités. Les hétérogénéités intra unités des modèles hydrogéologiques sont définies à partir des don- nées de K scalaire à une échelle locale. Cette échelle n’est généralement pas représentative de la taille des éléments d’un modèle hydrogéologique régional. La caractérisation stochastique de la K équivalente (tenseur-K ) à l’échelle de l’élément est un défi en soi. Elle nécessite plu- sieurs étapes incluant des contraintes à respecter dont la préservation de la structure spatiale et la contrainte d’inégalité Kverticale ≤ Khorizontale des tenseurs ainsi que les corrélations inter composantes. À cela s’ajoute, la mise à l’échelle de la K qui tient compte de l’effet des hétérogénéités locale régionalisée, la régionalisation des composantes corrélées des tenseurs-K des différentes unités ainsi que leur calage avec les variables d’état (p. ex. charge hydraulique) dans un contexte multivariable contraint. Une méthodologie complète devrait également fournir des réalisations significativement différentes afin de quantifier adéquatement l’incertitude associée aux tenseurs-K. L’objectif général de la thèse consiste à développer une démarche complète de caractérisation stochastique de l’incertitude des hétérogénéités inter et intra unités hydrostratigraphiques appliquées à la modélisation hydrogéologique régionale. Plus spécifiquement, la thèse vise à développer une méthodologie de simulation géostatistique des unités hydrostratigraphiques (UHS) afin de fournir des modèles équivalents d’aspect réaliste pour un contexte de sédimentation directionnelle (asymétrie) et à développer une méthodologie de caractérisation stochastique, de mise à l’échelle et de calage des tenseurs-K équivalents à partir de la K scalaire tout en préservant la structure spatiale des tenseurs, la contrainte d’inégalité ainsi que les corrélations inter composantes. La méthodologie globale s’oriente autour de trois étapes principales : i) la compilation et la préparation des données, ii) la simulation de modèles hydrostratigraphiques directionnels et iii) la caractérisation stochastique des tenseurs-K. La première étape a porté sur la compilation et le traitement des données (hydrostratigraphie et K) pour les développements méthodologiques et les différents tests méthodologiques appliqués sur le sous-bassin Innisfil Creek (comté de Simcoe, Ontario). La deuxième étape est basée sur le développement d’une démarche de simulation par la méthode MCP afin de générer des modèles équivalents dans un contexte hydrostratigraphique 3D complexe dont la séquence des unités est ordonnée verticalement. La méthodologie proposée inclut l’estimation des probabilités bivariables par transformée de Fourier à partir du modèle déterministe et des critères de qualité de simulation. La méthode est testée pour montrer sa capacité à quantifier la variabilité locale de la solution et son incertitude. La troisième étape fournit une méthodologie complète de caractérisation stochastique des hétérogénéités intra UHS qui inclut : la mise à l’échelle d’hétérogénéités locales (K) en tenseur-K équivalent, l’utilisation du modèle linéaire de corégionalisation pour la simulation des composantes des tenseurs-K équiprobables, la préservation des corrélations entre les composantes des tenseurs-K par nouvelle transformation gaussienne bivariée et l’utilisation de la méthode de déformation graduelle (GDM) pour caler les tenseurs-K en fonction des observations de charge hydraulique. Les résultats montrent la capacité de la méthode MCP à simuler efficacement un contexte hydrostratigraphique régional 3D complexe avec de multiples unités. De plus, les résultats montrent la capacité de la méthode à reproduire l’ordination verticale des unités et tout ça malgré la forte hypothèse d’indépendance conditionnelle. Les résultats des simulations conditionnelles et non conditionnelles pour le cas synthétique montrent clairement l’absence de tout biais substantiel. Dans l’exemple plus complexe du contexte réel, le biais est moins important pour le cas non conditionnel et pour l’échantillonnage aléatoire des données que pour le biais observé dans les données dû à l’échantillonnage préférentiel. La surreprésentation de données conditionnantes n’a pas d’impact significatif sur les proportions moyennes simulées mais l’ajout de données supplémentaires réduit la variabilité. Soulignons que la méthode MCP fournit de bons résultats malgré une image d’entraînement fortement non stationnaire. Avec une autre méthode (PGS ou MPS), rien n’aurait garanti un tel succès dans le contexte étudié. Les résultats de la méthodologie de caractérisation de la K montrent que la nouvelle méthode proposée fournit des modèles équivalents calés des tenseurs-K. La méthode de mise à l’échelle de la K scalaire en tenseur-K équivalent tient compte de la structure des hétérogénéités locales et des propriétés statistiques de chacune des UHS. La régionalisation des composantes des tenseurs-K effectuée par simulation géostatistique à partir de modèles linéaires de corégionalisation (LCM) rend le temps de calcul raisonnable. La simulation reproduit aussi, pour chacune des unités, les caractéristiques statistiques des différentes composantes de K. Pour rétablir les corrélations non linéaires des tenseurs-K non-gaussiens ainsi que les histogrammes des composantes souhaités, une nouvelle transformation bivariée est appliquée sur les réalisations. Les résultats du calage des tenseurs-K montrent que la méthode des déformations graduelles permet le calage multivariable contraint tout en préservant la covariance spatiale et croisée de UHS. Il a été démontré que la méthode conserve les relations non linéaires entre les composantes du tenseur-K grâce à l’intégration de la transformation bivariée. L’approche proposée a ainsi permis de préserver les propriétés spatiales des tenseurs ainsi que les relations entre les composantes. Les résultats démontrent que la GDM s’applique à un contexte multivariable contraint et pas seulement au contexte univariable. Soulignons que les inversions convergent rapidement en termes du nombre d’itérations grâce à l’optimisation du paramètre de déformation de la GDM. Il a aussi été montré que la méthode permet de caler la recharge séparément ou simultanément avec les tenseurs-K. Les exemples de calage montrent qu’il est avantageux d’inclure la recharge stochastique dans l’inversion. Par ailleurs, il a été observé qu’un nombre élevé d’observations de charge n’améliore pas significativement les résultats. Une comparaison, portant sur la définition de zones de captage autour de deux puits de pompage, entre la méthode GDM proposée et la méthode d’inversion PEST montrent des tenseurs-K géologiquement plus réalistes avec GDM que ceux obtenus avec PEST. L’ensemble des résultats permet de conclure que la méthode MCP simule efficacement un contexte régional complexe d’une succession stratigraphique directionnelle d’un bassin sédimentaire glaciaire. La méthode prend efficacement en charge la forte non-stationnarité spatiale des unités ainsi que d’imposantes contraintes comme la déposition verticale ordonnée en tenant compte de l’asymétrie stratigraphique à partir des probabilités de transition entre les unités. La méthode MCP s’est révélée sans biais dans le cas non conditionnel. La méthodologie proposée de caractérisation des tenseurs-K permet le calage individuel ou simultané des tenseurs-K de différentes UHS et/ou de recharge. La méthodologie inclut une mise à l’échelle des conductivités hydrauliques quasi ponctuelles vers un tenseur équivalent en tenant compte de la structure spatiale à l’échelle locale. Bien que la mise à l’échelle soit intensive en temps de calcul, il a été démontré qu’il est possible de limiter au minimum le nombre d’éléments à mettre à l’échelle grâce à la régionalisation des tenseurs-K par simulation géostatistique. La thèse a apporté quelques contributions originales. D’abord, il s’agit d’une première appli- cation de la méthode MCP pour la simulation 3D. À ma connaissance, aucune autre méthode de simulation stochastique ne semble être en mesure de solutionner avec autant de satisfaction un contexte hydrostratigraphique complexe 3D comme celui du comté de Simcoe. La méthode des déformations graduelles permet le calage efficace des composantes des tenseurs- K avec préservation des corrélations non linéaires des composantes et de leurs fonctions de covariance. À ma connaissance, il s’agit de la première approche utilisant la GDM pour un contexte multivariable avec contraintes. Enfin, cette thèse contribue à l’avancement des connaissances et des techniques de la caractérisation stochastique des hétérogénéités inter et intra hydrostratigraphiques. Dans une plus grande perspective, la thèse contribue à fournir des outils adaptés pour mieux quantifier l’incertitude associée aux hétérogénéités présente dans les modèles hydrogéologiques régionaux. ---------- Abstract One of the main sources of uncertainty in hydrogeological models relates to hydrostrati- graphic heterogeneities. They play a major role in quantification of groundwater dynamics and mass transport. Yet hydrogeological characterization rarely includes the assessment of the effect of uncertainty and spatial variability of hydrostratigraphic units and hydrogeologic parameters. This is partly explained by the difficulty of reproducing these physical charac- teristics satisfactorily in equivalent models. There is at present little or no comprehensive and efficient methodology that allows regional hydrogeological modelling to be carried out in a stochastic framework while quantifying the uncertainty of inter- and intra-unit hetero- geneities. The difficulty to effectively simulate inter-unit heterogeneities under such condi- tions is a hindrance to stochastic regional hydrogeological modelling. As well, the difficulty associated with the regional stochastic intra-unit parameterization (hydraulic conductivity or K ) of these models further aggravates the problem. These difficulties constitute a major issue when it comes to the need to adequately quantify model uncertainty. Modern geostatistical simulation methods, such as plurigaussian and multipoint simulations, make it possible to produce equivalent models of inter-unit heterogeneities, referred to as realizations. Despite their recent developments, these methods retain various shortcomings such as the lack of geological rigour in certain settings. For example, the plurigaussian method cannot easily incorporate asymmetrical transitions between units, whereas the multipoint method has dif- ficulties with non-stationarity and strong thickness contrasts. On the other hand, the MCP (Markov-type categorical prediction) method has the capacity to overcome these difficulties thanks to the transiogram and its 0/1 probability forcing property. The intra-unit heterogeneities of the hydrogeological models are defined from scalar K data at a local scale. This scale is generally not representative of the size of the elements used in regional hydrogeological models. Stochastic characterization of the equivalent K (K -tensor) at the element scale is a challenge in itself. It requires several steps including constraints to be respected such as the preservation of the spatial structure, the natural constraint Kvertical ≤ Khorizontal, and the inter-component correlations. In addition, the upscaling of K takes into account the effect of local heterogeneities, the regionalization of the correlated components of the K -tensors of the different units and their calibration with state variables (e.g., hydraulic head) in a constrained multivariate context. A comprehensive methodology should also provide significantly different realizations in order to adequately quantify the uncertainty associated with the K -tensors. The general objective of this thesis is to develop a comprehensive approach to stochastic characterization of the uncertainty of inter- and intra-hydrostratigraphic unit heterogeneities applied to regional hydrogeological modelling. More specifically, the thesis aims to develop a methodology for geostatistical simulation of hydrostratigraphic units (HSU) in order to provide equivalent models of realistic appearance for a directional sedimentation context (asymmetry) and to develop a methodology for stochastic characterization, upscaling and calibration of equivalent K -tensors from the scalar K while preserving the spatial structure of the tensors, the inequality constraint as well as the inter-component correlations. The overall approach revolves around three main steps: (i) data compilation and preparation, (ii) simulation of directional hydrostratigraphic models, and (iii) stochastic characterization of K -tensors. The first step focuses on data compilation and processing (hydrostratigraphy and K ) for methodological developments and various methodological tests applied to the Innisfil Creek sub-watershed (Simcoe County, Ontario). The second stage is based on the development of a simulation approach using the MCP method to generate equivalent models in a complex 3D hydrostratigraphic context with a vertically ordered sequence of units. The proposed methodology includes the estimation of bivariate probabilities by Fourier transform from the deterministic model and the simulation quality criteria. The method is tested to show its ability to quantify the local variability of the solution and its uncertainty. The third step provides a complete methodology for stochastic characterization of intra-HSU heterogeneities which includes: upscaling of local heterogeneities (K ) into equivalent K - tensors, use of the linear co-regionalization model for simulation of equiprobable K -tensor components, preservation of correlations between K -tensor components by new bivariate Gaussian transformation and use of the Gradual Deformation Method (GDM) to calibrate the K -tensors according to hydraulic head observations. The results show the ability of the MCP method to effectively simulate a complex 3D regional hydrostratigraphic with multiple units. Furthermore, the results show the ability of the method to reproduce the vertical ordering of the units, despite the strong assumption of conditional independence. The results of the conditional and unconditional simulations for the synthetic case clearly show the absence of any substantial bias. In the more complex example of field hydrogeological settings, the bias is less important in the unconditional case and in the random sampling case than the bias observed in the data due to preferential sampling. Over-representation of conditioned data does not have a significant impact on the mean simulated proportions, but the addition of supplementary data reduces variability. The MCP method provides good results despite a highly non-stationary training image. With another method (PGS or MPS), nothing would have guaranteed such success in the studied context. The results of the K characterization methodology show that the presented original method provides equivalent calibrated models of the K -tensors. The scalar K upscaling method to equivalent K -tensor takes into account the structure of local heterogeneities and the sta- tistical properties of each of the HSU. The regionalization of the K -tensor components by geostatistical simulation using Linear Co-Regionalization Models (LCMs) keeps the compu- tation time reasonable. The simulation also reproduces, for each of the units, the statistical characteristics of the different K components. To restore the non-linear correlations of the non-Gaussian K -tensors as well as the histograms of the desired components, a new bivariate transformation is applied to the realizations. The results of the calibration of the K -tensors show that the method of gradual deformations allows the multivariate constrained calibration while preserving the spatial and cross-covariance of HSU. The method has been shown to preserve the non-linear relationships between the K -tensor components through the integra- tion of the bivariate transformation. The proposed approach has thus preserved the spatial properties of the tensors and the relationships between the components. The results show that the GDM applies to a constrained multivariate context and not only to the univariate context. It should be noted that the inversions converge rapidly in terms of the number of iterations thanks to the optimization of the deformation parameter of the GDM. It has also been shown that the method allows to set the recharge separately or simultaneously with the K -tensors. The calibration examples show that it is advantageous to include stochastic recharge in the inversion. It has also been observed that a high number of head observations does not significantly improve the results. A comparison of the definition of capture zones around two pumping wells between the proposed GDM method and the PEST inversion method shows geologically more realistic K -tensors with GDM than those obtained with PEST. The overall results suggest that the LCM method effectively simulates a complex regional context of a directional stratigraphic succession in a glacial sedimentary basin. The method effectively supports the strong spatial non-stationarity of the units as well as imposing con- straints such as ordered vertical deposition by taking into account stratigraphic asymmetry from the transition probabilities between the units. The MCP method proved to be unbiased in the unconditional case. The proposed K -tensor characterization methodology allows in- dividual or simultaneous calibration of K -tensors from different HSU and/or recharge units. The methodology includes an upscaling of quasi-point hydraulic conductivities to an equiva- lent tensor taking into account the spatial structure at the local scale. Although the upscaling is computationally intensive, it has been shown that it is possible to minimize the number of elements to be scaled by regionalizing the K -tensors by geostatistical simulation. The thesis provided substantial original contributions. This is a first application of the MCP method for 3D simulations. To the author’s knowledge, no other stochastic simulation method seems to be able to solve as satisfactorily a 3D hydrostratigraphic context as complex as that of Simcoe County. The gradual deformation method allows for efficient calibration of the K - tensor components while preserving the non-linear correlations of the components and their covariance functions. To the author’s knowledge, this is also the first time where GDM has been used in a multivariate context with constraints. Finally, this thesis contributes to the advancement of knowledge and techniques for stochastic characterization of inter- and intra- hydrostratigraphic heterogeneities. In a broader perspective, the thesis contributes original methods and tools adapted to better quantify the uncertainty associated with heterogeneities present in regional hydrogeological models

Principled methods for mixtures processing

This document is my thesis for getting the habilitation à diriger des recherches, which is the french diploma that is required to fully supervise Ph.D. students. It summarizes the research I did in the last 15 years and also provides the short­term research directions and applications I want to investigate. Regarding my past research, I first describe the work I did on probabilistic audio modeling, including the separation of Gaussian and α­stable stochastic processes. Then, I mention my work on deep learning applied to audio, which rapidly turned into a large effort for community service. Finally, I present my contributions in machine learning, with some works on hardware compressed sensing and probabilistic generative models.My research programme involves a theoretical part that revolves around probabilistic machine learning, and an applied part that concerns the processing of time series arising in both audio and life sciences