Evaluation of reanalysis precipitable water vapor under typhoon conditions using multi-source observations

Precipitable Water Vapor (PWV) exhibits rapid and complex variations during typhoons, and its evaluation under typhoon conditions remains challenging due to sparse observations over the oceans. This study systematically evaluates PWV estimates from three state-of-the-art reanalyses during 113 typhoons between 2020 and 2024 over the Northwest Pacific and East Asia. The fifth-generation European Centre for Medium-Range Weather Forecasts Reanalysis (ERA5), the Modern-Era Retrospective Analysis for Research and Applications Version 2 (MERRA-2), and the Japanese Reanalysis for Three Quarters of a Century (JRA-3Q) are compared with ground-based Global Navigation Satellite System (GNSS), radiosonde, and Radio Occultation (RO) observations. Results based on ground-based GNSS, which serve as the primary reference in this study, show that ERA5 achieves the best overall accuracy, with a bias changing from −1.48 mm in non-typhoon periods to −0.18 mm during typhoons and root mean square error (RMSE) decreasing from 3.47 to 3.30 mm during typhoons. JRA-3Q also exhibits smaller errors during typhoons than during non-typhoon periods, although its bias and RMSE remain relatively large. Conversely, MERRA-2 shows higher error during typhoons, with the bias shifting from −0.63 mm in non-typhoon periods to 0.76 mm during typhoons, indicating a transition from modest underestimation to modest overestimation. Comparisons with radiosonde and RO observations show consistently high correlations for all three reanalyses, supporting their general reliability under typhoon conditions. Overall, the results provide an accuracy reference during typhoons and support the applicability of reanalysis PWV for typhoon-related studies, with ERA5 appearing the most robust among the datasets evaluated.

1 Introduction

Water vapor is primarily distributed within the troposphere and is one of the most important greenhouse gases. It plays a critical role in energy exchange within weather systems, the hydrological cycle, and climate change (e.g., Schneider et al., 2010; Sherwood et al., 2010). Precipitable Water Vapor (PWV), defined as the total amount of water vapor contained in a vertical column of the atmosphere per unit area, is a key variable for characterizing atmospheric water vapor. It is an important input parameter for Numerical Weather Prediction (NWP) models and is widely used in meteorological and climatological monitoring and forecasting (e.g., Zhao et al., 2020; Zhang et al., 2022). The spatio-temporal variation and distribution of PWV (also referred to as column water vapor) are related to the vertical structure of atmospheric moisture, tropical deep convection, and precipitation (Muller et al., 2009; Holloway and Neelin, 2009).

Tropical Cyclones (TCs) are among the most destructive types of extreme weather. They occur frequently, exhibit high intensity, and cause widespread impacts (Emanuel, 2005; Walsh et al., 2012; Chan et al., 2018; Wang et al., 2020a; Shi et al., 2021; Xi et al., 2023). TCs are often accompanied by heavy rainfall and secondary disasters such as flooding, landslides, and debris flows (Woodruff et al., 2013; Cogan et al., 2018; Utsumi and Kim, 2022). Under a warming climate, the average translation speed of TCs has been reported to decrease globally by about 10 % over 1949–2016. Meanwhile, TC precipitation rates are projected to increase, with a global median increase of approximately 14 % across models, which may contribute to more prolonged and damaging impacts in affected regions (Kossin, 2018; Elsner, 2020; Knutson et al., 2020) . In the Northwest Pacific and the South China Sea, where TCs are referred to as typhoons, typhoons have frequently affected the Asia–Pacific region in recent years, including China, South Korea, and Japan, with southeastern China being particularly vulnerable and experiencing substantial casualties and economic losses (Esteban and Longarte-Galnares, 2010; Jung et al., 2024; Wang et al., 2024). Some intense typhoons and their residual circulations have even penetrated deep into inland China. For example, Typhoon In-Fa (2106) in 2021 and the double typhoons Doksuri (2305) and Khanun (2306) in 2023 caused record-breaking rainfall and devastating flooding in Henan province and Beijing–Tianjin–Hebei region in China (Shi et al., 2022; Zhao et al., 2024). Typhoons are moisture-sensitive systems whose development depends strongly on the availability and transport of atmospheric water vapor. PWV therefore provides an effective metric for depicting the pathways and strength of moisture transport during typhoons and for tracking moisture-related variations associated with their evolution and movement. Accordingly, high-accuracy PWV estimates are crucial for understanding the mechanisms underlying extreme precipitation during TCs, improving TC monitoring and forecasting, and supporting disaster risk assessment and mitigation efforts.

Various measurement techniques have been employed to retrieve PWV, including radiosondes, water vapor radiometers, satellite-based microwave/infrared remote sensing, and sun photometers (Ichoku et al., 2002; King et al., 2003; Li et al., 2003; Turner et al., 2007). In recent decades, advances in Global Navigation Satellite System (GNSS) technologies and Low Earth Orbit (LEO) satellites have enabled the widespread application of ground-based GNSS and space-based GNSS Radio Occultation (RO) for atmospheric observations (e.g., Melbourne et al., 1994; Kursinski et al., 1996; Li et al., 2017). However, these techniques have limitations when applied to fast-moving, ocean-based, and moisture-complex weather systems such as TCs, and obtaining PWV with both high resolution and temporal continuity remains challenging. In this context, gridded global reanalysis datasets provide atmospheric fields with high spatial and temporal resolution and no gaps, making it possible to obtain PWV at any time and location through interpolation. Reanalysis data therefore offer valuable resources for retrospectively investigating moisture transport and evolution during typhoons and for characterizing the spatio-temporal variability of atmospheric water vapor. Currently, several research institutions provide atmospheric reanalysis datasets that are widely used, including the fifth-generation European Centre for Medium-Range Weather Forecasts (ECMWF) Reanalysis (ERA5) (Hersbach et al., 2020), the Modern-Era Retrospective Analysis for Research and Applications Version 2 (MERRA-2) (Gelaro et al., 2017), and the Japanese Reanalysis for Three Quarters of a Century (JRA-3Q) (Kosaka et al., 2024). Despite their broad use, the accuracy of reanalysis-based PWV can be affected by differences in data assimilation strategies and the underlying numerical models, as well as by the uneven spatio-temporal distribution of assimilated observations. Therefore, reanalysis-based PWV products should be systematically evaluated for accuracy and applicability before practical application.

