Next Article in Journal
Correlations of Stormwater Runoff and Quality: Urban Pavement and Property Value by Land Use at the Parcel Level in a Small Sized American City
Next Article in Special Issue
Spatial and Temporal Investigation of Dew Potential based on Long-Term Model Simulations in Iran
Previous Article in Journal
Assessment of the Effects of Wastewater Treatment Plant Modernization by Means of the Field Olfactometry Method
Previous Article in Special Issue
Characteristics and Drivers of Reference Evapotranspiration in Hilly Regions in Southern China
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Space–Time Kriging of Precipitation: Modeling the Large-Scale Variation with Model GAMLSS

by
Elias Silva de Medeiros
1,*,
Renato Ribeiro de Lima
2,
Ricardo Alves de Olinda
3,
Leydson G. Dantas
4 and
Carlos Antonio Costa dos Santos
4,5
1
FACET, Federal University of Grande Dourados, Dourados 79.825-070, Brazil
2
Department of Statistics, Federal University of Lavras, Lavras 37.200-000, Brazil
3
Department of Statistics, Paraíba State University, Campina Grande-PB 58.429-500, Brazil
4
Academic Unity of Atmospheric Sciences, Center for Technology and Natural Resources, Federal University of Campina Grande-PB, Campina Grande 58.109-970, Brazil
5
Daugherty Water for Food Global Institute, University of Nebraska, Lincoln, NE 68588-6203, USA
*
Author to whom correspondence should be addressed.
Water 2019, 11(11), 2368; https://doi.org/10.3390/w11112368
Submission received: 16 October 2019 / Revised: 1 November 2019 / Accepted: 7 November 2019 / Published: 12 November 2019
(This article belongs to the Special Issue Assessment of Spatial and Temporal Variability of Water Resources)

Abstract

:
Knowing the dynamics of spatial–temporal precipitation distribution is of vital significance for the management of water resources, in highlight, in the northeast region of Brazil (NEB). Several models of large-scale precipitation variability are based on the normal distribution, not taking into consideration the excess of null observations that are prevalent in the daily or even monthly precipitation information of the region under study. This research proposes a novel way of modeling the trend component by using an inflated gamma distribution of zeros. The residuals of this regression are generally space–time dependent and have been modeled by a space–time covariance function. The findings show that the new techniques have provided reliable and precise precipitation estimates, exceeding the techniques used previously. The modeling provided estimates of precipitation in nonsampled locations and unobserved periods, thus serving as a tool to assist the government in improving water management, anticipating society’s needs and preventing water crises.

1. Introduction

An effort has emerged in several areas of knowledge, especially in the climate sciences, to analyze a phenomenon by both space and time [1,2,3]. It is already established that the spatial–temporal behavior of rainfall is of great importance for the water resources of a region and that it has a direct influence on human activities, such as agriculture and commerce [4].
The northeast region of Brazil (NEB) is characterized by irregular distributions of rainfall, resulting in high spatial and temporal variability of precipitation [4]. In this region, it is common to observe high rainfall indices in a particular location and no record of rain in its surroundings [5]. The state of Paraíba, located in this region, presents the same characteristics, with periods of drought in rainy seasons and behavior different from the precipitation between the mesoregions that constitute the region [6]. Consequently, Paraíba presents a problem concerning the spatiotemporal variability of rainfall. The high variability of recorded rainfall in this region may be attributed to the characteristics of the meteorological systems and their dates of installation in different regions of the state. The periodicity is related to the general behavior of sea surface temperature (SST) in the tropical Pacific region, contributing to the events of El Niño and La Niña, and most notably, by the behavior of the SST in the tropical Atlantic region due to its proximity to the study area. In the literature, several studies have evaluated precipitation in this region using spatial statistics techniques [7,8]. However, these previous works did not take spatial and temporal dependence into account simultaneously.
Several previous works used multiple linear regression to model the deterministic component, such as the Gaussian model [9,10,11]. These works assumed that the variable response to be studied presented a behavior that could be correctly modeled by a Gaussian approach.
However, recent studies have evaluated non-Gaussian models, including generalized linear models (GLM), with variable response Gamma [12,13]. Once the Gamma model is assumed, the variable domain has to be strictly positive, and it does not include the value zero. For example, Menezes et al. (2016) [12] and Monteiro et al. (2017) [13] analyzed the concentration of NO2, which had asymmetric behavior with excesses of zero. However, when using GLM with Gamma distribution, it is not possible to model the frequency of zeroes in the variable under study. To circumvent this limitation, some models captured the frequency of zeroes in continuous variables. Stauffer et al. (2017) [3] analyzed the spatial–temporal distribution of daily precipitation at 117 stations located in Tirol, Austria, over 42 years. These authors pointed to the high frequency of zero observations in the data, suggesting the use of the normal zero-censored model. However, only the trend component was modeled, not taking into account the spatial–temporal dependence of precipitation in this region.
In this study, we analyzed the total monthly rainfall measured in the State of Paraíba, Brazil, which presents an asymmetrical behavior and has a high incidence of zeros, which is a climate characteristic of the region under study. Thus, to model this behavior, a new approach in the regression adjustment is proposed by applying generalized additive models for location, scale, and shape (GAMLSS) [14]. These models have several advantages, such as all distribution parameters can be modeled as covariate functions, cases of over-dispersion and excess of zeroes in the data can be modeled, and any distribution to the response variable can be adjusted.

2. Materials

The study region was Paraíba State, Brazil (Figure 1), located in the NEB, which is positioned between parallels 6° and 8° S and meridians 34° and 39° W, with an area of approximately 56,500 km2. Paraíba is included in the Tropical Region, and its area is divided into four geographical mesoregions: Zona da Mata, Agreste, Borborema, and Sertão.
The total monthly precipitation data measured at 54 pluviometric stations were considered during the period from 1994 to 2014. Table 1 shows the identification number, name of pluviometric station, latitude, longitude, altitude above sea level, and quantiles of precipitation (mm) at 2.5%, 50%, and 97.5%.

3. Methods

In space−time geostatistics, observations are modeled as a stochastic process using a random function { Z ( s , t ) : ( s , t ) D R d × R } . We assumed that Z has the first and second moments. Thus, Z(s, t) was decomposed by the sum of the trend components and the stochastic residue:
Z ( s , t ) = m ( s , t ) + ε ( s , t ) .
In Equation (1), the trend component m(s, t) = E[Z(s, t)] is not constant in space and time, representing the large-scale variation. The residual component ε(s, t) describes random fluctuations on a small scale and encompasses the three components: spatial, temporal, and interaction.

3.1. Trend Component

