A Strategy for Merging Objective Estimates of Global Daily Precipitation from Gauge Observations, Satellite Estimates, and Numerical Predictions

Precipitation is one of the most vital components of the global climate system. Gridded precipitation at the global scale is important for advancing various applications of numerical weather prediction, climate model simulation, water resource management, agricultural science, and disaster risk management (Huffman and Klepp, 2011). Gauge observations, satellite estimates (SEs), and numerical simulations are three existing principal sources of precipitation at the global scale. Gauge observations can provide accurate intensity measurements and long-term records from station locations, but they suffer from incomplete global coverage, especially over vast ocean areas (Schneider, 1993). In addition, the lack of a uniform time standard of reporting between different areas increases the difficulty of incorporating daily precipitation into global analyses (Huffman et al., 2001). Based on indirect relationships between the measurement quantities of satellite observation and actual precipitation intensity (PI), infrared (IR) and passive microwave (PMW) SEs have the capacity to retrieve spatial distributions of precipitation over regions that are not covered by dense gauge observational networks (Behrangi et al., 2009). IR sensors can provide cloud-top temperature measurements, which are converted to rainfall rates based on empirical algorithms. However, these algorithms occasionally perform poorly when the indirect relationships are weak. PMW sensors can estimate precipitation more directly by sensing the thermal emission of raindrops and the scattering of upwelling radiation from the earth to space. However, because these PMW sensors are restricted to only polar-orbiting satellites, the spatial and temporal samplings associated with these products are significantly limited (Joyce et al., 2004). The numerical simulation of precipitation performs relatively well at middle and high latitudes, but performs poorly at low latitudes (Arpe, 1991) because great uncertainties exist in the intensities and locations of precipitation distributions owing to the accumulation of model errors caused by various physical processes (Silva et al., 2011).

The deficiencies in individual precipitation sources have naturally led to attempts to combine their information to exploit the advantages of each for producing optimal estimations of global precipitation with complete global coverage. Using the threshold-matched index (TMPI) and daily-rescale method calibrated with monthly gauge data, the Global Precipitation Climatology Project (GPCP) produced a higher resolution dataset, named GPCP one-degree daily (GPCP-1DD), by merging daily sounding data from low-earth polar-orbit satellites (Huffman et al., 2001). Capturing the spatial and temporal variability of daily precipitation observations to a high degree with complete global coverage, the GPCP-1DD product has been widely used in many fields such as weather process diagnostic analyses, hydrologic streamflow model validation, and climate change applications (Crow, 2007; Ploshay and Lau, 2010, Medvigy and Beaulieu, 2012).

GPCP-1DD was the first product developed for estimating daily precipitation with complete global coverage from different data sources, but only daily precipitation information from SEs and monthly precipitation information from gauge observations were merged. Because daily information from gauge observations was not included, biases still existed in the intra-month-scale precipitation estimates of GPCP-1DD compared to actual gauge observations (Bolvin et al., 2009). Furthermore, because the error estimation technique in GPCP-1DD relied heavily on several empirical relationships, systemic errors of daily precipitation estimates could be found when comparing them to ground truth data in terms of their magnitude and spatial structures (McPhee and Margulis, 2005). Moreover, because the GPCP-1DD dataset required a monthly scale gauge-based precipitation product to remove its monthly bias, it was mainly delivered in post-real-time (with an approximate three-month lag). Therefore, the latest information from the available real-time global daily precipitation datasets [from SEs, model predictions (MPs) and gauge measurements] is still missing in the GPCP-1DD product. All of these factors limit the application of GPCP-1DD data in global near-real-time climate monitoring and surface hydrology forecasting.

Recently, merging studies of daily precipitation (Nie et al., 2015) and hourly precipitation (Pan et al., 2012) data from different sources have been carried out in China. However, because of their focus on China, these studies do not provide reliable global precipitation products for operational application in global climate monitoring and global climate model forecasting. In this paper, a strategy is presented for merging global daily precipitation information from various sources to develop a new method for constructing a gridded global daily precipitation product. Relying on a bias correction procedure and an objective merging algorithm, daily-scale gauge observations, advanced SEs combined from IR and PMW sensors, and predictions from a numerical model were used to take advantage of each source and produce a better gridded daily precipitation product with complete global coverage. Further, because all precipitation sources used in this study can be achieved in real time, the global merged daily precipitation output analysis has great potential for climate applications by various organizations, such as the National Climate Center of the China Meteorological Administration (NCC/CMA).

This study is organized into five sections: section 2 presents the three individual data sources and independent validation data used in the study; section 3 describes the global daily precipitation merging scheme in detail; section 4 provides the validation results and comparison among different global datasets; and section 5 provides conclusions and discussion.

The original source of gauge observations was the Global Telecommunication System (GTS) synoptic weather reports, which are archived at the National Meteorological Information Center of the CMA (Zhang et al., 2004). This dataset contains hourly (e.g., 0000, 0100, 0200 UTC), three-hourly (e.g., 0000, 0300, 0600 UTC), and daily observations at stations in different countries with all major meteorological variables in coding formats. A new global gauge daily precipitation dataset (NGDP) was obtained by applying multi-step automatic quality control procedures, including a duplication station check, internal consistency check, extreme value check, temporal consistency check, and spatial consistency check, on the original GTS precipitation dataset (Nie et al., 2012). The NGDP dataset contains approximately 15 000 land and island stations across the globe (Fig. 1), with substantial improvements in accuracy and completeness among the daily precipitation reports compared with the original GTS dataset. In this study, the NGDP dataset was used as a gauge observation source for the merging processes.

