In this work, we reconstruct the continuous function of LD from the mock GW or EM observational data with the Gaussian process. The advantages of the Gaussian process method are nonparametric and cosmological model-independent smoothing technique used to reconstruct a continuous function from the observed data. The reconstructed function, at each redshift point \(z\) , is related to a Gaussian distribution with a mean value and Gaussian error bands. The results from observational data points at any two redshifts \(z_i\)  and \(z_j\)  are correlated through a covariance function \(\kappa (z_i,z_j)\)  due to their nearness to each other. There is a wide range of possible candidates for this function, and it depends on a set of hyperparameters. The different choices of the covariance function may affect the reconstruction to some
extent. The detailed description and analysis of the Gaussian process can be found in Refs. [1]}, [2]}. This method has been widely used to reconstruct the equation of state of dark energy [1]}, [4]}, and test cosmography [5]} and the CDDR [6]}, [7]}. The simplest covariance function usually takes the squared exponential as
\(\kappa (z_i,z_j)=\sigma ^2_f \exp \bigg [-{(z_i-z_j)^2\over {2 l^2}}\bigg ],\) 
