Linear Model Of Coregionalization Research Articles

Multivariate simulation using a locally varying coregionalization model

Understanding the response of materials in downstream processes of mining operations relies heavily on proper multivariate spatial modeling of relevant properties of such materials. Ore recovery and the behavior of tailings and waste are examples where capturing the mineralogical composition is a key component: in the first case, to ensure reliable revenues, and in the second one, to avoid environmental risks involved in their disposal. However, multivariate spatial modeling can be difficult when variables exhibit intricate relationships, such as non-linear correlation, heteroscedastic behavior, or spatial trends. This work demonstrates that the complex multivariate behavior among variables can be reproduced by disaggregating the global non-linear behavior through the spatial domain and looking instead at the local correlations between Gaussianized variables. Local linear dependencies are first inferred from a local neighborhood and then interpolated through the domain using Riemannian geometry tools that allow us to handle correlation matrices and their spatial interpolation. By employing a non-stationary modification of the linear model of coregionalization, it is possible to independently simulate variables and then combine them as a linear mixture that locally varies according to the inferred correlation, reproducing the global multivariate behavior seen on input variables. A real case study is presented, showing the reproduction of the reference multivariate distributions, as well as direct and cross semi-variograms.

Multivariate Geostatistics for Mapping of Transmissivity and Uncertainty in Karst Aquifers

Due to their complex morphology, karst terrains are particularly more fragile and vulnerable to environmental damage compared to most natural systems. Their hydraulic properties, such as their transmissivity (T) and spatial variability, can be relevant for understanding groundwater flow and, consequently, for the sustainable management of water resources. The application of geostatistical methods allows for spatial interpolation and mapping based on observations combined with uncertainty quantification. Direct measurements of T are typically scarce, while those of the specific capacity (Sc) are more frequent. We established a linear and spatial relationship between the logarithms of T and Sc measured in 174 wells in a semi-arid karst region in northeastern Brazil. These relationships were used to construct a cross-variogram, whose Linear Model of Coregionalization proved valid. The values and the cross-variogram of logT and logSc were used to generate interpolations over 2554 values of logSc, which did not spatially coincide with logT. We used ordinary co-kriging (CO-OK) and conditional sequential Gaussian co-simulation (CO-SGS) to generate the interpolations. The cross-variogram of logT and logSc, when considering 174 wells, was isotropic with an exponential structure, a nugget effect of approximately 20% of the sill, and a range of 5 km. Cross-validation indicated an optimal number of 10 neighboring wells used in CO-OK, and we used 500 stochastic realizations in CO-SGS, which were then used to generate maps of logT estimates, deviations derived from the interpolations, and probabilistic scenarios. The resulting transmissivity maps are relevant for the design of groundwater management strategies, including stochastic approaches where the transmissivity realizations can be used to parameterize multiple executions of numerical flow models.

Fast matrix completion in epigenetic methylation studies with informative covariates.

DNA methylation is an important epigenetic mark that modulates gene expression through the inhibition of transcriptional proteins binding to DNA. As in many other omics experiments, the issue of missing values is an important one, and appropriate imputation techniques are important in avoiding an unnecessary sample size reduction as well as to optimally leverage the information collected. We consider the case where relatively few samples are processed via an expensive high-density whole genome bisulfite sequencing (WGBS) strategy and a larger number of samples is processed using more affordable low-density, array-based technologies. In such cases, one can impute the low-coverage (array-based) methylation data using the high-density information provided by the WGBS samples. In this paper, we propose an efficient Linear Model of Coregionalisation with informative Covariates (LMCC) to predict missing values based on observed values and covariates. Our model assumes that at each site, the methylation vector of all samples is linked to the set of fixed factors (covariates) and a set of latent factors. Furthermore, we exploit the functional nature of the data and the spatial correlation across sites by assuming some Gaussian processes on the fixed and latent coefficient vectors, respectively. Our simulations show that the use of covariates can significantly improve the accuracy of imputed values, especially in cases where missing data contain some relevant information about the explanatory variable. We also showed that our proposed model is particularly efficient when the number of columns is much greater than the number of rows-which is usually the case in methylation data analysis. Finally, we apply and compare our proposed method with alternative approaches on two real methylation datasets, showing how covariates such as cell type, tissue type or age can enhance the accuracy of imputed values.

Bayesian estimation of large-scale simulation models with Gaussian process regression surrogates