Previous studies have evaluated the accuracy of PWV estimates from reanalysis datasets across different regions. These evaluations commonly rely on independent observations, including ground-based GNSS, radiosondes, and GNSS RO, as reference data. In particular, ground-based GNSS-PWV typically achieves high accuracy (within 1–2 mm) and is not assimilated into any of the three reanalysis datasets examined in this study, making it an independent and reliable reference for validation (Wang et al., 2020b; Li et al., 2025). However, because GNSS stations are sparsely distributed over the oceans, radiosonde and GNSS RO observations are often used as complementary validation sources. At the regional scale, ERA5 generally achieves lower PWV errors (< 1 mm) over China and outperforms its predecessor ERA-Interim, with cross-validation using radiosonde data further supporting its reliability (Zhang et al., 2019a, b). In India, ERA5 also clearly outperforms MERRA-2 in PWV monitoring (Rani and Singh, 2025). Over the southern Tibetan Plateau, multiple reanalysis products exhibit systematic positive biases in the seasonal PWV cycle, likely linked to the persistent wet bias in regional models (Wang et al., 2017). In the Arctic, the Copernicus Arctic Regional Reanalysis (CARRA) shows good agreement with radiosonde observations, but with evident seasonal variability (Zhang et al., 2025). At the global scale, reanalysis PWV agrees well with GNSS, radiosonde, and RO observations (Zhang et al., 2018). However, considerable uncertainties remain in tropical and Southern Hemisphere regions, particularly for PWV estimates from the National Center for Environmental Prediction/Department of Energy (NCEP/DOE) dataset (Vey et al., 2010). Geographic and climatic factors have also been shown to influence the consistency between reanalysis and GNSS-derived PWV (Bock and Parracho, 2019). Among the various products, ERA5 generally performs better than others, while a newly released Chinese global reanalysis dataset, China Meteorological Administration-40 (CRA40), exhibits comparable performance to ERA5 in PWV estimation (Wang et al., 2020b; Li et al., 2025). However, as CRA40 data are not fully publicly available, it is not included in this study. GNSS RO data have also been used to compare ERA5 and MERRA-2 PWV estimates in tropical and subtropical regions (Johnston et al., 2021). Additionally, humidity modeling discontinuities in ERA5 at 09:00 and 21:00 UTC have been reported, which introduce diurnal jumps in Zenith Tropospheric Delay (ZTD) and subsequently affect PWV estimations (Yuan et al., 2025).

Existing evaluation studies have primarily focused on long-term and large-scale averages, with limited systematic evaluation of PWV accuracy from reanalysis datasets under extreme weather conditions such as typhoons. Moreover, comparative analyses of PWV estimation accuracy between typhoon and non-typhoon periods remain scarce. Driven by this research gap, we provide a systematic evaluation of PWV estimates from ERA5, JRA-3Q, and MERRA-2 during typhoon events using ground-based GNSS, radiosonde, and RO data from January 2020 to December 2024 over the Northwest Pacific and East Asia region. The data and methods used in this study are introduced in Sect. 2. Section 3 presents the evaluation results, and Sect. 4 provides the discussion and conclusions.

2 Data and methodology

In this section, the typhoon datasets, the three reanalyses, and the observational data from ground-based GNSS, radiosondes, and RO used in this study are introduced. The PWV retrieval method, data quality control procedures, and the spatio-temporal collocation strategy are also described. In addition, the evaluation metrics and their formulations are provided.

2.1 Typhoon data

The typhoon dataset used in this study is obtained from the International Best Track Archive for Climate Stewardship (IBTrACS), which provides a unified, publicly available best-track dataset by integrating historical and recent records from multiple operational agencies worldwide (Knapp et al., 2010; Gahtan et al., 2024). IBTrACS includes key parameters such as TC-center position, time, maximum sustained wind speed, and wind radius, typically reported at 3 and 6 h intervals. According to IBTrACS, the reporting agency currently responsible for the Northwest Pacific is the Japan Meteorological Agency (JMA), therefore, the JMA-reported typhoon information is used in this study. Because most GNSS stations used in this study are located in China and the China Meteorological Administration (CMA) and JMA intensity classifications are broadly comparable, typhoon categories are assigned following the CMA classification scheme, which defines six categories based on wind speed: tropical depression (10.8–17.1 m s⁻¹), tropical storm (17.2–24.4 m s⁻¹), severe tropical storm (24.5–32.6 m s⁻¹), typhoon (32.7–41.4 m s⁻¹), severe typhoon (41.5–50.9 m s⁻¹), and super typhoon (≥ 51.0 m s⁻¹). For simplicity, these are denoted as L1 to L6, respectively. In this study, L1 typhoons are not included since no stations are collocated with their tracks.

There are 113 typhoons recorded from January 2020 to December 2024, and the spatial distribution of their tracks is shown in Fig. 1a, together with the locations of the GNSS and radiosonde stations. Figure 1b illustrates the number of typhoons by category over the five-year period, summed by month. Typhoons occur most frequently from June to November, with relatively higher counts in August, September, and October. The largest number of tropical storms (L2) occurs in August, and the largest number of super typhoons (L6) occurs in September.

Figure 1 (a) Spatial distribution of typhoon tracks, ground-based GNSS stations from IGS (marked by circles) and CMONOC (marked by triangles), and IGRA radiosonde stations (marked by diamonds) with colorbars indicating wind speed (in m s⁻¹, left) and station elevation (in m, right). (b) Number of typhoons in each month and category from January 2020 to December 2024.

2.2 Reanalysis data

2.2.1 ERA5

ERA5 is the fifth generation of ECMWF reanalysis for global climate and weather, providing a consistent record since 1940 and representing a state-of-the-art product developed by ECMWF (Hersbach et al., 2020; Soci et al., 2024). ERA5 is produced using an ensemble-based four-dimensional variational (4D-Var) data assimilation system and model forecasts in CY41R2 of the ECMWF Integrated Forecast System (IFS), with 137 hybrid pressure levels in the vertical and the top level at 0.01 hPa. Atmospheric data are available on interpolated 37 pressure levels. The horizontal resolution is 0.25° × 0.25°, and the finest temporal resolution is hourly.

2.2.2 MERRA-2

MERRA-2, developed by the Global Modeling and Assimilation Office (GMAO) at the National Aeronautics and Space Administration (NASA), provides global atmospheric reanalysis data starting from 1980 (Gelaro et al., 2017). Compared to its predecessor MERRA, MERRA-2 incorporates improvements in the Goddard Earth Observing System (GEOS) model and the Gridpoint Statistical Interpolation (GSI) assimilation system, enabling the assimilation of modern satellite observations such as hyperspectral radiance, microwave radiances, RO, and NASA ozone profiles. MERRA-2 has a horizontal resolution of 0.5° × 0.625° and 72 vertical levels up to 0.01 hPa, and the 3-hourly instantaneous data are used in this study.

2.2.3 JRA-3Q

JRA-3Q is produced by the JMA using an advanced global NWP system to improve the quality and temporal coverage of long-term reanalysis (Kosaka et al., 2024). It builds on developments since JRA-55 and extending the record back to September 1947 and including many notable typhoon events. JRA-3Q assimilates a wide range of reprocessed observational datasets, including rescued historical observations and satellite data provided by global meteorological and satellite centers. It employs a 4D-Var data assimilation system and addresses many of the limitations found in JRA-55, resulting in a high-quality and consistent dataset spanning over 75 years. The vertical structure includes 45 pressure levels, with a horizontal resolution of 0.375° × 0.375° and a temporal resolution of 6 h. JRA-3Q dataset was officially released in 2022, and related evaluation studies are still relatively limited. The evaluation conducted in this study may serve as a helpful reference for its future application.

Table 1 summarizes the data providers, spatial and temporal resolutions, start time, update frequency, and assimilation strategies of the three reanalysis datasets. Further details can be found in the corresponding documentation (Gelaro et al., 2017; Hersbach et al., 2020; Bell et al., 2021; Kosaka et al., 2024).

