Overview

Seemingly, extreme weather events are becoming more severe and frequent, with their respective devastating impacts; in addition, perhaps they are caused by climate change. These two statements require that statistical analyses oriented to their prediction and forecast be based on quality meteorological and climatic data, that is, that they be accurate and homogeneous. In fact, obtaining perfectly homogeneous records or large series of climatic data is almost impossible, due to the inevitable change of the area surrounding the site of the climatological station, which affects the observed data. The loss of homogeneity occurs when there are changes in the record generated by non-climatic causes (^{Guijarro, 2014}; ^{Yozgatligil & Yazici, 2016}).

In summary, all records or data series that are measured at climatological stations often undergo alterations that are not related to climate variations. Such modifications arise from changes in the location of the station, by replacing its instruments or changing the measurement technique and by altering the physical conditions surrounding the site. Homogenization is a procedure to detect and correct the aforementioned artificial alterations. Regarding homogenization techniques, ^{Guijarro (2014)} briefly describes them and focuses on the available computing packages. A series of climatological data is considered homogeneous, when the measurement conditions at the station have not changed over time (^{Beaulieu et al., 2008}). ^{Guijarro (2014)} highlights the importance of the Meteorological Services in avoiding by all means, changes that affect the climatological measurements, since this makes the use of homogenization methods unnecessary.

The volcanic eruptions, being a natural phenomenon that causes an alteration of climate records, define change points that do not imply loss of homogeneity. On the contrary, break points are change points caused by non-climatic alterations and therefore imply loss of homogeneity (^{Toreti et al., 2011}).

The homogenization of the climatic series is a difficult task that must be carried out with extreme care, especially when additional data or historical data (metadata) of the climatological station are not available; common case. The objective of homogenization is to eliminate, or at least reduce, non-climatic alterations, while the climatic signal is preserved (^{Toreti et al., 2011}). Many methods have been proposed to homogenize the climatological series, a complete review that includes suggestions for their use can be found in ^{Peterson et al. (1998)}, and in ^{Beaulieu, Ouarda and Seidou (2007)}.

A first classification of the procedures and/or statistical tests that detect and correct the loss of homogeneity, establishes two groups: the first is based on the availability of additional data on the historical evolution of the climatological station, known as metadata and the second in the absence of such data. ^{Rhoades and Salinger (1993)} present techniques for the first group, under two approaches: (1) using neighboring stations and (2) in stations isolated or from old or primigenial records. In general, the historical data of the climatological station are essential to validate the changes detected by the tests and unfortunately, they are almost always not available.

Another classification of homogenization techniques divides them into: absolute tests, which exclusively use the series or record under analysis and relative tests, which use close or correlated records to form the reference series. Both approaches are valid and useful, and both also exhibit disadvantages. ^{Yozgatligil and Yazici (2016)} compare seven tests of each type and find that the relative tests have better performance, among them they cite first the standard normal homogeneity test, designated SNHT.

^{Syrakova and Stefanova (2009)} indicate that due to the difficulty in establishing a priori which stations are homogeneous, they choose to follow an iterative scheme of gradual homogeneity improvement, which was proposed by ^{Hanssen-Bauer and Forland (1994)}, and has been also applied by ^{Moberg and Alexandersson (1997)}; ^{González-Rouco, Jiménez, Quesada, and Valero (2001)}, and ^{Tuomenvirta (2001)}.

Instead, ^{Wijngaard, Klein Tank and Können (2003)} apply and cite two advantages of the absolute tests: (1) easy application in sparse networks and (2) superior to relative tests in networks with simultaneous changes, since these do not detect them. ^{Dhorde and Zarenistanak (2013)} use a mixed approach in the study of homogeneity of 20 climatological stations, applying the Pettitt and SNHT tests to each series and the SNHT test for gradual changes in the mean with reference series.