Large scale, computationally expensive simulation models pose a particular challenge when it comes to estimating their parameters from empirical data. Most simulation models do not possess closed-form expressions for their likelihood function, requiring the use of simulation-based inference, such as simulated method of moments, indirect inference, likelihood-free inference or approximate Bayesian computation. However, given the high computational requirements of large-scale models, it is often difficult to run these estimation methods, as they require more simulated runs that can feasibly be carried out. The aim is to address the problem by providing a full Bayesian estimation framework where the true but intractable likelihood function of the simulation model is replaced by one generated by a surrogate model trained on the limited simulated data. This is provided by a Linear Model of Coregionalization, where each latent variable is a sparse variational Gaussian process, chosen for its desirable convergence and consistency properties. The effectiveness of the approach is tested using both a simulated Bayesian computing analysis on a known data generating process, and an empirical application in which the free parameters of a computationally demanding agent-based model are estimated on US macroeconomic data.

Co-simulation of hydrofacies and piezometric data in the West Thessaly basin, Greece: A geostatistical application using the GeoSim R package

In the present study, we co-simulate hydrofacies and piezometric data in order to construct geostatistical realizations of underground geology in an area of the West Thessaly basin. This basin is of great importance in terms of sustainable water management and environmental perspective in Greece. Through Plurigaussian modeling, the hydrofacies are first transformed into Gaussian Random Fields. Then, a Linear Coregionalization Model is established to account for the dependencies between hydrofacies and the Normal scores of piezometric data. The effect of co-simulation shows an improvement of the facies transition probabilities in comparison with those of Plurigaussian simulation. For the purpose of this study, we use the GeoSim package in R developed by our team for the implementation of Plurigaussian simulation and co-simulation.

Learning Multitask Gaussian Process Over Heterogeneous Input Domains

Multitask Gaussian process (MTGP) is a well-known nonparametric Bayesian model for learning correlated tasks effectively by transferring knowledge across tasks. But current MTGPs are usually limited to the multitask scenario defined in the same input domain, leaving no space for tackling the heterogeneous case, i.e., the features of input domains vary over tasks. To this end, this article presents a novel heterogeneous stochastic variational linear model of coregionalization (HSVLMC) model for simultaneously learning the tasks with varied input domains. Particularly, we develop the stochastic variational framework with Bayesian calibration that: 1) infers posterior domain mappings to consider the effect of dimensionality reduction raised by domain mappings for achieving effective input alignment and 2) employs a residual modeling strategy to leverage the inductive bias brought by prior domain mappings for better-model inference. Finally, the superiority of the proposed model against existing heterogeneous LMC models has been extensively verified on diverse heterogeneous multitask cases and a practical multifidelity steam turbine exhaust case.

Spatial relationships between the van Genuchten soil-water retention curve parameters and photosynthetic gas exchange variables in a wheat field

Photosynthetic gas exchange indices can be used to assess spatial variation of crop growth and yield on the field. Spatial variation of soil parameters can be an important tool to assess and manage crop growth status and product quality in the field, given the impact of soil on photosynthetic function. The purpose of this study was to evaluate the spatial relationships between the van Genuchten (vG) soil water retention curve parameters, soil properties and the photosynthetic gas exchange indices in a wheat field with an extension of 200 ha. Different variables of soil and plant were determined at the 100 nodes of a nearly regular grid (150 m × 150 m) over the whole field. The range of variograms showed high variability for photosynthetic parameters and vG parameters, and greater continuity for soil parameters. Photosynthesis was significantly correlated with saturated soil water content (θs) (r = 0.26*), θ1500kPa (permanent wilting point) (0.21*), porosity (r = 0.31**), bulk density (r = −0.32**) and clay content (r = 0.24*). Cross correlogram and linear model of coregionalization of sub-stomatal CO2 concentration (Ci) and stomatal conductance (gs) with photosynthetic rate (A) and transpiration (E) showed opposite relation with soil and hydraulic properties. Regressing of maps showed that a large part of the total variation in A (66%), Ci (79%) and gs (37%) was explained by soil variables. Also, from the visual inspection of the interpolated maps the field can be divided in zones characterized by different soil properties and these maps could be effectively used by farmers to improve wheat yield. On the other hand, results of this study may enhance understanding of the soil properties and soil hydraulic parameters effects on gas exchange parameters of the plant and delineation of the field in homogeneous zones that could be advantageously used by farmers to make informed decisions in site-specific management.

