Detecting the Causal Nexus between Particulate Matter ( PM 10) and Rainfall in the Caribbean Area

: In this study, we investigate the interactions between particulate matter that have an aerodynamic diameter less than 10 µ m diameter ( PM 10) and rainfall ( RR ) in entropy framework. Our results showed there is a bidirectional causality between PM 10 concentrations and RR values. This means that PM 10 concentrations inﬂuence RR values while RR induces the wet scavenging process. Rainfall seasonality has a signiﬁcant impact on the wet scavenging process while African dust seasonality strongly inﬂuence RR behavior. Indeed, the wet scavenging process is 5 times higher during the wet season while PM 10 impact on RR is 2.5 times higher during the ﬁrst part of the high dust season. These results revealed two types of causality: a direct causality ( RR to PM 10) and an indirect causality ( PM 10 to RR ). All these elements showed that entropy is an efﬁcient way to quantify the behavior of atmospheric processes using ground-based measurements.


Introduction
In the literature, it is well known that aerosols have an impact on climate [1][2][3]. These effects can be direct by scattering and partly absorbing solar radiation [4,5], or indirect by modifying the microphysical properties of clouds [6,7]. Classically, clouds form when air is cooled and becomes supersaturated with water or ice [8]. The excess vapor that results from this process condenses on aerosol particles that will serve as cloud condensation nuclei (CCN) or ice nuclei (IN). The characteristics of aerosols are therefore a key parameter which define the properties of clouds, precipitation and the radiative effects of clouds. Thus, the aerosol-cloud interactions (ACIs) strongly depend on cloud types which are mainly controlled by atmospheric dynamics and thermodynamics.
There are several types of clouds in the atmosphere: shallow cumuli and stratocumuli (warm clouds), mixed-phase stratiform clouds, deep convective clouds and cirrus clouds [8]. Among this whole family of clouds, the impact of aerosols on warm clouds is the least complicated in ACIs as only the liquid phase is involved [8,9]. The "Twomey" effect is the process that comes into play for warm clouds [10]. This effect, relatively well understood in the scientific community, induces a reduction in droplets size and an increase in clouds reflectance due to the increase in the number of droplets for a constant liquid water path [8]. Hence, other aerosols' indirect effects have been introduced as increased cloud lifetime and cloudiness [6] and suppressed rainfall [11], which are both led by reduced droplet size and narrower droplet spectrum. When CCN concentrations are weak, rain forms faster without necessarily involving an ice phase, even for deep convective clouds with warm bases [7]. In the tropics, these types of clouds predominate.
Particulate matter with an aerodynamic diameter 10 µm or less (PM10) is one of the main aerosols in the Caribbean basin. These PM10 can be from natural (marine aerosols, mineral dust) [12][13][14][15] or anthropogenic origin [16] . Marine aerosols mainly come from the bursting of rising air bubbles that have been injected below the sea surface by breaking waves and the pulling of droplets from wave crests [17,18]. Wind speed is a key parameter to generate particles from these two mechanisms. Indeed, the production of spume droplets from wave tops only occurs for wind speeds >9 m·s −1 [7]. Every year, a large amount of mineral dust coming from the African deserts crosses the Atlantic to reach the Caribbean during the boreal summer [19]. Once these particles are emitted by their sources, they can be lifted up to 6 km of altitude and travel in the Saharan air layer (SAL), i.e., a hot and dry dust-laden layer confined between 1 and 5 km height [20,21], at a speed of 10 m·s −1 [22,23]. Depending on dust particle size and the atmospheric conditions during their transport, dry and wet depositions are the two main mechanisms that will enable them to reach the atmospheric boundary layer [24]. In the SAL, dust particles in suspension are a mixture of fine and coarse particles [25]. After their emission, coarse particles (PM > 70 µm) should settle in less than one day [24,26] while PM10 can reach the Caribbean area after a transport of 5 to 7 days [27]. As regards anthropogenic activity, PM10 may be mainly related to transport and industrial activity [28,29].
Over the past 20 years, many authors have studied ACIs using remote sensing measurements or aircraft observations [3,[30][31][32][33][34]. According to author, no study has yet investigated ACIs in entropy framework using ground-based measurements. Consequently, the feedback between PM10 concentrations and rainfall (RR) values was assessed for the first time with this approach. Contrary to correlation analyses, the causality methods enabled us to simultaneously determine the intensity and sign of interactions [35][36][37].

Experimental Data
To study the relationship between PM10 and RR, this work uses data from Guadeloupe archipelago (16.25 • N-61.58 • W). This French overseas region located in the central Caribbean Basin experiences a tropical rainforest climate [38,39].
Due to its insular context, there are many microclimates in this small territory (∼1800 km 2 ; 390,250 inhabitants) [40]. PM10 data were performed by Gwad'Air Agency (http://www.gwadair.fr/, accessed on 1 December 2021), which manages the Guadeloupe air quality network while RR data were provided by Météo France (https://meteofrance. gp/fr, accessed on 1 December 2021). Both measurements are collected simultaneously and continuously in the insular continental regime of the island, i.e., under the same atmospheric conditions [41]. To carry out this study, eleven years of PM10/RR data, respectively from January 2005 to December 2012 then from January 2015 to December 2017, are available at an hourly basis. To ensure data quality, both time series were previously preprocessed by Gwad'Air and Météo France. In order to investigate the impact of large scale dust event on PM10/RR interactions, PM10 data were converted into daily average values while RR to daily sum values. In the literature, several studies have already highlighted the cumulative effect of RR on atmospheric processes [42,43]. Furthermore, a recent study made in Guadeloupe has shown that Pearson correlation coefficient between the daily average PM10 and the daily RR (sum then average) exhibited the same result [44]. This is the reason why the sum was preferred to the average for the stochastic analysis. A sequence of synchronous PM10/RR time series is illustrated in Figure 1. Between both variables, one can notice that the groups of peaks seem to be present during the same period.

Coherence Function
In order to analyze a potential interaction between two random processes, the secondorder statistic is frequently performed [45]. Here, the coherence function is applied for estimating the correlations in the field of frequency between two processes x(t) and y(t) . The coherence function Λ xy is computed as follows [46,47]: where E xy ( f ) is the Fourier co-spectrum while E x ( f ) and E y ( f ) are respectively the Fourier spectra of x(t) and y(t) time series. According to Λ xy sign , x and y may be uncorrelated (Λ xy 0), perfectly correlated (Λ xy 1) or negatively correlated (Λ xy −1) [47].

Convergent Cross Mapping
To investigate the causality in complex environments, the convergent cross mapping (CCM) method was introduced by Sugihara et al. [35]. This approach is based on the state space reconstruction (SSR) technique which is an advanced nonparametric technique [48]. To detect the interactions between nonlinear time series, the SSR technique became an efficient framework after the introduction of the embedding theory by Takens [49]. Hence, based on the simplex projection proposed by Sugihara and May [50], the CCM was developed and applied for nonlinear dynamic systems to study the relationship between two variables by using time-delayed embedding to reconstruct two shadow manifolds of X and Y [48].
Considering two time series X and Y of length L, which are presented as follows [48]: By forming the lagged coordinate vectors of X and Y, the shadow attractor manifold M X and M Y can be reconstructed by the following [48]: where t = 1 + (E − 1)τ to t = L, E is the embedding dimension, which represents the size of the time window, and the time lag τ is positive. There are E-dimensional points between the first lagged coordinate X(t) and the last, i.e., X(t − (E − 1)τ). From the results of Sugihara and May [50], the smallest bounding simplex is formed from the E + 1 closest neighbors, so that the simplex can contain E-dimensional points. Therefore, an essential numerical limit on the potential embedding dimension E is needed. For n observations, the analysis requires n − 1 observations to characterize the historical dynamics and at least one observation to quantify the estimated values, hence E ≤ n − 1 [51].
From the cross mapping, a library with L points from M X may be employed to supply estimates of L points for the original time series Y. Then, the ability of the L cross-mapped estimates from X to describe the Ltrue value from Y can be quantified by the Pearson correlation coefficient ρ YŶt which is defined as follows: ρ YŶt represents the predictive skills, i.e., the strength of influences from one variable to another variable, and ranges from 0 to 1. Hence, if Y has a causal impact on X, thenŶ(t)M X will converge to Y(t), andX(t)M Y should converge to X(t). In theory, the CCM correlation should increase to L until infinity but in reality, the forecast skill of the cross-map estimates from X to Y can only reach a plateau [48].

Information Transfer
Another way to detect causality in stochastic dynamical systems is the information transfer method proposed by San Liang [52]. According to San Liang [53], for a two-dimensional system as follows: w i (i = 1, 2) represents the white noises, b ij and F i are functions of X and t. Liang and Kleeman demonstrated that the rate of information flowing from X 2 to X 1 is the following according to Shannon entropy [53]: where ρ 1 represents the marginal density of X 1 while ε is the mathematical expectation (units). Equation (6) naturally incorporates the rigorous principle of causality. Thereafter, this was improved by San Liang [36] with a concise formula. Thus, for two time series X 1 and X 2 , the estimate of the maximum likelihood of the information flow rate from X 2 to X 1 is expressed as follows [36]: with C ij the sample covariance between X i and X j , C i,dj the covariance between X i andẊ j , andẊ j the difference approximation of dX j dt using the Euler forward schemeẊ j,n = X j,n+k −X j,n k∆t with k = 2 [36]. T 2→1 is in nats per unit time (e.g., nats per day) with "nat" for natural unit of information. When T 2→1 = 0, X 2 does not cause X 1 ; if not, it is causal. If T 2→1 > 0 means that X 2 functions to make X 1 more uncertain, while T 2→1 < 0 indicates that X 2 tends to stabilize X 1 [36]. A recent study has shown that T 2→1 , arbitrarily termed Liang's causality, is a good way to estimate the sign of the information flow [54].
To better quantify the strength of the causality between two variables, a new approach to normalize the information flow formula proposed by San Liang [36] is introduced by Bai et al. [37]: where ε is the mathematical expectation, ρ 1 the marginal density of X 1 and g ij = ∑ k b ik b jk [52]. According to Bai et al. [37], the stochastic effects H noise 1 is equal to: | as the normalizer. Thus, T 2→1 , defined in Equation (7), becomes [37]: τ 2→1 , arbitrarily termed Bai's causality, quantifies the importance of the information flow from X 2 to X 1 in comparison to other stochastic processes [37]. Contrary to T 2→1 , τ 2→1 will always be positive due the absolute value applied during the normalization. For more details on information transfer theory and normalization computations, refer to San Liang [36,55] and Bai et al. [37].

Coherence Function Analysis
To assess an eventual relationship between PM10 and RR data in a dynamic way, the coherence function is firstly used. Figure 2 shows the coherence function Λ PM10/RR ( f ) for a frequency range 10 −9 f 10 −5 Hz. Overall, one can notice there is a correlation between PM10 and RR in the frequency domain. However, this correlation decreases as the frequency increases. Many peaks exceed the significant correlation threshold (Λ PM10/RR ( f ) ≥ 0.70) introduced by Calif and Schmitt [47] at 9.02 × 10 −9 f 6.92 × 10 −8 Hz corresponding to 7.0 T 53.5 months (2.2 years). At 9.77 × 10 −8 f 3.55 × 10 −7 Hz corresponding to 1.4 T 4.9 months, Λ PM10/RR ( f ) is between 50% and 70%. Λ PM10/RR ( f ) goes below 40% at 3.98 × 10 −7 f 5.57 × 10 −6 Hz corresponding to 2.1 T 29.1 days. These results show that for small time-scales, the interactions between PM10 and RR seem weak.
The coherence function analysis is an efficient method to highlight the dynamic correlation between PM10 and RR. Nevertheless, as correlation does not imply causality but causation implies correlation, this method may not be suitable to exhibit causality [52]. In the correlation analysis, it lacks the needed asymmetry or directedness between dynamic events [37] and the mirage correlation might not be detected [56]. Consequently, other approaches which aim to detect the causal relationship between PM10 concentrations and RR values will be used.

Overall Analysis
In order to detect the causal nexus between PM10 concentrations and RR values, the CCM and the information transfer methods are performed. To apply the CCM frame, one must first estimate the optimal embedding dimension (E) and the time lag (τ). As E depends on many factors, i.e., system complexity, time series length L and noise [57], it is therefore crucial to compute it in the most rigorous way. This is the reason why the algorithm proposed by Wallot and Mønster [58] is used. This latter is based on the average mutual information and the false nearest neighbors to determine respectively τ and E values. For illustration purposes, Figure 3 shows the results achieved for PM10 data. In Figure 3a, one can notice that the average mutual information drops below the default threshold value at τ = 1 while in Figure 3b the percentage of false nearest neighbors remains almost constant from E = 3. The same pattern was observed for RR data. Thus, τ = 1 and E = 3 were retained to compute the CCM. As regards the information transfer frame, this approach is a data-driven procedure. Figures 4 and 5 illustrate the causality computation between PM10 and RR respectively with CCM and information transfer methods. One can notice that Figure 4 exhibits an unidirectional causality (RR to PM10), while Figure 5  For RR to PM10, both Figures show the impact of the wet scavenging process [44]. In view of the studied data , it is difficult to define which scavenging process prevails between rainout or washout [59,60]. Everything will depend on PM10 location in the atmosphere. Indeed, the level of aerosols in and under the clouds during a rainy event determines whether the two phenomena occur simultaneously or not [61]. Contrary to the CCM, the information transfer shows simultaneously the strength (τ B ) and the sign of the causality (T L ). T L RR→PM10 < 0 indicates that RR values tend to stabilize PM10 concentrations. This is in agreement with the literature because during and after a rainy episode, the PM10 levels will decrease in the atmosphere making their concentrations more homogeneous [62]. A recent study by Plocoste and Calif [54] between the PM10 and the atmospheric temperature (T) found that τ B T→PM10 = 24.0% for the same period. This value is 8 times higher than the 3.0% obtained for τ B RR→PM10 . Thus, T values have more influence on PM10 concentrations than RR values. Such a difference may be due to the fact that wet deposition can only remove 30% of aerosols from the troposphere [63,64], while dry deposition can remove up to 70% [65]. In addition, during a rainy event, not all precipitation reaches the ground, as a certain fraction of the raindrops evaporate before reaching it.
As regards PM10 to RR, Figure 5 shows that PM10 also has an influence on RR behavior. Contrary to T L RR→PM10 < 0, T L PM10→RR > 0 which means that PM10 concentrations make RR values more uncertain. This may be due to PM10 nature which has several origins, i.e., mineral dust, marine aerosols or anthropogenic activity [7]. The impact of PM10 on T seems greater than the one of PM10 on RR as τ B PM10→T = 13.8% [54] and τ B PM10→RR = 1.1%. Several studies have shown that PM10 concentrations induce a greenhouse effect which increases T values [54,[66][67][68]. Conversely, ACIs are more complex mechanisms that involve several steps [8,9]. This may justify the weak information transfer from PM10 to RR. The author assumes that this weak exchange coupled with the fact that it is necessary to determine several parameters to compute the CCM can explain the non-detection of this interaction by this approach (see Figure 4). Given the drawbacks of the CCM method, only the information transfer frame would be carried out for the remainder of the study.

Seasonal Analysis
In the literature, many studies have demonstrated the impact of aerosol properties on ACIs [7,8,69,70]. In the Caribbean area, African dust seasonality has a significant impact on PM10 concentrations and air quality [19,28,71]. The low dust season is from October to April, while the high dust season is from May to September [12,19,27,72,73]. The Caribbean also experiences a rainfall cycle [74,75]. Classically, the dry season occurs from January to March, while the wet season takes place from July to November [40]. A recent study by Plocoste et al. [44] showed that the dry season can be extended from December to June. To discuss the causal results between the African dust seasonality and the rainy seasonality, the descriptive statistics for PM10 and RR time series are presented in Table 1. Table 1. The arithmetic mean (M), standard deviation (σ) and Kurtosis (K) of PM10 and RR over 11 years, for low (October to April) and high (May to September) dust seasons, for dry (December to June) and wet (July to November) seasons, and for high dust dry (May to June) and high dust wet (July to September) seasons. N is the sample size.

Study PeriodM (µg/m 3 ) σ (µg/m 3 ) KM (mm) σ (mm) K
Overall (N = 3847) 26. As benchmarks, PM10 and RR data over 11 years are also computed. In statistical analysis, the mean, the standard deviation and the kurtosis are parameters which respectively highlight the trend, fluctuation and intermittency of the data. The more intermittent the time series, the higher their kurtosis (K) [76]. The results of Table 1 show that the seasonality is well marked with average values 1.5 and 1.7 times higher respectively for the high dust season and the wet season. Figure 6a,b shows the results achieved according to African dust seasonality. In Figure 6a,b, one can notice that τ B PM10→RR and T L PM10→RR are almost the same for both seasons. During the low dust season, PM10 aerosol is mainly composed of sea spray [70] and anthropogenic pollution [77] as dust storms are weak [15]. This is the reason why in Table 1 K PM10 High < K PM10 Low . During high dust season, mineral dust coming from the dust belt in Africa is persistently added to this set of aerosols [78]. Despite this addition, PM10 influence on RR seems to remain the same during high dust season. One can observe that the whiskers are wider for the high season. This can be attributed to the inter-annual variability of African dust storms [79]. At first glance, it seems that there is no difference in PM10 impact on RR between low and high dust seasons. Regarding RR to PM10 (Figure 6a,b), the behavior is different as τ B RR→PM10 and T L RR→PM10 are 4 times higher in the low season. During the high season, RR average is higher (see Table 1) but the wet scavenging process efficiency seems less significant. Two complementary approaches can explain the difference between both seasons. At the synoptic scale, mineral dusts are known to reduce precipitation once they are far away from their emissions sources [7,[79][80][81]. From May to September, i.e., the high dust season, the large quantities of mineral dust from Africa both act on cloud microphysics and also reduce the convection process which enables their formation due to the extinction of solar radiation [7,8,80]. During this period, which coincides with the hurricane season [82,83], the rainfall will be mainly of synoptic origin as the African easterly waves always precede and follow a dusty event [84,85]. Furthermore, during the low season there are October and November, which are among the top five rainy months during the study period. October is even the rainiest month with a total of 2303 mm over 11 years [54]. All of these can explain the fact that the wet scavenging process is more significant in the low season than in the high one. Figure 7a,b illustrates the results achieved according to the rainy seasonality. τ B PM10→RR /T L PM10→RR and τ B RR→PM10 /T L RR→PM10 are 5 times higher respectively in Figure 7a,b. Although June is known to be the dustiest month [78,86], the dry season is dominated by the low dust season asM PM10 Dry <M PM10 Wet and K PM10 Dry > K PM10 Wet (see Table 1). Thus, it is rather the marine aerosols which are predominant at this period. The rainfall is mainly of convective nature from December to June, i.e., mesoscale and local scale origin, allowing a longer residence time of PM10 in the atmosphere due to RR strong intermittency (K RR Wet < K RR Dry in Table 1). Thus, τ B PM10→RR and T L

PM10→RR
show the impact of aerosol properties on RR behavior between both seasons. As expected, the wet scavenging process, i.e., τ B RR→PM10 and T L RR→PM10 , is more significant for the wet season due to higher and more persistent rainfall (see Table 1). On the other hand, PM10 influence on RR (τ B PM10→RR /T L PM10→RR ) is the lowest during that period. The more it rains, the greater the transfer of aerosols from the atmosphere to the ground. As a result, there will be less PM10 available to impact RR because their residence time in the atmosphere will be less. Furthermore, the strong atmospheric mixing generated by hurricane activity can also increase PM10 dispersion. It is important to underline that the whiskers are wider during the wet season. As this period mainly coincides with the high dust season and the hurricane activity, the author assumes that the inter-annual variability of African dust outbreaks coupled with PM10 dispersion can explain this behavior. For the wet scavenging process (RR → PM10), this phenomenon can be explained by the great variability of rainfall events in intensity and duration. Based on the results obtained for African dust seasonality and rains seasonality, the high dust season is divided into two parts, i.e., the high dust dry (May to June) and high dust wet (July to September) seasons, in order to better assess the impact of PM10 concentrations on RR values. Figure 8a,b shows that PM10 influence on RR (τ B PM10→RR ) is 3.5 times more significant during the high dust dry season, while the wet scavenging process (τ B RR→PM10 ) is slightly more important during the high dust wet season. In Figures 8a,b, one can notice that the PM10 (RR) make RR (PM10) more fluctuating in the high dust dry season (T L > 0), while they homogenize them during the high dust wet season (T L < 0 ). The author assumes that the residence time of PM10 in the atmosphere coupled with RR intermittency may explain this behavior. Indeed, Table 1 shows that PM10 and RR statistical parameters values between both periods are of the same order of magnitude, except for K RR , which is 6 times higher during the high dust dry period. In Figure 8a, the impact of PM10 on RR is greater than in Figure 6b. Thus, these results show that during the same season, the impact of mineral dust on RR can evolve. In addition to their residence time in the atmosphere aforementioned, the chemical properties of mineral dust may also affect clouds characteristics. During the high dust season, dust storms can come from different parts of Africa [15]. Mineral dust chemical characteristics will strongly depend on the soils' nature from which they originate [87][88][89]. A recent study showed that PM10 and PM2.5 broadly follow the same temporal pattern in the Caribbean area [78]. The highest PM2.5 values correspond to the highest PM10 values. Hence, dust outbreaks will bring PM2.5 of mineral origin. It is well known that PM2.5 is a better indicator of CCN activity than PM10 [90,91]. The CCN activity is determined by the aerosol number concentrations, which determine the cloud drop number concentrations. There is a good correlation between PM2.5 and CCN when all particles are small (near 0.1 µm), but once some of them are large (>1 µm) this correlation breaks. This is due to the fact that the mass of 1000 0.1 µm particles equals the mass of one 1 µm particle. As the high dust dry season corresponds to the highest PM10 concentrations (see Table 1), PM2.5 will also be significantly present. Even if the available data are not conducive to investigate this phenomenon, the author assumes that CCN activity could be more important during this period due to the strong presence of PM2.5 related to African dust.

Conclusions
The aim of this study was to assess the interactions between PM10 and rainfall (RR) using ground-based measurements. To carry out this analysis, PM10 concentrations and RR values from Guadeloupe archipelago were used. Our results highlighted the feedback between PM10 and RR data. In other words, PM10 influence RR (PM10 to RR), while RR induce the wet scavenging process (RR to PM10). We found that rainfall seasonality has a huge impact on wet scavenging process. During the wet season, this phenomenon is 5 times more significant, making PM10 concentrations more homogeneous. African dust seasonality also acts on RR behavior. During the first part of the high dust season, PM10 influence on RR is 2.5 times higher than in the low dust season, making RR values more heterogeneous. The highest wet scavenging process and the greatest impact of PM10 on RR correspond to the maximum average of RR and PM10 respectively. According to the author, these results highlight two types of causality: a direct causality (RR to PM10) and an indirect causality (PM10 to RR). Indeed, the intensity and duration of rainfall events will have a direct impact on PM10 concentrations. On the other hand, the seasonal variation of PM10 concentrations will be an indicator of the presence of other smaller mineral particles conducive to acting as CCN.
To conclude, the results of this study showed the natural balance that exists to harmonize the climate. Thus, the entropy framework is a robust approach to quantify the behavior of atmospheric processes using ground-based measurements. Nothing is done at random to maintain the earth's balance. Indeed, if mineral dusts were predominant during the dry season, it would accentuate the drought, while if the marine aerosols were predominant during the wet season, it would accentuate the rainfall, which would lead to more landslides due to soil saturation with water. In future research, we will examine the impact of climate change on CCN activity using PM2.5 number concentration.

Funding:
The present study has no external funding.

Institutional Review Board Statement: Not applicable.
Informed Consent Statement: Not applicable.

Data Availability Statement:
The data presented are available on request from the corresponding author. The data are not publicly available due to privacy or ethical reasons.