## 1 Introduction  
Soil hydraulic conductivity is essential information for soil water flow and solute transport simulations to represent hydrological processes such as water storage, release, and transfer in the soil, and is therefore a prerequisite and key link to ensure model effectiveness. Saturated hydraulic conductivity is often spatially heterogeneous and scale-dependent [1, 2]. Uncertainty in ground permeability can result in considerable alterations to the response of the porous medium [3, 4]. And the simulation process is often time-consuming and labor-intensive to obtain directly in field or laboratory experiments subject to scale constraints [5, 6]. In recent years, numerical inversion methods based on field observations have been increasingly used to estimate the hydraulic conductivity of soils [7, 8, 9, 10, 11, 12, 13, 14, 15, 16]. We note that most of the existing studies focus on the inversion of hydraulic parameters in homogeneous or stratified soils using observations at the same scale. However, strong spatial variability often requires us to obtain the spatial distribution of soil hydraulic parameters at a high resolution. Particularly, as the observations are obtained from multiple scales with different means of observation, it is necessary to integrate the dynamic observations at different scales to obtain the hydraulic parameters at the required scales. Thus, the problem boils down to integrating observations at different scales to invert the spatial distribution of soil hydraulic parameters at high resolution. This problem faces two difficulties: first, how to integrate observations on different scales for the inversion of hydraulic parameters on the scale of interest. Here, it is necessary to establish a kinetic relationship between observations on multiple scales and the hydraulic parameters on the scale of interest. The other is how to solve the problem of efficiency and uncertainty of the inversion simulation due to the huge parameter dimensions and the consequent large number of forward simulation calculations. This study uses Bayesian optimization \(\text{DREAM}_{\text{(ZS)}}\) as an inverse framework, and incorporates an efficient upscaling method and Karhunen-Loeve (KL) expansion to infer finer-scale saturated hydraulic conductivity distribution conditioned on coarse-scale measurements. First, the \(\text{DREAM}_{\text{(ZS)}}\) algorithm is used as a inverse framework to construct the posterior distribution of saturated hydraulic conductivity field in high-dimension soils conditioned on coarse-scale observations. The \(\text{DREAM}_{\text{(ZS)}}\) is an improved Markov chain Monte Carlo (MCMC) method proposed by [PERSON] et al. [17; 18; 19; 20], which greatly improves the convergence speed of the inverse model compared to the traditional MCMC method. Compared to the basic DREAM (DiffeRential Evolution Adaptive Metropolis) algorithm, it refers to the previous sampling state to adjust the current process in generating each chain candidate, which aims to reduce the number of parallel chains and accelerate the convergence of high-dimensional problems, and has been successfully applied to research in many fields [21; 22; 23; 24; 25]. Then, the upscaling method proposed by [PERSON] et al. [26] is used to link observations and hydraulic parameters, which may come from different scales. This method is based on the Kirchhoff transformation and the two-convergence theory to propose the macroscopic description of the Richards equation and the effective constitutive relations, rather than equivalent or upscaled flow parameters, as is the case with many parameter averaging approaches such as the stream tube approach [e.g., [27; 28; 29]] and the scaleway approach [e.g., [30; 31]].