^{Ducré-Robitaille, Vincent and Boulet (2003)} evaluated 8 homogenization techniques of temperature series, using homogeneous simulated series and containing one or several changes in the mean (steps). They found that most have a good performance, but two of them tend to be more efficient: the SNHT test and the multiple linear regression technique.

^{Beaulieu et al. (2008)} also compared eight statistical tests to detect lack of homogeneity in annual precipitation series, using several thousand synthetic series both homogeneous and non-homogeneous and found that no method is efficient in all types of loss of homogeneity, but three of them have a much better perform: the bivariate test, the Jaruskova method and the SNHT test.

In Mexico, and specifically for the state of Veracruz, ^{Guajardo-Panes, Granados-Ramírez, Sánchez-Cohen, Díaz-Padilla and Barbosa-Moreno (2017)}, have applied the absolute approach with the SNHT, Pettitt and Buishand tests, under a spatial verification scheme, to select neighboring climatological stations to the station being tested, within the same climatic unit, defined this based on the variation interval of the climate element being analyzed. The above is equivalent to processing climate information by sub-regions or climatic zones.

Quality control concept

^{González-Rouco et al. (2001)} have defined that the quality control processes of the climatic data series involve three stages: (1) detection and correction of outliers; (2) interpolation techniques of missing data and (3) homogenization of the series. In fact, the estimation of missing data occurs at monthly level, before and after homogenization, when trying to complete series to have a common information period. Due to the fact that all climatic data series are extremely sensitive to the presence of erroneous values and scattered data, ^{Eischeid, Baker, Karl and Diaz (1995)} address the detection and correction of scattered (outlier) data according to an objective approach of temporary verification in the analyzed and spatial record based on six interpolation techniques of monthly data.

Objective

A practical approach or strategy for verifying homogeneity in annual precipitation records has been formulated, which has the following relevant characteristics: (1) work by geographic regions or sub-regions, using monthly data; (2) missing monthly data are estimated using a statistical technique; (3) details how extreme maximum monthly values are verified and corrected, by truncation; (4) missing annual data are deduced, with a regional weighting technique; (5) an absolute contrast approach is used, using four statistical tests (Von Neumann, Pettitt, Buishand and SNHT) and (6) each record that was not homogeneous or reliable is checked graphically and the loss of homogeneity is corrected when there is a breaking point. In this study we worked on the semi-arid climatic region of the Potosino Plateau of Mexico, processing 16 long-term records, whose number of annual data varied from 42 to 53.

Correction of outliers

Outliers are values that are far from the trend shown in a set of data, which correspond to erroneous measurements or extreme meteorological events. There are several approaches to contrast and verify the temporal and/or spatial variability of the climatic variable studied, in order to identify the outliers and diagnose when they are erroneous data or feasible values to occur. When the outliers are undoubtedly wrong observations, such data is deleted and then there is a problem of missing data. When they are feasible values, it must be decided whether they are corrected and how it should be done (^{Eischeid et al., 1995}; ^{González-Rouco et al., 2001}).

Due to the practical difficulties involved in the verification of outliers in the original or field records, an approach to reducing the scattering in the right tail of the probability distribution was adopted, to seek a better performance of homogenization techniques that are not resistant (^{Lanzante, 1996}). Such approach is to limit the outliers to an extreme threshold magnitude (P _vd ), which is (^{Eischeid et al., 1995}):

where, P ₅₀ is the median or the percentile of the 50%, P ₇₅ and P ₂₅ are the percentiles of the 75% and 25% and FM is the interquartile range multiplication factor (RIC), with an empirical value adopted of 2.75 for average temperature data and 4.00 for precipitation, both monthly. ^{González-Rouco et al. (2001)} use a similar expression for their monthly precipitation data:

The outliers superior to P _vd are replaced by such limit. This truncation approach of the outliers reduces the bias that they originate and also maintains information relative to the extreme values. ^{Eischeid et al. (1995)}) indicate that RIC is used in quality control processes because it is resistant (^{Lanzante, 1996}) to outliers. In this study the criterion defined by equation 1 was adopted, having a greater theoretical and empirical justification.

