**3. Results**

#### *3.1. Seasonal Variations of GPS PWV in AACR and LIBE*

During the study period, the HYSPLIT trajectory analyses identified that air masses reaching Costa Rica predominantly came from the southeastern Caribbean Sea with a less frequent contribution from the Pacific Ocean (Figure 3). Overall, the wind direction and speed in Costa Rica are mostly influenced by the seasonal migration of the ITCZ. Thus, during the dry season (December–April) when the ITCZ is located south of Costa Rica, air mass trajectories were associated with the influence of the NE trade winds. During the wet season (May–November), NE trade winds were weaker due to the passage of the ITCZ over Costa Rica, and cross-equatorial winds from the southern hemisphere transported moisture from the Pacific Ocean to the Central American Isthmus. This moisture transport pattern controlled the precipitation regimes observed at the Central Valley (AACR, Figure 4A) and the northern Pacific region of Costa Rica (LIBE, Figure 4A). During the study period, approx. 23% (N = 12) of the air masses arrived from the Pacific Ocean and the rest (approx. 77%, N = 40) came from the Caribbean Sea. Air masses arriving from the Pacific Ocean and the Caribbean Sea predominantly traveled over the eastern Pacific Ocean and the central and southern Caribbean Sea basins, respectively. No significant di fferences were found between the mean sea levels of the air masses that reached the Central Valley in the dry and wet seasons, with typical mean sea levels of 1500 m to 2000 m.

Seasonal GPS PWV variations were clearly defined at AACR and LIBE (Figure 5A). During the dry season, GPS PWV values varied from 14.8 mm to 40.9 mm in AACR (mean value: 27.6 ± 6.3 mm), whereas in LIBE, the variation was in the range 20.2 mm–55.5 mm (mean value: 36.9 ± 7.6 mm). We observed an increment in the GPS PWV values at both sites at the end of April and at the beginning of May that coincides with the onset of the wet season in Costa Rica, namely the passage of the ITCZ. During the wet season, at AACR, the GPS PWV estimates ranged from 24.3 mm to 46.2 (mean value: 39.7 ± 3.7 mm), and at LIBE, these GPS PWV values fluctuated from 31.5 mm to 62.6 mm (mean value: 54.1 ± 5.2 mm). At the end of November, when the transition wet-to-dry season began, we registered a decrease in the GPS PWV estimations related to the beginning of the dry season and the influence of the NE trade winds. Overall, the GPS PWV values were greater at LIBE than at AACR due to the elevation di fference between the GPS stations ( Δ1,027 m a.s.l.). For example, the mean di fferences between the estimations for AACR and LIBE were −9.5 ± 4.5 mm in the dry season and −14.4 ± 3.0 mm in the wet season. As shown in Figure 5A, these observed di fferences in the GPS PWV measurements at the GPS stations were more evident during the wet season when the ITCZ predominantly influenced the air circulation over Costa Rica. During the dry season, on the other hand, the di fferences between the GPS PWV values for AACR and LIBE were relatively more di fficult to separate. However, although the GPS stations are situated approx. 160 km from each other (one in the Central Valley, AACR, and the other one in the northern Pacific region of Costa Rica, LIBE), we found a good Spearman's correlation (r = 0.929, *p* > 0.001) between the GPS PWV values estimated for AACR and the corresponding estimations calculated for LIBE (Figure 5B). The best performing linear regression model shown in Figure 5B overall explained 86.9% of the variance for the GPS PWV estimates calculated for LIBE using the GPS PWV values at AACR. Overall, this finding confirms that the PWV variations at both sites are controlled by the climatic conditions of the Pacific slope which is also reflected in the precipitation patterns and air temperature/relative humidity variations shown in Figure 4A–C. This is an important result that demonstrates the applicability of PWV to monitor changes in the hydrometeorological conditions at regions that share similar climatic conditions. Additionally, our HYSPLIT analysis is able to identify the seasonal PWV variations at AACR. For instance, air masses arriving from the Pacific Ocean between May and October are associated with high PWV estimations with values between 39 and 44 mm/day. These values are practically equal to or greater than the 75th percentile of our data set (41 mm/day). In turn, air masses coming from the Caribbean Sea were associated with greater variations in the PWV estimations registered between November and April but also to smaller PWV estimations (up to 14 mm/day, Figure 5A).

#### *3.2. GPS PWV Comparison to Other Estimations Methods and MRL Analysis*

