An Improved One-Dimension Variational Method for a Ground-Based Microwave Radiometer

1. Introduction

The Atmospheric temperature and relative humidity profiles are important parameters of the atmospheric state, especially in the field of numerical weather prediction (NWP) [1,2]. Radiosonde observation method is one of the useful methods for detecting atmospheric temperature and humidity. This method involves launching sounding balloons equipped with sophisticated sensors, enabling highly accurate measurements of various atmospheric parameters. However, the inherent constraints of launch procedures and the substantial costs involved preclude continuous, real-time monitoring of atmospheric conditions [3]. To improve the temporal coverage and resolution, microwave remote sensing techniques have been proposed [4,5]. Compared to RAOBs, microwave remote sensing obtains the atmospheric parameters continuously and data can be collected autonomously.

A microwave radiometer (MWR) is a passive remote sensing instrument with 22 observation channels. It obtains a vertical distribution of atmospheric temperature and relative humidity using the V and K band from the surface to an altitude of 10 km. In recent years, MWRs have achieved excellent performance when measuring atmospheric parameters, especially temperature and precipitable water vapor[6,7,8].

The retrieval of the temperature and relative humidity profile is an ill-posed problem. Various retrieval methods have been proposed which include physical algorithms, linear regression retrieval method, the neural network(NN) method and the one-dimensional variational (1D-VAR) method[9]. Currently, most MWR instruments use the BPNN algorithm to retrieve atmospheric temperature and water vapor profiles[10]. However, in some regions , there are not enough historical atmospheric datasets or the quality of the datasets does not satisfythe neural network training requirements. Expanding the sample set by adding data from other sensors is one approach to solve these issues[11]

The one-dimensional variational (1D-VAR) method is a optimal estimation method, that combines a precise forward atmosphere transport model with a priori information. Several studies have shown that 1D-VAR method have good performance on the retrieval of atmospheric parameters[9,12,13,14,15]. Clive D. Rodgers introduced the basic principles of BPNN and 1D-VAR retrieval methods in detail, and compared the two retrieval methods, elaborating on their advantages and limitations[16].

Recent research has investigated the integration of Backpropagation Neural Network (BPNN) algorithms with MWR systems[17]. However, the efficacy of this approach continues to be bounded by constraints related to dataset dimensions and quality. In the absence of viable Neural Network (NN) methods, an improved 1D-VAR method is introduced with an innovative iteration factor, which significantly improves the convergence efficiency of the algorithm. By incorporating this iteration factor, the method optimizes the iterative process, potentially leading to more accurate and computationally efficient retrievals of atmospheric temperature profiles from MWR observations.

The paper is organized as follows. The MWR and dataset used in this study are introduced in Section 2. the basic concept of 1D-VAR retrieval method is presented in Section 3. Section 4 introduces the experimental setup, The results of the 1D-VAR retrieval method are described in Section 5. Section 6 provides a brief summary and conclusion.

2. Instruments and Datasets

2.1. Microwave Radiometer

For this study, a 22-channel MP-3000A MWR manufactured by Radiometrics Corporation Boulder, CO, USA, was deployed in Harbin, at $45.{46}^{o}N$, $126.{41}^{o}E$[18]. The MWR is a passive instrument; it measures the downwelling radiance for two spectral ranges. The 22 to 30 GHz (K-band) range is used to retrieve water vapor density and relative humidity profiles with 8 channels, with 22.234, 22.500, 23.034, 25.000, 26.234, 28.000, and 30.000 GHz center frequencies. The 51 to 60 GHz (V-band) range is used to retrieve temperature profiles with the remaining 14 channels, with 51.248, 51.760, 52.280, 52.804, 53.336, 53.848, 54.400, 54.940, 55.500, 56.020, 56.660, 57.288, 57.964, and 58.800 GHz center frequencies. The MP-3000A MWR equipped with surface temperature,humidity and pressure sensors, as well as an infrared sensor for measuring cloud base temperature. A BPNN retrieval model has been built in MWR, and the MWR has been calibrated to ensure accuracy.

2.2. Datasets