Multivariate Modeling for Spatio-Temporal Radon Flux Predictions.

Nowadays, various fields in environmental sciences require the availability of appropriate techniques to exploit the information given by multivariate spatial or spatio-temporal observations. In particular, radon flux data which are of high interest to monitor greenhouse gas emissions and to assess human exposure to indoor radon are determined by the deposit of uranium and radio (precursor elements). Furthermore, they are also affected by various atmospheric variables, such as humidity, temperature, precipitation and evapotranspiration. To this aim, a significant role can be recognized to the tools of multivariate geostatistics which supports the modeling and prediction of variables under study. In this paper, the spatio-temporal distribution of radon flux densities over the Veneto Region (Italy) and its estimation at unsampled points in space and time are discussed. In particular, the spatio-temporal linear coregionalization model is identified on the basis of the joint diagonalization of the empirical covariance matrices evaluated at different spatio-temporal lags and is used to produce predicted radon flux maps for different months. Probability maps, that the radon flux density in the upcoming months is greater than three historical statistics, are then built. This might be of interest especially in summer months when the risk of radon exhalation is higher. Moreover, a comparison with respect to alternative models in the univariate and multivariate context is provided.

GeoSim: An R-package for plurigaussian simulation and Co-simulation between categorical and continuous variables

Plurigaussian simulation is widely used to model geological facies in geosciences and is predominantly applied in mineral deposits and petroleum reservoirs exploration. GeoSim package builds geostatistical models of categorical regionalized variables via conditional or unconditional Plurigaussian simulation and co-simulation. Co-simulation between Gaussian Random Fields representing the geological facies and other numerical variables accounting for auxiliary hydrological or geophysical quantities, is also available in this package with the definition of a linear coregionalization model. The algorithm is not restricted by the number of simulated facies and the number of truncated Gaussians, while parts of the code requiring heavy computations are compiled in C++ taking benefits of the integration between R and C++. In this work, we introduce the GeoSim package and demonstrate its capabilities. We present a 3D application focused on a lignite mine in Greece, where we investigate the Plurigaussian simulation and co-simulation of five geological facies (categorical variables) and the lower calorific value (continuous variable). The findings of our study highlight the significant benefits of Plurigaussian and co-simulation to capture the geological spatial heterogeneity.

A multivariate approach for mapping a soil quality index and its uncertainty in southern France

AbstractPedometricians have spent a lot of effort on mapping soil types and basic soil properties. However, end‐users typically need a more elaborate soil quality index for land management. Soil quality indices are typically derived from multiple individual soil properties by evaluating whether specific criteria are met. If this is based on individually mapped soil properties, then an important consequence is that cross‐correlations between soil properties are ignored. This makes it impossible to quantify the uncertainties associated with the mapped indices. The objective of this study was to map a soil potential multifunctionality index for agriculture (Agri‐SPMI) over a 12,125 km2 study region located along the French Mediterranean coast to help urban planners preserve soils of the highest quality. The index considered the ability of soils to fulfil four functions under five land use scenarios. A binary map represented each soil function fulfilment for a given scenario. The final soil quality index map was the sum of the 20 binary maps. A regression cokriging model was developed to map the basic soil properties first individually from legacy soil data and spatial soil covariates using a random forest algorithm, and next, interpolate the residuals using cokriging and the linear model of coregionalisation. The mapping uncertainties of soil properties were propagated by calculating the soil quality index over 300 stochastic simulations of soil properties derived from the linear models of coregionalisation. Results showed a poor prediction accuracy of the quality index, mainly because some soil properties were poorly predicted (notably available water capacity and coarse fragments) and used in combination with extreme thresholds that defined land suitability. Overall, the uncertainty was correctly quantified because the stochastic simulations reproduced the width of the observed distribution well, but the shapes of the distributions differed considerably from those of the observations. We envisage some ways for improvement, such as creating probability maps instead of the mean from simulations, and changing the prediction support from point to area.Highlights The study mapped a soil quality index (Agri‐SPMI) to help preserve soils of highest quality along the French Mediterranean coast. The Agri‐SPMI considered the ability of soils to fulfill four functions under five land use scenarios, derived from multiple individual soil properties. A regression cokriging model was developed to map the basic soil properties, followed by interpolating the residuals using cokriging and the linear model of coregionalisation. The study accurately quantified mapping uncertainties using stochastic simulations, but the soil quality index prediction accuracy was poor, with suggestions for improvement.