Table 1 Summary of three atmospheric reanalysis datasets used in this study.

2.3 GNSS data

2.3.1 Ground-based GNSS data

Ground-based GNSS provides continuous, all-weather, high-precision observations of atmospheric variables with high temporal resolution, making it well-suited for evaluating reanalysis data under typhoon conditions. The International GNSS Service (IGS) provides Zenith Path Delay (ZPD, also known as ZTD) products at 5 min intervals. However, most IGS stations in the Asia–Pacific region are concentrated in Japan and South Korea, with sparse coverage along China's southeastern coast. To improve spatial coverage, GNSS data from the Crustal Movement Observation Network of China (CMONOC) are also used. These data are processed using the Position and Navigation Data Analyst (PANDA) software developed by Wuhan University (Shi et al., 2008), based on the Precise Point Positioning (PPP) technique (Zumberge et al., 1997), to generate ZTD estimates at the same temporal resolution as the IGS products. A quality control procedure is applied to exclude loosely constrained ZTD estimates, defined as those deviating from the station's monthly mean by more than four standard deviations (Zhang et al., 2017). The locations and elevations of IGS and CMONOC stations are shown in Fig. 1a. To minimize errors in the vertical interpolation of PWV, stations with elevations greater than 500 m are excluded. In total, 32 IGS and 29 CMONOC stations are used in this study.

2.3.2 GNSS RO data

GNSS RO data of the Constellation Observing System for Meteorology, Ionosphere, and Climate-2 (COSMIC-2) mission, as the successor to COSMIC, are provided via the COSMIC Data Analysis and Archive Center (CDAAC) in Boulder, USA (Schreiner et al., 2020). The orbital inclination of the COSMIC-2 constellation is designed to increase RO sampling in tropical and subtropical regions, and most profiles are distributed within a latitude range of 45° N to 45° S. PWV estimation requires profiles of specific humidity or water vapor partial pressure profiles. This study uses near-real-time wet profile (hereafter “wetPrf”) from the Level 2 products. The wetPrfs provide atmospheric parameters with a vertical sampling of 50 m below 20 km and 100 m between 20 and 60 km (the upper limit of the profiles). These data are retrieved using a one-dimensional variational (1D-Var) technique, and the lowermost height varies among profiles (Wee et al., 2022). In tropical and subtropical regions, super-refraction can prevent signals from penetrating to the surface, resulting in variations in the lowermost height across COSMIC-2 profiles (Schreiner et al., 2020; Wang et al., 2022). To reduce vertical interpolation errors, only wetPrfs that pass the CDAAC quality control procedures and reach below 500 m are used in this study.

2.4 Radiosonde data

Radiosonde data used in this study are obtained from the Integrated Global Radiosonde Archive (IGRA), with routine observations conducted twice daily at approximately 00:00 and 12:00 UTC. PWV derived from radiosonde profiles typically has an uncertainty of 5 % to 8 % (Pérez-Ramírez et al., 2014; Turner et al., 2003). Despite certain limitations, radiosonde observations remain a standard reference for evaluating the PWV retrieved from other techniques (Gui et al., 2017). Quality control follows procedures described in previous studies (Wang and Zhang, 2008; Zhang et al., 2017), with additional criteria applied to ensure profile completeness and temporal coverage: (1) Profiles must reach at least 300 hPa, and only profiles with at least 30 vertical levels are used for PWV calculation; (2) Profiles containing large data gaps, defined as pressure intervals exceeding 200 hPa between successive humidity measurements, are discarded; (3) Stations must operate continuously from January 2020 to December 2024, with at least 200 observations per year; (4) Stations with elevations exceeding 500 m are excluded (Shi et al., 2025). Based on these criteria, a total of 25 radiosonde stations are retained for evaluation, and their locations are shown in Fig. 1.

2.5 PWV estimation methods

PWV can be derived using two approaches. The first integrates specific humidity or water vapor partial pressure and is applied to reanalysis, radiosonde, and RO data. The second converts the Zenith Wet Delay (ZWD) estimated from GNSS PPP into PWV using a conversion factor.

2.5.1 Reanalysis-PWV, radiosonde-PWV, and RO-PWV

The vertical integration for PWV can be expressed as follows:

$$\begin{array}{}\text{(1)}& \mathrm{PWV}={\int }_{p_{\mathrm{1}}}^{p_{\mathrm{2}}}{\displaystyle \frac{q}{{\mathit{\rho }}_{\mathrm{w}}\cdot g_{\mathrm{s}}}}\mathrm{d}p\end{array}$$

where p₁ and p₂ (in hPa) are the upper and lower pressure limits of the integration, respectively, q is the specific humidity (in g kg⁻¹), ρ_w is the density of liquid water, whose value is 1000 kg m⁻³, and g_s is the mean gravitational acceleration at the station. It is defined as:

$$\begin{array}{}\text{(2)}& \begin{array}{rl}{g}_{\mathrm{s}}\left(\mathit{\phi },h\right)& =g_{\mathrm{n}}(\mathrm{1}-\mathrm{0.0026373}\cos \left(\mathrm{2}\mathit{\phi }\right)\\ & +\mathrm{5.9}\cdot {\mathrm{10}}^{-\mathrm{6}}{\cos }^{\mathrm{2}}\left(\mathrm{2}\mathit{\phi }\right))\\ & \cdot \left(\mathrm{1}-\mathrm{3.14}\cdot {\mathrm{10}}^{-\mathrm{7}}\cdot h\right)\end{array}\end{array}$$

where g_n=9.80665 m s⁻² is the value of normal gravity, φ is latitude (in rad), h is elevation (in m). Specific humidity q can also be replaced by water vapor partial pressure e, with the conversion relationship as follows:

$$\begin{array}{}\text{(3)}& e={\displaystyle \frac{qp}{\mathrm{0.622}+\mathrm{0.378}q}}\end{array}$$

where p is pressure (in hPa).

2.5.2 GNSS-PWV

The ZTD, comprising Zenith Hydrostatic Delay (ZHD) and ZWD, can be accurately estimated from GNSS observations using PPP, with ZHD precisely computed using the Saastamoinen model (Saastamoinen, 1972; Elgered et al., 1991):

$$\begin{array}{}\text{(4)}& \mathrm{ZHD}={\displaystyle \frac{\mathrm{0.002277}\cdot p_{\mathrm{s}}}{\mathrm{1}-\mathrm{0.00266}\cdot \cos \left(\mathrm{2}\mathit{\phi }\right)-\mathrm{0.00028}\cdot H}}\end{array}$$

where φ is latitude (in rad) of the station, H is the ellipsoidal height of the station (in km). p_s represents the surface air pressure (in hPa). The ZWD is obtained by subtracting the ZHD from the ZTD:

$$\begin{array}{}\text{(5)}& \mathrm{ZWD}=\mathrm{ZTD}-\mathrm{ZHD}\end{array}$$

the PWV can be obtained by multiplying ZWD by the water vapor conversion factor:

