codes have been widely used to solve mathematical models that represent real-world processes in virtually every field of science and engineering. They are often referred to as *simulators* in Uncertainty Quantification (UQ) and computer experiments [@OHagan2006; @Santner2018]. In practice, computer codes can be too timing-consuming to be used for adequately addressing UQ tasks. To overcome this bottleneck, Gaussian processes have been widely used as surrogate models to approximate simulators due to its computational advantages and attractive theoretical properties [@Sacks1989].

In real applications, computer codes can be run at different levels of accuracy due to sophistication of physics incorporated in mathematical models, accuracy of numerical solvers and resolutions of meshes; see [@Peherstorfer2018] for formal definition of multifidelity models. Several works have been proposed to combine output from computer codes at different fidelity levels based on a well-known geostatistical method called *cokriging*; see Chapter 3 of [@Cressie1993]. The approach of cokriging to synthesizing multiple computer model outputs is originated in [@Kennedy2000], which is developed based upon a first order Markov assumption that given output from a low-fidelity code run at an input, no more information can be learnt about the high-fidelity code with output from the low-fidelity code at any other input. 