Silver origins and spatial distribution in soils of southwestern France: A first assessment with cokriging

Ag is a rare, non-essential metal with many applications in different forms, including nanoparticles (Ag-NPs). Its release into the soil environment can represent a risk, due to the toxicity of silver. This study is based on the analysis of samples from the two sampling campaigns of the Soil Quality Measurement Network, RMQS1 and RMQS2, to establish the spatial distribution and origin of Ag in the soils of southwestern France (area of 90,293 km2). Given the progress of the second campaign (completion in 2027), we chose to do a cokriging to estimate the spatial distribution, with Pb, Cd or Cu as potential secondary variables. The results of the analyzes of the 35 samples of the RMQS2 indicate a range of Ag concentration in the soils going from 0.015 to 0.222 mg kg−1 and a median of 0.081 mg kg−1. Correlations were found between the concentrations of Ag Cd, Cu and Pb, and also between Ag concentrations and distance to Pb mines and deposits (r = 0.62). This suggests a common origin for Ag and Pb linked to historical exploitations of argentiferous galena. No link could be demonstrated between the concentration of Ag and the high density of anti-hail ground generators in the region. The analysis of the (auto-)variograms showed that Pb is the most relevant variable compared to Cd and Cu to achieve a cokriging. Thus, based on Ag concentrations from RMQS2 and Pb concentrations from RMQS1 (356 samples), a linear co-regionalization model was established and the co-kriging of Ag concentrations in soils calculated. Compared to ordinary kriging, cokriging improved the quality of the interpolation (decrease in the root mean square error).

Use of Geostatistics for Multi-Scale Spatial Modeling of Xylella fastidiosa subsp. pauca (Xfp) Infection with Unmanned Aerial Vehicle Image

In recent years, the use of Unmanned Aerial Vehicles (UAVs) has been spreading widely, as in plant pest control. The collection of huge amounts of spatial data raises various issues including that of scale. Data from UAVs generally explore multiple scales, so the problem arises in determining which one(s) may be relevant for a given application. The objective of this work was to investigate the potential of UAV images in the fight against the Xylella pest for olive trees. The data were a multiband UAV image collected on one date in an olive grove affected by Xylella. A multivariate geostatistics approach was applied, consisting firstly of estimating the linear coregionalization model to detect the scales from the data; and secondly, of using multiple factor kriging to extract the sets of scale-dependent regionalized factors. One factor was retained for each of the two selected scales. The short-range factor could be used in controlling the bacterium infection while the longer-range factor could be used in partitioning the field into three management zones. The work has shown the UAV data potential in Xylella control, but many problems still need to be solved for the automatic detection of infected plants in the early stages.

Blind recovery of sources for multivariate space-time random fields

With advances in modern worlds technology, huge datasets that show dependencies in space as well as in time occur frequently in practice. As an example, several monitoring stations at different geographical locations track hourly concentration measurements of a number of air pollutants for several years. Such a dataset contains thousands of multivariate observations, thus, proper statistical analysis needs to account for dependencies in space and time between and among the different monitored variables. To simplify the consequent multivariate spatio-temporal statistical analysis it might be of interest to detect linear transformations of the original observations that result in straightforward interpretative, spatio-temporally uncorrelated processes that are also highly likely to have a real physical meaning. Blind source separation (BSS) represents a statistical methodology which has the aim to recover so-called latent processes, that exactly meet the former requirements. BSS was already successfully used in sole temporal and sole spatial applications with great success, but, it was not yet introduced for the spatio-temporal case. In this contribution, a reasonable and innovative generalization of BSS for multivariate space-time random fields (stBSS), under second-order stationarity, is proposed, together with two space-time extensions of the well-known algorithms for multiple unknown signals extraction (stAMUSE) and the second-order blind identification (stSOBI) which solve the formulated problem. Furthermore, symmetry and separability properties of the model are elaborated and connections to the space-time linear model of coregionalization and to the classical principal component analysis are drawn. Finally, the usefulness of the new methods is shown in a thorough simulation study and on a real environmental application.

A Fully Natural Gradient Scheme for Improving Inference of the Heterogeneous Multioutput Gaussian Process Model.