The observed BTs of MWR used in this study were collected in April 2017 at Harbin, at $45.{46}^{o}N$, $126.{41}^{o}E$. The retrieval products from MWR are profiles of temperature(K) and water vapor density($g/m^3$). The atmospheric vertical profile from surface to 10 km is divided into 58 layers with varying vertical resolution at different altitude ranges. The resolution of these profiles is measured as initial and final distances from the surface, e.g., from 0.05 km to 0.5 km, 0.1 km to 2 km, and 0.25 km to 10 km. Five-year (February 2012 to August 2016) RAOB data from Harbin station were used to train BPNN retrieval models (temperature and water vapor density) and one-year (August 2016 to September 2017) used for validation. The simulated BTs were used as inputs for the BPNN retrieval model, which were obtained via calculations using the Monochromatic radiative transfer model (MonoRTM)[19,20] . NCEP-2 is an improved version of the NCEP and has updated data assimilation techniques. In this study, the five-year(January 2012- December 2017) NCEP-2 data with a spatial resolution of $2.5^{\circ }\times 2.5^{\circ }$ were used as the background state of 1D-VAR method. The red hollow circle in Figure 1 signifies the installed position of the MWR (microwave radiometer), while the blue dot designates the location surrounding the MWR utilized for reanalyzing the data. The reanalysis data situated at the lower right corner position was used as the background state. The NCEP-2 reanalysis data encompasses global atmospheric profile information from 1979 to the present, with four daily atmospheric profiles containing temperature and humidity data from the surface up to 17 pressure levels: 1000, 925, 850, 700, 600, 500, 400, 300, 250, 200, 150, 100, 70, 50, 30, 20 and 10 hPa. The 58 layers are interpolated and serve as the atmospheric background field $x_b$. The radiosonde data used as the truth value to validate the retrieval results.

3. Methodology

3.1. 1D-VAR Retrieval Method

In this study, a 1D-VAR method is used to retrieve temperature and humidity profiles. The method is based on Bayesian theory to obtain the maximum probability of atmospheric parameter values given the BTs, using a prior knowledge,

$$P\left(x\right|y)={\displaystyle \frac{P\left(y\right|x\left)P\right(x)}{P\left(y\right)}}$$

(1)

where $P\left(x\right)$ and $P\left(y\right)$ refer to the PDF of state vector and measurement, respectively. $P\left(y\right|x)$ represents the conditional PDF of y for a given x.

To estimate the optimal state vector x, the prior information and the errors in state and measurement spaces are assumed to be independent and follow Gaussian distributions. The form of the atmosphere transfer model can be expressed by:

$$y=F\left(x\right)+ϵ$$

(2)

$F(·)$ and $ϵ$ denote the atmosphere transfer model and system error, respectively.

The maximum a posteriori probability of the state vector x can be obtained using an iterative approach, which is tantamount to minimizing the cost function J(x),

$$J\left(x\right)={\left(x-x_b\right)}^T{B}^{-1}(x-x_b)+{\left[y-F\left(x\right)\right]}^T{R}^{-1}\left[y-F\left(x\right)\right]$$

(3)

where $x_b$ is the background state vector, B is the error covariance matrix of the background state vector, R is the observation error covariance matrix, T and −1 represent matrix transpose and inverse, respectively.

3.2. Improevd 1D-VAR Retrieval Method

As described in Section 3.1, there are two key components in the 1D-VAR retrieval algorithm; one is Jacobian matrix of the atmospheric radiative transfer model and the parameters of J(x), which include $x_b$, B, R, $F\left(x\right)$ and $k_i$ (Jacobian matrix). The other is the method that derives the optimal state vector by minimizing the cost function.

the Newton-iteration method is widely employed to find the optimal solution for objective functions. In ideal circumstances, this method ensures that each iteration yields a value of the objection function smaller than the previous one.

$$x_{i+1}=x_i-{\left[{\nabla }_{x}g\left(x_i\right)\right]}^{-1}g\left(x_i\right)$$

(4)

where $x_i$ is the ith estimate of the atmospheric state vector; $g\left(x_i\right)$ represents the derivative of the cost function $J\left(x\right)$, and the $g\left(x\right)$ can be written as,

$$g\left(x\right)={\nabla }_{x}J\left(x\right)=B^{-1}(x-x_b)-{\left[{\nabla }_{x}F\left(x\right)\right]}^T{R}^{-1}\left[y-F\left(x\right)\right]$$

(5)

${\nabla }_{x}g\left(x\right)$ is the Hessian matrix of cost function $J\left(x\right)$,

$${\nabla }_{x}g\left(x\right)=B^{-1}+{\left[K\left(x\right)\right]}^T{R}^{-1}K\left(x\right)-{\left[{\nabla }_{x}K\left(x\right)\right]}^T{R}^{-1}[y-F\left(x\right)]$$