Figure 1.  Spatial distribution of the gauge stations included in the NGDP data. The marked boxes indicate six monsoon regions (M1-M6 represent the Indian monsoon region, the East Asian monsoon region, the Australian monsoon region, the north American monsoon region, the South American monsoon region, and the West African monsoon region, respectively) and six ocean regions (O1-O6 represent the Indian Ocean region, the Equatorial Pacific region, the South Pacific region, the North Pacific region, the South Atlantic region, and the North Atlantic region, respectively) for regional comparisons of daily precipitation in section 4.

The satellite-based estimates from the CPC's morphing technique (CMORPH) (Joyce et al., 2004) combine precipitation estimates from multiple IR and PMW sensors. The original CMORPH data were available at a 0.25^°× 0.25^° and 3-h resolution between 60^°N and 60^°S from December 2002. The daily-scale data used in this study were obtained based on an area average onto a 1^°× 1^° grid box and time accumulation on each calendar day.

The NCEP-NCAR reanalysis (hereafter, NCEP) (Kalnay et al., 1996) precipitation data included short-range forecast accumulations from the model based on physics and parameterizations, which were not adjusted by precipitation observations. The original NCEP precipitation data had four 6-h forecast times (0000, 0600, 1200, 1800 UTC) and a triangular 62 wave truncated ("T62") horizontal resolution (approximately 2.5^°× 2.5^°). The daily total precipitation used in this study was produced by summing four forecasts each day and by interpolating them onto a 1^°× 1^° grid box over the entire globe using the bilinear interpolation method.

The GPCP-1DD dataset provides global 1^°× 1^° gridded fields of daily precipitation from October 1996 to the present (Huffman et al., 2001). These data are based on a combination of satellite and gauge data. In the 40^°N-40^°S bands, the TMPI calculated from geostationary-satellite IR sensors was used to compute three-hourly estimates for each 1^°× 1^° grid box. Then, daily precipitation estimates were produced by summing the three-hourly TMPI estimates for each day. Outside the 40^°N-40^°S bands, precipitation estimates were computed based on TOVS data. Finally, the 1DD data were rescaled to match the monthly accumulation of the GPCP satellite-gauge combination precipitation estimates (Adler et al., 2003) on a monthly time scale for both datasets. Because the daily precipitation information from the CMORPH estimates, NCEP forecasts, and GTS gauge observations was not included in the GPCP-1DD product, it is reasonable to use them as the independent validation data to evaluate the quality of the merged daily precipitation analyses produced in this study.

The Tropical Rainfall Measuring Mission (TRMM) Multi-satellite Precipitation Analysis (TMPA) 3B42 near-real-time (RT) version 7 product is a popular and widely used product that was developed by the GSFC (Huffman et al., 2007). The TMPA 3B42RT product (hereafter, TMPA) has provided a three-hourly temporal resolution and 0.25^°× 0.25^° spatial resolution product since 2000. It was used as a comparison dataset for the multi-satellite sources in the multi-source merged BMEP-d analysis in this study.

Systemic bias and random error of all input precipitation sources are two major reasons for uncertainties in merging processes. The systemic bias is the error source that can be characterized statistically, whereas the random error is the stochastic component of the error that does not exhibit regularity (AghaKouchak et al., 2012). A multi-step merging technique was designed first, to remove systemic bias and then to objectively reduce random error from all individual sources to produce a global daily precipitation analysis that is superior to each individual source. The study period was from 1 January 2011 to 31 December 2014.

Figure 2.  The (a) G-none area-mean and (b) G-located area-mean ratios of monthly NCEP (solid lines) and CMORPH (dashed lines) data compared to the NGDP data, and (c) comparisons between the ratios of the G-none area and G-located area during the period 2011-14.

There was an assumption that the NGDP is unbiased at the observation points because gauge data yield relatively accurate measurements of precipitation at the observation locations. This assumption can be considered reasonable, as the NGDP data exhibited significantly improved performance against the GTS data after comprehensive quality control procedures (Nie et al., 2012). Based on this assumption, a two-step strategy was adopted to reduce biases in the CMORPH and NCEP data across the globe.

Over global land areas with gauge stations (hereafter, G-located areas), the bias correction procedure was similar to that employed by (Nie et al., 2015), who matched the cumulative distribution function (CDF) of the CMORPH and NCEP data with that of the NGDP data. A short summary of the CDF methodology is given here. In this study, the NGDP and CMORPH (or NCEP) data are defined as the "reference" and "target" data, respectively. According to the spatial scale of precipitation observations (Nie et al., 2008), the reference and target precipitation data pairs were matched over each 5^°× 5^° grid box with at least one gauge station located therein, and in each calendar month to consider the seasonal variation factors of the biases. To eliminate dry days from the statistics, only precipitation data larger than 0.1 mm d^-1 were included as valid records in the matching procedure.

Over oceans and remote land areas (hereafter, G-none areas), it is difficult to implement the CDF strategy because not enough reference data can be found in the target grid. A rescaled strategy was employed to correct the overall bias from the G-none areas in the target data. This strategy is based on an additional assumption that the monthly precipitation ratios between bias-corrected target data and original target data in the G-none areas were approximately equal to the ratios in the G-located areas. The purpose of this assumption was to make the bias-corrected NCEP and CMORPH monthly data over ocean areas have a somewhat low overall bias compared to the original NCEP and CMORPH data based on the existing gauge precipitation observation network. This assumption is not optimal because the intensity and patterns of daily precipitation are dependent on the regions of precipitation occurrences. However, because correcting bias in satellite and model precipitation data is exceedingly difficult (Rudolf and Rubel, 2005), and because the gauge data are usually not reported with sufficient density over most ocean areas, this assumption is considered appropriate for correcting the overall bias of satellite and model precipitation over ocean areas, to a certain extent.

According to this assumption, the spatial average monthly precipitation of the bias-corrected target data in the G-none areas can be scaled as follows: \begin{equation} P_{{bc, O},g,m}\approx P_{{ori, O},g,m}\times (P_{{bc, L},g,m}/P_{{ori, L},g,m}) ,(1) \end{equation} where the subscripts O and L denote the G-none and G-located areas, respectively; the subscripts g and m denote the spatial average monthly precipitation for month m; and the subscripts ori and bc represent precipitation before and after the bias correction processes. The overall bias over ocean areas can be reduced, and a typically smaller bias of monthly precipitation from the target data can be obtained with this equation. Then, the bias-corrected monthly precipitation P _{bc, O,m} in grid box i is scaled by the corresponding ratio of bias-corrected target data to the original target data over the G-none areas using the following formula: \begin{equation} P_{{bc, O},i,m}=P_{{ori, O},i,m}\times (P_{{bc, O},g,m}/P_{{ori, O},g,m}) , (2)\end{equation} where P _{bc, O,i,m} and P _{ori, O,i,m} are the bias-corrected and original monthly precipitation data in grid box i over the G-none areas for month m, respectively. This equation was used to transfer the spatial inhomogeneity of the bias correction effects to each 1^°× 1^° grid box over ocean areas. Subsequently, the bias-corrected daily precipitation in grid box i was computed using the daily/monthly ratios to scale each daily field during the month as follows: \begin{equation} P_{{bc, O},i,m,d}=P_{{ori, O},i,m,d}\times (P_{{bc, O},i,m}/P_{{ori, O},i,m}) , (3)\end{equation} where P _{bc, O,i,m,d} and P _{ori, O,i,m,d} are the bias-corrected and original monthly precipitation over the G-none areas at day d in the month m, respectively. This equation was used to transfer the monthly bias correction effects from each grid box to each day in each calendar month.

Figure 2 shows the ratios of the monthly spatial average NCEP (CMORPH) target data to the NGDP reference data over the G-located areas and G-none areas during the study period. The ratios of both the NCEP and CMORPH data show similar values over global land areas (Fig. 2a) and global ocean areas (Fig. 2b). The comparisons between the ratios in the G-none areas and the G-located areas for NCEP and CMORPH data in each month are shown in Fig. 2c. The comparison curves are relatively stable at approximately one (with a variability ≈ 10%) during the entire period for both the NCEP and CMORPH data. This result indicates that an approximate 10% uncertainty bias might still exist in the bias-corrected monthly NCEP and CMORPH data over all of the ocean areas. This uncertainty of overall-ocean bias in the target data exhibits a similar magnitude to that of the TMPA data over ocean areas for the 30^°S-30^°N latitude band (Huffman et al., 2007). Therefore, the additional assumption of Eq. (1) can be considered reasonable and acceptable for both the NCEP and CMORPH target data during the study period.

Figure 3.  A flow chart of the HL-OI gauge-satellite-model merging strategy. The bias correction data flows are shown with dashed lines, and the combination data flows are shown with solid lines. The gridded data with complete global coverage are shaded. The final output analysis fields are shown as boxes with bold borders.

The bias-corrected NCEP and CMORPH data are combined with NGDP data using an algorithm based on optimal interpolation (OI) (Gandin, 1965). The daily precipitation analysis x _a,i at the i-th grid box was calculated by modifying the background x _b,i using the observation vectors x _o around the gird box as follows: \begin{equation} x_{{a},i}=x_{{b},i}+{w}_i^{T}({x}_{o}-{x}_{b}) , (4)\end{equation} where w_i is the weighting coefficient for the analysis increment in the grid box, which is determined by the background error σ _b and observation error σ _o, and is estimated as follows: \begin{equation} {w}_i=({B}+{O})^{-1}{b}_i , (5)\end{equation} where B and O are the background and observation error covariance matrices and b_i is the background error covariance vector. In this study, the NGDP data are used as the observation to modify the background (CMORPH and NCEP data) across the globe.

The quantification of magnitudes for σ _b and σ _o is the key to successfully carrying out the OI-based combining algorithm. An objective error estimation method, which was developed by (Hollingsworth and Lönnberg, 1986) and applied to the OI method (referred to as the HL-OI scheme) to merge high-resolution precipitation over China (Nie et al., 2015), was employed here to quantity the observation and background errors of daily precipitation at the global scale. The HL-OI scheme has been described in detail by (Nie et al., 2015), but a short summary of the scheme is nevertheless also provided here: Under the assumption that no correlation exists between the observation and background errors, the background departure covariance at any two points l and k can be calculated as follows: \begin{equation} {cov}_{l,k}\approx{H}_l{BH}_k^{T}+{O}_{l,k}=\left\{ \begin{array}{l@{\quad}l} \sigma_{{b},l}^2+\sigma_{{o},l}^2 & k=l\\[2mm] \sigma_{{b},l}^2\rho_{{b},l}(r_{l,k}) & k\neq l \end{array} \right. , (6)\end{equation} where H is the observation operator, ρ _b is the spatial correlation of the background error covariance, and r is the spatial distance. By calculating a histogram of the background departure covariance against spatial distance, the background and observation errors can be separated objectively. Then, using Eqs. (4) and (5), a daily precipitation analysis was produced as a weighted average of the background and observations.

The strategy used to merge gauge, satellite and model data in this study is outlined in Fig. 3. The SEs and MPs were adjusted by the CDF matching technique to reduce the biases for each grid box in the G-located areas. Alternately, the ratio of spatially averaged monthly precipitation estimates between the global G-none areas and the G-located areas was used to provide consistent and unbiased global SE and MP data in the G-none areas. This eliminates the bias in the SE and MP data on the global scale while maintaining the local details and temporal variation of the data. The unbiased SEs (used as observation fields) were then combined with unbiased MPs (used as background fields) with the HL-OI scheme, thus producing merged SE-MP daily precipitation data with complete global coverage. Finally, the gauge-satellite-model merged analysis, named the Beijing Climate Center Merged Estimation of Precipitation with daily resolution (BMEP-d), on a 1^°× 1^° grid box covering 90^°S-90^°N was merged from a combination of NGDP data and SE-MP data using the HL-OI scheme. Benefiting from this strategy, a daily precipitation analysis at the global scale was produced by using daily-scale gauge observation information and applying the HL-OI scheme for the first time. Although the NGDP global gauge observations, CMORPH SEs and NCEP predictions were chosen from individual daily precipitation sources in this study, this HL-OI merging technique is not only suitable for these datasets. Because the background errors and observation error estimates in the analysis were objectively calculated by the HL-OI scheme without additional ad hoc assumptions, this daily precipitation merging strategy can also be conveniently applied to other datasets from different daily precipitation data sources.

Inferences about the quality of the merged BMEP-d analysis are made with respect to the existing popular satellite-gauge merged GPCP-1DD dataset with complete global coverage and the TMPA multi-satellite merged dataset in this section. In addition, daily precipitation from NCEP, CMORPH and NGDP data are used as comparisons with the BMEP-d merged analysis for each data source. The NGDP gauge data were interpolated onto the same regular 1^°× 1^° grid boxes as the other datasets using the inverse distance weighted interpolation method from (Cressman, 1959). Because the NGDP gauge stations were much sparser than the grid boxes in other gridded datasets, a unified spatial interpolation radius of 400 km was used to make the gauge observations influence more grid boxes. A comparison of the gridded NGDP data and the BMEP-d data can provide a clear understanding of the advantages of the HL-OI technique compared to a simple spatial interpolation method for reducing interpolation-based representative errors in the grid boxes around the gauge station.

Figure 4 shows the spatial distributions of biases for NCEP, CMORPH, NGDP, TMPA and BMEP-d against GPCP-1DD. In contrast to the GPCP-1DD data, the NCEP data exhibit wet biases in global tropical regions, whereas the CMORPH data exhibit obvious dry biases in the latitude bands between 40^°N (^°S)-60^°N (^°S). The NGDP gridded data exhibit significant wet biases compared to the GPCP data over areas where the gridded NGDP data have values (hereafter, G-value areas). Because gauge data always exhibit higher PI than satellite data and model data at the station locations, significant positive biases in the gridded NGDP data indicate that the simple spatial interpolation method might have led to a large amount of interpolation-based representative errors in the grid boxes around gauge stations. As a comparison, biases in both the annual and seasonal means of BMEP-d data remain less than 1.0 mm d^-1 in most parts of the globe, which are significantly less than with other datasets. This result indicates that the wet biases in NCEP and NGDP data, and the dry biases in CMORPH data, have been successfully reduced by the HL-OI merging scheme.

Comparisons of global-mean bias and absolute bias in different datasets against the GPCP-1DD product are shown in Table 1. BMEP-d had the smallest magnitude of biases compared to NCEP, NGDP and CMORPH, regardless of the annual and seasonal means (DJF for December-February; MAM for March-May; JJA for June-August; SON for September-November). Similarly, the absolute biases of BMEP-d in all seasons were also greatly reduced compared with those of other precipitation source data, which indicated substantial increases in the global consistency of the magnitude of daily precipitation between BMEP-d and GPCP-1DD. Compared to TMPA, BMEP-d also had much better bias performance against GPCP-1DD in terms of both spatial pattern (Fig. 4) and annual and seasonal mean intensity (Table 1). Because TMPA is a multi-satellite merged product, whereas the GPCP-1DD merged monthly gauge precipitation information with its estimate, this result indicates that more daily precipitation bias, compared to GPCP-1DD, can be reduced when precipitation information from both a gauge source and model source is merged, compared to from multiple satellite sources only.

Figure 4.  Time-averaged biases (mm d$^-1$) against the GPCP-1DD data of the NCEP, CMORPH, NGDP, TMPA and BMEP-d data for January-December (ANN) means (left column), December-February (DJF) means (middle column) and June-August (JJA) means (right column) during the period 2011-14. Gray coloring represents missing values.

Figure 5.  Time-averaged correlations (left) and RMSEs (right; units: mm d$^-1$) of NCEP, CMORPH, NGDP, TMPA and BMEP-d compared with GPCP-1DD during the period 2011-14.

The time correlation and RMSE of the different datasets versus GPCP-1DD are shown in Fig. 5. NCEP had the smallest correlation with GPCP-1DD, which was below 0.5 for most parts of the globe. Because CMORPH and TMPA used similar satellite sensors, these two datasets had comparable performance for both time correlation and RMSE against GPCP-1DD. The CMORPH and TMPA also had significantly higher correlation than NCEP, especially in tropical regions, where the correlations were almost consistently above 0.7. The gridded NGDP data had a small correlation and large RMSE with GPCP-1DD, especially in low-latitude (30^°S-30^°N) regions, where the RMSE could reach 10 mm d^-1. Generally, the gauge data always had higher PI than both the satellite data and model data. When the NGDP gauge data were interpolated onto each 1^°× 1^° grid box around the gauge station with a large interpolation radius (400 km), the gridded NGDP data in most of the grid boxes were actually obtained by spatial interpolation from the nearby gauge observations. Therefore, the large RMSE values of the gridded NGDP data in these regions were primarily due to large interpolation-based representative errors from the relatively simple spatial interpolation algorithm. The spatial distribution of the time correlation for BMEP-d was an approximate combination of NCEP, NGDP and CMORPH. Between the 60^°S-60^°N bands, where the NGDP and CMORPH data were available, BMEP-d had a similar correlation pattern to CMORPH in oceanic areas and to NGDP over most of land areas. Outside the 60^°S (60^°N) bands, the magnitudes and spatial structures of BMEP-d were close to those of NCEP because much less precipitation information from the gauge and satellite sources was available in those regions. Table 2 summarizes the annual and seasonal mean time correlations and RMSEs of different datasets against the GPCP-1DD. There was no significant seasonal variability in either the time correlations or RMSEs in all of the datasets. The global-mean RMSEs of BMEP-d were obviously smaller than those of NCEP, CMORPH and NGDP in all periods. These results indicate that the HL-OI technique can effectively reduce the daily PI errors from each input data source, leading to the RMSEs of BMEP-d being smaller than those of each input dataset across the globe.

Figure 6.  Zonal-average (a-c) annual-mean, (d-f) DJF-mean, and (g-i) JJA-mean precipitation (mm d$^-1$) over all global areas where the NGDP gridded data have values (left column), the global land areas (middle column), and the global ocean areas (right column) for GPCP-1DD (black lines), BMEP-d (red lines), NCEP (green lines), CMORPH (purple lines), TMPA (blue lines), and gridded NGDP (yellow lines).

Figure 6 provides the zonal-averaged latitudinal profiles of the BMEP-d merged analysis and the other five gridded daily precipitation datasets. Over the G-value areas, the profiles of NGDP data were always higher than those of other data. The profiles of GPCP-1DD and BMEP-d exhibited reliable consistency over G-valued areas for the annual mean and the DJF (JJA) means. Over global land and ocean areas, the general shapes of the annual, DJF and JJA means for GPCP-1DD, BMEP-d, TMPA and CMORPH were relatively close between the 40^°N-40^°S bands, where the PMW-derived precipitation estimates were used in the GPCP-1DD, TMPA and CMORPH products. Outside the 40^°N-40^°S bands, both CMORPH and TMPA seriously underestimated the zonal-averaged precipitation compared to the other data, especially over global ocean areas (Figs. 6c, f and i). In contrast to the underestimated profiles of CMORPH and TMPA between the 40^°N (^°S)-60^°N (^°S) bands, the BMEP-d and NCEP profiles exhibited high consistency with GPCP-1DD. Outside the 60^°N (60^°S) areas, where no SEs and fairly rare gauge observations exist, the model data dominated the merged analysis, which caused the profiles of BMEP-d to almost coincide with those of NCEP. Overall, depending on the daily precipitation information from satellite and gauge data in low- and midlatitude regions and the information from model and gauge data in high-latitude regions, the BMEP-d data provided highly consistent zonal-averaged precipitation with GPCP-1DD, not only for annual means but also for the DJF and JJA means during the study period.

The time series of spatial correlations, spatial biases, and spatial RMSEs of BMEP-d, NCEP, CMORPH, TMPA, and NGDP data against GPCP-1DD products are shown in Fig. 7. The clear gaps of BMEP-d and CMORPH around 540 days are due to an abnormally large number of missing data in the original CMORPH data (Figs. 7a-f). The spatial correlations of CMORPH, TMPA and BMEP-d are consistently above 0.7 and are significantly higher than those of NCEP and NGDP. TMPA has the highest spatial correlations with GPCP-1DD over the G-value areas and the global land areas. The spatial correlations of BMEP-d and TMPA show almost the same performance for the global oceanic areas and better performance than those of CMORPH and NCEP. Over the G-value areas, significantly positive (negative) biases are presented in NGDP (CMORPH) compared to GPCP-1DD. The spatial bias curves of BMEP-d are among the curves of the three input data sources and are closest to the zero line, especially over the G-value areas and the global ocean areas. Over global ocean areas, BMEP-d has the smallest spatial bias, at about -0.05 mm d^-1, compared to GPCP-1DD, whereas NCEP (CMORPH and TMPA) overestimated (underestimated) the ocean-mean precipitation by approximately 0.14 (-0.39 and -0.46) mm d^-1 (Fig. 7f). The spatial RMSEs of BMEP-d were much smaller than other datasets in both the G-value areas and the global land (ocean) areas (Figs. 7g-i). These results indicate that a combination of precipitation information from gauge observations, SEs and model simulations can reduce the random errors in each source dataset and improve the accuracy of the merged analysis.

Table 3 summarizes a comparison of the time-mean statistics of the NCEP, CMORPH, TMPA, and BMEP-d datasets. Because CMORPH and TMPA have no values outside the 60^°N-60^°S bands, the comparison mainly reflects the overall performance of each dataset between the bands. BMEP-d and GPCP-1DD exhibited the most consistency among these statistics. The global-mean spatial RMSEs between BMEP-d and GPCP-1DD decreased by 0.78 and 0.98 mm d^-1 from CMORPH and NCEP, respectively. The spatial correlation of TMPA was higher than for the other datasets because similar calibration methods were used in the input satellite data for TMPA and GPCP-1DD (Huffman et al., 2001, 2007). The spatial correlations between BMEP-d and GPCP-1DD were 0.70 (0.73) over the global land (ocean) regions, which were larger than for NCEP and almost equal to that of CMORPH. The spatial bias of CMOPRH (NCEP) was significantly negative-estimated (positive-estimated) by -0.38 (0.18) mm d^-1 compared to GPCP-1DD for the global average. After the merging processes, the spatial bias of BMEP-d significantly decreased to about -0.07 mm d^-1.

Figure 7.  Time series of spatial correlations, spatial biases, and spatial RMSEs for BMEP-d (red lines), NCEP (green lines), CMORPH (purple lines), TMPA (blue lines), and NGDP (yellow lines) versus GPCP-1DD over all global areas where the NGDP gridded data have values (left column), the global land areas (middle column), and the global ocean areas (right column) during the period 2011-14. All lines are obtained by carrying out a 21-day moving average onto the original time series.

Figure 8.  Seasonal cycles of monthly precipitation during the period 2011-14 at the global scale [(a) G-value areas; (b) global land areas; (c) global ocean areas] and the (d-i) regional monsoon regions and (j-o) oceanic regions. The monsoon regions in the G-value areas [(d) Indian monsoon region; (e) East Asian monsoon region; (f) Australian monsoon region; (g) North American monsoon region; (h) South American monsoon region; (i) West African monsoon region] and the oceanic regions [(j) Indian Ocean region; (k) Equatorial Pacific region; (l) South Pacific region; (m) North Pacific region; (n) South Atlantic region; (o) North Atlantic region] are displayed in the marked boxes of Fig. 1.

Figure 9.  The PDF (%) of daily PI (units: mm d$^-1$) for GPCP-1DD (black), BMEP-d (red), NCEP (green), CMORPH (purple) and TMPA (blue) between the 60$^\circ$N-60$^\circ$S bands during the period 2011-14.

Figure 10.  Statistical indices of BMEP-d, NCEP, CMORPH and TMPA versus GPCP-1DD for (a-c) POD, (d-f) FAR, and (g-i) ETS computed at different precipitation thresholds.

Figure 8 shows the seasonal variability of the different datasets at both the global and regional scales. At the global scale, the curves of the BMEP-d data were closest to that of GPCP-1DD (Figs. 8a-c). All of the datasets exhibited obvious seasonal cycles with seasonal peaks around July and annual minima around January in the global land areas (Fig. 8b). The temporal variabilities of all of the datasets in ocean areas (Fig. 8c) were significantly smaller than on land areas. Over the G-value areas, the large representative errors from spatial interpolation made the curves of the gridded NGDP data much higher than for the other datasets (Fig. 8a).

Several monsoon and oceanic regions (labeled boxes in Fig. 1) were selected to facilitate regional comparisons of seasonal variability (Figs. 8d-o). The seasonal variability amplitudes of NCEP were overestimated in the East Asian monsoon region, the South American monsoon region, and the West African monsoon region with respect to the amplitudes of GPCP-1DD, whereas they were underestimated in the North American monsoon region (Figs. 8d-i). The seasonal cycles of CMORPH were exaggerated in the North American region and the West African monsoon region (Figs. 8g and i) but were shifted downward compared to other data in the Indian and East Asian monsoon regions (Figs. 8d and e). Over the oceanic regions, seasonal cycles of BMEP-d represented more features of daily precipitation from CMORPH and NCEP because the NGDP stations are quite scarce in these regions. The cycles of BMEP-d, GPCP-1DD, TMPA and CMORPH exhibited similar seasonal variations over the Indian Ocean region, the Equatorial Pacific region, and the South Pacific region (Figs. 8j-l), whereas BMEP-d, GPCP-1DD and NCEP had consistent seasonal cycles over the North Pacific region and the North Atlantic region (Figs. 8m and o). This result shows that the satellite and model data have a greater capacity for characterizing the seasonal variations of daily precipitation in low- and high-latitude regions, respectively. By combining the advantages of all of these input precipitation data, the BMEP-d dataset exhibited the most consistency in seasonal curves with GPCP-1DD in almost all of the monsoon and oceanic regions. The smallest differences between the BMEP-d and GPCP-1DD data can likely be attributed to the location-dependent objective error reductions using the HL-OI scheme during the merging processes.

The PDF is another important statistic of the mean amount and spatiotemporal variation patterns of global daily precipitation. The PDFs of daily PI computed from each dataset in all of the 1^°× 1^° grid boxes between the 60^°N-60^°S bands are shown in Fig. 9. The PDFs were computed in six PI classes as follows: 0≤ PI<5 mm d^-1, 5 mm d\(^-1\le{PI}<10\) mm d^-1, 10 mm d^-1≤ PI<15 mm d^-1, 15 mm d^-1≤ PI<20 mm d^-1, 20 mm d\(^-1\le{PI}\le 25\) mm d^-1, and PI>25 mm d^-1. Because the PDF values of the class 0≤ PI<5 mm d^-1 of all of the datasets were approximately 80% and much higher than those of other classes, the probabilities of this class are displayed in Fig. 9a separately. The frequencies of BMEP-d in the light precipitation class (0≤ PI<5 mm d^-1) were closest to that of GPCP-1DD, especially in the global land areas (Fig. 9a). BMEP-d, TMPA and CMORPH tended to overestimate the frequency of light precipitation across the entire global area and the global land area, whereas NCEP tended to underestimate compared to GPCP-1DD. In the entire global area (Fig. 9b), the occurrence of moderate precipitation (5 mm d^-1≤ PI<15 mm d^-1) in NCEP (CMORPH) was overestimated (underestimated) by approximately 41% (-24%) compared to GPCP-1DD, whereas the occurrence of heavy precipitation ( PI>25 mm d^-1) in NCEP (CMORPH) was underestimated (overestimated) by about -60% (106%) against GPCP-1DD, respectively. BMEP-d performed better for the PDFs of moderate precipitation and heavy precipitation than NCEP, CMORPH, and TMPA. The occurrence of moderate precipitation and heavy precipitation estimated by BMEP-d was generally equivalent to that of GPCP-1DD in both global land and oceanic areas (Figs. 9c and d).

To further elucidate the accuracy of the merged BMEP-d analysis, a quantitative evaluation of it versus the GPCP-1DD product, compared with other datasets, was also performed. Three statistical indices for evaluating daily precipitation events, including the probability of detection (POD), false alarm ratio (FAR), and equitable threat score (ETS), were computed at different precipitation thresholds in all grid boxes across the globe. The definitions and calculations for these indices were described in detail by (Ebert et al., 2007).

Figure 10 shows the results of the POD, FAR and ETS scores for BMEP-d, NCEP, CMORPH and TMPA at different precipitation intensity thresholds. The PODs of all datasets had similar trends, which decreased with an increasing threshold (Figs. 10a-c). For thresholds of 0.1-5 mm d^-1, the POD values of BMEP-d and NCEP were larger than those of CMORPH (in the entire global area and the global ocean area) and TMPA; then, the POD values of NCEP decreased rapidly and were much smaller than those of BMEP-d, CMORPH and TMPA for thresholds larger than 10 mm d^-1. The FAR results indicate that NCEP had obviously higher false alarm scores than BMEP-d, TMPA and CMORPH for precipitation thresholds larger than 5 mm d^-1 (Figs. 10d-f). TMPA and CMOPRH had smaller FAR scores than BMEP-d for moderate precipitation (0.1-5 mm). For precipitation thresholds larger than 10 mm d^-1, BMEP-d exhibited better performance than other datasets over global land and oceanic areas. The overall capacity for detecting precipitation events was evaluated in terms of ETS scores. BMEP-d, TMPA and CMORPH exhibited increasing ETS scores for precipitation thresholds up to 5 mm d^-1 for the global land area and 10 mm d^-1 for the entire globe and global oceanic areas. Then, the ETS scores declined slowly at higher thresholds (Figs. 10g-i). The ETS scores of BMEP-d were higher than those of CMORPH, TMPA and NCEP at most precipitation thresholds larger than 10 mm d^-1. Summing up the results from POD, FAR and ETS, BMEP-d has better capability to capture intense precipitation events than CMORPH and NCEP when compared with GPCP-1DD.

In this study, a merging strategy using a bias-correction procedure and the HL-OI scheme was adopted to combine three individual kinds of global daily precipitation data from NGDP gauge observations (Nie et al., 2012), CMORPH SEs (Joyce et al., 2004), and NCEP reanalysis (Kalnay et al., 1996). This strategy was designed to produce a daily precipitation analysis with complete global coverage by three steps. First, a CDF matching procedure was performed to remove systemic biases over G-located areas. Then, based on the assumption that the monthly ratios between bias-corrected target data and original target data in the G-none areas were approximately equal to those in the G-located areas, the biases of SEs and MPs over ocean areas were successfully eliminated. Third, the HL-OI scheme was used to combine unbiased gauge observations, SEs, and MPs to reduce random error from each data source and to produce the gauge——satellite-model sources merged daily precipitation analysis, called BMEP-d, on 1^°× 1^° grid boxes with complete global coverage. The reliability of BMEP-d was evaluated against independent daily precipitation verification data across the globe during a four-year period from 2011 to 2014. The results of this study reveal the following:

(1) The strategy developed in this study can merge daily-scale precipitation information from gauge observations, SEs, and MPs to successfully construct a gridded global daily precipitation product. The merged BMEP-d data provide improved daily precipitation distributions with complete global coverage, and are of improved quality compared to each input data source from NCEP, NGDP and CMORPH. Compared with GPCP-1DD, the magnitudes of the global-averaged bias, absolute bias, and RMSE of BMEP-d are significantly less than those of NCEP, NGDP and CMORPH. Moreover, BMEP-d can reproduce the spatial distribution and temporal variability of global daily precipitation with higher quality. The spatial bias, spatial RMSE and spatial correlation of the merged analysis consistently improved compared to each source data.

(2) The general shapes of the zonal-mean latitudinal profiles for GPCP-1DD and BMEP-d were relatively close at the middle and low latitudes. In high-latitude areas, where the satellite data (CMORPH and TMPA) seriously underestimate the zonal-averaged precipitation, the profiles of BMEP-d exhibited reliable consistency with those of NCEP. Using the precipitation information from satellite and gauge data in low- and midlatitude areas, and the information from model and gauge data in high-latitude areas, the BMEP-d data exhibited highly consistent zonal-averaged precipitation with GPCP-1DD. The HL-OI technique clearly functioned as intended, as better source data can obtain more weights and the zonal-mean errors of BMEP-d were smaller than each input data source.

(3) The HL-OI technique also resulted in significant improvements in seasonal variability of daily precipitation at both the global and regional scales. The magnitudes and phases of the BMEP-d seasonal cycles matched those of GPCP-1DD more closely than not only the source data from NCEP, NGDP, and CMORPH, but also the TMPA multi-SEs in most of the selected monsoon regions and oceanic regions. The smallest differences between BMEP-d and GPCP-1DD can likely be attributed to the location-dependent objective error reductions in the HL-OI scheme.

(4) In terms of statistical probability indices, BMEP-d demonstrated the best performance in estimating the frequency of daily precipitation, as well as with respect to its POD, FAR, and ETS scores. The PDFs of BMEP-d were reliably consistent with GPCP-1DD over most of the daily PI thresholds. The resulting analyses of POD, FAR, and ETS scores over both the entire globe and global land (oceanic) areas indicated an improved capability (over CMORPH and NCEP) of BMEP-d to capture intense precipitation events, based on comparison with GPCP-1DD.

This study demonstrates the capacity of a multi-source strategy to provide merged global daily precipitation analysis of substantially improved quality compared to data from each input data source. Because this strategy only requires daily-scale precipitation information to achieve bias correction, error estimation, and objective merging processes, the BMEP-d data have strong potential to provide an RT global daily precipitation product in the future. The plan is to produce the final BMEP-d data for the period 2003 to the present. Currently, a preliminary product of 13 years, from 1 January 2003 to 31 December 2015, has been constructed. A current version of the data is now available from NCC/CMA, China.

Although the BMEP-d data show reliable consistency with the GPCP-1DD data, errors and inconsistencies remain, especially over mid- and high-latitude areas. Because the quality of the merged analysis is determined primarily by the accuracy of the merged source data, better source data are essential for further improving the performance of BMEP-d in future work. Because gauge observations play an important role in reducing magnitude errors of daily precipitation, efforts to merge more gauge observations worldwide are worthwhile. The Global Historical Climatology Network daily dataset (Menne et al., 2012) contains high-density daily data from over 80 000 stations to be used for this purpose in the future. Meanwhile, gauge observations from TAO/TRITON (Hayes et al., 1991) are expected to provide more accurate magnitude information over ocean areas. Furthermore, because satellite data are important for determining spatial patterns of daily precipitation, especially over oceanic areas and sparsely populated land areas with very few gauge observations, the latest high-quality satellite precipitation estimates, such as the GPM mission products (Huffman et al., 2014), are expected to be merged to improve the accuracy of the next version of BMEP-d. Additionally, accurate definitions of the error structures for each data source are another key point for improving the quality of the merged precipitation. The development of a more objective and advanced merging algorithm is also under consideration for the future.