A recent novel extension of multioutput Gaussian processes (GPs) handles heterogeneous outputs, assuming that each output has its own likelihood function. It uses a vector-valued GP prior to jointly model all likelihoods' parameters as latent functions drawn from a GP with a linear model of coregionalization (LMC) covariance. By means of an inducing points' framework, the model is able to obtain tractable variational bounds amenable to stochastic variational inference (SVI). Nonetheless, the strong conditioning between the variational parameters and the hyperparameters burdens the adaptive gradient optimization methods used in the original approach. To overcome this issue, we borrow ideas from variational optimization introducing an exploratory distribution over the hyperparameters, allowing inference together with the posterior's variational parameters through a fully natural gradient (NG) optimization scheme. Furthermore, in this work, we introduce an extension of the heterogeneous multioutput model, where its latent functions are drawn from convolution processes. We show that our optimization scheme can achieve better local optima solutions with higher test performance rates than adaptive gradient methods for both the LMC and the convolution process model. We also show how to make the convolutional model scalable by means of SVI and how to optimize it through a fully NG scheme. We compare the performance of the different methods over the toy and real databases.

Correlated Chained Gaussian Processes for Modelling Citizens Mobility Using a Zero-Inflated Poisson Likelihood

Modelling the mobility of people in a city depends on counting data with inherent problems of overdispersion. Such dispersion issues are caused by massive amounts of data with zero values. Though traditional machine learning models have been used to overcome said problems, they lack the ability to appropriately model the spatio-temporal correlations in data. To improve the modelling of such spatio-temporal correlations, in this work we propose to model the citizens mobility, for the Chinese city of Guangzhou, by means of a Zero-inflated Poisson likelihood in conjunction with Gaussian process priors generated from convolution processes. We follow the idea of chaining the likelihood’s parameters to latent functions drawn from Gaussian process priors; this way allowing a higher flexibility to model heteroscedasticity. Additionally, we derive a stochastic variational inference framework that allow us to use two types of convolution process models in the context of large datasets: correlated chained Gaussian processes with a convolution processes model, and correlated chained Gaussian processes with variational inducing kernels. We present quantitative and qualitative results comparing the performance between Negative Binomial and Zero-inflated Poisson likelihoods, both in combination with three types of Gaussian process priors: a linear model of coregionalisation, and our two proposed methods based on a convolution processes model and variational inducing kernels.

Scalable Marginalization of Correlated Latent Variables with Applications to Learning Particle Interaction Kernels

Marginalization of latent variables or nuisance parameters is a fundamental aspect of Bayesian inference and uncertainty quantification. In this work, we focus on scalable marginalization of latent variables in modeling correlated data, such as spatio-temporal or functional observations. We first introduce Gaussian processes (GPs) for modeling correlated data and highlight the computational challenge, where the computational complexity increases cubically fast along with the number of observations. We then review the connection between the state space model and GPs with Matérn covariance for temporal inputs. The Kalman filter and Rauch-Tung-Striebel smoother were introduced as a scalable marginalization technique for computing the likelihood and making predictions of GPs without approximation. We introduce recent efforts on extending the scalable marginalization idea to the linear model of coregionalization for multivariate correlated output and spatio-temporal observations. In the final part of this work, we introduce a novel marginalization technique to estimate interaction kernels and forecast particle trajectories. The computational progress lies in the sparse representation of the inverse covariance matrix of the latent variables, then applying conjugate gradient for improving predictive accuracy with large data sets. The computational advances achieved in this work outline a wide range of applications in molecular dynamic simulation, cellular migration, and agent-based models.

Bayesian Optimization With Transfer Learning: A Study on Spatial Variability of Rock Properties Using NMR Relaxometry

AbstractNuclear magnetic resonance measurements of sedimentary rocks are used to extract various transport properties including hydraulic conductivity and water retention curves. These estimates are controlled by intrinsic physical quantities like surface relaxivities, and effective relaxation time and restricted self‐diffusion coefficient of water in clay. Sampling these properties on a set of core plugs presents a series of inverse problems where some of the extracted parameters are expected to be similar. To leverage such valuable information, we extend a previously developed single‐task inverse solution workflow (ISW) to the multi‐task case, transferring the knowledge gained from previous optimization tasks. Two multi‐task kernels: intrinsic model of coregionalization (ICM) and linear model of coregionalization are compared to capture the underlying correlations. We consider three micro‐CT images of Bentheimer sandstone from two different blocks imaged at different resolutions, following different segmentation pathways, and demonstrate our approach for the case of low and high task similarity. In both scenarios the multi‐task ISW finds lower fitting residuals and uses only one‐third to one‐half of the function evaluations required by the single‐task ISW. The scalability of the multi‐task ISW is demonstrated by transferring knowledge of two completed optimization tasks to a third task, which outperforms the single‐task ISW, with ICM showing faster convergence. The observed 4% difference for the values identified for samples from the same block and around 28% difference across blocks indicates significant spatial variability in surface relaxivity of the main mineral component, while effective clay parameters show a significantly higher variability.