The trend component is responsible for the large-scale variation of the stochastic process. In spatiotemporal datasets, this component is not a deterministic function. It was, therefore, necessary to specify stochastic patterns that represented the observed variability. The specification of the trend component included the geographic coordinates as covariables, as well as time. The purpose of this study was to adjust this component using GAMLSS models. In the GAMLSS models, it is assumed that the observations z i are independent, with i = 1 , 2 , , n and conditioned to θ k = ( μ , σ , ν , τ ) , with k = 1 , 2 , ,   n . So, Z | θ k G ( μ , σ , ν , τ ) , where G represents Z distribution. The parameters μ , σ , ν , and τ are classified as the parameters of the location, scale, asymmetry, and kurtosis, respectively, where the latter two represent the distribution shape [14]. Thus, according to these authors, the GAMLSS models were expressed in the form,
g k ( θ k ) = η k = X k β k + j = 1 J k Y j k γ j k ,
where g k ( ) is a monotonic link function, and θ k is related to explanatory variables and random effects; θ k and η k are vectors of size n; β k is a parameter vector of size J k ; X k and Y j k are the fixed (covariate) matrices of size n × J k and n × q i j , respectively; and γ j k is a random variable of dimension q i j . In a purely parametric linear model, we have that j = 1 J k Y j k γ j k = 0 , that is, there are no additive terms associated with the distribution parameters in the model. In the space–time case, we have η ( s , t ) = g ( m ( s , t ) ) .
In the modeling of the trend, the component was considered Y j k = I n , where I n is the identity matrix of order n × n and γ j k = f j k = f j k ( x j k ) . Thus, Equation (2) is rewritten as:
g k ( θ k ) = η k = X k β k + j = 1 J k f j k ( x j k ) .
The model given in Equation (3) was fitted to the data of this study, where f j k is an unknown function of the explanatory variable x j k and f j k ( x j k ) is a vector that evaluates the function f j k in x j k . This model is called the semiparametric additive GAMLSS and may contain parametric, nonparametric and random effects terms [14].
We considered the geographic coordinates (latitude and longitude) and the temporal index (time) as covariables to model the trend. This time index was constructed as follows: an index equal to 1 denoted Jan/1994, the index equal to 2 denoted Feb./1994, and so on. In all covariates, cubic splines (cs) were inserted, including the covariate time that was designed to model the seasonal effect of precipitation. Thus, the terms of the trend component were expressed as:
g 1 ( μ ) = ln ( μ ) = β 10 + β 11 x 11 + β 12 x 12 + β 13 x 13 g 2 ( σ ) = β 20 g 3 ( ν ) = ln ( ν 1 ν ) = β 30 x 13 ,
where x 11 , x 12 , and x 13 represent the covariates of longitude, latitude, and time, respectively (Equation (4)). The parameters β 10 ,   β 11 ,   β 12 , and β 13 are related to the lease vector of the GAMLSS model, β 20 is the parameter related to the intercept for the scale vector, and β 30 is the parameter associated with covariable time and that is responsible for modeling the occurrence of no zeroes in the response variable.
In order to obtain the parameter estimates for this regression model, the penalized maximum likelihood method was used, which differs from conventional methods because it assigns the value zero to the coefficients of nonsignificant variables. However, this method does not take into account the assumption that errors are uncorrelated. Thus, it was necessary to evaluate the residuals of this model, taking into account the correlation structure. For this, we proposed to fit a valid nonseparable space–time variogram model.

3.2. Spatiotemporal Variogram

