Precipitation forecast estimation applying the change point method and ARIMA

Abstract

In this study, a time series is analyzed with daily historical data from 1989 to December 2021 of the precipitation variable with a total of 12,053 observations, these data are obtained from the Tunguativa weather station located in Paipa a small town in Colombia. For this research, the records of the ‘precipitation’ variable were taken into account. The objective was to analyze the trends, use the data until 31 December 2020, to estimate a forecast for the year 2021. The observed data is provided by the weather station, and the data of the year 2021 are compared with the data forecast by the Arima method and the change point. Finally, statistical tests are carried out to contrast the degree of similarity of the data obtained from the forecasts with the real data produced by the station. This study can be used to estimate precipitation forecasts and avoid environmental damage.

Keywords:

• ARIMA
• change point
• forecast
• precipitation
• time series

Reviewing Editor:

• Qingsong Ai, Wuhan University of Technology, China

Subjects:

• meteorology
• environmental studies
• research methods in environmental studies
• statistics & computing
• general engineering education

1. Introduction

The small town of Paipa is located in the Sogamoso valley, one of the most important internal valleys of the Andean region, in the eastern central part of Colombia at 2525 m above sea level, approximately 184 km from Bogotá (Colombia’s capital). The economy of this town largely depends on tourism. One of its main tourist attractions is Lake Sochagota. The Lake Sochagota it is an artificial body of water built in 1956 by the government as a means of promoting tourism and promoting nautical sports in the region. It has an area of ​​163 Hectares, a volume of 4,557,000 cubic meters, an average depth of 2.80 m, and a maximum depth of 3.20 m. Containing in its first meter of depth, approximately 1,623,000 cubic meters, which this means that for every centimeter that the water level drops in this first meter, the reservoir loses 16,230 cubic meters of water (https://www.corpoboyaca.gov.co/sochagota/).

The lake has a lentic condition, which implies resting or accumulated water, this facilitates the precipitation of solids, the greater penetration of light and, on occasions, its heating, which favors greater algal and biological activity in general. The water quality of Lake Sochagota depends largely on the flows of the Salitre stream and the characteristics of its waters. The collapse of the perimeter sewer system that causes the overflow of wastewater by the municipality’s tourism sector, which eutrophicates (mass of water naturally or artificially enriched with nutrients, such as phosphates and nitrates) from the reservoir and from the thermo-mineral water deposits in the lake, which provide a considerable amount of salts that can be harmful to biota (https://www.corpoboyaca.gov.co/sochagota/).

All of the above leads to the enrichment of the aquatic environment with nutrients, favoring the proliferation of producing organisms, such as algae. This constitute phytoplankton and which, being present in abundance, are their photosynthetic activity, abundantly consuming oxygen from the water. With their photosynthetic activity, consume oxygen from the water at night, being able to exhaust it and then appear in states of anoxia (a condition characterized by the absence of oxygen). Which causes suffocation and death of animals (such as fish), and chemical reduction reactions that release substances with bad odors, such as hydrogen sulfide, which indicates that the environmental capacity of the aquatic ecosystem has been exceeded, due to the excessive presence of nutrients and phytoplankton, caused by anthropic pressure and its repetitive actions, which produce cumulative impacts on the ecosystem.

All of the above happens in seasons where precipitation is high. It must be taken into account that Lake Sochagota is a constructed system that acts as an oxidation pond. Its level is regulated by gates, it has an outcrop of thermo-mineral water and receives runoff from acid- sulfated soils and a portion of fresh water from the Chicamocha River. The gates help move the water and oxygenate the lake water. It is for this reason that Lake Sochagota has been the subject of many proposals, studies and plans for its management, recovery and protection. This research can contribute to predict in which seasons there are precipitations to open the gates of the lake and avoid environmental problems.

The observed data provided by the Tunguavita weather station go up to the year 2021, for this reason we use the data up to the year 2020 to forecast the year 2021 and compare the observed data with the results of the forecasts for that year.

2. Theoretical framework

This article covers the theoretical framework to contextualize the methods used; the methodology used to achieve the objectives; the analysis and interpretation of the results to identify the behavior of the models calculating the precipitation forecast; and finally, the conclusions are presented, summarizing the most relevant aspects, such as the scope and limitations of the research.

Time series analysis has been widely used to predict precipitation events. Predicted future weather conditions play an important role in the decision-making process of many organizations. In order to predict future events of weather variables, it is necessary to rely on information obtained from past events. Forecasting weather variables resulting from numerous interactions that create complex systems is a difficult but important task. Modeling time series to generate data and forecast precipitation is an important step in the design and analysis of water resources.

Forecasting precipitation variables is a very useful tool in the management of natural resources. The importance of forecasting in precipitation issues makes them use more accurate statistical methods to study weather and climate change, the main objective of this study is to investigate the use of additive and multiplicative forms of ARIMA time series model and change point to predict precipitation variables such as precipitation.

Precipitation is directly produced by cloud microphysics, but the occurrence of precipitation is associated with precipitation dynamics and thermodynamics of climatic events (Sui et al., 2007). Forecasts used in weather variables benefit the agricultural sector, since, in many regions, interannual climate variability can affect agricultural production, decrease farmers’ income and increase market prices (Esquivel et al., 2018).

Limited water resources for water supply in the agricultural and industrial sector have caused major problems, particularly in arid and semi-arid areas. Therefore, forecasting and estimation of precipitation is one of the most important climatic parameters for efficient use of water resources (Feidas et al., 2007). In general, measurement and prediction of this variable are required to study its behavior. Moreover, the prediction of precipitation for any type of territory is considered one of the most important climatic parameters for the optimal use of water resources (Nastos et al., 2014).

Climate prediction is the process of predicting the likely trend of future climate development based on the changing laws of past climate. In recent years, the progressive role of climate prediction in disaster prevention and mitigation has been increasingly recognized, with the needs of social and economic development emerging in the world, climate prediction research should be urgently improved (Li et al., 2022).

Next, forecasts are studied using the Bayesian model together with the change point theory versus the traditional ARIMA method. The study uses a quantitative approach and compares the results of the change point model with the ARIMA model used to forecast the behavior of precipitation during one year (2020), as a result, the values obtained in the forecasts are compared with the observed data of the original series.

As described by Phan et al. (2018) the algorithms commonly used to perform forecasting studies of meteorological variables are: Single Exponential Smoothing (SES), Seasonal Naive (Snaive), Seasonal ARIMA (SARIMA), Feedforward Neural Network (FFNN), Dynamic Time Warping Based Imputation (DTWBI) and Bayesian Structural Time Series (BST). There are other methods that are used by airport weather stations as stated by Jacobs and Maat (2005) where numerical weather prediction (NWP) data are used which provide an acceptable quality of forecasts, or as mentioned by Arroyo and Maté (2009) where by applying the nearest neighbor algorithm (KNN) has the ability to yield promising results in forecasting meteorological variables, and according to Heydari et al. (2020) the ARIMA model is an effective method for prediction and modeling of weather parameters. ARIMA is most commonly used time series model for short-term load forecasting. The ARIMA model forecasts a linear time series data and uses ARIMA to transform the stationary series. ARIMA method is most commonly used, as it predicts the load purely based on the historical loads and no other assumptions are considered.

The change point model is also a tool that has been used to analyze time series data in areas such as the stock market (Zhitlukhin & Ziemba, 2016; Lenardon & Amirdjanova, 2006), geological data (Barry & Hartigan, 1993; Cmejla et al., 2013), weather data (Zhitlukhin & Ziemba, 2016; Lenardon & Amirdjanova, 2006), genetic data analysis (Chen & Gupta, 2012) and to predict catastrophes, such as earthquakes (Alawadhi & Alhulail, 2016) among others.

Those looking for a change point mainly investigate a problem in two different ways: 1. Detect a failure that does not affect the mean of the observations, but changes their covariance structure. When this happens, they recommend considering procedures that take into account smooth or deep perturbations, perform a non-parametric type test, such as the Kolmogorov–Smirnov test. And 2. Detect a failure and estimate the parameters of change for a failure occurring at the same time in the mean and covariance of a bounded-order autoregressive process. It is recommended to apply a likelihood ratio test and make appropriate normalizations (Picard, 2013).

Statistical inference on the change points has two aspects, the first is to detect if there is any change in the sequence of observed random variables. And the second is to estimate the number of changes and their corresponding locations. The methods for CP estimation are mainly likelihood ratio, nonparametric and Bayesian. Some authors also considered the change point problem in other model configurations, such as gamma and exponential (Chen & Gupta, 2012).

The likelihood model is used for the location of the change point, as well as a developed information criterion, for the case of known variance (Plummer, 2012). ‘This model is a hypothesis test that helps to choose the “best” model between two nested models. “Nested models” means that one is a special case of the other. The best model is the one that makes the data more likely, or maximizes the likelihood function. However, although the concept is relatively easy to understand, the calculations to find the inputs to the procedure are not’ (Lehmann, 2006). Likelihood ratio tests use log-likelihood functions, which are difficult and time-consuming to calculate by hand; most statistical software packages have built-in functions to handle them (Wang et al., 2020).

While nonparametric methods are designed without a particular parametric assumption (Zhou et al., 2017), the nonparametric method is a flexible tool that mechanizes the construction of prediction models, selecting relevant variables, transforming predictor variables, treating missing values and preventing overfitting by means of a self-test. It also allows prediction taking into account structural factors that could influence the response variable, generating hypothetical models (Vanegas & Vásquez, 2017). Nonparametric methods have traditionally been recommended for smaller data sets and/or non-normally distributed data (Grech and Calleja, 2018). However, these methods tend to be less accurate in high dimensional problems due to the so-called dimensionality effect (Liu et al., 2013).

Another CP method is the Bayesian method (BM) for change-point problems, introduced by Barry and Hartigan in 1993. With this model-based method, an unknown partition of blocks is estimated for a data set. Therefore, the probability of a partition is calculated using the so-called cohesion (Herrmann et al., 2014). ‘A sequence of observations undergoes sudden changes at unknown times. We model the process by assuming that there is an underlying sequence of parameters partitioned into contiguous blocks of equal parameter values; the beginning of each block is said to be a change point. The observations are then assumed to be independent in different blocks given the parameter sequence. In a Bayesian analysis, it is necessary to give probability distributions to both the change points and the parameters’ (Barry and Hartigan, 1993).

The MB has become particularly widespread because it is useful for the solution of problems in decision making. The usefulness of this method consists basically in the use in situations where there is limited information about a large number of variables or when the information comes from different sources. The method primarily incorporates prior knowledge to be able to estimate useful models within a master space and thus be able to estimate parameters that come from experience or from a probabilistic theory (Mesa et al., 2011).

Forecasting is a tool that provides a quantitative estimate of the probability of future events (Contreras et al., 2016). The objective of a forecast is to enable decisions about the future and to provide an estimate of the risk involved in a decision (Masterton, 2014).

Forecasting is divided into three standardized methods: qualitative, historical projection and causal. Qualitative methods focus on determining which intangible and subjective factors affect the forecast based on judgment, intuition, surveys or comparative techniques to generate quantitative estimates about the future (Ballou, 2004). Historical projection methods are developed from past data and serve to predict future behavior (Mas-Machuca et al., 2014). Finally, causal methods ‘assume that the factor to be predicted exhibits a cause-effect relationship with one or more independent variables. The purpose of causal models is to describe the form of the relationship between the variables and use it to predict future values of the dependent variable. Among the most commonly used causal methods are regression techniques and econometric techniques’ (Contreras et al., 2016).

The quantitative methods normally used to perform forecasts are: linear regression, moving averages and exponential smoothing, among others. As for the linear regression method, it is an important and useful tool in many statistical analyses to study the relationship between variables, which is mainly used to ‘predict values of the response variable to interesting values of the predictor variables, discover the predictors associated with the response variable, and estimate how changes in the predictor variables affect the response variable’ (Eck, 2018).

To model time series, one can work with traditional statistical models that include moving average, exponential smoothing and ARIMA. These models are linear in that future values are linear functions of past data. Over the past few decades, researchers have focused heavily on linear models as they have demonstrated their simplicity in understanding and application. Time series forecasting models are mainly used to predict demand (Ferbar et al., 2016).

The accuracy of the mentioned models can be determined by comparing the observed values used with the predicted values, likewise they assure (Taylor, 2011) ‘the forecast error tells us how well the model performs when compared to itself using historical data’. Therefore, analyzing the prediction errors is a task that allows a comparison of the forecast models and the selection of the most appropriate model (De Oliveira & Oliveira, 2018). Within the literature there are several procedures to calculate the overall error of a forecast which are used to measure the performance, compare several forecast models and monitor the forecasts. The most common are the mean absolute deviation (MAD), mean absolute percent error (MAPE) and mean squared error (MSE), among others (Heizer & Render, 2004).

Precipitation is a fundamental component of the global water cycle and a key hydrological variable of the water cycle for meteorology, climatology and hydrology. Accurate observations of precipitation and its regional and global distributions have long been scientific challenges. Precipitation forecasts have matured over four decades of development (Frías-Paredes et al., 2018). From a physical point of view, precipitation is defined as ‘all aqueous particles in liquid or solid phase that originate in the atmosphere and fall to the Earth’s surface’ (Shcherbakov et al., 2013).

A fundamental component of the global water cycle, precipitation is the moisture fluxes from the atmosphere to the Earth’s surface, and is a key hydrological variable of the water cycle for meteorology, climatology and hydrology, the atmosphere derives approximately three-quarters of its heat energy from the release of latent heat by precipitation (Jindi et al. 2020). Precipitation varies from year to year and over decades, and changes in amount, intensity, frequency and type affect the environment and society. Steady moderate rainfall penetrates the soil and benefits plants, while the same amount of rainfall in a short period of time can cause flooding, leaving soils much drier (Michaelides et al., 2009).

3. Methods

The methodology for this research is to take the daily data from 1 January 1989 to 31 December 2020, these data are obtained from the Tunguativa weather station located in Paipa a small town in Colombia. Make the forecast for one year (2021) using the ARIMA method and point of change. Take the predicted results for the year 2021 by the two methods and compare them with the observed data produced by the weather station.

The procedure performed for the calculation of the forecast includes a series of six consecutive steps (Figure 1). Once the precipitation data for the time series under study has been cleaned, the first step consists of plotting the data in order to illustrate its behavior in the time window; then, in the second step of the methodology, components of the series, such as trend, seasonality and randomness are identified; the third step is dedicated to the calculation of the forecast by means of ARIMA; In the fourth step, the change points of the series in the mean and variance dimensions are calculated; the fifth step of the methodology uses the results obtained in the previous step to forecast precipitation using the change point method; finally, in the sixth step, the forecasts obtained by the two methods are compared, applying inferential statistical methods, such as dependence tests and residual analysis. In the following section, each of these steps of the methodology used for this research is developed.

Figure 1. Methodology flow chart.

The data used were supplied by a meteorological station in charge of monitoring the study area, which guarantees both scientific rigor and methodology for data collection. The data for this meteorological variable of precipitation is a daily time series compiled from 1989 to 31 December 2020. With this historical daily data, forecasts are made for the year 2021 considering two calculation methods, such as ARIMA and change point.

Precipitation according to Kummerow et al. (1998) is classified depending on the intensity, when it is less than 2 mm/h it is weak, between 2 mm/h and 15 mm/h it is moderate, between 15 mm/h and 30 mm/h it is strong, very strong if it is between 30 mm/h and 60 mm/h and when it is greater than 60 mm/h it is classified as torrential. Taking into account the above, the data according to the data of this research are grouped as illustrated in Figure 2, 69.78% of the data are weak precipitation, 28.31% to moderate precipitation, 1.76% and 0.14% to strong and very strong precipitation, respectively. Since the change point model has proved to be useful for the prediction of stock market variables, hydrology or natural phenomena as in the studies conducted by Michaelides et al. (2009), Kummerow et al. (1998), Trenberth (2011) and Bourgouin (2000), one of the contributions of this work is to verify that the model is robust for the prediction of atmospheric variables.

Figure 2. Grouping of data according to precipitation intensity.

4. Results and discussion

Step 1. Plotting the time series. Plotting a time series allowed visualizing trends of numerical values as a function of time, each point on the graph corresponds to both a time and a quantity that is being measured for this study, precipitation given in millimeters (mm). In Figure 3, it can be established that most of the time where precipitation obtained high values is in the decade from 1995 to 2005, confirming the information provided by the weather station.

Figure 3. Time series precipitation data from 1989 to December 2020.

Step 2 of the methodology identifies the components of the series, such as trend, seasonality and randomness.

The series presents fluctuations movements that occur year after year with more or less the same intensity, this phenomenon allows affirming that the time series of the precipitation variable of the data collected by the meteorological station from 1989 to 2020 is seasonal (Figure 4). It is evident that the trend allows predicting the future states of the series over time where there are versatile trends (increase, decrease and stability of the precipitation coefficient), the behavior of precipitation undergoes changes that are not cyclical, i.e. its behavior is irregular and causes a high degree of randomness.

Figure 4. Components of the time series 1989–2020.

In Step 3, the forecast is calculated using the ARIMA model. Figure 5 illustrates the forecast using the ARIMA model since it uses a weighted average of past and current values (1989–2020) to provide predictions (2021). For this calculation more weight was given to more recent data to increase the importance of this data compared to older data. The time series used contains seasonal variation, and using the ARIMA model is a good decision for this type of series. The ARIMA model can be used to make short-term predictions because most of them place more emphasis on the recent past than on the distant past.

Figure 5. Forecast results for the year 2021 ARIMA method.

In order to calculate the forecasts using the change point model, the change points must first be identified and estimated as proposed in Step 4 of the methodology established for this research.

Most of the change points calculated in the mean of the original time series are reflected in the days where the precipitation variable had values less than 10 mm (Figure 6), while the change points located in the variance take values greater than 10 mm (Figure 7).

Figure 6. Estimated points of change in the mean.

Figure 7. Variance change points estimation.

Once the change points have been estimated, we proceed to Step 5 of the methodology, which consists of calculating the forecast using the change point model. The change point model shows a trend to the mean, that is, it considers the days without rain for its projection, but in the final series it represents the trend, which does not show the days without rain or days of extreme events, which are of great importance for decision making in the sectors that may be affected or benefit from the fluctuations of the variable under study (Figure 8).

Figure 8. Forecast results for the year 2021 using the change point method.

Finally, in Step 6, a comparison of the forecast results obtained with the application of the two models is made. The series with the highest significance were plotted according to Figure 9, where it is evident that the two series represent the same bimodal regime trend, but it can be seen that the change point model does not represent the extreme events, that is it does not represent the dry days or high precipitation, and the ARIMA model takes the extreme values, this information can offer a viable and reasonably accurate solution to the studies of precipitation behavior.

Figure 9. Comparison of observed data 2021 with forecasted data 2021 from January 1 to December 31.

To obtain these results of the forecasts for the year 2021, the data of the time series up to 31 December 2021 were used to analyze the calculated data with the observed data obtained for the meteorological station during the year 2021. To estimate the validity of the results obtained by applying the ARIMA model and the change point model, they were compared with the observed data using statistical tests to calculate the degree of similarity of the forecasts with the real data, as shown in Table 1.

Table 1. Tests of statistical dependence of observed data on predicted data.

Test	Statistic value	Reliability (%)	p Value	Hypothesis statement	Hypothesis not accepted
Fisher	0.27	95	2.2E-16	H0: the variance of the CP forecast is homogeneous to the variance of the observed data.	H0
				H1: The variance of the CP forecast is not homogeneous to the variance of the observed data.
Wilcoxon	51,339	95	2.2E-16	H0: The ARIMA forecast data are dependent on the observed data.	H0
				H1: ARIMA forecast data are independent of observed data.
Wilcoxon	6554	95	2.2E-16	H0: The CP forecast data are dependent on the observed data.	H0
				H1: The CP forecast data are independent of the observed data.
T	4443	95	1.11E-05	H0: the variance of the ARIMA forecast is homogeneous to the variance of the observed data.	H0
				H1: The variance of the ARIMA forecast is not homogeneous to the variance of the observed data.
T	36,602	95	2.20E-16	H0: the variance of the CP forecast is homogeneous to the variance of the observed data.	H0
				H1: The variance of the CP forecast is not homogeneous to the variance of the observed data.

H₀: Null hypothesis ARIMA: ARIMA.

H_1: alternative hypothesis; CP: change point

According to (Frías-Paredes et al., 2018) the prediction of meteorological variables generally implies that the probability will be higher than normal. Such a statement requires knowledge of ‘normal behavior’, it requires a null hypothesis. The hypotheses calculated in Table 1 were tested in three ways: comparing the observed data with the predicted data, comparing the statistical independence of the observed data with the distribution of the predicted data, and comparing the likelihood ratio with the formulation of the null hypotheses.

To obtain a comparison index between the forecasts of the two models, the mean absolute error (MAE) was calculated, which is a robust measure when considering the high volume of data used for forecasting, the MAE for the forecasts obtained with the ARIMA model is 7.92 while the MAE for the results of the change point model was 2.56, which allows us to evidence that the forecast yielded by the latter model is better adjusted to the observed data, corroborating what is illustrated in Figure 9.

Another way to verify the results of the tests in Table 1, is the analysis of residuals as shown in Figures 10 and 11, these graphs show that the Naive method analyzes which forecasts seem to take into account all the available information. The residuals were useful to check that the prediction of precipitation for the year 2021 using the ARIMA and change point models are a good forecast, taking into account that the results comply with the following properties: the mean of the residuals is close to zero and there is no significant correlation in the series. The time plot of the residuals shows that the variance of the residuals remains almost the same throughout the historical data, apart from the outliers and therefore the residual variance can be treated as constant. This can also be seen in the histogram of the residuals, the histogram indicates that the residuals may not be normal, considering that the right tail in Figures 10 and 11 appear too long for the two cases, even when outliers are ignored.

Figure 10. Residuals of data obtained by the ARIMA method.

Figure 11. Residuals of data obtained by the change point method.

5. Conclusions

The results obtained from this research can be a tool that can be taken into account to prevent the impact that rain may have on Lake Sochagota, so that in the rainy seasons the right decisions are made to preserve the integrity of the lake and avoid environmental damage.

In this study, the data obtained from the change-point model show higher accuracy and fit the observed data relatively well. However, it may appear that the optimization procedure plays an important role in the convergence of the model.

The data obtained from the change point model show higher precision and fit the observed data relatively well. However, it may seem that the optimization procedure has an important role in the convergence of the model. Since algorithms are used in this study, it is difficult to state the accuracy of these procedures. Therefore, it is advisable to use several algorithms focused on forecasting weather variables in which better results can be obtained.

When comparing the forecasts obtained with the ARIMA model versus the results of the change point model, it can be verified that the two series represent the same trend of bimodal regime, considering that the study has a high volume of data, a parameter to conclude which of the forecasts is better was the calculation of the MAE, obtaining a value of 2.56 and 7.92 for the change point and ARIMA forecasts, respectively, therefore, it is evident that the forecast calculated through the change point estimation is better adjusted to the observed data.

By applying the change point model to calculate the precipitation forecast for the year 2021, the forecast obtained shows a trend to the mean, that is it considers the days with weak precipitation for its projection, but in the final series it represents the trend, in which the data with weak precipitation and data with strong or very strong precipitation are not visible.

The results of the forecasts obtained through the application of the ARIMA model are adjusted to a lesser extent to the observed data; however, according to the residual analysis, it is considered a good forecast because it was found that there is no relationship between the observed data and the calculated data.

This study is considered preliminary and for the results to be considered conclusive they must be applied to a significant number of time series of meteorological variables. It is also advisable to test the use of algorithms focused on the forecast of weather variables in order to verify that better results are obtained.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Additional information

Notes on contributors

Juan Camilo Valderrama Balaguera

Juan Camilo Valderrama Balaguera, a professional in aeronautical engineering with a postgraduate degree in production and operations management and a master’s degree in industrial engineering. Currently director of the industrial engineering program and member of the LOGyCA research group at the Universidad de Boyacá, Tunja, Colombia. Research areas Logistics, simulation, production and operations planning and scheduling.