Potential of GPR data fusion with hyperspectral data for precision agriculture of the future

Precision Agriculture (PA) requires accurate spatial and temporal information of soil properties at a very fine scale. Traditional soil characterization methods are time consuming, laborious and invasive and do not allow long-term repeatability of measurements. Ground Penetrating Radar (GPR) appears to be a particularly suitable methodology for characterizing soil and subsurface from a physical property point of view. Visible-near infrared-shortwave infrared (VIS-NIR-SWIR) reflectance spectroscopy has now become a widespread technique in soil analysis. Information on soil variability can be improved by the integration of data from multiple sensors. The overall objective of this paper was to examine the potential of fusing GPR data with hyperspectral data using multivariate geostatistics for delineating the management zones in the soil of an olive grove of centuries-old trees in Italy.A linear model of coregionalization (LMC) was individually fitted for the raw hyperspectral data and for GPR data including for each case a nugget effect and two spherical models at short scale and at longer scale. After that, one data set was obtained from the fusion of the two sensor data sets and a LMC was fitted for the combined data to be then used in factor cokriging. The application of this technique produced a delineation of the field into homogeneous zones, highlighting a wide southern-central zone, characterized by different granulometric and chemical properties. The proposed approach was then effective to discriminate areas with different properties by using multi-sensor data. It then has the potential to be used in PA.

A data science approach for analysis and reconstruction of spinodal-like composition fields in irradiated FeCrAl alloys

A statistical method for the analysis of continuously distributed data representative of composition fluctuations in irradiated FeCrAl alloys acquired using Energy Dispersive X-ray Spectroscopy (EDS) method is presented. Using probability distribution functions, direct and cross-covariances between the elemental compositions, the effects of alloy composition and irradiation dose were investigated on the spatial distribution and length scale of composition fluctuations at the nanoscale. It was observed that, for neutron-irradiated FeCrAl alloys, the distribution of Fe and Cr followed a left-skewed and right-skewed distribution, respectively for all (average) alloy compositions and irradiation doses. The analysis also revealed enhanced spatial gradients in the elemental compositions at higher irradiation dose. Direct and cross-covariance estimates of the experimental data were also utilized for reconstruction of composition data through fitting it to a parametric form of the covariance functions. Linear Model of Coregionalization was used to determine the parameters of the covariance functions. Subsequently, a spectral method was utilized for simulating a realization of the alloy compositions. Close correspondence was observed between the experimental and the reconstructed data which was analyzed using probability distribution functions and covariance functions. Composition space of the experimental and reconstructed data and dislocation velocities as a function of applied stress and line directions over the entire composition maps were also examined.

Scalable multi-task Gaussian processes with neural embedding of coregionalization

Multi-task regression attempts to exploit the task similarity in order to achieve knowledge transfer across related tasks for improving the quality of prediction and alleviating the demand of big data. The application of Gaussian process (GP) in this scenario yields the non-parametric yet informative Bayesian multi-task regression paradigm. Multi-task GP (MTGP) provides not only the prediction mean but also the associated prediction variance to quantify uncertainty, thus gaining popularity in various scenarios. The linear model of coregionalization (LMC) is a well-known MTGP paradigm which exploits the dependency of tasks through linear combination of several independent and diverse GPs. The LMC however suffers from high model complexity and limited model capability when handling complicated multi-task cases. To this end, we develop the neural embedding of coregionalization that transforms the latent GPs into a high-dimensional latent space to induce rich yet diverse behaviors. Furthermore, we use advanced variational inference as well as sparse approximation to devise a tight and compact evidence lower bound (ELBO) for higher quality of scalable model inference. Extensive numerical experiments have been conducted to verify the higher prediction quality and better generalization of our model, named NSVLMC, on various real-world multi-task datasets and the cross-fluid modeling of unsteady fluidized bed.