$$\begin{array}{}\text{(6)}& \mathrm{PWV}=\mathrm{\Pi }\times \mathrm{ZWD}\end{array}$$

where Π is the water vapor conversion factor, expressed as:

$$\begin{array}{}\text{(7)}& \mathrm{\Pi }={\displaystyle \frac{\mathrm{1}}{{\mathrm{10}}^{-\mathrm{6}}{\mathit{\rho }}_{\mathrm{w}}{R}_{\mathrm{v}}\left[\left(k_{\mathrm{3}}/T_{\mathrm{m}}\right)+k_{\mathrm{2}}^{\prime }\right]}}\end{array}$$

where ρ_w=1000 kg m⁻³ is the density of liquid water, R_v=461.51 J K⁻¹ kg⁻¹ is the specific gas constant for water vapor, $k_{\mathrm{2}}^{\prime }=\mathrm{17}$ ± 10 K hPa⁻¹ and k₃=3.776 ± 0.004 × 10⁵ K² hPa⁻¹ are atmospheric refractivity constants, respectively. T_m denotes the weighted mean temperature, which can be determined using either an empirical linear model (Bevis et al., 1994) or by integration based on meteorological reanalysis data. It has been indicated in previous studies that calculating T_m using the integration method provides greater accuracy (Wang et al., 2005). The formula is presented as Eq. (8) (Davis et al., 1985; Bevis et al., 1992):

$$\begin{array}{}\text{(8)}& T_{\mathrm{m}}={\displaystyle \frac{{\int }_{h_{\mathrm{s}}}^{\mathrm{\infty }}\left(e/T\right)\mathrm{d}h}{{\int }_{h_{\mathrm{s}}}^{\mathrm{\infty }}\left(e/T^{\mathrm{2}}\right)\mathrm{d}h}}={\displaystyle \frac{{\sum }_{\mathrm{1}}^n\stackrel{\mathrm{‾}}{\left(\frac{e_i}{T_i}\right)}\left(h_i-h_{i-\mathrm{1}}\right)}{{\sum }_{\mathrm{1}}^n\stackrel{\mathrm{‾}}{\left(\frac{e_i}{T_i^{\mathrm{2}}}\right)}\left(h_i-h_{i-\mathrm{1}}\right)}}\end{array}$$

where e is the water vapor pressure at the station's zenith (in hPa), T is the temperature (in K), h_s is the height of the station (in m), and h is the height (in m).

2.6 Data collocation schemes and PWV horizontal interpolation method

Wind radius can be used to characterize the spatial extent of a typhoon's influence. For the Northwest Pacific, recent IBTrACS records are primarily based on information from the JMA, which includes wind radius reported for different wind speed thresholds. In this study, a dynamic collocation scheme is adopted using the 30 kt wind radius to associate GNSS and radiosonde stations as well as RO profiles with the typhoon track. Specifically, the 30 kt wind radius includes the longest and shortest radius of 30 kt winds, and their mean is taken as the effective influence radius at each time. For GNSS and radiosonde observations, a station is considered collocated with the typhoon center when the distance between the station and the typhoon center is smaller than the time dependent effective radius. Temporally, the value at the target time is obtained by linear interpolation using the two observations closest in time. For RO profiles, a profile is regarded as collocated when the distance between the mean tangent point location and the typhoon center is smaller than the corresponding time dependent effective radius and the time difference does not exceed 30 min.

For PWV horizontal interpolation, bilinear interpolation is used. The target location is given by (lon,lat), where $\mathrm{lon}\in \left[{\mathrm{lon}}_{\mathrm{1}},{\mathrm{lon}}_{\mathrm{2}}\right]$, and $\mathrm{lat}\in \left[{\mathrm{lat}}_{\mathrm{1}},{\mathrm{lat}}_{\mathrm{2}}\right]$. The bilinear interpolation is expressed as:

$$\begin{array}{}\text{(9)}& \left\{\begin{array}{l}{w}_{\mathrm{lon}}=\frac{\mathrm{lon}-{\mathrm{lon}}_{\mathrm{1}}}{{\mathrm{lon}}_{\mathrm{2}}-{\mathrm{lon}}_{\mathrm{1}}},w_{\mathrm{lat}}=\frac{\mathrm{lat}-{\mathrm{lat}}_{\mathrm{1}}}{{\mathrm{lat}}_{\mathrm{2}}-{\mathrm{lat}}_{\mathrm{1}}}\\ \mathrm{PWV}(\mathrm{lon},\mathrm{lat})=\left(\mathrm{1}-w_{\mathrm{lon}}\right)\left(\mathrm{1}-w_{\mathrm{lat}}\right){\mathrm{PWV}}_{\mathrm{11}}\\ +w_{\mathrm{lon}}\left(\mathrm{1}-w_{\mathrm{lat}}\right){\mathrm{PWV}}_{\mathrm{21}}\\ +\left(\mathrm{1}-w_{\mathrm{lon}}\right)w_{\mathrm{lat}}{\mathrm{PWV}}_{\mathrm{12}}+w_{\mathrm{lon}}{w}_{\mathrm{lat}}{\mathrm{PWV}}_{\mathrm{22}}\end{array}\right.\end{array}$$

where a w_lon and w_lat are the weights along the longitudinal and latitudinal directions, respectively. PWV₁₁ denotes the PWV value at (lon₁,lat₁), and PWV₂₁, PWV₁₂, and PWV₂₂ are defined analogously at (lon₂,lat₁), (lon₁,lat₂), and (lon₂,lat₂).

2.7 Statistical metrics

The statistical metrics used in this study include the systematic deviation of PWV from reanalysis with respect to observations, denoted as bias, Relative Bias (RB), Root Mean Square Error (RMSE), and de-biased RMSE (dRMSE). The dRMSE represents the residual error after bias removal and is used to isolate the random component of PWV errors. They are calculated by the following formulas:

$$\begin{array}{}\text{(10)}& \left\{\begin{array}{l}\mathrm{Bias}={\mathrm{PWV}}_{\mathrm{reanalysis}}-{\mathrm{PWV}}_{\mathrm{reference}}\\ \mathrm{RB}=\frac{{\mathrm{PWV}}_{\mathrm{reanalysis}}-{\mathrm{PWV}}_{\mathrm{reference}}}{{\mathrm{PWV}}_{\mathrm{reference}}}\times \mathrm{100}\mathit{\%}\\ \mathrm{RMSE}=\sqrt{\frac{\mathrm{1}}{N}\sum _{i=\mathrm{1}}^N{\left({\mathrm{PWV}}_{\mathrm{reanalysis},i}-{\mathrm{PWV}}_{\mathrm{reference},i}\right)}^{\mathrm{2}}}\\ \mathrm{dRMSE}=\sqrt{\frac{\mathrm{1}}{N}\sum _{i=\mathrm{1}}^N{\left[\left({\mathrm{PWV}}_{\mathrm{reanalysis},i}-{\mathrm{PWV}}_{\mathrm{reference},i}\right)-\stackrel{\mathrm{‾}}{\mathrm{Bias}}\right]}^{\mathrm{2}}}\end{array}\right.\end{array}$$

where PWV_reanalysis represents the PWV from different reanalysis data, and PWV_reference is the reference PWV, which can be GNSS-PWV, radiosonde-PWV, or RO-PWV, $\stackrel{\mathrm{‾}}{\mathrm{Bias}}$ is the mean value of bias.

In addition, to characterize the local spatial variability of reanalysis data, the Neighborhood Standard Deviation (NSD) is adopted as one of the evaluation metrics (Wei et al., 2013). NSD is used as a descriptive metric of subgrid scale spatial variability in reanalysis PWV around the target location. Larger NSD indicates stronger spatial heterogeneity within neighboring grid cells and thus a higher potential for representativeness mismatch when gridded fields are compared with point observations (Bock and Parracho, 2019). Because TCs are asymmetric and can exhibit sharp moisture gradients, high NSD does not necessarily imply a large reanalysis error, but provides context for interpreting the evaluation results under typhoon conditions. As illustrated in Fig. 2, NSD is computed at two spatial scales using the standard deviation of PWV from the 4 nearest grid points (NSD-4) and the 16 nearest grid points (NSD-16) around the typhoon center.

Figure 2 Schematic of NSD calculation. The orange line shows the typhoon track; the red marker indicates the typhoon center. The 16 nearest grid points are shown in dark blue and light blue, with the 4 closest points in light blue.

NSD-4 and NSD-16 can be calculated according to the following equation:

$$\begin{array}{}\text{(11)}& \left\{\begin{array}{l}\stackrel{\mathrm{‾}}{\mathrm{PWV}}=\frac{\mathrm{1}}{M}\sum _{j=\mathrm{1}}^M{\mathrm{PWV}}_j\\ \mathrm{NSD}=\sqrt{\frac{\mathrm{1}}{M}\sum _{j=\mathrm{1}}^M{\left({\mathrm{PWV}}_j-\stackrel{\mathrm{‾}}{\mathrm{PWV}}\right)}^{\mathrm{2}}}\end{array}\right.\end{array}$$

where $\stackrel{\mathrm{‾}}{\mathrm{PWV}}$ denotes the mean PWV of the neighboring points, and M is set to 4 or 16.

3 Results

This section presents the evaluation results of PWV estimates from the three reanalysis datasets. For clarity, prefixes are used to distinguish PWV derived from different data sources and their corresponding statistical metrics. Specifically, GNSS-PWV, radiosonde-PWV (RS-PWV), and RO-PWV denote PWV retrieved from ground-based GNSS, radiosonde, and RO data, respectively. For the reanalysis products, the initial letter of each dataset is used as a prefix. Accordingly, E-PWV, M-PWV, and J-PWV refer to PWV estimated from ERA5, MERRA-2, and JRA-3Q, respectively. The prefix “REA-” is used as a collective term for PWV estimates and associated statistical metrics derived from the three reanalysis datasets.

3.1 Evaluation using ground-based GNSS data

In this subsection, subscripts are used to distinguish results from different GNSS networks when needed. Specifically, GNSS-PWV_I and GNSS-PWV_C denote PWV estimated from IGS and CMONOC stations, respectively.

3.1.1 Monthly evaluation

To provide a comprehensive understanding of the PWV accuracy of the three reanalysis datasets, evaluations are conducted using data from 64 GNSS stations spanning January 2020 to December 2024. Monthly mean PWV from the reanalyses is compared to GNSS-PWV, with results shown in Fig. 3.

Figure 3 Monthly mean evaluation results. Panels (a)–(c) show REA-PWV (mm); red and blue bars denote REA-PWV sampled at the locations of CMONOC (REA-PWV_C) and IGS (REA-PWV_I) stations, respectively. Panels (d)–(f) show REA-bias and REA-RMSE (mm); red and yellow bars denote bias and RMSE at CMONOC stations (REA-bias_C and REA-RMSE_C), respectively, and blue and green bars denote bias and RMSE at IGS stations (REA-bias_I and REA-RMSE_I), respectively. Panels (g)–(i) show REA-RB (%); red and blue bars denote RB at CMONOC (REA-RB_C) and IGS (REA-RB_I) stations, respectively. Subscripts C and I indicate that reanalysis PWV is sampled at CMONOC and IGS station locations, respectively, and that the corresponding statistical metrics (bias, RMSE, and RB) are computed against GNSS observations from the same network.

Monthly mean PWV from the three reanalyses (E-PWV, M-PWV, and J-PWV) show high consistency and pronounced seasonality in Fig. 3a–c, peaking at about 55 mm in July–August and dropping to about 20 mm in boreal winter. Notably, the PWV maximum slightly precedes the main typhoon season (August–October). REA-PWV_C is generally higher than REA-PWV_I from February to November but lower in the remaining months, mainly reflecting the different spatial distribution of GNSS stations from these two networks. Among the reanalyses, E-PWV and M-PWV are very similar, whereas J-PWV remains consistently lower than both E-PWV and M-PWV across all months. Biases in Fig. 3d–f indicate that ERA5 and JRA-3Q have negative biases in every month (red and blue bars), with the absolute value of the bias positively correlated with PWV, consistent with previous studies indicating that mean REA-PWV is negatively biased in low-latitude regions (Wang et al., 2020b). The largest E-bias_C and E-bias_I occur in October and July (−1.70 and −1.66 mm), while the largest J-bias_C and J-bias_I both occur in August (−4.59 and −3.69 mm). In contrast, M-PWV agrees well with GNSS-PWV, with M-bias_C and M-bias_I remaining below 1 mm and approaching near-zero from May to September. RMSE (yellow and green bars in Fig. 3d–f) follows a PWV-like seasonal cycle, with ERA5 having the smallest RMSE, followed by MERRA-2, and JRA-3Q the largest. RB in Fig. 3g–i is smaller in summer and larger in winter, as expected because it represents bias normalized by PWV, so its seasonality is largely governed by the PWV distribution and aligns with the bias pattern.

For the months with more than ten typhoons within the five-year period (June–November), network-weighted mean bias and RMSE are obtained by combining the CMONOC- and IGS-based results, with weights determined by each network's proportion of the total number of used GNSS stations. The resulting E-bias, M-bias, and J-bias are −1.26, −0.73, and −3.16 mm, respectively, while the E-RMSE, M-RMSE, and J-RMSE are 3.18, 4.25, and 4.55 mm, respectively (not shown). Among the three datasets, MERRA-2 shows the smallest systematic deviation, whereas ERA5 yields the lowest RMSE, suggesting greater stability for E-PWV. For JRA-3Q, both the absolute value of J-bias and J-RMSE are the largest among the three reanalyses in most months. The next section presents a more detailed evaluation of REA-PWV under typhoon conditions.

3.1.2 Composite evaluation considering adjacent periods of typhoons

Typhoon monitoring agencies release typhoon data based on wind speed. Data recording begins or ends when the wind speed reaches or falls below a specified threshold. It has been found that water vapor plays a critical role in TC formation, high column water vapor appears near the pouch center and starts to increase about 42 h prior to genesis, while a substantial increase in precipitation occurs within 24 h before genesis (Wang, 2014; Wang and Hankes, 2016). Moreover, even after typhoon dissipation or passage, residual circulation can continue to exert influence (Duan et al., 2014). Therefore, in addition to the recorded typhoon period (hereafter “r-typhoon”), REA-PWV is also evaluated during a 3 d Adjacent Period (AP), defined as the 3 d before and the 3 d after r-typhoon. A schematic of r-typhoon and AP is provided in the top panel of Fig. 4.

Figure 4 The top panel shows the schematic timeline of r-typhoon and the AP. Panels (a)–(l) present biases and RMSEs under different AP settings across typhoon categories – ERA5 (a–d); MERRA-2 (e–h); JRA-3Q (i–l). Within each bar cluster, color from light to dark correspond to AP ± 3, AP ± 2, AP ± 1 d, and r-typhoon, respectively. “L2-6” and “Non” denote the results for all typhoons and for non-typhoon periods, respectively. Panels (c) and (d) report the numbers of CMONOC and IGS station–typhoon pairs. Panels (m) and (n) summarize the mean bias and RMSE averaged over all typhoons as a function of AP.

Figure 4a–l shows bar charts of REA-biases and REA-RMSEs under different AP settings and typhoon categories. Here, “Non” denotes the mean results for the same stations and calendar periods in non-typhoon years during 2020–2024; for example, for typhoon Doksuri in 2023, “Non” refers to the same stations and dates in 2020, 2021, 2022, and 2024, when no typhoon occurred. The magnitude of the REA-bias generally decreases as AP shortens from ± 3 d to r-typhoon period, meaning that biases tend to move closer to zero or shift from negative to positive. For ERA5 and JRA-3Q, both E-bias and J-bias decrease overall with a shorter AP, and E-bias at CMONOC stations becomes slightly positive during r-typhoon period (Fig. 4a, b, i, j). In contrast, M-bias increases slightly as AP shortens and turns more clearly positive during r-typhoon period (Fig. 4e, f), with a larger absolute value than under the ± 3 d AP. Under non-typhoon conditions, the “Non” group in Fig. 4a, b, e, f, i, j shows consistently negative REA-biases and no pronounced dependence on AP length. For RMSE, E-RMSE and J-RMSE exhibit a slight overall decrease as AP shortens (Fig. 4c, d, k, l), although this tendency is weak for E-RMSE computed at IGS stations. In addition, J-RMSE decreases more markedly with shorter AP for weaker typhoons (L2–L4; Fig. 4k). By contrast, M-RMSE during r-typhoon period is slightly higher than during AP (Fig. 4g, h), indicating slightly larger errors when only r-typhoon period is considered. Similar to the bias, RMSE under non-typhoon conditions shows no pronounced dependence on AP (the “Non” clusters in Fig. 4c, d, g, h, k, l).

When aggregated across all typhoons, E-bias and J-bias are smaller than those during non-typhoon periods regardless of whether AP is included (the two rightmost bar groups in Fig. 4a, b, i, j). The corresponding RMSEs are also slightly lower than their non-typhoon counterparts (the two rightmost bar groups in Fig. 4c, d, k, l). This indicates that E-PWV and J-PWV exhibit higher overall accuracy during typhoon periods than during non-typhoon periods. Moreover, as AP shortens, E-bias, J-bias, E-RMSE, and J-RMSE further decrease, indicating that r-typhoon-only evaluation shows smaller errors than evaluations that include AP. For MERRA-2, M-bias is negative and close to zero when AP spans ± 1 to ± 3 d, whereas it becomes positive during r-typhoon period, with a larger absolute value than that during non-typhoon periods (Fig. 4e, f, two rightmost bar groups). Meanwhile, M-RMSE is slightly higher during typhoon periods than during non-typhoon periods and is consistently larger than E-RMSE, while remaining comparable to J-RMSE (comparison of the two rightmost bar groups in Fig. 4g, h with those in Fig. 4c, d, k, l). The summaries in Fig. 4m, n quantify these differences: during r-typhoon period, the mean E-bias, M-bias, and J-bias are −0.18, 0.76, and −1.86 mm, respectively, and the mean E-RMSE, M-RMSE, and J-RMSE are 3.30, 5.13, and 4.44 mm, respectively. Overall, E-PWV shows the best performance during r-typhoon periods. M-PWV shows a positive shift in bias, but the mean bias remains within 1 mm. J-PWV improves relative to non-typhoon periods, although a negative bias remains.

3.1.3 Evaluation during recorded typhoon periods

This section focuses on the accuracy of PWVs during r-typhoon from reanalyses, without considering AP. Table 2 presents mean REA-biases, REA-RMSEs and dRMSEs of all typhoons during both r-typhoon and non-typhoon.

Table 2 Mean REA-bias, REA-RMSE, and REA-dRMSE (mm) during r-typhoon period and non-typhoon periods at all GNSS stations, aggregated over all typhoons.

Based on the bias results in Table 2, E-bias and J-bias during r-typhoon period are both closer to zero than during non-typhoon periods, indicating reduced underestimation in E-PWV and J-PWV under typhoon conditions. Among the three reanalyses, E-bias is closest to zero. In contrast, M-bias is larger during r-typhoon period than during non-typhoon periods and changes sign from negative to positive. Even so, r-typhoon M-bias (0.76 mm) is smaller in absolute value than r-typhoon J-bias (−1.86 mm). During non-typhoon periods, all three reanalyses underestimate PWV (E-bias = −1.48 mm, M-bias = −0.63 mm, and J-bias = −3.74 mm), with MERRA-2 showing the smallest bias magnitude. During r-typhoon period, the underestimation in ERA5 and JRA-3Q is alleviated, whereas MERRA-2 shifts to a moderate overestimation.

The RMSE results further show that, relative to non-typhoon periods, E-RMSE and J-RMSE during r-typhoon period decrease by 0.17 and 0.69 mm, respectively, whereas M-RMSE increases by 0.37 mm. E-RMSE is also the smallest among the three reanalyses, indicating the lowest overall error for E-PWV during r-typhoon period, followed by J-PWV, while M-PWV shows comparatively larger uncertainty. Since RMSE reflects both systematic and random components of error, Table 2 also reports dRMSE to characterize the random component after removing the bias. The dRMSE results show that all three REA-dRMSEs are higher during r-typhoon period than during non-typhoon periods, indicating an overall increase in the random error component under typhoon conditions. Notably, for ERA5 and JRA-3Q in Table 2, RMSE decreases from non-typhoon periods to r-typhoon period while dRMSE increases, indicating that the lower RMSE during r-typhoon is mainly driven by a reduction in systematic bias, which offsets part of the increase in the random error component. Conversely, the larger RMSEs during non-typhoon periods are primarily associated with larger biases rather than larger random errors. Overall, Table 2 suggests improved performance for E-PWV and J-PWV during r-typhoon period, whereas MERRA-2 shows a sign change in bias and a modest increase in RMSE, although its bias magnitude remains smaller than that of JRA-3Q.

Figure 5 shows the spatial distributions of REA-biases, REA-RMSEs, and REA-dRMSEs at GNSS stations during r-typhoon period for all typhoons. For the biases (Fig. 5a, d, g), E-bias exhibits both positive and negative values across stations and is overall closer to zero. In contrast, M-bias is predominantly positive, whereas J-bias is dominated by negative values (Fig. 5d, g), consistent with the earlier finding that M-bias shifts from negative to positive during r-typhoon period. Stations with positive E-bias and M-bias show some spatial similarity and are relatively concentrated in southern and southeastern China, while positive J-bias occurs at only a few stations. For the RMSEs (Fig. 5b, e, h), E-RMSE is generally the smallest and shows a comparatively homogeneous spatial distribution (Fig. 5b). M-RMSE is relatively low over southern China but exceeds 6 mm at a few stations in southeastern China (Fig. 5e). J-RMSE falls between E-RMSE and M-RMSE overall (Fig. 5h), with isolated stations showing comparatively large values (e.g., ∼ 7.91 mm at GXHC), which may be related to the limited number of available stations–typhoon matchups, making the statistics more sensitive to individual events. For dRMSE (Fig. 5c, f, i), E-dRMSE remains the lowest, followed by J-dRMSE, whereas M-dRMSE is the largest and shows a latitudinal gradient, with larger values at higher latitudes. Notably, the difference between J-RMSE and J-dRMSE is comparatively pronounced (Fig. 5h, i), indicating a larger contribution of systematic bias to J-RMSE.

Figure 5 Spatial distribution of mean REA-biases, REA-RMSEs, and REA-dRMSEs (in mm) at GNSS stations during r-typhoon period, aggregated over all typhoons.

3.2 Evaluation using radiosonde observations

Although all three reanalyses assimilate radiosonde observations, RS-PWV is still used as a reference because of its high accuracy, to further evaluate REA-PWVs. Figure 6 compares REA-PWVs against RS-PWV across different typhoon stages and categories, including the 3 d before r-typhoon (pre-typhoon, Fig. 6a), r-typhoon (r-typhoon, Fig. 6b), and the 3 d after r-typhoon (post-typhoon, Fig. 6c), as well as r-typhoon results for individual categories (Fig. 6d–h). Figure 6i provides the non-typhoon counterpart for comparison, and Fig. 6j–k summarizes the corresponding biases and RMSEs using bar charts. For convenience, E-r, M-r, and J-r denote the correlation coefficients (r) between RS-PWV and E-PWV, M-PWV, and J-PWV, respectively.

Figure 6 Correlation statistics of REA-PWVs and RS-PWVs for (a) pre-typhoon period (3 d before r-typhoon), (b) r-typhoon period, (c) post-typhoon period (3 d after r-typhoon), (d–h) r-typhoon periods for different categories, and (i) non-typhoon periods. Panels (j) and (k) summarize the corresponding REA-biases and REA-RMSEs for panels (a)–(i). Correlation coefficient (r) and the number of samples are shown in each panel (a)–(i).

Comparing the three typhoon stages (Fig. 6a–c), E-r, M-r, and J-r are identical for pre-typhoon and r-typhoon periods, with E-r and J-r exceeding 0.9 and M-r equal to 0.88, consistent with the relatively larger M-RMSE during r-typhoon noted in Sect. 3.1. Correlations increase modestly in the post-typhoon period. For r-typhoon results stratified by typhoon category (Fig. 6d–h), all three correlations are higher for L2, L5, and L6 and lower for L3 and L4. Overall, E-r is the highest across categories and remains above 0.91; J-r is slightly lower and matches E-r during L6; M-r is generally the lowest but peaks during L6, reaching 0.93. Both J-r and M-r reach their minima during L4, at 0.82 and 0.88, respectively. When comparing the aggregated r-typhoon and non-typhoon results, E-r, M-r, and J-r are essentially the same, with values of 0.94, 0.88, and 0.92. To further quantify accuracy differences between r-typhoon and non-typhoon conditions, Table 3 reports the corresponding summary statistics, which align with the “R” and “Non” bar groups in Fig. 6j–k.

The conclusions drawn using RS-PWV as the reference are generally consistent with those based on GNSS-PWV. As shown in Table 3, during r-typhoon period, E-bias and J-bias are closer to zero than during non-typhoon period, decreasing by 0.30 and 0.78 mm, respectively. E-RMSE and J-RMSE are also smaller during r-typhoon period, with reductions of 0.21 and 0.51 mm. By contrast, M-RMSE is nearly unchanged between r-typhoon and non-typhoon periods, differing by only 0.02 mm. For dRMSE, M-dRMSE and J-dRMSE are slightly higher during r-typhoon period than during non-typhoon period, whereas E-dRMSE is slightly lower during r-typhoon period. It should be noted that the RS-PWV- and GNSS-PWV-based statistics do not always rank the three reanalyses in the same way. For example, during r-typhoon period, M-bias is the smallest among the three REA-biases. Such differences may be related to the relatively limited number of radiosonde observations and a smaller set of matched samples available for evaluation, which makes the statistics more sensitive to individual samples. Therefore, the results in this section are best interpreted as complementary that supports the GNSS-referenced evaluation.

Table 3 Mean REA-bias, REA-RMSE, and REA-dRMSE (mm) during r-typhoon and non-typhoon periods at all radiosonde stations, aggregated over all typhoons.

3.3 Comparison with COSMIC-2 RO data

Evaluation using GNSS-PWV and RS-PWV as references is constrained by station locations, whereas RO-PWV can be used to compare REA-PWVs over the ocean. Based on the collocation scheme described in Sect. 2.6, 160 and 1299 COSMIC-2 RO profiles are collocated with typhoon centers during r-typhoon and non-typhoon periods, respectively, and these RO profiles are mainly distributed over the ocean. Given that only 160 samples are available during r-typhoon period, results stratified by typhoon category are not further reported in this section. Figure 7 shows the comparison results for (a) r-typhoon and (b) non-typhoon periods. During r-typhoon period, E-r, M-r, and J-r are 0.81, 0.77, and 0.80, respectively, which are lower than those during non-typhoon period by 0.09, 0.12, and 0.10. Table 4 further summarizes biases, RMSEs, and dRMSEs for r-typhoon and non-typhoon periods.

Figure 7 Correlation statistics of REA-PWVs and RO-PWVs for (a) r-typhoon and (b) non-typhoon periods. Correlation coefficients (r) and the number of samples are shown in each panel.

Table 4 Mean REA-bias, REA-RMSE, and REA-dRMSE (mm) during r-typhoon and non-typhoon periods when compared with RO-PWV, aggregated over all typhoons.

The error characteristics derived from comparisons with RO-PWV differ from those in the previous two sections, where GNSS-PWV and RS-PWV were used as references. During r-typhoon period, all three REA-biases are positive and are generally larger than those during non-typhoon period, indicating that the reanalyses tend to overestimate PWV relative to RO-PWV under typhoon conditions. Specifically, during r-typhoon period, E-bias and M-bias are 1.27 and 1.12 mm, respectively, are comparable, whereas J-bias exceeds 3 mm. The REA-RMSEs are also larger during r-typhoon period than during non-typhoon period, suggesting increased overall uncertainty. In particular, M-RMSE and J-RMSE are 5.17 and 5.25 mm, respectively, exceeding E-RMSE by 1.05 and 1.13 mm. For dRMSE, M-dRMSE and J-dRMSE increase during r-typhoon period relative to non-typhoon period, whereas E-dRMSE decreases slightly. Overall, across both periods, MERRA-2 shows the smallest bias, while ERA5 yields the lowest RMSE. It should be noted that although more than 4500 COSMIC-2 RO profiles are available per day in tropical and subtropical regions, the number of RO profiles used for comparison remains limited under the current collocation scheme. The differences between the RO-based results and those based on GNSS-PWV and RS-PWV may be related to the limited number of collocated samples and differences in land–ocean sampling, because most RO profiles used here are located over the ocean, whereas GNSS and radiosonde stations are mainly distributed over land. In addition, under typhoon influence, the water vapor field over the ocean tends to be more variable and spatially inhomogeneous, which may also contribute to larger discrepancies in this comparison.

3.4 Local spatial heterogeneity of REA-PWVs

NSD is used to characterize the local spatial heterogeneity of gridded REA-PWVs and thus reflects the consistency among neighboring grid points. Comparing NSD between typhoon and non-typhoon periods helps quantify the complexity of PWV spatial structure under typhoon conditions. Table 5 summarizes NSD-16 and NSD-4 for the three REA-PWVs during r-typhoon and non-typhoon periods across all typhoon events.

Table 5 Mean NSD-16 and NSD-4 (mm) of three REA-PWVs during r-typhoon and non-typhoon periods over all typhoons.

According to Table 5, NSD-16 is generally larger than NSD-4 because it spans a broader neighborhood and includes more grid points, indicating that spatial variability is larger when evaluated over a broader neighborhood. Relative to non-typhoon periods, NSD increases for all three reanalyses during r-typhoon period, suggesting a more complex local PWV structure and stronger spatial inhomogeneity under typhoon conditions. Specifically, E-NSD-16, M-NSD-16, and J-NSD-16 increase by 15.5 %, 50.9 %, and 38.1 %, respectively, while the corresponding increases for NSD-4 are larger, reaching 30.4 %, 69.3 %, and 47.1 %. The relative increase is larger for NSD-4 than for NSD-16, indicating that NSD changes more strongly when evaluated at a smaller neighborhood scale. Systematic differences are also evident among the reanalyses. MERRA-2 exhibits the largest NSD at both scales and the strongest increase from non-typhoon to r-typhoon, whereas ERA5 and JRA-3Q show lower NSD values and more moderate changes. It should be noted that NSD describes spatial variability of the local PWV field rather than PWV error itself. Therefore, higher NSD during r-typhoon period can be interpreted indicating stronger local spatial variability when comparing REA-PWVs.

4 Discussion and conclusions

This study presents a systematic evaluation of PWV during typhoons from the ERA5, MERRA-2, and JRA-3Q reanalyses using three independent observational references: ground-based GNSS, radiosonde, and COSMIC-2 RO data. A total of 113 typhoon events from January 2020 to December 2024 are analyzed. The evaluation focuses on contrasts between typhoon and non-typhoon periods in bias, RMSE, and dRMSE, and further examines dependence on typhoon category when sample sizes permit.

In the ground-based evaluations, GNSS provides higher temporal resolution and a larger sample size and is therefore used as the primary reference in this study. The radiosonde-based evaluation serves as an independent external check. Using GNSS-PWV as the reference, ERA5 shows the best overall performance under typhoon conditions, with the bias closest to zero and the lowest RMSE (−0.18 and 3.30 mm). M-bias lies between E-bias and J-bias, while M-RMSE is the largest among the three reanalyses. J-PWV exhibits a more pronounced negative bias, and J-RMSE falls between E-RMSE and M-RMSE. Relative to non-typhoon period, the absolute values of E-bias and J-bias decrease, indicating reduced mean systematic bias in E-PWV and J-PWV during typhoons. M-bias increases slightly but remains below 1 mm, suggesting that MERRA-2 maintains a bias level comparable to that in non-typhoon period. In terms of RMSE, E-RMSE and J-RMSE decrease by 0.17 and 0.69 mm during typhoon period, respectively. Together with the bias changes, this indicates an overall reduction in total error for E-PWV and J-PWV, whereas M-RMSE is higher during typhoon period than during non-typhoon period. By construction, dRMSE is no greater than RMSE, however, all three REA-dRMSEs are higher during typhoon period than during non-typhoon period, which may reflect more pronounced short-term variability of water vapor and a stronger random error component under typhoon conditions. Therefore, changes in REA-PWV performance during typhoons should be interpreted by jointly considering both systematic and random error components. When RS-PWV is used as the reference, the three REA-biases during typhoon period are overall closer to zero than those during non-typhoon period, with M-bias being the smallest, followed by E-bias, while J-bias remains relatively pronounced. For RMSE, E-RMSE and J-RMSE are lower during typhoon period than during non-typhoon period, whereas M-RMSE is broadly comparable between the two periods.

Using RO-PWV as the reference yields patterns that are not fully consistent with the ground-based evaluations. Under typhoon conditions, the three REA-biases are positive overall, and both their absolute values and RMSEs are larger than those under non-typhoon conditions. This indicates that, for the RO collocations available in this study, REA-PWVs tend to overestimate PWV relative to RO-PWV. This contrast does not contradict the GNSS-based and radiosonde-based conclusions because the RO collocations are concentrated over oceanic regions, introducing a clear land–ocean sampling difference relative to the samples used in the ground-based evaluations. Moreover, under the relatively strict collocation scheme, only 160 collocated samples are available under typhoon conditions, making the statistics more sensitive to sample composition. In terms of inter-product comparison, E-bias and M-bias are comparable under typhoon conditions, with M-bias lower by 0.15 mm, whereas J-bias is larger. For RMSE, E-RMSE is the lowest and is 1.05 mm lower than M-RMSE, while J-RMSE is the largest. Overall, these results continue to indicate the relatively robust performance of ERA5. In addition, the NSD results show enhanced local spatial variability in REA-PWV under typhoon conditions, reflecting a more complex spatial structure of the water vapor field during typhoons.

Consistent evidence from three independent references suggests that the reanalyses can provide continuous PWV information that is practically useful and reasonably reliable under typhoon conditions. Overall, ERA5 shows more consistent performance and lower error levels. JRA-3Q improves relative to non-typhoon conditions in most cases, although a systematic bias can remain in some situations. MERRA-2 tends to maintain a bias close to zero, but its uncertainty becomes more pronounced under typhoon conditions, consistent with its higher NSD and its larger increases in RMSE and NSD from non-typhoon to typhoon conditions.

This study supports the reliability and applicability of reanalysis-derived PWV under typhoon conditions and provides continuous water vapor information for typhoon monitoring and process-based analyses. It also provides a performance reference for reanalysis PWV under typhoon weather. Over oceanic regions where ground-based observations are sparse, reanalysis PWV can serve as an important reference for PWV retrieved from other technologies, facilitating product evaluation and bias diagnosis. This reference can also provide contextual information for interpreting discrepancies among different observations under typhoon conditions, thereby improving the interpretability and utility of oceanic PWV estimates in both scientific research and operational applications related to typhoons.