(6)

where $K_i$ donates the Jacobian matrix. Equation (15) represents the sensitivity of the observed state y to perturbations in the state vector x at different altitudes (surface-10km) and frequencies (22.234-58.8GHz).

the Newton-Gauss iterative formula is derived:

$$x_{i+1}=x_i+{\left[B^{-1}+K_i^T{R}^{-1}{K}_i\right]}^{-1}\left[K_i^T{R}^{-1}(y-F\left(x_i\right))+B^{-1}(x_i-x_b)\right]$$

(7)

However,the Newton iteration method occasionally exhibits convergence issues when applied to optimize the objective function J(x). This convergence issue arises due to problems with the iteration step size. To address this, we introduce a step size factor $\alpha$ into the original Newton iteration method.

$$x_{i+1}=x_i-\alpha *{\left[{\nabla }_{x}g\left(x_i\right)\right]}^{-1}g\left(x_i\right)$$

(8)

Substituting equations (5) and (6) into (8), the Newton iteration method is derived as follows:

$$x_{i+1}=(1-\alpha )x_i+\alpha \left[x_b+B{K}_i^T{(K_{i}B{K}_i^T+R)}^{-1}[y-F\left(x_i\right)+K_i(x_i-x_b)]\right]$$

(9)

in this paper, $\alpha =0.5$

3.3. Optimal Configuration of 1D-VAR Method

The background state $x_b$ is part of the a priori knowledge and represents the initial atmospheric state. As shown in Section 2.2, a background state vector derived from NCEP-2 reanalysis data.

The background error covariance matrix B represents the correlation between the true atmospheric state and the background at varying altitudes. In this study, the background error covariance matrix B can be obtained via the NCEP-2 reanalysis data and RAOB data, and written as:

$$B=E\left[(x_{NCEP}-x_{RAOB}){(x_{NCEP}-x_{RAOB})}^T\right]$$

(10)

where $x_{NCEP}$ represents the NCEP reanalysis data, and $x_{RAOB}$ represents the radiosonde data.

Figure 2 shows the background error covariance of temperature and water vapor density. Both the abscissa and ordinate axis represents altitude from the surface up to 10 km. The background error covariance matrix B comprises a $58\times 58$ dimensional matrix, The elements within the matrix represent the error covariances between the NCEP reanalysis data and the radiosonde data (true atmospheric state) at various altitudes. As shown in Figure 2, the temperature error covariance matrix presents a large covariance near the diagonal, especially within the 0 to 5 km altitude range, while the water vapor density background error covariance matrix exhibits large covariance at the lower altitude levels. Above 2 km, the low water vapor content results in smaller covariances for the water vapor density.

Given these assumptions, there is no correlation between the measurements of the 22 channels in the MWR. The observation error covariance matrix R is comprised of the microwave radiometer noise N and the forward model error F, and can be written as follow:

$$R=N+F$$

(11)

where the diagonal matrix N represents the instrument noise, and the diagonal matrix F represents the error of MonoRTM.

4. Experimental Setup

To eliminate an obstruction above the antenna, the MWR was positioned on the rooftop. An antenna shroud was placed below the antenna to shield it from background radiation reflecting from the ground surface.

However, the detection by ground-based MWR sensor was susceptible to interference from solar radiation and the surface materials present at the observation location.

For the experimental validation, the impact of ground radiation background on the MWR requires consideration. Comparative analysis is conducted between the inversion outcomes obtained during daytime and nighttime.

As shown in Figure 4, significant discrepancies were observed between ground sensor and radiosonde data at 08:00 and 20:00 during April 2017. At 08:00, the surface temperatures measured by the MWR were notably higher than those recorded by radiosondes, with an average error of 8.2 K and a maximum error of 16.2 K. Conversely, at 20:00, the MWR-measured surface temperatures were lower than the radiosonde data, with an average error of 4.1 K and a maximum error of 12.2 K. The bias in ground sensor measurements significantly impacts the retrieval results of the MWR. This issue can be effectively corrected using the improved 1D-VAR method.

5. Results

In this section, the experiments utilize radiosonde data as the truth value reference for validation analysis. The output resulting from the 1D-VAR retrieval algorithm are evaluated quantitatively based on the mean error(ME) and root mean square error(RMSE).

