and variance parameters.

Inference about model parameters has been approached in several different ways in autoregressive cokriging models with focus on empirical Bayesian approaches. A particular challenge in autoregressive cokriging models is to estimate correlation parameters. To tackle this issue, [@Kennedy2000] assume independent noninformative priors for all the parameters, and then carry out numerical maximization for the marginal likelihood functions after integrating out regression parameters with respect to scale discrepancy parameters, variance parameters and correlation parameters [@Kennedy2000]. In [@Qian2008], conjugate priors are assumed for regression parameters, scale discrepancy parameters and variance parameters. For correlation parameters, proper gamma priors are assumed. The correlation parameters are then estimated by maximizing the corresponding marginal posteriors. For the code at the first level, standard nonlinear optimization is performed; while for the code at the second level, the corresponding posterior does not have a closed form and its evaluation requires numerical integration. This optimization procedure is then recast into a stochastic programming problem. To alleviate computational difficulties in [@Kennedy2000; @Qian2008], [@Gratiet2013] develop an efficient joint Bayesian estimation approach with either non-information priors or informative priors for all the model parameters except the correlation parameters. Without further assuming prior distributions for correlation parameters, [@Gratiet2013] maximizes a 