Estimation of seasonal boundaries using temperature data: a case of northwest part of Bangladesh

Seasons are the divisions of the year into months or days according to the changes in weather, ecology and the intensity of sunlight in a given region. The temperature cycle plays a major role in de ning the meteorological seasons of the year. This study aims at investigating seasonal boundaries applying harmonic analysis in daily temperature for the duration of 30 years, recorded at six stations from 1988 to 2017, in northwest part of Bangladesh. Year by year harmonic analyses of daily temperature data in each station have been carried out to observe temporal and spatial variations in seasonal lengths. Periodic nature of daily temperature has been investigated employing spectral analysis, and it has been found that the estimated periodicities have higher power densities of the frequencies at 0.0027 and 0.0053 cycles/day. Some other minor periodic natures have also been observed in the analyses. Using the frequencies between 0.0027 to 0.0278 cycles/day, the observed periodicities in spectral analysis, harmonic analyses of minimum and maximum temperatures have found four seasonal boundaries every year in each of the stations. The estimated seasonal boundaries for the region fall between 19-25 February, 19-23 May, 18-20 August and 17-22 November. Since seasonal variability results in imbalance in water, moisture and heat, it has the potential to signi cantly a ect agricultural production. Hence, the seasons and seasonal lengths presented in this research may help the concerned authorities take measures to reduce the risks for crop productivity to face the challenges arise from changing climate. Moreover, the results obtained are likely to contribute in introducing local climate calendar.


Introduction
Bangladesh experiences diverse threats including oods, cyclones, hurricanes, waves, sea level rise, and landslides invariably in various locations, and accordingly, has been underlined as one of the most vulnerable countries in the world towards climate change [30,31,64]. The northwest part, focus of this research, experiences extreme weather, recurrent drought and irregular precipitation [2,5,23], which is also in line with the evidenced a rmation made by the Intergovernmental Panel on Climate Change (IPCC) about global warming and climate change [21,22]. Since extreme dry and drought weather regulates agricultural production, which plays pivotal role in food security and national economy, impact of temperature on weather is much required investigation to explore local climate [1,60]. Climate of the northwest part of Bangladesh and the neighboring regions have undergone for various investigations using temperature and rainfall analyses [3,29,38,40,44,46,51,52,65]. However, seasons in this region in terms of boundary and duration in a year, and the respective variations are underexplored as of yet. This research thus is intended to assess and comprehend the seasons and seasonal variability in temperatures.
Bangladesh is located in the tropical monsoon region and characterized by high temperature, heavy rainfall, often excessive humidity, and fairly marked seasonal variations. While about 68% of the total land areas are used for agriculture and village settlement, approximately 80% of the total populations are directly or indirectly involved in agriculture and related activities [52], although a constant decrease in agricultural employment is evidenced during last 30 years [64]. Agricultural yield in Bangladesh largely depends on climate variability [45]. Along with other environmental and climate parameters, temperature and rainfall are the major visible climate variables that play pivotal roles in agriculture [17,61,62], and hence, require detail investigation in order to assess climatic pattern of the region.
Bangladesh is known as the land of six seasons [8,24], in contrast to the usual four-season concept, having a temperate climate. Within these six seasons, there are three distinct spells: the pre-monsoon hot season from March to May, the rainy monsoon season from June to October, and a cool dry winter season from November to February (Table 1). There are di erent views and multiple ways in separating and characterizing the seasons over the calendar year in Bangladesh. Table 1 shows six, three, four and two di erent seasons and timeframes of the seasons over a calendar year. Seasonal trend analysis is considered as one of the most vital analyses in climate change research [41,55]. Although trend analysis along with estimation of climate variables is frequently found in the literatures studying the climate of Bangladesh, seasons in Bangladesh are yet to be well-de ned both in terms of duration and boundary. Therefore, this research intends to investigate the seasons through identifying the beginning, ending and duration of a particular spell within a calendar year. Although perception and practices of farmers, various adaptation measures, impacts of climate change on agricultural yield are among the top researched areas (For instance, [12,19,25]), empirical evidences about the changes in season and respective duration in local scale require more investigations to provide the agriculture department of any government with loads of information about how the seasons are being shifted.
Estimation of annual seasonal cycle varies over land and ocean due to anthropogenic forcing along with the time derived from equinoxes [15,49,58]. Time series records of climate variables available in meteorological stations are of di erent types. However, the seasons in di erent time series in a yearly cycle should be unique. It is thus rational to expect the seasons in time series of di erent climate variables would be similar. Hence, this work estimates the seasonal boundaries using harmonic analysis, which has already been tested in various applications in the eld of geoscience [6,39,50,54,67]. The extension of harmonic analysis has also been implemented in meteorology and climatology; for instance, VanLoon et al. (1973) [59] used harmonic analysis to describe zonal standing waves in the atmosphere-pressure waves that described the ridges and troughs exhibited by the isobars on weather maps. Legates and Willmott (1990b) [34] employed harmonic analysis to explain seasonal trends in global surface air temperature and showed areas that exhibited a strong seasonal cycle in air temperature at the middle and upper latitudes from areas where air temperature had a strong biannual cycle. Successful applications of harmonic analysis have also been observed in describing rainfall pattern [32,33,66].
Fourier or harmonic analysis can decompose a time-dependent periodic function into a series of oscillatory components and can rebuild the function using the oscillatory components de ned as unique amplitude and phase values [7]. The proportion of variance in the original time series data set accounted for by each term of the harmonics can also be calculated [26]. Justino et al. (2010) [27] have shown the meaningful application of harmonic analysis in seasonal variability, and the authors have estimated the largest seasonal variability over the Weddell and Ross seas along with the amplitude of the seasonal cycle over the West Antarctic ice sheet. In most of the contributions the analysis has been re ected as a viable tool for investigating seasonal cycles of temperature and the uctuation of annual and daily rainfall [11,31,68]. Accordingly, this research implemented harmonic analysis over daily temperature time series in order to show the boundaries and dura- Hossain, 2019 [20] Apr to Oct (spring, summer and fall) Nov to Mar