When inserting structures in the trend component, one advantage is that it describes the large-scale variation of the phenomenon being studied, causing the residuals to present small-scale variations, consequently reducing the prediction error in the kriging. After obtaining the estimates of the trend component parameters, the stochastic residue is obtained by difference, ε ^ ( s , t ) = Z ( s , t ) m ^ ( s , t ) . Therefore, the next step of the analysis was the empirical estimation of the covariance function or spatiotemporal variogram to model this residue [2].
Spatiotemporal covariates are usually described using a space–time variogram ( γ s t ) , which measures the mean difference between separated data in the space–time domain using the distance vector. The estimate of this variogram, obtained by the moment’s method, was given by:
γ ^ s t ( r s , r t ) = 1 2 | L ( r s , r t ) | L ( r s , r t ) [ ε ^ ( s + r s , t + r t ) ε ^ ( s , t ) ] 2 ,
where | L ( r s , r t ) | is the cardinality for the set L ( r s , r t ) . It is noted in Equation (5) that the empirical variogram depends only on spatial and temporal distances, and this may not result in a valid variogram. Thus, to obtain estimates of γ ^ s t in any arbitrary lag, the empirical variogram must be smoothed. As a solution to this limitation found in the empirical variogram, one should use a theoretical model, γ s t ( h s , h t , ω ) , which fits, as best as possible, the space–time dependence structure of the residues [15]. The vector ω contains all the unknown parameters to be estimated.
The typical approach is to fit a theoretical model to the empirical variogram. The least-squares estimation method is the most common approach. However, there are other methods, such as maximum likelihood, generalized estimation equations, and composite likelihood, among others [16]. In this study, the least-squares method was used.
There are no assumptions about the probability distribution of the empirical variogram’s generated values in the least-squares method. The approach considers the statistical model of the form γ ^ ( h s , h t ) = γ ( h s , h t ) + ε ( h s , h t ) . Assume that ε ( h s , h t ) has null mean and covariance matrix R = R ( ω ) , which depends on ω . We suggested that the diagonal of R ( ω ) can be approximated to:
V a r [ γ ^ ( h s i , h t j ) ] 2 γ ( h s i , h t j , ω ) 2 | N ( h s i , h t j ) | ,
where N ( h s i , h t j ) is the number of pairs at each space–time distance. Thus, the weights assigned in the weighting matrix are proportional to the number of pairs in a given space–time distance.
The adjustment of the theoretical space–time variogram model was made by using the weighted least-squares (WLS) method, which replaces R ( ω ) with the diagonal matrix W ( ω ) . The elements of W ( ω ) were obtained from Equation (6). Thus, using the weights given in this Equation (7), the WLS method is given by:
( γ ^ ( h s , h t ) γ ( h s , h t , ω ) ) W ( ω ) 1 ( γ ^ ( h s , h t ) γ ( h s , h t , ω ) ) = i = 1 I j = 1 J | N ( h s i , h t j ) | γ ( h s i , h t j , ω ) 2 { γ ^ ( h s i , h t j ) γ ( h s i , h t j , ω ) } 2 .
Generally, statistical assumptions are considered, which comprise a space–time covariance function and ensure that this function is defined as positive [17,18,19]. Several researchers have studied the construction of valid space–time covariance functions [15,20,21]. Although this type of construction is a difficult task, the complexity becomes more significant when it is desired to determine valid models of space–time covariance. In the space–time approach, there are two types of general model classes: the separable and the nonseparable.
The separable models assume that spatial and temporal processes are uncorrelated and that variograms are composed of purely spatial and purely temporal models. As an example, the space–time covariance function can be expressed as the product between purely spatial and purely temporal components. The major drawback of this type of model is the nonincorporation of the spatial–temporal interaction component [18]. However, the nonseparable models assume that the space–time phenomenon is correlated in space–time [22]. We used the generalized sum–product model, which belongs to the nonseparable covariance model class.
The stationary covariance model product–sum [17,23] is given by C ( h s , h t ) = k 1 C s ( h s ) C t ( h t ) + k 2 C s ( h s ) + k 3 C t ( h t ) . Assuming that ε ( s , t ) is second-order stationary, we have that: γ s ( h s ) =   C s ( 0 ) C s ( h s ) C s ( h s ) = C s ( 0 ) γ s ( h s ) and γ t ( h t ) = C t ( 0 ) C t ( h t ) C t ( h t ) = C t ( 0 ) γ t ( h t ) , in which, based on the assumption of second-order stationarity, we have that the product–sum covariance function is defined as:
C ( h s , h t ) = γ s ( h s ) [ k 1 C t ( 0 ) + k 2 ] γ t ( h t ) [ k 1 C t ( 0 ) + k 3 ] + k 1 γ s ( h s ) γ t ( h t ) + C ( 0 , 0 ) ,
where k 1 > 0 , k 2 0 e are constants that ensure that the covariance function is positive definite (Equation (8)). C s ( 0 ) e C t ( 0 ) are the spatial and temporal thresholds, respectively. Applying the relation, γ ( h s , h t ) = C ( 0 , 0 ) C ( h s , h t ) , we have the product–sum variogram. This model gives origin to the following relations:
γ ( h s , 0 ) = γ s ( h s ) ( k 2 + k 1 C t ( 0 ) ) = k s γ s ( h s )   and γ ( 0 , h t ) = γ t ( h t ) ( k 3 + k 1 C s ( 0 ) ) = k t γ t ( h t ) ,
where k s and k t are the proportionality coefficients between the space–time variograms, γ s ( h s , 0 ) and γ t ( 0 , h t ) , and the spatial and temporal variograms, γ s ( h s ) and γ t ( h t ) , respectively. Therefore, modeling γ ( h s , 0 ) is equivalent to modeling γ s ( h s ) . Likewise, estimating and modeling γ ( 0 , h t ) is equal to estimating and modeling γ t ( h t ) . The two relations established in Equation (9) can be simplified by imposing three constraints: k 2 + k 1 C t ( 0 ) = 1 , k 3 + k 1 C s ( 0 ) = 1 , and k 1 + k 2 + k 3 = 1
These constraints facilitate the estimation of the model parameters γ ( h s , 0 ) and γ ( 0 , h t ) using γ s ( h s ) and γ t ( h t ) , for this it is necessary to determine k 1 , k 2 , and k 3 . These three coefficients of the model can be written in terms of the thresholds and parameters k s and k t [17]. In this particular case, this leads to modeling γ s ( h s ) and γ t ( h t ) , using the models of γ ( h s , 0 ) and γ ( 0 , h t ) , respectively. In addition, these two parameters were combined into a single parameter k, resulting in the generalized product–sum model (Equation (10)).
By combining the parameters k s and k t to form a single parameter k, the space–time variogram expression, γ ( h s , h t ) , was simplified to,
γ ( h s , h t ) = ( k 2 + k 1 C t ( 0 ) ) γ s ( h s ) + ( k 3 + k 1 C s ( 0 ) ) γ t ( h t ) k 1 γ s ( h s ) γ t ( h t ) = k s γ s ( h s ) + k t γ t ( h t ) k 1 γ s ( h s ) γ t ( h t ) = γ s ( h s , 0 ) + γ t ( 0 , h t ) k 1 γ s ( h s ) k s γ t ( h t ) k t = γ s ( h s , 0 ) + γ t ( 0 , h t ) k γ s ( h s ) γ t ( h t ) ,
where k = k 1 k s k t . In Equation (11), the parameter k is expressed as:
k = k t C t ( 0 ) + k s C s ( 0 ) C ( 0 , 0 ) C s ( 0 ) C t ( 0 ) k s k t = k t C t ( 0 ) + k s C s ( 0 ) C ( 0 , 0 ) [ k s C s ( 0 ) ] [ k t C t ( 0 ) ]
We used the covariance generalized product–sum model, and the variogram function of this model is expressed by:
γ s t ( h s , h t ) = ( k C t ( 0 ) + 1 ) γ s t ( h s , 0 ) + ( k C s ( 0 ) + 1 ) γ s t ( 0 , h t ) k γ s t ( h s , 0 ) γ s t ( 0 , h t ) ,
where γ s t ( h s , 0 ) and γ s t ( 0 , h t ) are the marginal spatial and temporal variograms, respectively (Equation (12)). The great advantage of this type of model is that sample-based marginal adjustments are used, and only one global threshold parameter is incorporated for space–time interaction [24]. In Equation (13) the parameter k is positive and has an identity involving the global threshold, C s t ( 0 , 0 ) together with the spatial and temporal threshold, C s ( 0 ) and C t ( 0 ) , respectively, given by:
k = C s t ( 0 , 0 ) C s ( 0 ) C t ( 0 ) C s ( 0 ) C t ( 0 ) .

3.3. Geostatistical Prediction