Absolute Approach Tests

The tests that are applied to the series studied are those of the ^{Wijngaard et al. (2003)} who uses four methods to test the loss of homogeneity, these are: the Von Neumann test, the Pettitt range test, the Buishand range test and the SNHT test. These four tests assume that under the null hypothesis the annual values of the series (Y _i ) of the variable Y that are studied, are independent and are identically distributed. Under the alternative hypothesis, the first test assumes that the series is not randomly distributed and the remaining three consider that there is a change in the mean or break point. Under the alternative hypothesis -the first test-, it is assumed that the series is not randomly distributed and the remaining three consider that there is a change in the mean or break point. Logically, the homogeneity or consistency of the series implies that their data belong to a single population that has an invariant mean over time (^{Machiwal & Jha, 2012}).

The Pettitt, Buishand and SNHT tests can locate the year where the break point probably occurs and, in general, they are different in their approach or conception. It has been found that the first two are more sensitive to the break points that occur around the middle of the series; on the contrary, the SNHT test detects changes more easily at the start and end of series.

Von Neumann, Buishand and SNHT tests assume that the Y _i data have normal distribution; on the other hand, Pettitt's test does not require such consideration, since it is based on the ranges and not on the values of the series, which also makes it less sensitive to outliers. The Von Neumann test is considered complementary to the other three, because it detects loss of homogeneity due to causes other than the change in the mean (^{Wijngaard et al., 2003}).

Classification of the series tested

The classification depends on the number of tests that rejected the null hypothesis, defining three categories whose designations according to ^{Wijngaard et al. (2003)}, and ^{Dhorde and Zarenistanak (2013)} are: useful, doubtful and suspect and according to ^{Guajardo-Panes et al. (2017)} are reliable, moderately reliable and unreliable. In this study the following designation is adopted: in class 1 or “reliable” the series that were homogeneously detected by the four tests are grouped; or, only one test detected loss of homogeneity; class 2 or “less reliable” includes series with two tests that detect a lack of homogeneity and class 3 or “unreliable” are series with three or four tests that reject the null hypothesis.

With the class 1 series, trend analysis and climate variability studies can be performed; with class 2, the results of the aforementioned analyses should be kept in reserve, but they are feasible if there are no other reliable series available and the class 3 series should not be used for the aforementioned estimates.

There is a notable difference in establishing the previous classification, since ^{Wijngaard et al. (2003)} define a level of significance (α) in the tests of 1% and instead, in the other two works mentioned 5% is used. In this study, α = 5% was accepted, to coincide with the recommendations of the SNHT test (^{Alexandersson, 1986}).

Selection of climatological stations

This study used the monthly climatological information from the Excel file provided by the San Luis Potosí Local Agency of the National Water Commission (CONAGUA), which is grouped into the three geographical areas of the state: Potosino Plateau, Middle Zone and Huasteca Region. Processed data records correspond to the climatic variable annual precipitation in millimeters, each data integrated by the sum of the twelve monthly rainfall values.

The Potosino Plateau archive covers more than 60 climatological stations, but some are repeated because they belong to the extinct Ministry of Water Resources and the National Meteorological Service. From such availability, stations with less than 40 years of registration and lapses of missing years were eliminated. With such restrictions a sample of 16 stations was obtained whose initial year of their records is indicated with AI in Table 2 and Table 3 of data to be processed. All records cover until 2016, with the exception of the Los Pilares station that was suspended in 2009.

To establish the longest common period of data, the year 1964 was selected for its start and with it, 13 records processed have 53 or 52 values and three: La Presa, El Grito and Los Pilares, cover 42, 48 and 45 years, respectively. In the Moctezuma station, the seven missing years from 1964 to 1970 were taken from Climate Bulletin No. 3 (^{SARH, 1980}), with data on the beginning of the stations in 1978. These values have been indicated at the beginning of the Table 3 with key parenthesis. In the aforementioned bulletin, the beginnings of the La Presa and El Grito stations were verified.

Figure 1 shows the geographical location of the 16 processed climatological stations of the Potosino Plateau, Mexico, high lighting that the first four stations form a north-extreme subgroup, the next three make up a northern subgroup and then the five that continue define the subgroup of the center. Finally, the last four establish the southern subgroup.

Figure 1 Geographic location of the 16 processed climatological stations of the Potosino Plateau, Mexico. 

Estimation of missing monthly data

In this study, the respective monthly mode, estimated based on the fitting of the Mixed Gamma distribution of two parameters (^{Campos-Aranda, 2005}), was adopted as the missing monthly value of the incomplete years, under three specific conditions: (1) it is estimated only in the rainy semester, that is, from May to October; (2) a maximum of three months is completed in each incomplete year and (3) when the mode does not exist, because the distribution is of the type J inverted, monthly precipitation is used with a probability of not exceeding 25% (P ₂₅); October case at Charcas station (see Table 4). In Table 2 and Table 3, the annual values completed with the previous procedure have been indicated in circular parenthesis.

Table 4 shows results of the monthly fitting, only in two months of each incomplete record, due to space limitation. From Table 4 it follows that only four pluviometric stations had complete records: Vanegas, La Maroma, Los Pilares and Los Filtros.

Correction of maximum monthly values

The Mixed Gamma distribution fitting allowed estimating the elements of equation 1, to calculate the maximum monthly dispersed value (P _vd ), which is cited in Table 4. The maximum values that should be truncated to the P _vd value due to exceeding it are also indicated in circular parenthesis. The respective annual values corrected for the truncation described stand out among asterisks in Table 2 and Table 3. Due to the adopted process of deduction of missing monthly data, only extreme maximum values were corrected in the rainy semester from May to October, but erroneous values were sought in the rest of the months, to eliminate them.

Estimation of missing annual data

The pluviometric stations: Matehuala, Charcas, Palo Blanco and El Mezquite presented missing annual data; the first three, the last two and the rest a value. By observing in Table 2 and Table 3 that the annual values available at the surrounding stations of each incomplete record were quite similar, it was decided to deduct each missing annual value based on the three values of the nearest stations, according to a weighted factor with the inverse of the distance in kilometers, between each neighboring and incomplete stations. This procedure has been described and applied by ^{Eischeid et al. (1995)}, ^{Mishra and Singh (2008)}, and by ^{Campos-Aranda (2013)}, it is suitable in areas of plain topography and for few missing data. In Table 2 and Table 3, the seven annual values estimated with the described approach have been cited in rectangular parenthesis.

When several missing annual data are deducted at climatological stations in mountainous areas, the normal-ratio method should be used (^{Paulhus & Kohler, 1952}), applying defined weighting factors based on the correlation coefficients between annual values of the incomplete station and neighboring stations, as defined by ^{Young (1992)} and applied by ^{Eischeid et al. (1995)}.

The numerical concepts and techniques included in the proposed strategy are cited in Figure 2, for the verification of the homogeneity in annual precipitation records of a given geographical area.

Figure 2 Block diagram of statistical concepts and techniques that includes the suggested strategy to test and achieve homogeneity of climatological records. 

Statistical Parameters of records

In the last ten lines of Table 2 and Table 3 the statistical indicators relative to each completed record are shown. With AI the year of beginning of each annual precipitation series is indicated whose number of data, within the established common period, is n. After the median (Med) of the record, the coefficients of variation (Cv), of asymmetry (Cs) and kurtosis (Ck), are cited as well as the serial correlation of order one (r ₁).

Verification of normality

The classic test of ^{Shapiro and Wilk (1965)} was applied to the annual precipitation records with less than 50 data and that of ^{Shapiro (1998)} in those with more than 50 values. It was found that the records of the six following stations: Santa Ma. del Refugio, Palo Blanco, Reforma, El Mezquite, El Peaje and Villa de Arriaga do not come from a normal distribution.

To these six records the Geary Quotient test was applied (^{Machiwal & Jha, 2012}; ^{Campos-Aranda, 2015}) and with this it was confirmed that the first four records do not come from the Normal distribution. It is observed in Table 2 and Table 3, that the asymmetry and kurtosis coefficients of these four records are the ones that differ most from the normal values of zero and three. Such records will be processed using the natural logarithm of their data, to make them normal.

Critical values of statistical tests

Table 1 shows that the critical values of the N, P _k , R/(n s) and T _k statistics of each test are a function of the size of the record (n) and how this varies from 42 to 53, in Table 5 their respective interpolated magnitudes are exposed.

Results of statistical tests

In Table 6 the maximum quantities that were obtained in each test for their indicator have been concentrated (N, P _k , R/(ns) and T _k ) and whose contrast against the values in Table 5, allowed us to indicate whether each record resulted homogeneous (H) or non-homogeneous (NH). It was found that 11 records are reliable, that is, 68.8% of the annual precipitation series analyzed. Three are less reliable (Vanegas, Palo Blanco and Villa de Arriaga) and two are unreliable (Moctezuma and El Mezquite). For these last two stations the Pettitt test provides its maximum values in the 16 records processed, 236 and 234.

Analysis of non-homogeneous records

With respect to the less reliable records of the Vanegas, Palo Blanco and Villa de Arriaga stations, it is likely that the Von Neumann test detects non-randomness due to the persistence shown by such series, whose serial correlation coefficients of order one (r ₁) were 0.241, 0.347 and 0.489; that is, the largest ones of Table 2 and Table 3, together with the largest of all (r ₁ = 0.587) of the El Mezquite station. Anderson's test (^{WMO, 1971}; ^{Linsley, Kohler, & Paulhus, 1988}) was applied to the records of the three stations mentioned, with the result that their coefficients r ₁ are significant.

Figure 3 shows the graphs of the statistics of Pettitt, Buishand and SNHT tests, applied to the Moctezuma station record. It is observed that the first two define maximum positive and negative extreme values and that the SNHT test does define a positive maximum at k = 45 (year 2008). Therefore, there is no single break point, according to the three tests; but two of them define it in k = 45. In fact, in 2008, an extremely dry period of seven years begins, whose presence was ratified in such an area, for a regional drought that affected the entire Mexican Highlands. Therefore, there is no application of equation 10 in k = 45 to correct the prior record.

Figure 3 Graphs of the indicated statistical tests applied to the annual precipitation record of the Moctezuma station of the Potosino Plateau, Mexico 

On the other hand, Figure 4 shows the graphs of the Pettitt, Buishand and SNHT tests, applied to the record of the El Mezquite station. The first two define maximum positive and negative extreme values and the third a single positive maximum, in k = 37 (year 2000). Accepting such year as a break point, a period from 2000 to 2016 is defined with 17 values and a mean of 579.365 mm after the change and 36 data and a mean of 327.811 mm before the change. With equation 10 a value of 1.7674 is obtained as the correction factor (f).

Figure 4 Graphs of the indicated statistical tests applied to the annual precipitation record of the El Mezquite station of the Potosino Plateau, Mexico.  

In the corrected record of the El Mezquite station it has the following statistical parameters: minimum = 59.9 mm; maximum = 1344.1 mm; X- = 579.4 mm; median = 567.3 mm; Cv = 0.480; Cs = 0.467, Ck = 3.486, and r ₁ = 0.503. In addition, according to the Shapiro and Geary tests, it turned out that this may come from a normal distribution. The new maximum values of the Von Neumann, Pettitt, Buishand and SNHT test statistics were: N = 1.016; P _k = 182.0; R/(ns) = 2.03, and T _k = 8.02. Therefore, it is not homogeneous according to the first and third tests, becoming a less reliable record, with the correction made.