3
Pre-monsoon season Monsoon season Dry winter season Rakib, 2018 [46] Mar to May Jun to Oct Nov to Feb tion of the season in days within a year, and the corresponding variability of seasons in Rajshahi and Rangpur, the two northern divisions in Bangladesh.

Materials and methods . Study area
The study area is located in the northwestern part of Bangladesh that extends from 23

. Data
Daily minimum and maximum temperature data recorded at six stations for 30 years from 1988 to 2017 were obtained from Bangladesh Meteorological Department. Obtained data shows a very minimal record of overall missing data, less than 1%, which indicates an acceptable level of data quality. However, the missing data obtained for Ishwardi station were slightly more. Missing records were estimated using the smoothing spline technique [14]. It is to note that daily minimum and maximum surface temperature have been used for analyzing annual cycle of surface temperature [36].

. Harmonic analysis
Harmonic analysis decomposes a signal into an in nite series of harmonic components, and comprises initially of a sine wave and a cosine wave of equal frequency. These two waves are then combined into a single cosine wave, which has characteristic amplitude and phase angle. Harmonic analysis lets a complex curve to be expressed as the sum of a Fourier series of cosine waves and an additive term [13,47]. Successive harmonic terms are added to produce a complex curve, and each component curve or term accounts for a percentage of the total variance in the original time series data set.
A continuous function f (x) in the interval [ , l] according to Fourier series can be written as [26]: where, the term a , an and bn are Fourier coe cients. The coe cients a , an and bn can be de ned as followed.
Assuming the ith harmonic to the ith term in the Fourier series for i ≥ , the terms can be written as: The ith harmonic can be converted to a single cosine term and represent it as follows: where, c i = a i + b i is the modulus or length of the vector a i , b i in the two dimensional plane, is the amplitude and ϕ i = tan − b i a i is the angle between the vector a i , b i . The reason is that the inverse tangent function only returns values in the interval − π , π , whenever a i < the modi ed de nition ϕ i = tan − b i a i + π must be used to obtain a true phase angle. The phase angle follows that ϕ i ∈ − π , π . If the Fourier coe cient c = a , then the continuous function f (x) takes the following form: where, the parameter cn is the amplitude, and ϕn is the phase angle of the nth harmonic.
Since can compare the quality of models, outputs from the model are justi ed with Akaike Information Criterion (AIC), formulated by the Japanese statistician Akaike Hirotsugu in the early 1970 [4]. It is one of the best suitable techniques to select better model using AIC: AIC may over t, particularly when sample size becomes small. It is extended further with nite correction and is said AICc, as formulated below: Where, the parameter n and K denote the sample size and the number of free parameters in the model, respectively. When n → ∞, the extra penalty term converges to 0, and thus AICc converges to AIC. Bayesian Information Criterion (BIC) is formally de ned as-BIC = − (ln (likelihood)) + ln (n) K

Temperature data analysis
Harmonic analysis of daily temperature data and determination of seasonal boundaries in the analysis for a year are explained in this section, while the detail results from harmonic analyses using data from di erent stations have been provided in the following section. Before the analysis, the periodicity or the periodic nature of data is examined using spectral analysis. Figure 2(a) shows the one-year time series of daily minimum temperature. The power spectrum of the daily minimum temperature time series for 30 years from 1988-2017 has been estimated using Fourier transformation for Rajshahi station as shown in Figure 2(b), which describes the periodic nature of the data. There are only two higher peaks observed in the analysis at 0.002732 cycles/day and 0.005464 cycles/day, and the peaks are interpreted as the period of 366 days and 183 days, respectively. However, there are many small peaks of frequencies all over the ranges from greater than 0 to more than 1 cycles/day; however peaks in the range from 0 to 0.05 cycles/day only are shown in Figure 2(b). This work likes to visualize the period greater than a month over the data. The periods equal to 36 days or more, equivalent frequency 0.0278 cycles/day or less, are considered in this work. Harmonic analysis is conducted for the periods from 366 days to 36 days or alternatively the frequencies from 0.0027 cycles/day to 0.0278 cycles/day. Used frequency range in harmonic analysis lies between the judged frequencies as shown in Figure 2(b) Harmonics are the waves or signals of individual frequencies in the time series. It is assumed that all successive frequencies together with appropriate amplitude and phase can reproduce the real time series. Frequency of the harmonics can be de ned, as shown in Equation (6), n/l, where, l is the total number of days in the year and n = 1, 2, 3. . . 365 or 366. This way the total harmonics can be employed for the frequencies 1/366 (=0.002732) to 366/366 (=1) cycles/day with respect to the number of days in the year. On other hand, estimated spectral analysis of the climate time series has gured that the power densities are high at adjacent values to the frequencies of 0.002732 cycles/day and 0.005464 cycles/day (Figure 2(b)). The remaining frequencies are of very low power densities. Therefore, setting the frequency range from 0.0027 cycles/day to 0.0278 cycles/day for the harmonic analysis of climate time series appears sensible to reconstruct the basic time series as shown in Figure 3(a) with the error, and to visualize the seasons over few harmonics in a year. Figures 3(a-c) show the harmonic analysis over daily minimum temperature data recorded at Rajshahi station in the year 1988 for studying the seasonal boundaries. Figure 3(a) shows the reconstructed daily temperature along with original records using 366 harmonics. Estimations of total of 10 successive harmonics (n = 1,2,3. . . 10) are to be applied to de ne the seasonal boundaries. Among 10, rst ve individual estimated The lower values of AIC and modified AICc or BIC indicate better model fit (the estimation using 5H is better than 1H; mean square error, MSE is also indicating the same having smaller at 5H). According to estimated MSE, the uses of the estimated harmonics, 1H to 5H are reasonable to define seasonal boundaries as estimated four boundaries keeping common 1H with others.
harmonics are shown in Figure 3(b). Cumulative or successive ve harmonics are shown in Figure 3(c) and marked as 1H, 2H, 3H, 4H and 5H. In Figure 3(c), the days in which successive harmonics are same or close to the same value of amplitude are seasonal boundaries, and marked as B1, B2, B3 and B4. Table 2 shows a detail understanding of harmonic analysis to help determine the boundaries in the year. AICs were estimated for each successive harmonics, and it has found that MSE and AICc are lower with the increasing number of harmonics (Table 2). Figure 3(c) shows the boundaries indicated with an arrow mark, termed as B1, B2, B3 and B4.

Results and discussions
First ten, instead of ve as shown in Figure 3(c) above, successive harmonics of temperature time series recorded at six stations were estimated year by year using the equations (1-6) as described earlier. The boundary of a season is the time when the amplitude of the harmonics becomes identical, and in the yearly analyses four boundaries are obtained. Estimated yearly boundaries termed as B1, B2, B3 and B4 in minimum and maximum temperature with harmonic analysis are basically in days. Estimations of such seasonal boundaries and the respective conversion to calendar date in minimum and maximum temperature over time are shown in Figures 4-5. The rst, second, third and fourth boundaries vary between 49 th and 54 th day, 137 th and 146 th day, 240 th and 246 th day, and 325 th and 338 th day, respectively, throughout the time series. Variation in estimated boundaries within a year is trivial, however boundaries along the year have little ± trends.
For Rajshahi station, rst, second, third and fourth boundaries, marked B1, B2, B3 and B4 in the analysis using minimum temperature, are obtained at 57 th , 143 rd , 243 rd and 332 nd day, respectively. The same boundaries are obtained at 54 th , 144 th , 241 st and 336 th day in the analysis using maximum temperature. Estimations of seasonal boundaries using both minimum and maximum temperature are found to be very close.
For Ishwardi station, rst, second, third and fourth boundaries are obtained at 56 th , 142 nd , 241 st and 329 th day, and at 54 th , 144 th , 241 st and 337 th day using minimum and maximum temperature, respectively.
For Bogura station, rst, second, third and fourth boundaries are obtained at 51 st , 137 th , 240 th and 331 st day, and at 54 th , 146 th , 246 th and 338 th day using minimum and maximum temperature, respectively.
For Rangpur station, rst, second, third and fourth boundaries are obtained at 52 nd , 140 th , 240 th and 330 th , and at 50 th , 142 nd , 241 st and 335 th day using minimum and maximum temperature, respectively.
For Syedpur station, rst, second, third and fourth boundaries are obtained at 49 th , 138 th , 240 th and 329 th , and at 50 th , 143 rd , 241 st and 335 th day using minimum and maximum temperature respectively.
For Dinajpur station, rst, second, third and fourth boundaries are obtained at 49 th , 137 th , 238 th and 325 th , and at 51 st , 142 nd , 243 rd and 338 th day using minimum and maximum temperature, respectively.
Identi ed seasons at six stations over 30 years appear to be almost unvarying. Thus, it may be argued that the study locations belong to an identical climate zone that has comparable four di erent seasons.
The boundaries estimated for di erent stations appear to start little earlier in the northern stations than that of the south. Overall boundaries estimated in calendar date using harmonic analysis in daily minimum and maximum temperature in di erent stations are shown in Figure 6 and Table 3. Spatial variation in the boundaries of the seasons in the study area is found to be very small and varying from 1 to 4 days.

N-number, C-calendar
Estimated duration of the four seasons are mostly March-May, June-August, September-November and December-February, and the ndings are in line with the seasons as argued by Shahid (2010a, b) [51,52] and Syed and Amin (2016) [53]. Seasonality may not be uniquely de ned due to di erent methods employed to isolate seasonal components [42]. The traditional approach assumes a periodic response. Franses (1996) [16] and Lund et al. (1995) [35] used periodic correlation to di erentiate the seasonal variation. This work is also   the approach of periodicity analysis of climate time series that has used the phases with the same amplitude of the estimated harmonics, and identi ed four seasons in the northwest region of Bangladesh. The estimated  [9,57]. Along with summer monsoon four meteorological seasons in Bangladesh are pre-monsoon (March, April and May), monsoon (June to September), post-monsoon (October and November) and winter (December, January and February) [24,28]. This work has also visualized similar four seasons with an average seasonal length of 90.5 days (from 82 to 101 days) as shown in Table 3. However, the seasonal lengths found are not uniform. Length of the fourth season appears 84 days, while the second season appears 100 days and the remaining lengths fall between 84 to 100 days. Little early-start of season and variation in seasonal length coupled with observed temperature pattern provide the evidences of change in local climate and weather.
Higher temperatures along with changing rainfall, frequent ood and drought, rising salinity are likely to reduce crop yields and crop production, and eventually income of the farmers [18,37]. During the last decades, researchers around the world have shown massive interests about adaptation interventions in agriculture in the context of climate change. IPCC (2014) [22] estimated that, by 2050, rice and wheat production in Bangladesh could decline by 8 and 32%, respectively (against a base year of 1990). Since adaptation is a priority to address the adverse impacts of climate change, agriculture-based economy like Bangladesh has started to intervene various adaptation measures including crop diversi cation, change in harvesting time, salinetolerant and drought-resistance varieties, as well as combined crop-livestock farming [10,43,56]. Among these interventions, change in harvesting time requires special attention, since core agriculture is primarily season-based. Harvesting should be done according to the shift in seasons.
National policies and practices related to farming and agriculture may be ne-tuned over long-term based on analyses this study attempted. For example, if seed delivery time arranged by the agricultural extensions services was the rst week of March 30 years back, change in seasonality may require the service to be shifted according to the actual change in season after 30 years. Besides, if specialized nancial institutions provide any facility to the farmers, the time for planning the fund requirement would change with the shift in season. For many of the developing countries, agriculture is the core contributor to gross domestic product. If it is not practiced considering seasonality, agricultural output may decline which would result in reduced national income.

Conclusion
Investigation about seasonal variability is very common in climate change research. Understanding seasons along with the variations is important to comprehend meaningful climate of a region. This work decomposed the temperature time series into successive harmonics, and harmonics had also reconstructed the time series. Using decomposed temperature data harmonics seasonal boundaries were estimated. Four boundaries were estimated in each of the year from 1988 to 2017 and presented in this article. First, second, third and fourth seasonal boundaries were visualized within 19-25 February, 19-23 May, 18-20 August and 17-22 November, respectively. Hence, apparent seasons found were from March to May, from June to August, from September to November and from December to February. Results have showed that the seasonal boundaries in a year are varying from 2 days to 6 days in a 30-year time frame. Similar research using data from all stations in the country would help in developing a national climate calendar to assist policymaker.
Further observations can also be noted that there were little trends close to a slope of ±0.05 of the estimated boundary in a year; however, 30-year long time series is not su cient enough to make an error-proof trend over the boundary. For reasonable trend analysis, the length of the time series should be of 100 years or more. Nonetheless, in absence of data for a longer time span in Bangladesh context, this study is likely to contribute to the formation of local climate calendar, a near-future necessity of agricultural extensions services provided by the government. The estimation of realistic season consequently may assist in increased agricultural and industrial productivity.