GPS PWV observations at the Central Valley of Costa Rica (AACR) compared well to the atmospheric sounding measurements during the dry and wet season but only to the MODIS Terra estimations during the dry season. As shown in Figure 6A, GPS PWV observations followed the seasonal variations registered using the radiosonde data. The best performing satellite-based estimations were those retrieved from the MODIS Terra, which also followed the seasonal variations in the GPS and radiosonde PWV observations. Unlike MODIS Terra, MODIS Aqua PWV estimations showed a systematic positive bias with respect the GPS PWV values and the radiosonde data. To better identify the seasonal di fferences found after applying these PWV estimation methods, we split our data set into two groups: dry season and wet season estimations (Figure 6B,C). For the dry season, the GPS PWV median value (26.5 mm) was not significantly di fferent from the radiosonde PWV median value and the MODIS Terra PWV median estimation (27.0 mm and 25.8 mm, respectively; *p* > 0.05). However, it was significantly di fferent from the median value estimated using the MODIS Aqua PWV values (29.7 mm, *p* > 0.001). In turn, for the wet season, the GPS PWV median value (40.3 mm) was significantly di fferent from the MODIS Terra and MODIS Aqua PWV median estimations (36.0 mm and 51.4 mm, respectively; *p* < 0.001) but not significantly different from the radiosonde PWV median value (41.4 mm, *p* > 0.05). The mean relative biases for MODIS Aqua PWV and MODIS Terra PWV were also calculated using the GPS PWV as a reference. During the dry season, these values corresponded to 0.16 ± 0.24 mm and 0.02 ± 0.30 mm, respectively, and were equivalent to RMSE values of 7.43 mm and 7.21 mm, in that order. During the wet season, the mean relative biases were 0.30 ± 0.24 mm and −0.06 ± 0.19 mm, respectively, corresponding to 15.2 mm and 8.05 mm, respectively.

**Figure 5.** (**A**) The time series of GPS precipitable water vapor (PWV) (mm/day) estimated for AACR (blue circles) and LIBE (red circles): The wet season period 2017 (May–November) is delimited in the graph. (**B**) The graph shows the relationship between the GPS PWV estimated for ACCR and LIBE (*p* < 0.001, N = 365). The 95% confidence limits are also shown (dashed lines).

**Figure 6.** (**A**) The time series of PWV (mm/day) estimated for AACR using GPS (blue circles), atmospheric sounding (red squares), MODIS Aqua (green inverted triangles), and MODIS Terra (blue triangles). (**B**,**C**) Box plots of the PWV (mm/day) estimated using GPS, atmospheric sounding, MODIS Aqua, and MODIS Terra for the dry season and wet season in AACR, respectively: The grey box indicates the 25th and 75th percentiles with the median in middle. The error bars indicate the minimum and maximum values. The black circles indicate outliers (1.5 times the central box).

Using the available surface meteorological data at AACR, we conducted a multiple linear regression (MLR) analysis to identify the major drivers controlling the seasonal variability of GPS PWV measurements made at the Central Valley of Costa Rica (AACR, Figure 7A). The best performing MLR model was calculated as

$$\text{GPS PWN} = 4.257(\text{T}) + 0.355(\text{RH}) - 0.0486(\text{FLUX}) - 0.257(\text{P}) - 999.125, \text{R}^2 = 0.597 \tag{7}$$

where T is the mean daily air temperature (K), RH is the mean daily relative humidity (%), FLUX is the mean daily downward solar radiation flux (W/m2), and P is the mean daily atmospheric pressure (hPa). The mean relative bias associated with the estimations of GPS PWV values using this model during the dry season was 0.10 ± 0.25 mm, whereas for the wet season the mean relative bias was −0.04 ± 0.09 mm. The corresponding RMSE estimated for the dry and wet seasons were 6.09 mm and 4.02 mm, respectively. When this model was applied to estimate the sounding PWV measurements, the mean relative bias during the dry season and wet season were 0.12 ± 0.27 mm and −0.06 ± 0.10 mm, respectively. The RMSE values calculated for the dry and wet season were 6.83 mm and 4.72 mm, respectively. As shown in Figure 7B, the correlation between the GPS PWV data at MROC and the corresponding PWV values estimated from the MLR model was better for the values between 30 and 45 mm. As these values were mostly registered during the wet season, it seems that our MLR model performs better when the atmospheric conditions in the Central Valley are controlled by the seasonal migration of the ITCZ and worse during the less stable atmospheric conditions linked to the influence of NE trade winds.

**Figure 7.** (**A**) The simulated PWV (mm/day) time series (blue triangles) in relation to the estimated GPS PWV (mm/day) at AACR (red circles). (**B**) The GPS PWV (mm/day) vs. simulated PWV at AACR shows the goodness-of-fit of the multiple linear regression (MLR) model.