$$RMSE=\sqrt{{\displaystyle \frac{1}{n}}\sum _{i=1}^n{(R_i-O_i)}^2}$$

(12)

$$MAE={\displaystyle \frac{1}{n}}\sum _{i=1}^n\left|(R_i-O_i)\right|$$

(13)

After excluding anomalous and cloudy instances, 52 clear sky temperature and humidity soundings remain, which are separately assessed for daytime and nighttime to account for varying effects of background radiation on the ground sensor.The resulting statistics are presented in Figure 5Figure 6.

5.1. Statistics in Daytime

As shown in Figure 5a,d, The statistical errors of the retrieval results for 52 cases during daytime are summarized in Figure 5. The Figure 5a,b show the rms errors and the mean errors for the improved 1D-VAR, 1D-VAR and BPNN method. The improved 1D-VAR method achieves a better performance than the BPNN and 1D-VAR method, the temperature rms errors of the BPNN have a great deviation from surface to 10km height with a maximum rmse greater than 8.9K. The Figure 5a shows the rms errors for the improved 1D-VAR, 1D-VAR and BPNN, the rms errors of improved 1D-VAR method are reduced 0.5-1.4K randomly within the range of 1-3km compare to 1D-VAR method. The rms errors and mean errors for water vapor density are shown in Figure 5c,d. Similar to the temperature profile retrievals, the BPNN retrieval method has the worst performance above 0.7km. the maximum rms errors of the improved 1D-VAR method and 1D-VAR method compared to BPNN method are 0.2$g/m^3$ and 0.4$g/m^3$, respectively. The retrieval errors of the improved 1D-VAR method are significantly improved compared with the 1D-VAR method and BPNN method at almost height levels.

5.2. Statistics in Nighttime

Figure 6a–d show the temperature and water vapor density errors profile retrieved via the improved 1D-VAR, 1D-VAR and BPNN approach during nighttime. As shown in Figure 6a,b, the improved 1D-VAR method achieves a best performance below 6km. In Figure 6c, the rms errors of improved 1D-VAR are smaller than the 1D-VAR method below 1km. In the range of 1-6km, the results of the improved 1D-VAR method and 1D-VAR are better than BPNN method.

Table 1 presents the comparative temperature and water vapor density retrieval outcomes from 0 to 10 km for the BPNN ,1D-VAR and improved 1D-VAR methods. Evidently, the improved 1D-VAR approach demonstrates superior performance over the BPNN and 1D-VAR method for temperature and water vapor density profile retrieval.

6. Summary and Conclusions

The optimization and validation of numerical weather prediction models need accurate atmospheric vertical profiles with high temporal resolution. This requirement underscores the need for advanced measurement techniques capable of providing detailed atmospheric parameters across various altitudes. Ground-based MWRs have emerged as promising instruments for meeting the aforementioned requirements, exhibiting remarkable capabilities in providing high-temporal resolution atmospheric vertical profiles within the range of 0-10km. However, the MWR built-in BPNN retrieval method has large retrieval errors both in temperature and humidity profiles in some cases. To meet the requirements of accuracy, we improved the 1D-VAR algorithm to retrieve the vertical profile of temperature and water vapor density by combining with NCEP2 reanalysis profiles. The improved 1D-VAR retrieval method proposed in this paper improves the iteration efficiency by adding an iteration factor $\alpha$, and performs more accurate atmospheric temperature and water vapor density profiles. Experimental validation of the proposed technique is also demonstrated.

Through statistical comparisons between the three retrieval methods(BPNN, 1D-VAR and improved 1D-VAR), we conclude that:

• The results of the BPNN retrieval method show large deviations when the ground measurement sensors are faulty, and the improved 1D-VAR method can effectively correct these errors;
• Compared to the BPNN and 1D-VAR methods, the RMSE of temperature of improved 1D-VAR method are reduced by $60.8\%$ and $29.4\%$ during daytime, and by $54.2\%$ and $49.7\%$ during nighttime,respectively.
• Compared to the BPNN and 1D-VAR methods, the accuracy of the improved 1D-VAR method regarding water vapor density is not significant. Although the rms error during surface to 10km is reduced, the rms errors during surface to 1km are higher than BPNN method. In the next step, the improved 1D-VAR method and the BPNN method can be fused to achieve a high retrieval results of water vapor density.

Moving forward, we will improve the parameter tuning methods of the 1D-VAR inversion model to address the problem of adaptive tuning of stepping factors.