Once the trend component was determined and estimated and a valid space–time variogram model fitted, the next step consisted of the space–time kriging of the residuals using the ordinary kriging method. The variogram model is crucial in space–time kriging to calculate the best invariant linear predictor. Thus, it was necessary to obtain the space–time kriging equations, which provide the weights of the observations in such a predictor and the prediction variation that indicates the prediction accuracy [1].
Ordinary space–time kriging [15] consists of predicting the value ε ( s 0 , t 0 ) , at the nonsampled space–time point ( s 0 , t 0 ) , from the stochastic component ε ( s , t 0 ) . For this, the linear predictor, ε * ( s 0 , t 0 ) = i = 1 n λ i ε ( s i , t i ) , is used, where λ i are the kriging weights obtained by imposing the condition that the prediction error expectation is zero and that it has minimum variance; that is, it is the best linear unbiased predictor (BLUP). To obtain the weights, we used the ordinary space–time kriging equations that were calculated in terms of the variogram:
{ j = 1 n λ j γ s t ( s i s j , t i t j ) + α = γ s t ( s i s 0 , t i t 0 ) , i = 1 , 2 , , n i = 1 n λ i = 1
In the system of Equations (14), α is the Lagrange multiplier required to minimize the prediction variance. The prediction variance (Equation (15)) is given by:
V ( Z * ( s 0 , t 0 ) Z ( s 0 , t 0 ) ) = i = 1 n λ i C ( s i s 0 , t i t 0 ) + α
in which the final predictor ( z ^ ( s 0 , t 0 ) ) for the precipitation variable ( Z ), at the location ( s 0 , t 0 ) , is defined as:
z ^ ( s 0 , t 0 ) = m ^ ( s 0 , t 0 ) + ε ^ ( s 0 , t 0 )
where m ^ ( s 0 , t 0 ) is the estimated value for the location ( s 0 , t 0 ) obtained by adjusting GAMLSS models (Equation (16)). Figure 2 shows the flowchart of regression kriging with GAMLSS model.
The interpolation performance in space–time kriging was evaluated using the leave-one-out method, which consisted of removing an observed point from the precipitation in space–time ( p ( s 0 , t ) ) and then predicting its value ( p ^ ( s 0 , t ) ) . This process was then repeated for all remaining points, and the residue (E) of this procedure was then obtained by the difference between the observed and predicted values at each location ( E = p ^ ( s 0 , t ) p ( s 0 , t ) ) . After that, the RMSE (root mean squared error), MAE (mean absolute error), and R2 (coefficient of determination) were extracted as selection criteria.
All statistical analyses were performed in software R, using the libraries: gamlss [14], gstat [25], and sp [26].

4. Results

4.1. Descriptive Analysis

The precipitation variable, measured monthly in each of the 54 locations, in the period 1994–2014, presented a median value of 29.05 mm and an interquartile range of 3.60–90.20 mm. In the whole historical series, approximately 20% of the observations presented 0.00 mm of total monthly precipitation. In order to evaluate the precipitation distribution of the studied region, a histogram was initially developed for the entire dataset (Figure 3).
In Figure 3, it is possible to observe the excess of zeroes values in the dataset and, consequently, the presence of asymmetric behavior, indicating that the Gaussian model is not the most suitable for modeling.
Box plots were built to explore the behavior of precipitation at the 54 pluviometric stations (Figure 4).
Due to the asymmetric behavior and the high variability of precipitation, the State of Paraíba suffers from irregular distribution of rainfall throughout the region. This distribution presented high variability among the mesoregions analyzed in this work. A positive asymmetry was observed in all the locations, as well as the presence of atypical points (Figure 4).

4.2. Trend Analysis

The distribution of precipitation relative to geographical coordinates (latitude and longitude) and temporal index, with the corresponding adjustments through cubic splines, is presented in Figure 5.
An exploratory analysis by adjusting the cubic smoothness function of the spline showed that precipitation occurred more frequently in the extremes of the studied region, corresponding to the mesoregions of Sertão and Zona da Mata, and that there was limited precipitation in the areas located near the meridian 36° W, corresponding to the Borborema mesoregion (Figure 5a). The spatial distribution of rainfall in the region presented a spatial tendency that varied with the latitudinal and longitudinal gradients (Figure 5a,b). For the analysis of the time trend, the graph of the time series of precipitation was calculated, based on the average of 54 pluviometric stations (Figure 5c). This series showed an evident seasonal behavior of 12 months. Thus, for this series, a GAMLSS model was fitted using cubic smoothing splines function.
To model the trend of precipitation data, the GAMLSS model was used, via zero-adjusted gamma (ZAGA) distribution, as indicated in the system of Equations (4). Table 2 presents the coefficient estimates, standard error, t value, p-value, and the coefficient of determination (R2) in the model adjustment.
The results presented in Table 2 indicate, at the level of 0.01 of significance, that the covariates longitude, latitude, and the temporal index are related to the precipitation. The R2 shows that 48% of the large-scale variation of rainfall was explained by the trend component, corroborating with the results found through the descriptive analysis.
After estimating the trend component, which was described in Equation (1), we then analyzed the stochastic residue, ε ^ ( s , t ) = Z ( s , t ) m ^ ( s , t ) . Next, the space–time sampling variogram of these residues was calculated, taking as reference the adjustment of the generalized product–sum theoretical model.

4.3. Geostatistical Analysis

The residuals, obtained from the regression via GAMLSS models, presented an evident pattern of spatiotemporal correlation (Figure 6a). There was a strong spatial correlation between nearby sites, and the spatial structure became weaker as time differences increased. Similarly, the same occurred with the temporal structure (Figure 6b). Thus, spatial–temporal kriging of the residue was required.
The behavior of the empirical space–time variogram estimates is presented in terms of the marginal effects (Figure 7).
The spatial variability exhibited by the marginal spatial variogram (Figure 7a) differed from the temporal variability shown in the marginal temporal variogram (Figure 7b), since a difference, concerning the threshold, was observed between the two marginal variograms. For this reason, the spatiotemporal variability of the residues could be modeled through the generalized product–sum model. One of the advantages of constructing marginal variograms is the fact that it is possible to identify some theoretical variogram models (for example, Gaussian, Exponential, and Spherical) that can compose the structure in the generalized product–sum model (Equation (12)).
As described above, the Normal, Gamma, and ZAGA models have been adjusted for large-scale variations. The different structures of the generalized product–sum model were adjusted for the residues obtained in each of these models, followed by the spatial–temporal kriging of that model. The cross-validation results are provided in Table 3. The results of the validation were calculated after the values adjusted by the trend component were added to the values interpolated in space–time kriging.
The results presented in Table 3 indicate that the Exponential model was the best fit for the spatial ( γ s t ( h s , 0 ) ) and temporal ( γ s t ( 0 , h t ) ) structures in all regression models, since it presented higher R2 and smaller RMSE and MAE. Note that the ZAGA model (EXP + EXP) presented the best results, with RMSE of 34.598 mm, MAE of 20.970 mm, and R2 equal to 82.8%. These results demonstrated the effectiveness of adjusting the precipitation data to an adequate model that contemplates the spatial–temporal characteristics presented by the phenomenon.
The components of the generalized product–sum model (Equation (12)) and their estimates are shown in Table 4. The components were chosen based on the results in Table 3.
The parameter estimates (Table 4) indicated that all components act on the spatial–temporal pattern of the ZAGA regression residuals. The parameter estimation related to the spatial extent was high, indicating that the residues were correlated in distances of up to 171 km. The temporal scale estimate indicated that sites with a delay of up to 47 days had residual autocorrelation.
The main advantage of the methodology proposed in this study is that with the space–time kriging technique, predictions can be made for unobserved locations and times [21,24]. As an example of this methodology applicability, the spatial–temporal kriging of precipitation was carried out in the period 2015 (Figure 8).
In Figure 8, it is possible to evaluate the precipitation behavior in nonsampled sites, as well as in times not observed in the sample, helping policymakers manage the allocation of water in the present and the future. It is noted that the spatial distribution of rainfall presented variability within each mesoregion and also among the four mesoregions. Precipitation in the region occurred more intensely during the first half of the year, but during the same period, there was an irregularity in the distribution of rainfall between the mesoregions. In the Sertão, the high spatial variability can be justified by meteorological systems in the area, such as the the Intertropical Convergence Zone (ITCZ) and the Upper Tropospheric Cyclonic Vortices (UTCV). It is observed that in this region, the rainfall distribution during February to May occurred irregularly. In years when El Niño and La Niña occurred with strong intensity, there was substantial evidence that rainfall occurred, respectively, above or below the expected value, resulting in high temporal variability of precipitation [27]. The rainy season in this area is concentrated between February and May, presenting high variation in space and time. The low volume of rainfall and the irregularity of its distribution have contributed, for example, to the occurrence of the desertification phenomenon [5,28].

5. Discussion

The state of Paraíba has abnormalities in the distribution of precipitation, with high spatial and temporal variability, resulting in heavy rains with short duration and prolonged periods of drought. Consequently, the scarcity of rain directly affects water resources, resulting in severe problems in the region’s economy, which focuses on family agriculture and livestock [29]. Also, drought has affected water distribution for human consumption, electricity supply, the occurrence of diseases and migration, among other issues, resulting in an unprecedented water crisis [30]. Thus, a space–time distribution modeling of rainfall becomes of paramount importance for the water resources management to reduce the impacts caused by drought in this region.
Irregular rainfall, resulting in periods of prolonged drought, is associated with meteorological phenomena. Strong El Niño years have caused intense droughts in Paraíba [27]. Rainfall periods are concentrated in a few months, occurring mainly in winter, which is associated with the ITCZ action and the passage of cold fronts from the South Atlantic Ocean [31]. Due to the high rainfall variability of the region, rainfall behavior can be classified into three periods, namely: pre-rainy season, rainy season and post-rainy season. During the pre-rainy season, which runs from December to January in the Sertão mesoregion, meteorological systems such as UTCV and instability associated with cold fronts occur. The rainy season itself persists from February to May, when the ITCZ acts as the primary meteorological system responsible for the recharge of water resources in the state. The post-rainy season usually occurs between May and July in the Zona da Mata mesoregion and some localities of the Southeast. The Eastern Wave Disturbances (EWDs) are the primary meteorological system acting in this period.
In this study, we used the gamma distribution, which shows asymmetric behavior, to model precipitation data. However, it is known that this probability distribution does not include values equal to zero, not having a parameter that can be used to model the excess of zero that is common in data of monthly total precipitation throughout the NEB. As an advantage, the methodology allowed modeling the occurrence of null values, resulting in better estimates when compared with models that do not contemplate the excess of zeroes in the dataset.
Several studies that deal with spatial–temporal modeling and which have null values in the database have used the regression model with Normal distribution [11,32] or the Gamma regression model [12,13]. However, the Normal regression does not solve the data asymmetry problem, and the Gama regression does not contemplate in its domain the occurrence of zeroes. Thus, we propose a new way of modeling the trend component through GAMLSS models, especially using the ZAGA distribution, allowing one of its three parameters to model the excess of zeroes that occurs with high frequency in data of monthly total precipitation.
The trend component alone does not model the spatiotemporal dependence that is present in the regression residuals. Thus, the need for space–time modeling of these residues arises. Previous studies did not take into account this spatial and/or temporal residue dependence [7,33]. Consequently, we modeled the spatiotemporal dependence structure of this stochastic component, using a covariance function. However, this approach has numerous parameters present in the space–time theoretical variogram, resulting in great computational difficulties to determine the parameter estimates.
A significant benefit of the suggested methodology is that it can predict in nonsampled places and in periods that have not been observed, thereby assisting water-related initiatives in the region. The Agreste and Borborema mesoregions, with a reduced amount of precipitation, were discovered to be an aggravating factor in the uneven rainfall distribution in the study region. The majority of the state of Paraíba is part of the Brazilian semi-arid region, which has a climate characterized by low humidity and rainfall. Studies in recent years have revealed the process of desertification in the semi-arid Northeast regions [30]. This research can, therefore, identify areas that will soon be affected by water supplies for human and animal consumption. These findings can serve as a supporting tool to help government agencies in water management to predict societal needs and avoid water shortages that are rapidly intensifying, not only in Paraíba State, but also in the NEB as a whole.

6. Conclusions

In this study, total monthly precipitation in the state of Paraíba was analyzed using spatiotemporal geostatistical tools. Rainfall data of 54 pluviometric stations from 1994 to 2014 were analyzed. Large-scale variation was modeled using spatial and temporal covariates to remove the spatiotemporal trend in the regression residuals. It was proposed to use the adjustment of GAMLSS models with different distributions for the response variable (Normal, Gamma, and ZAGA). The results of the cross-validation indicated that the ZAGA distribution best fits the data, with the exponential model for the spatial component and also for the time component of the generalized product–sum variogram model. Also, the results indicated the presence of the space–time correlation in the rainfall phenomenon, implying in this way the proper use of space–time kriging proposed in this study. Finally, space–time kriging in Paraíba state for the year 2015 evidenced the great irregularity in the distribution of rainfall throughout the region, providing a detailed visual analysis of sectors suffering from scarcity and excess of rain, favoring government policies that address water resources.

Author Contributions

Conceptualization, E.S.d.M., R.R.d.L., R.A.d.O., L.G.D. and C.A.C.d.S.; Methodology, E.S.d.M., R.R.d.L. and R.A.d.O; Software, E.S.d.M.; Validation, E.S.d.M. and R.R.d.L.; Formal analysis, E.S.d.M., R.R.d.L., R.A.d.O., L.G.D. and C.A.C.d.S.; Investigation, E.S.d.M., R.R.d.L., R.A.d.O., L.G.D. and C.A.C.d.S.; Resources, E.S.d.M.; Data curation, E.S.d.M. and C.A.C.d.S.; Writing—original draft preparation, E.S.d.M., R.R.d.L., R.A.d.O., L.G.D. and C.A.C.d.S.; Writing—review and editing, E.S.d.M., R.R.d.L., R.A.d.O., L.G.D. and C.A.C.d.S.; Visualization, E.S.d.M.; supervision, E.S.d.M., R.R.d.L., R.A.d.O. and C.A.C.d.S.

Funding

This research was funded by the Brazilian National Council for Scientific and Technological Development (CNPq) and the Brazilian Federal Agency for the Support and Evaluation of Graduate Education (CAPES).

Acknowledgments

The authors acknowledge the Brazilian National Council for Scientific and Technological Development (CNPq) and the Coordination for the Improvement of Higher Education Personnel (CAPES)-Finance Code 001 (Pró-Alertas Research Project-Grant No. 88887.091737/2014-01 and 88887.123949/2015-00). We thank Lacey Bodnar from the Robert B. Daugherty Water for Food Global Institute (DWFI) at the University of Nebraska for her invaluable assistance, and the support of DWFI is also appreciated.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Kilibarda, M.; Tadić, M.P.; Hengl, T.; Luković, J.; Bajat, B. Global geographic and feature space coverage of temperature data in the context of spatio-temporal interpolation. Spat. Stat. 2015, 14, 22–38. [Google Scholar] [CrossRef]
  2. Martínez, W.A.; Melo, C.E.; Melo, O.O. Median Polish Kriging for space-time analysis of precipitation. Spat. Stat. 2017, 19, 1–20. [Google Scholar] [CrossRef]
  3. Stauffer, R.; Mayr, G.J.; Messner, J.W.; Umlauf, N.; Zeileis, A. Spatio-temporal precipitation climatology over complex terrain using a censored additive regression model. Int. J. Climatol. 2017, 37, 3264–3275. [Google Scholar] [CrossRef] [PubMed]
  4. Marengo, J.A.; Torres, R.R.; Alves, L.M. Drought in Northeast Brazil—past, present, and future. Theor. Appl. Climatol. 2017, 129, 1189–1200. [Google Scholar] [CrossRef]
  5. Borges, B.C.; de Baptista, G.M.; Meneses, P.R. Identification of hydromorphic areas by means of spectral analysis of remote sensing data, as support for the preparation of forensic reports issued to evaluation of rural properties. Rev. Bras. Geogr. Física 2014, 7, 1062–1077. [Google Scholar] [CrossRef]
  6. De Almeida, H.A.; Medeiros, E.A. Variabilidade no regime pluvial em duas mesorregiõres da Paraíba e sua relação com o fenômeno EL Niño Oscilação Sul. J. Environ. Anal. Prog. 2017, 2, 177. [Google Scholar] [CrossRef]
  7. Dos Santos, C.A.C.; Gomes, O.M.; de Souza, F.A.S.; Paiva, W. De Análise Geoestatística da Precipitação Pluvial do Estado da Paraíba (Geostatistical Analysis of Precipitation in Paraiba State). Rev. Bras. Geogr. Física 2012, 4, 692. [Google Scholar] [CrossRef]
  8. Francisco, P.R.M.; da Mello, V.S.; Bandeira, M.M.; de Marcedo, F.L.; Santos, D. Discriminação de cenários pluviométricos do estado da Paraíba utilizando distribuição Gama Incompleta e Teste Kolmogorov-Smirnov. Rev. Bras. Geogr. Física 2016, 9, 41–61. [Google Scholar]
  9. Gräler, B.; Gerharz, L.; Pebesma, E. Spatio-temporal analysis and interpolation of PM10 measurements in Europe. ETC/ACM Tech. Pap. 2011/10 2012, 8, 37. [Google Scholar]
  10. Hengl, T.; Heuvelink, G.B.M.; Perčec Tadić, M.; Pebesma, E.J. Spatio-temporal prediction of daily temperatures using time-series of MODIS LST images. Theor. Appl. Climatol. 2012, 107, 265–277. [Google Scholar] [CrossRef]
  11. Hu, D.; Shu, H.; Hu, H.; Xu, J. Spatiotemporal regression Kriging to predict precipitation using time-series MODIS data. Cluster Comput. 2017, 20, 347–357. [Google Scholar] [CrossRef]
  12. Menezes, R.; Piairo, H.; García-Soidán, P.; Sousa, I. Spatial–temporal modellization of the NO2 concentration data through geostatistical tools. Stat. Methods Appt. 2016, 25, 107–124. [Google Scholar] [CrossRef]
  13. Monteiro, A.; Menezes, R.; Silva, M.E. Modelling spatio-temporal data with multiple seasonalities: The NO2 Portuguese case. Spat. Stat. 2017, 22, 371–387. [Google Scholar] [CrossRef]
  14. Stasinopoulos, M.; Rigby, R.A.; Heller, G.Z.; Voudouris, V.; De Bastiani, F. Flexible Regression and Smoothing Using GAMLSS in R; CRC Press: New York, NY, USA, 2017. [Google Scholar]
  15. Montero, J.-M.; Fernández-Avilés, G.; Mateu, J. Spatial and Spatio-Temporal Geostatistical Modeling and Kriging; John Wiley & Sons, Ltd: Chichester, UK, 2015. [Google Scholar]
  16. Schabenberger, O.; Gotway, C.A. Statistical Methods for Spatial Data Analysis; CRC Press: New York, NY, USA, 2004. [Google Scholar]
  17. Iaco, S.D.; Myers, D.E.; Posa, D. Space–time analysis using a general product–sum model. Stat. Probab. Lett. 2001, 52, 21–28. [Google Scholar] [CrossRef]
  18. Huang, H.-C.; Martinez, F.; Mateu, J.; Montes, F. Model comparison and selection for stationary space–time models. Comput. Stat. Data Anal. 2007, 51, 4577–4596. [Google Scholar] [CrossRef]
  19. Cressie, N.; Wikle, C.K. Statistics for Spatio-Temporal Data; John Wiley & Sons, Inc.: Hoboken, NJ, USA, 2015. [Google Scholar]
  20. Iaco, S.D.; Myers, D.E.; Posa, D. Nonseparable space-time covariance models: Some parametric families. Math. Geol. 2002, 34, 23–42. [Google Scholar] [CrossRef]
  21. Kilibarda, M.; Hengl, T.; Heuvelink, G.B.M.; Gräler, B.; Pebesma, E.; Perčec Tadić, M.; Bajat, B. Spatio-temporal interpolation of daily temperatures for global land areas at 1 km resolution. J. Geophys. Res. Atmos. 2014, 119, 2294–2313. [Google Scholar] [CrossRef]
  22. Sherman, M. Spatial Statistics and Spatio-Temporal Data; Wiley Series in Probability and Statistics; John Wiley & Sons, Ltd.: Chichester, UK, 2010. [Google Scholar]
  23. Cesare, L.D.; Myers, D.; Posa, D. Estimating and modeling space–time correlation structures. Stat. Probab. Lett. 2001, 51, 9–14. [Google Scholar] [CrossRef]
  24. Iaco, S.D.; Palma, M.; Posa, D. Spatio-temporal geostatistical modeling for French fertility predictions. Spat. Stat. 2015, 14, 546–562. [Google Scholar] [CrossRef]
  25. Graler, B.; Pebesma, E.; Heuvelink, G. Spatio-temporal interpolation using gstat. Wp 2016, 8, 1–20. [Google Scholar] [CrossRef]
  26. Bivand, R.S.; Pebesma, E.; Gómez-Rubio, V. Hello World: Introducing Spatial Data. In Applied Spatial Data Analysis with R; Springer: New York, NY, USA, 2013; pp. 1–16. [Google Scholar]
  27. Marengo, J.A.; Alves, L.M.; Alvala, R.C.; Cunha, A.P.; Brito, S.; Moraes, O.L.L. Climatic Characteristics of the 2010–2016 Drought in the Semiarid Northeast Brazil Region. 2017. Available online: http://www.scielo.br/scielo.php?script=sci_arttext&pid=S0001-37652018000501973 (accessed on 7 November 2019).
  28. Souza, B.I.; Menezes, R.; Artigas, R.C. Efeitos da desertificação na composição de espécies do bioma Caatinga, Paraíba/Brasil. Investig. Geográficas 2015, 45–59. [Google Scholar] [CrossRef]
  29. Macedo, M.J.; Guedes, R.; de Sousa, F.A.; Dantas, F. Analysis of the standardized precipitation index for the Paraíba state, Brazil. Ambient. e Agua-An Interdiscip. J. Appl. Sci. 2010, 5, 204–214. [Google Scholar] [CrossRef]
  30. Gutiérrez, A.P.A.; Engle, N.L.; De Nys, E.; Molejón, C.; Martins, E.S. Drought preparedness in Brazil. Weather Clim. Extrem. 2014, 3, 95–106. [Google Scholar] [CrossRef] [Green Version]
  31. Rodrigues, R.R.; McPhaden, M.J. Why did the 2011–2012 La Niña cause a severe drought in the Brazilian Northeast? Geophys. Res. Lett. 2014, 41, 1012–1018. [Google Scholar] [CrossRef]
  32. Guo, L.; Lei, L.; Zeng, Z.-C.; Zou, P.; Liu, D.; Zhang, B. Evaluation of spatio-temporal variogram models for mapping Xco2 using satellite observations: A case study in China. IEEE J. Sel. Top. Appl. Earth Obs. Remote Sens. 2015, 8, 376–385. [Google Scholar] [CrossRef]
  33. Gomes, O.M.; dos Santos, C.A.C.; de Souza, F.A.S.; de Paiva, W.; de Olinda, R.A. Análise comparativa da precipitação no estado da paraíba utilizando modelos de regressão polinomial. Rev. Bras. Meteorol. 2015, 30, 47–58. [Google Scholar] [CrossRef]
Figure 1. Pluviometric stations over the four mesoregions of Paraíba State, Brazil. The information for each respective number is presented in Table 1.
Figure 1. Pluviometric stations over the four mesoregions of Paraíba State, Brazil. The information for each respective number is presented in Table 1.
Water 11 02368 g001
Figure 2. Regression kriging flowchart with GAMLSS model. GAMLSS: generalized additive models for location, scale, and shape.
Figure 2. Regression kriging flowchart with GAMLSS model. GAMLSS: generalized additive models for location, scale, and shape.
Water 11 02368 g002
Figure 3. Histogram of the precipitation data of the 54 pluviometric stations sampled in the period from 1994 to 2014.
Figure 3. Histogram of the precipitation data of the 54 pluviometric stations sampled in the period from 1994 to 2014.
Water 11 02368 g003
Figure 4. Box plot of monthly rainfall for 54 pluviometric stations, listed in Table 1, in this study.
Figure 4. Box plot of monthly rainfall for 54 pluviometric stations, listed in Table 1, in this study.
Water 11 02368 g004
Figure 5. Distribution of precipitation data (mm) (blue dots) in the longitude direction (a) and in the latitude direction (b) and respective trend lines (gray line) using cubic smoothing splines function over the study region; time series of precipitation data (mm) (blue line) and adjusted temporal trend (gray line) (c).
Figure 5. Distribution of precipitation data (mm) (blue dots) in the longitude direction (a) and in the latitude direction (b) and respective trend lines (gray line) using cubic smoothing splines function over the study region; time series of precipitation data (mm) (blue line) and adjusted temporal trend (gray line) (c).
Water 11 02368 g005
Figure 6. Empirical space–time variogram (a) and adjusted variogram (b).
Figure 6. Empirical space–time variogram (a) and adjusted variogram (b).
Water 11 02368 g006
Figure 7. Marginal experimental variograms of the residuals: purely spatial (a) and purely temporal (b).
Figure 7. Marginal experimental variograms of the residuals: purely spatial (a) and purely temporal (b).
Water 11 02368 g007
Figure 8. Spatial–temporal prediction of total monthly rainfall (mm) in the study region for the year 2015. The crosses (+) represent the positions of each pluviometric station.
Figure 8. Spatial–temporal prediction of total monthly rainfall (mm) in the study region for the year 2015. The crosses (+) represent the positions of each pluviometric station.
Water 11 02368 g008
Table 1. Pluviometric stations with an identification number (ID), name of the mesoregion, name of rainfall station, latitude (Lat), longitude (Long), Altitude (Alt), and precipitation quantiles 2.5%, 50%, 97.5%.
Table 1. Pluviometric stations with an identification number (ID), name of the mesoregion, name of rainfall station, latitude (Lat), longitude (Long), Altitude (Alt), and precipitation quantiles 2.5%, 50%, 97.5%.
IDMesoregionStationLatLongAlt (m)Quantiles
2.5%50.0%97.5%
1Zona da MataAlhandra−7.43−34.9156.286.49119.85481.71
2Zona da MataJacaraú−6.61−35.29181.190.8855.90311.01
3Zona da MataMamanguape−6.84−35.1216.622.6066.70403.02
4Zona da MataPedras de Fogo−7.40−35.12175.364.7179.25399.52
5Zona da MataSapé−7.09−35.22116.841.8664.35324.14
6AgresteAlagoinha−6.96−35.55168.640.0070.50306.56
7AgresteAraçagi−6.83−35.39104.620.0056.95290.73
8AgresteAraruna−6.53−35.74575.410.0053.40235.32
9AgresteAreia−6.98−35.72572.284.3495.55331.73
10AgresteAreial−7.05−35.93693.090.0038.55182.72
11AgresteBoa Vista−7.26−36.24486.000.0020.30127.07
12AgresteCampina Grande−7.23−35.90544.370.9348.00244.28
13AgresteCasserengue−6.79−35.89396.630.0018.00137.88
14AgresteCuité−6.49−36.15668.960.0033.35190.79
15AgresteDona Inês−6.61−35.63421.720.0051.25242.04
16AgresteIngá−7.29−35.61155.560.0044.35191.81
17AgresteMogeiro−7.31−35.48108.280.0040.95194.68
18AgrestePocinhos−7.08−36.06650.320.0021.45140.05
19AgresteSoledade−7.06−36.36523.630.0017.75139.61
20BorboremaBarra de São Miguel−7.75−36.32488.630.0015.30142.21
21BorboremaBoqueirão−7.49−36.14355.080.0021.65148.11
22BorboremaCamalaú−7.89−36.83519.160.0011.15185.53
23BorboremaCaraúbas−7.73−36.49442.240.0010.60161.65
24BorboremaCongo−7.80−36.66491.990.0012.25164.82
25BorboremaJuazeirinho−7.07−36.58553.960.0020.10175.97
26BorboremaJunco do Seridó−7.00−36.71589.560.0022.40210.18
27BorboremaPedra Lavrada−6.76−36.46521.820.0013.55180.72
28BorboremaPrata−7.70−37.08584.000.0018.85220.19
29BorboremaSalgadinho−7.10−36.85430.450.0013.60240.85
30BorboremaSanta Luzia−6.87−36.92311.240.009.80246.99
31BorboremaSão João do Tigre−8.08−36.85572.840.0013.80161.03
32BorboremaSão José dos Cordeiros−7.39−36.81530.490.0014.30313.51
33BorboremaSão Seb. do Umbuzeiro−8.15−37.01595.790.0020.10200.60
34BorboremaSumé−7.67−36.90519.880.0015.70264.53
35BorboremaVárzea−6.77−36.99267.640.0010.65270.29
36SertãoAgua Branca−7.51−37.64732.800.0041.40274.52
37SertãoBom Sucesso−6.44−37.93289.120.0027.20296.80
38SertãoBrejo do Cruz−6.35−37.50200.670.0027.70320.60
39SertãoCajazeiras−6.89−38.54299.440.0038.80409.20
40SertãoCatolé do Rocha−6.34−37.75298.890.0034.55301.66
41SertãoConceição−7.56−38.50388.280.0027.90272.48
42SertãoCondado−6.92−37.59260.820.0023.75336.19
43SertãoLagoa−6.59−37.91275.560.0031.30314.85
44SertãoMãe D’Água−7.26−37.43411.100.0015.65277.71
45SertãoManaíra−7.71−38.15767.400.0028.95284.65
46SertãoNova Olinda−7.48−38.04321.000.0027.55326.82
47SertãoPassagem−7.14−37.05305.040.0014.65273.24
48SertãoPatos−7.00−37.31256.690.0026.25303.92
49SertãoPiancó−7.21−37.93261.160.0024.90321.85
50SertãoPombal−6.77−37.80191.560.0022.85353.35
51SertãoRiacho dos Cavalos−6.44−37.65206.930.0024.45269.48
52SertãoSanta Teresinha−7.08−37.45307.160.0022.55368.42
53SertãoSão J. do Rio do Peixe−6.73−38.45248.200.0031.90325.64
54SertãoSousa−6.77−38.22235.440.0030.30327.17
Table 2. Estimates, standard errors (SE), p-values, and the coefficient of determination (R2) of the GAMLSS regression model for the dataset used in this study.
Table 2. Estimates, standard errors (SE), p-values, and the coefficient of determination (R2) of the GAMLSS regression model for the dataset used in this study.
ParameterEstimativeSEt Valuep-ValueR2
μ β 10 −16.7600.517−32.399< 10 4 0.48
β 11 < 10 4 < 10 4 12.717< 10 4
β 12 0.002< 10 4 39.081< 10 4
β 13 < 10 4 < 10 4 −5.955< 10 4
σ β 20 0.8550.005163.700< 10 4
ν β 30 −0.0100.002−59.350< 10 4
Table 3. Space-time kriging results combining the different trend models (NORMAL, GAMMA and ZAGA) with the different variogram models: Gaussian (GAU), Exponential (EXP) and Spherical (SPH). The RMSE (root mean squared error), MAE (mean absolute error), and R2 (coefficient of determination) estimates were obtained by "leave-one-out" cross-validation.
Table 3. Space-time kriging results combining the different trend models (NORMAL, GAMMA and ZAGA) with the different variogram models: Gaussian (GAU), Exponential (EXP) and Spherical (SPH). The RMSE (root mean squared error), MAE (mean absolute error), and R2 (coefficient of determination) estimates were obtained by "leave-one-out" cross-validation.
MODELNORMALGAMMAZAGA
RMSEMAER2RMSEMAER2RMSEMAER2
GAU + EXP36.09921.8810.81336.65222.0650.80736.11621.9200.812
GAU + SPH36.88122.4140.80536.70822.0900.80737.59622.8950.797
EXP + EXP34.71021.0110.82735.38121.2530.82134.59820.9700.828
EXP + SPH34.99121.1810.82435.38821.2620.82135.34221.4300.820
Table 4. Parameter estimates of the generalized product–sum variogram model adjusted to precipitation residues. The final model is obtained by replacing these estimates in Equation (12).
Table 4. Parameter estimates of the generalized product–sum variogram model adjusted to precipitation residues. The final model is obtained by replacing these estimates in Equation (12).
ComponentVariogram ModelNuggetThresholdRangeK
Spatial ( γ s t ( h s , 0 ) )Exponential0.7574.350171 km15.689
Temporal ( γ s t ( 0 , h t ) )Exponential17.47952.05847 days

Share and Cite

MDPI and ACS Style

Medeiros, E.S.d.; de Lima, R.R.; Olinda, R.A.d.; Dantas, L.G.; Santos, C.A.C.d. Space–Time Kriging of Precipitation: Modeling the Large-Scale Variation with Model GAMLSS. Water 2019, 11, 2368. https://doi.org/10.3390/w11112368

AMA Style

Medeiros ESd, de Lima RR, Olinda RAd, Dantas LG, Santos CACd. Space–Time Kriging of Precipitation: Modeling the Large-Scale Variation with Model GAMLSS. Water. 2019; 11(11):2368. https://doi.org/10.3390/w11112368

Chicago/Turabian Style

Medeiros, Elias Silva de, Renato Ribeiro de Lima, Ricardo Alves de Olinda, Leydson G. Dantas, and Carlos Antonio Costa dos Santos. 2019. "Space–Time Kriging of Precipitation: Modeling the Large-Scale Variation with Model GAMLSS" Water 11, no. 11: 2368. https://doi.org/10.3390/w11112368

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop