A new rain detection method to complement high-resolution global precipitation products

The Global Satellite Mapping of Precipitation (GSMaP) project provides semi-global and hourly rain rate estimates based on multi-satellite measurements. However, the quality of GSMaP estimates degrades when no microwave radiometers are available or when a cloud moving vector is applied. This study proposes a new rain detection method for the interval between microwave measurements based on variations in surface soil moisture. The new method is tested with Tropical Rainfall Measuring Mission Microwave Imager data and evaluated with rain gauge measurements from the Oklahoma Mesonet. Generally, 10 GHz with horizontal polarization (10H) is the optimal channel to detect temporal increases in surface soil moisture induced by rainfall. Compared with GSMaP, the performance of hourly rain detection in the new method is not high. However, the new method is expected to perform better if measurement intervals are shortened by use of multiple microwave radiometers. The new method is good for less vegetated areas where GSMaP performance is reduced. Currently, the new method is mostly independent of GSMaP as 10H is not used for rain retrieval over land. Incorporation of the new method into GSMaP is expected to improve the latter.


INTRODUCTION
With recent progress in satellite remote sensing of precipitation, such as the successful long-term observations of the Tropical Rainfall Measuring Mission (TRMM), global precipitation products with high spatial and temporal resolutions have become available and have been used in various hydro-meteorological applications (e.g., Andermann et al., 2011;Ferguson et al., 2011;Su et al., 2011). The Global Satellite Mapping of Precipitation (GSMaP), one of the leading products, provides hourly rain rate estimates with a resolution of 0.1 degrees longitude by 0.1 degrees latitude for the entire world, excluding polar areas outside 60 degrees north and south (Okamoto et al., 2005;Kubota et al., 2007). The near real-time version of GSMaP (called JAXA Global Rainfall Watch) is published with a latency of less than four hours (Kachi et al., 2011). The reanalysis version of GSMaP is processed with additional measurements and by using the latest algorithm (Version 5.222) and is available for nearly 11 years (from March 2000to December 2010. Major inputs to GSMaP are observations from microwave radiometers including the TRMM Microwave Imager (TMI), the Advanced Microwave Scanning Radiometer for the Earth Observing System (AMSR-E), the Special Sensor Microwave Imager (SSMI), the Special Sensor Microwave Imager and Sounder (SSMIS), and the Advanced Microwave Sounding Unit (AMSU). Over land, scattering signals measured at higher frequency channels (mainly between 85 and 90 GHz) are used for rain retrieval. The rain area is detected by the method of Seto et al. (2005Seto et al. ( , 2008 and the rain rates are quantitatively estimated as described by Aonashi et al. (2009).
However, microwave radiometers do not cover the entire globe in one hour. To fill the gap in the microwave observations, cloud moving vectors are derived from infrared imagery of geostationary meteorological satellites. The horizontal movement of rain systems is then simulated with the cloud moving vector and the Kalman filter technique (called the MVK method; Ushio et al., 2009). The MVK method also simulates decaying and growing of existing rain systems by checking the change in brightness temperature of infrared imagery. However, if a rain system is not detected by any microwave sensors, MVK cannot simulate the rain system. In other products such as the TRMM Multisatellite Precipitation Analysis (Huffman et al., 2007), rain rates are related to cloud-top temperatures measured by infrared imagery. However, there are sometimes no strong relations between surface rain rate and cloud-top temperature.
In this study, a new method is proposed to detect rainfall for the interval of microwave measurements. The proposed method is based on the principle that rainfall temporally increases surface soil moisture. Conveniently, surface soil moisture can be estimated from lower frequency channels of microwave imagers such as TMI, AMSR-E, and SSMI (e.g., Koike et al., 2004;Dorigo et al., 2010). For simplicity, if we neglect the effects of vegetation and the atmosphere, brightness temperature can be approximated from the product of surface physical temperature and surface emissivity, the latter of which decreases when surface soil moisture increases. Therefore, by examining the temporal change in measured brightness temperature at a lower frequency channel, rainfall in the measurement intervals can be detected. However, in reality, temporal variations in vegetation, the atmosphere, and the surface physical temperature may hamper the signal of surface soil moisture. In this paper, the proposed method is tested with TMI measurements over Oklahoma for one year and rain detection performance is evaluated. Crow et al. (2009) developed a method to correct one by one degree and three-day rainfall accumulation by assimilating surface soil moisture estimates derived from AMSR-E. Their method is based on a similar principle as our method. However, our new method focuses on rain detection and is applied at much higher spatial and temporal resolutions.

DATA
The test region was between 34 and 36 degrees north and between 100 and 94 degrees west. This region largely belongs to the state of Oklahoma in the United States. The vegetation cover is lower in the western compared to the eastern part of the region according to the annual average of the normalized differential vegetation index (NDVI; Figure 1). Here, the Moderate Resolution Imaging Spectroradiometer (MODIS) Terra product (MOD13C2) was used for the year 2006. The region was divided into a 0.1 by 0.1 degrees latitude by longitude grid for the analysis. The test period was the whole year of 2006.
A matched-up dataset for Precipitation Radar (PR) and TMI was used. This dataset is basically the same as that used by Seto et al. (2005). The PR is on the same platform as the TMI and measures the middle band (220 km width; referred to as narrow swath) of the TMI swath (760 km width; referred to as wide swath). The matched-up dataset, rather than a TMI-only dataset, was required for later analysis. When the TRMM visits this region, it flies mostly from west to east along the 35 degrees north line. Hence, the region is well covered by the narrow swath. Even if we use the wide swath data, the number of samples does not increase substantially in this region. For each orbit, the brightness temperature was converted into 0.1-degree grid data. If multiple TMI pixels existed in a grid, a simple average of the brightness temperature was taken. Usually, the measurement interval was around 24 hours. In some cases, the interval was around 1.5 hours where two successive orbits overlap. In some other cases, the measurement interval was around 48 hours or longer.
For the validation, rain gauge data of the Oklahoma Mesonet (Brock et al., 1995;McPherson et al., 2007) was used. Data from 68 stations were available for this study. The temporal resolution was five minutes and the measurement unit of tipping buckets was 0.254 mm (0.01 inch). For each station, the data was converted to hourly and 0.1-degree grid data by including rain gauges in surrounding five by five grids using the square inverse distance method. For comparison with the proposed method, the reanalysis version of GSMaP was also evaluated. In addition to hourly rain rate estimates by GSMaP, a measurement time (0.01 hours) was given if a microwave radiometer was used for the current estimate. If no microwave radiometers were used for the current estimate, the time of the temporally nearest microwave observation was given.

METHOD
The TMI has nine channels: 10 GHz, 19 GHz, 37 GHz, and 85 GHz with both vertical and horizontal polarization and 21 GHz with vertical polarization only. Hereafter, the brightness temperature is denoted by the frequency and polarization (H for horizontal and V for vertical) of the channel. For example, 10H stands for the brightness temperature of the 10 GHz horizontal polarization channel. For each frequency with dual polarization, the ratio of the brightness temperature at the horizontal polarization to the brightness temperature at the vertical polarization was calculated. This value multiplied by 300 K is defined as the polarization ratio. For example, the polarization ratio at 10 GHz (denoted by 10P) is given by 10P = 10H/10V × 300 K. This definition makes it easy to handle the polarization ratio in the same way as the brightness temperature in the proposed method. By taking the ratio of dual polarization, the effects of the variation in surface physical temperature were expected to be reduced.
With each of 13 indices (10V, 19V, 21V, 37V, 85V, 10H, 19H, 37H, 85H, 10P, 19P, 37P, 85P), rainfall detection for TMI measurement intervals was tested. The change of index was calculated by subtracting the value at the first measurement from the value at the second measurement. When the change of index was smaller than a given threshold, rainfall was judged to occur in the interval of the two measurements (referred to as "detected"). Thresholds between −30 K and 20 K with a step of one K were tested for each index. Therefore, the number of tests was 13 × 51 = 663.
For evaluation, rain gauge data for the hours which are bounded by hours with TMI measurements were used. If the first measurement was taken at 00:30 and the second measurement at 23:15 on the same day, the rain amount was estimated between 01:00 and 23:00. If the rain amount is non-zero, rainfall is regarded as existing. The performance of rain detection was evaluated with "equitable threat score (ETS) weighted by the rain amount" (rETS) defined as follows. (1) (3) where H is the number of "hits" (rain exists and is detected), M is the number of "misses" (rain exists but is not detected), F is the number of "false alarms" (rain does not exist but is falsely detected), N is the number of "nones" (rain does not exist and is not detected). R(H) and R(M) are the sum of rain amounts for the cases of hits and misses, respectively. If H and M are used instead of H r and M r , respectively, rETS becomes ETS. rETS is used instead of ETS to evaluate the rain detection performance in terms of the rain amount rather than the number of cases. Like ETS, rETS can range from −1/3 to 1, with larger rETS meaning better performance.

EVALUATION
With the same index and threshold, rain detection was tested at 68 stations and rETS was calculated by summing data from all stations. Figure 2a shows the dependence of rETS on the index and threshold. rETS was largest when 10H was used with a threshold of −5 K. Generally, better performance was obtained with lower frequency and horizontal polarization channels. This result corresponds to the fact that the sensitivity of surface emissivity to surface soil moisture is higher for lower frequency and horizontal polarization channels. The polarization ratio does not give better results than the brightness temperature at horizontal polarization. This suggests that the variation in surface physical temperature is not large enough to disturb the signal of surface soil moisture change.
Some positive performance was shown for 85V and 85H, although they are less sensitive to surface soil moisture. It is possible that the changes in 85V and 85H are affected by rainfall at the time of TMI observations. If rain exists at the time of the second observation, the changes in 85V and 85H can be negative, and "rain is detected" can be returned by the proposed method. In this case, it could be expected that rain actually exists in the interval, as rain events are likely to last for several hours. Thus, unexpectedly, 85V and 85H showed positive rain detection performance. Figure 2b shows results by excluding cases when the matched-up PR detects rainfall at the time of one or both of the two TMI observations. With 85V and 85H, as well as with the other indices, rETS was closer to zero when the threshold was near the low and high edges. This suggests that the effects of rainfall at the time of TMI observations were successfully excluded, and the rain detection relies on the change in surface soil moisture change in Figure 2b. For lower frequency channels, rETS became slightly higher after exclusion of TMI observations affected by rainfall. This means that these channels are slightly affected by rainfall. For later analysis, all available TMI observations were used.
The evaluation is shown for each station. We used 10H with a threshold of −5 K, and rETS for each station is shown in Figure 3a. In the western part of this region, which is less vegetated, rETS was higher. However in the eastern  part of this region, which has much more vegetation, rETS was lower. This result corresponds to the fact that variation in surface soil moisture is difficult to detect when the surface is covered by vegetation. Next, the threshold of 10H was optimized for the best rETS at each station. The optimal threshold of 10H and the best rETS are shown in Figure 3b. In the southwestern part of the region, with the optimal threshold (which is smaller than −5 K) the best rETS was close to 1.0. In the eastern part, even with the optimal threshold, rETS was not substantially improved. Next, both the index and threshold were optimized for each station. The optimal index and the best rETS are shown in Figure  3c. At most stations, 10H was selected, although at some stations in the western part 10P was selected. Where vegetation is sparse, the variation in surface physical temperature is large and the polarization ratio may be effective. At some other stations, 19P and 10V were selected, but they did not increase rETS significantly.

COMPARISON WITH GSMAP
The rain detection of GSMaP was evaluated in the same way as in the previous section, but the evaluation was performed every hour as GSMaP has hourly resolution. The results for all stations are summarized below. rETS was 0.357 averaged over the whole study, but increased to 0.408 for hours with microwave measurements (45% of all times), and decreased to 0.323 for hours without microwave measurements (55% of all times). This result confirms that the accuracy of GSMaP degrades when no microwave measurements are available.
To allow for comparison with GSMaP, the proposed method was also evaluated every hour. As the proposed method does not have hourly resolution, if rain is (is not) detected for an interval of TMI measurements, we assumed that rain is (is not) detected for all hours in the interval. rETS for the hourly evaluation (Figure 4a) was much lower than rETS shown in Figure 2a, and also much lower than that of the GSMaP without microwave measurements. However, we can expect improvements in the proposed method. As only TMI is used in the present test, the interval is generally as long as 24 hours, but it is sometimes two hours. Figure 4b shows hourly evaluation only for the cases with intervals of two hours. rETS increased up to nearly 0.3 when 10H was used, and the proposed method was comparable to the GSMaP without microwave measurements. This result suggests that rETS is higher when the interval is shorter. Relatively high rETS was shown for 85V and 85H, but this is probably because scattering signal was detected as explained before. If multiple microwave sensors are used with proper calibration of brightness temperatures, the proposed method is expected to show higher rain detection performance.
GSMaP and the proposed method were evaluated for each station and the results compared in Figure 5. GSMaP was evaluated for hours without microwave measurements. For the proposed method, the index and threshold were optimized for each station. Except for two stations, GSMaP showed better performance than the proposed method. In Figure 5, different symbols are used for lower NDVI (NDVI < 0.45), medium NDVI (0.45 < NDVI < 0.6), and higher NDVI (0.6 < NDVI) areas. As explained before, the proposed method tends to give better performance in less vegetated areas. On the other hand, GSMaP shows slightly worse Figure 5. Comparison of GSMaP and the proposed method in terms of hourly rain detection. The horizontal axis is rETS of GSMaP when no microwave radiometers were available. The vertical axis is rETS of the proposed method with the optimal index and threshold. Different symbols are used for different categories of annual average NDVI. performance in less vegetated areas where rain detection by microwave sensors is degraded because of the larger variation in surface physical temperature and significant effects of soil particle scattering. In a sense, GSMaP and the proposed method are complementary in terms of the rain detection performance. Therefore, we expect that GSMaP could be improved by incorporating the proposed method.

SUMMARY
This paper describes a new rain detection method based on the variation in surface soil moisture. The proposed method was tested with TMI data and evaluated with rain gauge measurements of the Oklahoma Mesonet. Generally, rain can be taken to occur when 10H decreases by more than five K. The method was more effective in less vegetated areas. These results reflect the fact that lower frequency and sparser vegetation increase the sensitivity to surface soil moisture. The optimal index and threshold varied with location. For less vegetated areas, a lower threshold can be selected and 10P is sometimes better than 10H for reducing the effects of variation in surface physical temperature.
In terms of the performance of hourly rain detection, the proposed method is generally not as good as GSMaP with the MVK method as shown in Figure 5. This is because the interval of TMI observations is usually around 24 hours or longer, and the proposed method does not have hourly resolution. The use of multiple microwave imagers with proper calibration can shorten the measurement interval and better performance could be expected. GSMaP performs better in more vegetated areas, whereas the proposed method performs better in less vegetated areas. GSMaP and the proposed method use mostly independent input data over land, with the former using high frequency channels and the latter using lower frequency channels. Therefore, we can expect to improve the performance of GSMaP by incorporating the proposed method.