Assessing seismicity in Bangladesh: an application of Guttenberg-Richter relationship and spectral analysis

Abstract

Bangladesh has a high risk of earthquakes because the Dauki, Jamuna, and Chittagong-Myanmar faults are still active. However, the assessment of seismicity remains a big challenge due to the complex geologic setting of Bangladesh. This study employed the Guttenberg-Richter relationship and the spectral models to assess and analyze the earthquake conditions in Bangladesh. Besides, an instrumental earthquake catalogue, obtained from the Bangladesh Meteorological Department (BMD), covering 1985–2017, is established. The results revealed that the Guttenberg-Richter constants of a and b were 2.981 and 0.392, which propagated a strain release from 1992 to 2017. The spectral model analyses, e.g. wavelet transform (WT), short-time Fourier transformation (STFT), and multitaper model (MTM), demonstrated the magnitude and strain release anomalies of the same magnitude ranging from 4.8 to 5.7, indicating the probable precursor of an upcoming earthquake. Notably, magnitudes have been running around 4.5–5.8, which may act as a signal to major earthquakes that have not been evident before. The proposed models allowed for the completion of the Bangladesh earthquake catalogue and provided a platform for future seismicity assessment and earthquake probability analysis. These results should be considered in determining how likely earthquakes are to happen in an area or region.

Keywords:

• Seismicity
• strain release
• Guttenberg-Richter relationship
• wavelet transformation
• earthquake probability

1. Introduction

Bangladesh is located where two active Indo-Burmise tectonic plates meet. When these plates move or push against each other, they can cause earthquakes (Steckler et al. 2016; Rahman et al. 2020). Therefore, numerous destructive massive earthquakes have been well-recognized around Bangladesh (Alam 2019; Tabassum and Ansary 2020). Though the country has not witnessed any major earthquakes in recent decades, the probability of earthquakes cannot be underestimated (Haque et al. 2020; Kamal et al. 2021). Scientists anticipated possible earthquakes in and around Bangladesh by interpreting the global occurrence of earthquakes, particularly at the plate margins of India and Eurasia, next to Bangladesh (Steckler et al. 2008; Kayal et al. 2012; Rahman et al. 2018).

A large earthquake can devastate the country, considering its dense population and poor infrastructure. Therefore, the seismicity assessment is further needed to identify the possible vulnerability. However, seismicity assessment in Bangladesh faces some problems. The inadequate opportunity of earthquake data from the past centuries and in the new era has the same data incompleteness problem. More surveys and preliminary research must be needed to assess the entire seismicity and its possible threats. The obsolete data and related sources form a barrier to accessing the seismicity problem. Finally, insufficient drilling to identify seismicity is also a big issue for conducting research. Though some approaches have been made to drilling, completing a rigorous seismicity assessment is inadequate. In addition, seismic catalogues are the essential products that act as the primary source for most studies related to seismicity. Considering historical seismicity and the lack of instrumental seismicity measurements and knowledge of tectonically active features in the country, it is critical to assess the statistical likelihood of severe earthquakes occurring in and around Bangladesh using the most relevant techniques.

Enormous bodies of literature have been written based on observational relationships between earthquake frequency and magnitude, including Gutenberg and Richter (1944), Kanamori (1983), and Bilham (2009). Such relationships have been employed to estimate earthquake recurrence intervals and map earthquake-prone zones’ spatial patterns. The Guttenberg-Richter Law (G-R) is one such empirical association between the magnitude (M) and the number of events (N) with the same magnitude, also known as the Frequency-Magnitude Relationship (FMR) (Asfahani and Darawcheh 2017). The studies showed the G-R FMR holds for virtually all magnitude ranges at all sites and all times. Ogata (1988) investigated the G-R models for earthquakes and provided a clear view of their occurrences. The gamma, log-normal, Weibull, and exponential distributions have also been applied by Utsu (1984) to describe the probability distribution of the interference time of massive earthquakes using G-R and show its effectiveness. However, the G-R model consists of two seismicity parameters, a and b, which vary widely with space. It is essential to determine their values for a region for a reliable assessment of seismicity in the area.

The spatial definition of a particular area with various attributes is defined by spectral modeling analysis. Therefore, the earthquake magnitude and strain release condition are better identified through a scalogram with several units to identify the intensity of the observations and relevant values. The most common spectral models, including the wavelet transformation (WT), short-time Fourier transformation (STFT), and multitaper models, revealed an excellent spatial relationship between an earthquake’s magnitude and its strain release mechanism. For instance, the wavelets are natural functions with the same shape, concentrating on time and frequency. The signal used in WT is multiplied and transformed into a particular form calculated for diverse parts of the time-domain signal (Heidari and Salajegheh 2006). In the short-time Fourier transform (STFT), no resolution issue is noticed in the frequency domain, as it is known precisely what frequencies exist. Chakraborty and Okaya (1995) worked on and compared wavelet-based and Fourier-based models for performing frequency analysis on the seismic dataset. Sinha et al. (2005) proposed and demonstrated a noble approach for creating a time-frequency map by doing a Fourier transform on the inverse cross WT. On the contrary, the time resolution in the STFT and the frequency resolution in the time domain are both zero. The multitaper model of spectral analysis is adapted to the case of data irregularity or missing value, which illustrates not only the lack of data in past centuries but also in the twentieth and twenty-first centuries of sampling.

Several research scholars in Bangladesh have made numerous efforts to evaluate the seismic hazard (Ansary and Sharfuddin 2002; CDMP 2009; Al-Hussaini and Al-Noman 2010; Trianni et al. 2014; Al-Hussaini et al. 2015; Carlton et al. 2018; Rasel et al. 2019; Rahman et al. 2020). Many studies have been performed with the statistical analysis of earthquakes in Bangladesh (Saha 2005; Akhtar 2010; Paul and Bhuiyan 2010; Haque et al. 2020; Kamal et al. 2021). According to Akhtar (2010), the studies revealed that the megacity of Dhaka, Bangladesh’s rapidly expanding capital city with a large population (12.8 million as of 2008), presents an exceptionally high risk for an earthquake. Bangladesh’s northeastern cities are highly susceptible to earthquake threats. Most of the earlier works have possible inaccuracies due to data inaccuracy or the use of limited data. Seismicity analysis based on the Bangladesh earthquake is vital for assessing the potential earthquake risk. However, Bangladesh’s assessment of seismicity based on the proposed Guttenberg-Richter (G-R) relationship and spectral analysis is minimal. Besides, the WT, STFT, and multitaper approaches are relatively new in the seismicity assessment of Bangladesh. Such scalograms or spectrograms are rarely found in analyzing seismicity in Bangladesh.

From the above-mentioned problem and earlier research gaps, the following three primary research questions are raised: First, what is the value of a and b of the proposed G-R relationship regarding the frequency magnitude analysis in Bangladesh? Second, how does the strain release, and what is its probable amount of discharge from 1992 to 2017 based on the seismicity catalog of Bangladesh? Third, how can the actual earthquake magnitude and strain release condition be identified through a scalogram using wavelet transformation (WT), short-time Fourier (STFT), and multitaper models? The answers to these questions are unknown in Bangladesh. Therefore, the primary goals of this study are (i) to determine the constants a and b and to estimate the frequency and magnitude of the earthquake in Bangladesh using the G-R relationship; (ii) to estimate the strain release from the Bangladesh earthquake since 1985 to 2017; and (iii) to assess the spectrogram and related scalogram of earthquake magnitude and strain release mechanism using WT, STFT, and multitaper models. A generalized hypothesis has been developed to identify the probable earthquake condition and the total accumulation of strain release in Bangladesh since 1985. This empirical research presents the new facts and aspects of the Bangladesh earthquake by focusing on the proposed G-R relationship and spectral representation, including WT, STFT, and multitaper model analysis. The findings of this work are mainly crucial in seismic hazard studies in and around Bangladesh, positioned in the eastern portion of the Indian plate.

2. Data and methods

2.1. Study area description

Bangladesh was selected as the study area due to being one of the most vulnerable countries to natural disasters in the world. The epicenter locations of earthquakes are adjacent to the country. A relatively wider study area, from 20°34′2 N–28°38′2 N and 86°01′2 E–94°412 E, a site covering 722,904.41 sq. km, was selected to get a complete earthquake catalog from 1992 to 2017. The whole area includes the three major fault lines that run through Bangladesh from north to south (at the east): the Dauki fault, the Sylhet-Assam fault, and the Chittagong–Myanmar fault, which cause large earthquakes (Figure 1). A section of the world’s most significant river delta is at sea level, which increases the danger of tsunamis and floods following earthquakes (Alam et al. 2003). Bangladesh occupies most of the Bengal Basin, a geotectonic portion of the Assam-Himalayan area and the world’s most significant depositional feature (Kuehl et al. 1989). The Bengal Basin, the world’s biggest delta (60,000 km²), drains most of the Himalayan runoff (Johnson and Alam 1991). Three tectonic plates cover the nation: the Indian, Eurasian, and Burmese plates (Figure S1). Two active tectonic plates, the Indian and Eurasian, have created two subduction zones near Bangladesh (Alam 2019). An active fault in Haluaghat, Mymensingh, was recently established, increasing the vulnerability risk (Saha 2005). Table 1 summarizes the devastating earthquakes in the study region between 1663 and 2011, along with the areas impacted, magnitudes, geological risks, and death/injuries. At least 12 earthquake events occurred outside Bangladesh (in India and Myanmar), triggering significant internal damage. The 1663 Assam, 1869 Cachar, 1897 Assam, 1923 Meghalaya, 1930 Dhubri, 1934 Bihar-Nepal, 1950 Assam, 1954 north Myanmar, and 1977 Bangladesh earthquakes are only a few of these historical earthquakes (Alam and Dominey-Howes 2016). Other vital historical earthquakes with epicenters in Bangladesh, like those in 1762 and 1885, also triggered severe damage.

Figure 1. Seismotectonic map of Bangladesh and the surrounding areas displaying earthquake epicenters from 1762 to 2016 (declustered catalogue). the tectonic characteristics are taken from Kayal et al. (2012); Steckler et al. (2016) and Kamal et al. (2021). The global multi-resolution terrain elevation dataset 2010 (GMTED 2010) of the USGS served as the source for the backdrop digital elevation model (DEM). the basic geological map of the country is placed on the seismotectonic map, which was modified by Alam et al. (1990).

Table 1. Occurrence of large earthquakes in history around Bangladesh (Source: Akhtar 2010; Alam 2019; Kamal et al. 2021).

Date	Name of earthquake	Latitude/ longitude	Affected regions	Magnitude	Death/ injuries
19-02-1663	Assam earthquake	26.1/92.56	Assam, India, and Bangladesh	8	–
02-04-1762	Chittagong earthquake	22/92	India, Bangladesh, and Myanmar	7.8?	200
10-01-1869	Cachar earthquake	24.75/93.25	India and Bangladesh	7.5	–
14-07-1885	Bengal earthquake	24.8/89.5	Bangladesh and India	7.0	75
12-06-1897	Great Indian earthquake	26/91	India and Bangladesh	8.7	1626
1923	Meghalaya earthquake	24.3/91.7	Bangladesh and India	7.6	9
02-07-1930	Dhubri earthquake	25.8/90.2	India and Bangladesh	7.1	1
15-01-1934	Bihar-Nepal earthquake	26.5/86.5	Nepal, India, and Bangladesh	8.3	13,772
15-08-1950	Assam earthquake	28.12/94.05	Assam, India, and northern Bangladesh	8.3	–
22-03-1954	North Myanmar	24.5/95.3	Myanmar and Bangladesh	7.3	–
12-05-1977	Bangladesh-Myanmar	21.75/92.99	Myanmar and Bangladesh	5.7	–
06-08-1988	Myanmar	25.14/95.12	Myanmar and Bangladesh	7.3	–
21-08-1988	Bihar-Nepal	26.7/86.8	Nepal, India, and Bangladesh	7.8	998
18-09-2011	Sikkim Earthquake	27.73/88.15	India and Bangladesh	6.9	97

2.2. Geology and seismotectonic setting

The Bengal Basin is in the eastern part of the Indian subcontinent. It is the largest fluvio-deltaic to the shallow marine sedimentary basin. It comprises the riverine channel, floodplain, and delta plain ecosystems. Common indicators of deformation in the basin are tilting, warping, and folding of the strata with anomalous drainage (Reimann 1993). Plate-to-plate contact changes throughout time, affecting basin architecture and sedimentation style. The Bengal Basin and neighboring areas have withstood 23 severe earthquakes, indicating geotectonic activity and the reactivation of ancient faults (Kayal et al. 2012). The Pre-Cambrian Shillong Plateau surrounds the Bengal Basin, the Indian Platform, the Arakan-Yoma-Naga folded system, and the Bay of Bengal to the south. The Bengal Basin is an exo-geosyncline with thick detrital deposits from the craton uplift. The Bengal foredeep is one of the world’s biggest exo-geosynclines (Alam et al. 1990). Earlier investigations of the regional stratigraphic and tectonic scenarios (e.g. Evans 1932; Sengupta 1966) created the groundwork for understanding basin development and sediment-fill history. Bakhtine (1966) recognized Bangladesh’s tectonic characteristics. Alam (1972) described the basin’s geosynclinal evolution.

The active edge of eastern Bangladesh is a brittle-plastic transition zone characterized by ‘Pelagic Mud,’ whereas intra-fold decollement is ‘Shefal Mud’ (Khan 1980). This mud may dampen earthquakes. Bangladesh’s earthquake history is connected to tectonically active north and eastern seismic zones. The Bengal Basin is tectonically unstable since most geological activity occurs around plate boundaries. Choudhury (1993) described ‘active’ and ‘passive’ Bengal Basin margins. The ‘active margin’ is the deep basinal portion and its folded eastern edge. The active margin of the Bengal Basin moves because of subduction, while the passive margin is unstable because of changes in the crust. The NW-SE Padma fault, Karotoa-Banar fault, Tista-Old Brahmaputra fault, N-S Dubri-Jamuna-Madhupur fault, and E-W Dauki fault are in charge of controlling the tectonics of the passive margin (Figure S2). A shaking core in the earth causes an earthquake, and physical jolts on the earth’s surface cause new ruptures and reactivate old ruptures.

Based on these faults and ruptures, the researcher noted fracture propagation is limited to the northwest, northeast, and southeast zones, which coincide with seismic fault source zones (Reimann 1993). Fracture propagation and fault rupture are expected in the northeast and southeast Bangladesh. The faulting in the Dauki source zone is vertical to near vertical basement induced, whereas faulting in the southeast occurs along detachment zone slip planes. The Rajshahi Division features neotectonic faults in Bangladesh’s northwest. The northeast has seen moderate-to-large earthquakes. The Dauki Fault Zone is the largest crustal displacement. Habiganj’s western section includes Raghunanda Hill. Habiganj has a clay-and-clay-gall fault escarpment. The Lalmai Hills in Comilla is a fault-bound elevated terrace from the Neogene era. The fault escarpment, terracotta, pottery, and stone weapons indicate recent reactivation and uplift along the Lalmai fault. Borkol’s eastern margin is faulted in Rangamati. Active faults are in Chittagong, Bandarban, and Teknaf.

2.3. Data analysis and tabulation

Earthquake databases are obtained from Bangladesh Meteorological Department (BMD). Data processing after collecting earthquake data usually involves the manipulation of items, including locations, magnitude, and the tabulation of data with possible sorting. The tabulation of data involves an orderly arrangement of the data through various classifications. Data processing involves multiple steps as follows

• ✓ Validation—the valid seismic data were collected from BMD, Dhaka.

• ✓ Sorting—the items were arranged and sorted by magnitude and location.

• ✓ Summarization—a summarized table of Bangladesh earthquake data from 1908–2017

• ✓ Aggregation combining multiple types of datasets

• ✓ Analysis—the collection, organization, analysis, interpretation, and presentation of data

• ✓ Classification—separation of datasets into different types for applying techniques

After the compilation of the data, it is processed and then finally analyzed.

2.4. Instrumental earthquake catalog

Compiling earthquake data is crucial to understanding earthquake dynamics in any location. To describe its correctness, one must include the whole magnitude evaluation and the estimation errors for the parameters (locations, magnitudes, and depths). An instrumented seismic inventory for the research region was generated for 1992–2017, as presented in Figures S3–S5.

2.5. Frequency-magnitude analysis

The most extensively used empirical association between magnitude (M) and frequency of occurrence of earthquakes (N) in a specific area is given by the Gutenberg and Richter relationship (Gutenberg and Richter 1944) as follows in Equation (1) (1) $$\mathit{Log}\mathrm{ }N=\mathrm{ }{a-b}^M$$(1)

Where N is the number of earthquakes with M_s equal to or greater than M per year, and a and b are the constants for the region (Gutenberg and Richter 1944). Where a is a measure of the level of seismic activity and b (b value) is the rate at which the events occur within a given magnitude range. Generally, b takes values between 0.67 and 1. Higher values of b denote that smaller magnitude events are more abundant than the larger ones for that particular region. The b-value plays an essential role in the GR relationship. As it helps to structure a clear idea about the seismic pattern and the seismic stress level of an area (Ahmed et al. 2016; Ray et al. 2019), the calculation of the b-value is a fundamental step for the seismic hazard analysis of any zone (Tinti et al. 1987; Khalid et al. 2014).

2.6. Strain release accumulation

Magnitude is related to the energy that emits from an earthquake source as elastic waves. In this paper, strain accumulation and release patterns were obtained from the geotectonic zone of Bangladesh. The seismic energy was calculated based on Gutenberg-Richter equation (Equation 2(2) $$\mathit{Log}{J}_{10}{}^{1/2}=5.9+0.75\mathrm{\hspace{0.33em}}M$$(2) ). Benioff (1955) method and Gutenberg and Richter (1954) formula were used where the annual cumulative strain release is (2) $$\mathit{Log}{J}_{10}{}^{1/2}=5.9+0.75\mathrm{\hspace{0.33em}}M$$(2) where J is the energy released in ergs by an earthquake of surface wave magnitude M. The plots of cumulative values of seismic energy against a given period for Bangladesh revealed that this area is characterized by the amount of strain accumulation that can precipitate in an earthquake of body wave magnitude.

2.7. Wavelet transformation analysis

A robust data analysis technique often utilized in research is spectral analysis. Several approaches have been devised to analyze extended records of stationary processes of varying magnitude and data on strain release. This study uses three spectral analyses to estimate leakage-free spectra: wavelet Transformation (WT), short-time Fourier analysis, and Multitaper approach. The Fourier transformation (FT) divides a signal into a succession of sine waves of various frequencies, while the wavelet transformation (WT) divides a signal into its components. The so-called mother wavelet is scaled and shifted into wavelets by the WT, which divides the signal into these wavelets. WT permits exceptional frequency (scale) localization through dilations and wavelet translations in the time domain. The wavelets have the same structure and are real or complex functions concentrating on time and frequency. The signal is multiplied by the wavelet, and the transform is independently calculated for various time domain signal segments in the WT. In general, the WT of the signal, x(t), is defined as the following inner production Equation (3)(3) $$WT \left(\tau , b\right)={\int }_{-\infty }^{+\infty } x \left(t\right)\mathit{\text{ gt }}(ᴦ, b) dt\mathrm{ }$$(3) (3) $$WT \left(\tau , b\right)={\int }_{-\infty }^{+\infty } x \left(t\right)\mathit{\text{ gt }}(ᴦ, b) dt\mathrm{ }$$(3)

The family of continuously transformed wavelets is generated from the mother wavelet (Equation 4(4) $$g \left(\tau , b\right)=\frac{1}{\sqrt{b}} g \left(\frac{t-ᴦ}{b}\right)$$(4) ) (4) $$g \left(\tau , b\right)=\frac{1}{\sqrt{b}} g \left(\frac{t-ᴦ}{b}\right)$$(4) where t is the transform parameter, corresponding to the position of the wavelet as it is shifted through the signal, b is the scale dilation parameter defining the width of the wavelet. The scale b > 1 dilates (or stretches out) the signals, while scale b < 1 compresses the signal. The correlation between the wavelet and a particular localized area of the signal is represented by the wavelet coefficients, WT (t,b). The wavelet at this scale is near to the signal at the specific place if the signal includes a significant frequency component that corresponds to the provided scale, and the associated wavelet transform coefficient (Equation 5(5) $${g \left(t\right)=e}^{{-0.5t}^2} e^{\mathit{ict}}$$(5) ) obtained at this point has a comparatively high value. One of the most extensively used mother wavelets is Morlet’s wavelet (Morlet 1981) (5) $${g \left(t\right)=e}^{{-0.5t}^2} e^{\mathit{ict}}$$(5)

Using the representation Equation (5)(5) $${g \left(t\right)=e}^{{-0.5t}^2} e^{\mathit{ict}}$$(5) , the Morlet wavelet takes the form of Equation (6)(6) $$g (\tau ,b)\frac{1}{\surd b} e^{{-0.5\left(\frac{t-\tau }{b}\right)}^2} e^{ic\frac{t-\tau }{b}}$$(6) : (6) $$g (\tau ,b)\frac{1}{\surd b} e^{{-0.5\left(\frac{t-\tau }{b}\right)}^2} e^{ic\frac{t-\tau }{b}}$$(6)

Wavelet transform is a relatively new technique, and in recent years, there has been a lot of interest in using wavelets in many applications. On the other hand, research publications seldom use the wavelet transform for earthquake engineering.

2.8. Short-time Fourier transformation

The accuracy of the time localization is governed by the window width and the frequency localization despite the constant product of these two (Vetterli and Kovačević 1995). The window width should be carefully selected to detect all frequencies. The tradeoff between time and frequency localization is ideal for all frequencies in the study range, where f is the smallest frequency of interest, and the temporal frame should be bigger than 1/f. However, too large windows would degrade the time localization’s accuracy. In Equation 7(7) $$\text{STFT θ}x \left(\omega ,t\right)={\int }_{-\infty }^{\infty }w \left(\tau -t\right){\theta }_x \left(\tau \right)e^{i\omega \tau } d\tau$$(7) , the short-time Fourier transform (STFT), e.g. of $\theta$x(t), is by definition equal to (7) $$\text{STFT θ}x \left(\omega ,t\right)={\int }_{-\infty }^{\infty }w \left(\tau -t\right){\theta }_x \left(\tau \right)e^{i\omega \tau } d\tau$$(7)

Where w ($\tau -t$) = a real and usually used as a symmetric window centered at time t. It is a surface defined in the frequency-time plane, representing a Fourier transform that is localized at time t. The ripple effects of sharp window margins are avoided using w (t), which eventually disappear. A trapezoidal-shaped window with linear ramps at both ends and constant amplitude in the middle was utilized as w (t) in Udwadia and Trifunac (1973).

2.9. Multitaper model analysis

The multitaper spectrum analysis approach was first presented by Thomson (1982) and has been extensively used in seismogram analysis (Park 1987). In a multitaper analysis, the data are multiplied by some leakage-resistant tapers rather than just one. As a result, one record produces a very tapered time series. Each of these time series’ DFTs (discrete Fourier transforms) yields several ‘eigen spectra,’ which are then averaged to create a single spectral estimate. Many different multitapers have been suggested, such as Slepian tapers, discrete prolate spheroidal sequences, sinusoidal tapers, etc., to name a few of them. The fundamental concept of this multitaper technique is that, given moderate circumstances, the spectral estimations would be independent of each other at every frequency, provided the data tapers are constructed as suitably orthogonal functions. The design of the multitapers maximizes resistance to spectral leakage while allowing each taper to sample the time series differently. The second taper partly recovers the statistical data that the first taper abandoned; the third taper partially recovers the data that the first two tapers partially discarded, and so on. The higher-order tapers allow for an unacceptably large amount of spectral leakage; hence, only a few lower-order tapers are used. The trade-off between leakage and variation that plagues single-taper estimates may be avoided by using these tapers to construct an estimate. Group of orthonormal tapers, including harmonically connected sinusoidal tapers. These tapers are also known as minimal bias tapers or sinusoidal tapers. The continuous time minimum bias tapers are given in Equation (8)(8) $$V_k\mathrm{ }\left(t\right)=\sqrt{2}\sin \hspace{0.17em}\left(\frac{\pi kt}{\tau }\right)\mathrm{ }\left(k=1,\mathrm{ }2,\dots \right)$$(8) (8) $$V_k\mathrm{ }\left(t\right)=\sqrt{2}\sin \hspace{0.17em}\left(\frac{\pi kt}{\tau }\right)\mathrm{ }\left(k=1,\mathrm{ }2,\dots \right)$$(8)

and its Fourier Transform can be expressed by Equation (9)(9) $$V_k\left(\omega \right)=\sqrt{2j}\left[\frac{\sin (\omega +\frac{\pi k}{\tau })}{(\omega +\frac{\pi k}{\tau })}+\frac{\sin \hspace{0.17em}\left(\omega -\frac{\pi k}{\tau }\right)}{\left(\omega -\frac{\pi k}{\tau }\right)}\right]$$(9) (9) $$V_k\left(\omega \right)=\sqrt{2j}\left[\frac{\sin (\omega +\frac{\pi k}{\tau })}{(\omega +\frac{\pi k}{\tau })}+\frac{\sin \hspace{0.17em}\left(\omega -\frac{\pi k}{\tau }\right)}{\left(\omega -\frac{\pi k}{\tau }\right)}\right]$$(9)

The discrete analogs of the continuous time minimum bias tapers are known as sinusoidal tapers. The k-th sinusoidal taper is given in Equation (10)(10) $$n-1,\mathrm{ }2,\dots .N;\mathrm{ }k=1,\mathrm{ }2,\dots .\mathrm{ }K,$$(10) (10) $$n-1,\mathrm{ }2,\dots .N;\mathrm{ }k=1,\mathrm{ }2,\dots .\mathrm{ }K,$$(10)

Where the normalization factor in the amplitude term on the right assures the orthonormality of the tapers. These sine tapers feature side lobes that are significantly taller and a main lobe that is substantially smaller. Due to the main lobe’s smoothing, they can attain a reduced bias, albeit at the price of side lobe suppression. This performance is adequate if the spectrum changes gradually. The kth order sinusoidal taper has its spectral energy focused in the frequency bands, (11) $$\frac{\pi (k-1)}{N+1}\le \left|\omega \right|\le \frac{\pi \left(k+1\right)}{N+1}\mathrm{ }k=1,\mathrm{ }2,\dots K$$(11)

The sinusoidal taper of order factor k = 1, 2, and 3 is shown in time. Higher-order tapers significantly weigh the sample values of the data that the first taper only moderately weighs. The third and fourth tapers heavily weigh the data samples that are only slightly weighted by the first and second tapers. With the extra flexibility these tapers have, because they have a slope of fall for the weighting function, and so limit the leakage, the data is weighted equally at all places, like the rectangular taper. A multitaper technique can estimate the spectrum of a synthetic time series so that it looks like the actual spectrum. ArcGIS (version 10.2.1) was used for preparing maps. Other software, including SPSS, R software, PAST, and Microsoft Excel 2007, was used to compute the result, prepare data sets, and analyze relevant data statistically. The flowchart adopted to carry out this research is provided in Figure S6.

3. Calculation and results

3.1. Calculation of b-value using GR relationship

The Guttenberg and Richter frequency magnitude relationship method was used to compute the b-value for the entire Bangladesh Equation (12)(12) $$\mathit{Log}\mathrm{ }N={2.981-0.39225}^M$$(12) implied to obtain coefficients a and b, which establishes the following GR relationship.

The constants a and b are evaluated first for the 442 main recorded events (M_s ≥ 3.5) for Bangladesh and its vicinity from 1992 to the end of 2017 using the least squares technique for computing the recurrence periods of different earthquake magnitudes. The frequency-magnitude relationship for Bangladesh suggests that the frequency linearly decreases with increasing magnitude and follows Equation (12)(12) $$\mathit{Log}\mathrm{ }N={2.981-0.39225}^M$$(12) (12) $$\mathit{Log}\mathrm{ }N={2.981-0.39225}^M$$(12)

Four hundred forty-two seismic events were used in this case, and an estimated b-value of 0.39225 was obtained. Though the value found in the study appears to be satisfactory, generally, the b-value of 1 is more expected when the study area is a seismically active zone with numerous occurrences of events (Ray et al. 2019). Furthermore, as the foreshocks and aftershocks were excluded, it eventually decreased the value of b (Ahmed et al. 2016). The value of b may increase if the foreshocks and aftershocks are included in the calculation, increasing the number of seismic events. It is safe to choose the main events only to lessen the uncertainty (Ray et al. 2019). However, the regression coefficient of R² value is 0.69. This study is concerned with solving for the long-term average values of the and b parameters for the entire Bangladesh.

3.2. Frequency-magnitude relationship

The frequency-magnitude relation presents the number of earthquake with magnitude lying between a pair of values occurring in a given region during a certain period of time is expressed by Equation (13)(13) $$N\left(M\right)={10}^{(a-b\mathrm{ }M)}$$(13) : (13) $$N\left(M\right)={10}^{(a-b\mathrm{ }M)}$$(13)

Where N(M) is the number of shocks of magnitude M or greater, per unit time, a is a measure of the level of seismic activity, and b (b value) is the rate at which the events occur within a given magnitude range.

Under the assumption that the magnitude data are random samples from a population obeying the Gutenberg-Richter relation (Figure 2), the method of moment and method of maximum likelihood (Aki 1965) both yield solution in Equation (14)(14) $$b=\frac{{\mathit{log}}_{10}e}{M-M_z}$$(14) (14) $$b=\frac{{\mathit{log}}_{10}e}{M-M_z}$$(14)

Figure 2. Frequency-Magnitude distribution of the instrumental earthquakes recorded along the country in the study period; (a) shows the linear trend; (b) shows the future Predictive line using different variable values and (c) normal P-P plot of regression standardized residuals which depict the residual vs. leverage using cook’s distance.

Where M is the mean magnitude of earthquakes of M ≥ M_z and log e = 0.434294.

The generalized linear model also shows the same linear relation with two probable linear lines in Figure 2(b). The cumulative probable observation in Figure 2(c) also illustrates the linearity against the number of magnitudes. Scales M_s, M_b, and M_l were used in Figures 2(a,b).

Table S1 shows that the regression coefficients (R²) and standard error of the estimate are 0.686 and 0.29041, respectively. This study is concerned with solving for the long-term average values of these parameters for Bangladesh. Table S2 shows the Frequency-Magnitude relation analysis, which is significant at a 99% confidence level.

Figure 3 shows the residual-frequency curve where the sum of residuals will always be zero. The residual median point is 0, which touches the highest peak of the frequency, where the mean is −4.37, and the standard deviation is 0.986. The Cook’s distance is related to the residuals where magnitudes are the independent variables. The highest distance was found at 0.0 residual, as shown in Figure 3(a). A catalog of earthquakes with M_s ≥ 3.5 from 1992 to 2017 in Bangladesh was considered for this study. It demonstrated a strong log-linear relationship for events of M_s ≥ 3.5. Moreover, it includes 442 main earthquake events, whose locations and magnitude values are shown in Figure S3.

Figure 3. Estimation of residual-frequency curve and probability of earthquake magnitude in Bangladesh; (a) Residual-Frequency curve presents the highest peak of frequency; (b) presents the cook's distance using regression analysis to find influential outliers.

3.3. Strain release analysis

The highest stress release was recorded in 1987, with body magnitudes M_b of 5.5, 5.1, and 4.8 (Figure 4(a)). The highest value was found in 1987 when only five events were recorded. The second highest was recorded in 1992 when only 25 events were found. The other years revealed significant events also. Figure 4(b) presents the parallel linear line higher than in the previous 8 years. In this case, 1995 shows the higher stress release, where 12 events were recorded. The following years, 1996, 1997, 1998, 1999, and 2000, also significantly indicate a significant rise in pressure. Figure 4(c) continues the higher release of energy shown in 2004 when 12 events were recorded. The following five years also showed parallel high-strain release. The parallel fluctuation is shown in Figure 4(d), but 2011, 2012, and 2016 showed higher energy release for high magnitude.

Figure 4. Strain release accumulation pattern of Bangladesh from 1985 to 2017; (a) shows stress emits from 1985 to 1992; (b) shows stress emits from 1993 to 2001; (c) shows stress emits from 2002 to 2009 and (d) shows stress emits from 2010 to 2017.

Figure 5 represents the strain energy released over 32 years from 1985 to 2017 and offers a valuable way of comparing Bangladesh’s seismic activities. Figure 5(b) shows the highest strain released in 1995, where body magnitude M_b was found at 5.1, 5.2, 5.9, and 6.4. Figure 5(c) presents the highest strain release in 2004 because of moderate body magnitude (Mb) of 4.8, 4.9, and 5.5. Figure 5(d) shows the higher strain value in 2016 only when M_b 4.8, 5.5, 5.8, 6.6, and 7.2 are present. The strain accumulation and its release pattern in Bangladesh’s geotectonic zones showed a linear increasing strain buildup. Records from six places in Bangladesh showed how quickly the Bengal Delta moves and the resulting stress builds up (Table S3).

Figure 5. Strain accumulation and release pattern of the geotectonic zone in Bangladesh during the study period.

3.4. Wavelet analysis

The wavelet analysis was computed to present magnitude and strain release on different scales. The colorization of the shade presents magnitude and strain release on different scales where the reconstruction factor Cδ, decorrelation factor for time averaging γ, and factor for scale averaging δj₀ are 1.966, 1.37, and 0.97, respectively.

The average wavelet analysis is combined into a two-dimensional contour plot, as shown in Figure 6(a), which shows consistent fluctuation. A zoomed-in view to show the self-similar (fractal) magnitude pattern with more clarity. Figure 6(b) depicts a spectral presentation using the derivative of Gaussian (DOG) wavelet with m = 6. Linear, curved regions on either end indicate the ‘cone of influence,’ where edge effects become important. The shaded contour has a normalized variance of 8.0. The method surprisingly displays the cone of impact (COI), where the e-folding time for the autocorrelation of wavelet power at each scale was defined as the area of the wavelet spectrum in which edge effects became essential (Figure 6(c)). This e-folding time is selected to guarantee that the edge effects are minimal beyond this point by causing the wavelet power for a discontinuity at the edge to decrease by a factor of e-2. For a single spike in a time series, the magnitude of the COI at each scale also provides a measurement of the decorrelation time of 1.37 s. The red spots show that the scale is consistent from 64.00 to 22.63. It also shows that the magnitude scale goes from 4.8 to 5.5 (Figure 6(c)).

Figure 6. Wavelet analysis based on magnitude; (a) the anomaly depicts the magnitudes plotted against the overall observation; (b) cross correlation between magnitude and overall observation and (c) scalogram representation of earthquake magnitudes in Bangladesh during the study period.

Figure 7(a) shows the continuous fluctuation of energy release. It presents a zoomed view to observe a better stress pattern. The average wavelet analysis is combined into a two-dimensional contour plot, as shown in Figure 7(a). The thick contour encloses regions greater than 95% confidence for a red-noise process with a lag-0 significance level of 0.05. Linear, curved areas on either end indicate the ‘cone of influence,’ where edge effects become important. As in Figure 7(b), this energy is shown on the lag and correlation scale. Figure 7(c) illustrates the uniform pattern of red shade that refers to the highest energy release in a specific range. Figure 7(c) also shows a zone of red patches, illustrating a range of particular magnitude, 4.8, 4.9, 5.2, 5.5, 5.8, and 6, respectively.

Figure 7. Wavelet analysis based on propagated strain release; (a) the anomaly depicts the release of strain plotted against the overall observation; (b) cross correlation between magnitude and overall observation and (c) scalogram representation of strain release in Bangladesh earthquake.

3.5. Short-time Fourier transformation

Figure 8(a) shows the magnitude in different frequency bands. It demonstrates the magnitude at different frequencies for the window width of 32. In this paper, 442 earthquake datasets have been presented to illustrate the performance of Fourier transformation through STFT operation. STFT operated the strain release for the recorded magnitude. In Figure 8(b), the hammering style was adapted to perform the STFT of strain release, where 32 windows are used. As in Figure 8(a), the red patches were found at 0.05 levels, illustrating the higher magnitude, containing 4.9, 5.5, and 5.8, respectively. The Fourier and wavelet patterns are changed for their equation differences. The p-value of 0.095 illustrates the maximum significance of the spectrogram.

Figure 8. (a) Short-time Fourier analysis performed on magnitude; (b) short-time Fourier analysis performed on strain release in the study area.

The STFT has been applied to analyze, modify, and synthesize non-stationary or time-varying signals. Figure 8(b) illustrates the strain release in the Bangladesh earthquake from 1985 to 2017. It shows the red patches contain a vast area. The specific problem with the STFT is the possible occurrence of leakage (Proakis and Manolakis 1996), which has been reduced using Hammer style 32 windowing. It may demonstrate the higher strain release linked to the high magnitude range. The higher degrees of magnitude, 4.9, 5.2, 5.5, 5.6, 5.8, and 6.8, are responsible for higher strain release. The scale shows the dark red patches contain 0.10–0.20 and 0.20–0.30 for 5.8, 6.5, and 6.8 magnitude.

3.6. Multitaper model analysis

Figure 9(a) represents the multitaper, simple, and smooth periodograms using blue, green, and black lines. The frequency and power grid illustrated the magnitude, which shows a parallel fluctuation. Multiplying the time series by a taper before operating a discrete Fourier transform (DFT) to reduce spectral leakage is standard practice. The initial findings contain the fluctuation of magnitude in the overall change, indicating that the magnitude of the Bangladesh earthquake follows the same range of occurrence in specific periods (Figure 9). The study showed that at different frequencies, there are breaking points at 0.05, 0.20, 0.25, 0.32, 0.40, and 0.47, where the exact size of change keeps happening.

Figure 9. Multitaper spectral analysis of earthquake magnitude; (a) multitaper spectra, smooth and simple periodogram conjugated lines of earthquake magnitude in Bangladesh; (b) F value using F test, (c) degree of freedom line.

Figure 9(b) presents the F value and degree of freedom. The upper and bottom rows show the F value of the Multitaper. The adequate degree of freedom increases with the number of tapers used, and the F value necessary to reach a certain significance level decreases accordingly. In the case of Bangladesh, the significance levels of 95 and 99% are shown as dotted lines.

For statistical correctness, the F test was done. To overlay the spectrum, the F test scale shows the visual correlation in bands of high F value and power. Figure 9(b) shows that the F values were in the 5–10 modulation range. The adequate degree of freedom increases with the number of tappers. The degree of freedom is shown in Figure 9(c). The degree of freedom indicates the drop-off point of frequency in 0.20, 0.25, 0.30, 0.40, and 0.50. The drop points may be shown only using three types of tapers. Figure 10 illustrates the continuous fluctuation of the strain release in Bangladesh. It also shows the breaking points where the same oscillation pattern continues are 0.20, 0.30, 0.40, and 0.47, quite similar to the previous magnitude break points.

Figure 10. Multitaper spectral analysis of earthquake strain release; (a) multitaper spectra, smooth and simple periodogram conjugated lines of earthquake strain release in Bangladesh, (b) F value using F test, (c) degree of freedom line.

In this case, the average yearly stress from 1985 to 2017 follows the Guttenberg-Richter relationship formula. Figure 10 shows the Multitaper spectral with a smooth and straightforward periodogram. Figures 10(b,c) also indicate the F value and the degree of freedom line, where the f value shows a 99% significant level. The findings also show the F value of the Multitaper conducted by the F test, considering the continuous fluctuation of strain release. F is between 4 and 5; the degree of freedom is also shown in Figure 10(c).

4. Discussion

4.1. Frequency-magnitude association

The findings show that the frequency-magnitude relationship always follows a linear relationship. Gutenberg and Richter (1954) found that the distribution of earthquakes in a given location on the Earth typically followed the Gutenberg-Richter power law, which was approximately linear. The linearity of this law depends on the earthquake magnitude and the logarithmic number of magnitudes in a particular year (Asfahani and Darawcheh 2017). Crampin et al. (2013) said that some events bigger than any magnitude apply logarithmically and would be proportional to the magnitude. Therefore, the linearity is of moderate magnitude for many occurrences. The findings show that Bangladesh’s b value is 0.392, indicating the maximum shear stress condition prevailing in Bangladesh.

In Bangladesh, earthquakes with an average of around five on the Richter scale are frequent. Although there has been paleoseismic activity inside Bangladesh, as suggested in recent years, for the last 200 years (Kamal et al. 2021). The physiography of Bangladesh makes it more vulnerable to earthquakes. It also indicates recent tectonic activity or the propagation of fractures from the adjacent seismic zone of Bangladesh. The historical evidence of an earthquake in Bangladesh shows the magnitude averaging 5, 6, and 7 on the Richter scale. The Bengal Delta is situated on the converging plate margin, and tectonic activity is the major factor associated with the evolution of the delta, making it genetically vulnerable to earthquakes. The deposition of quaternary sediments was markedly influenced and controlled by tectonic activity. The presence of several faults in Bangladesh, namely the NW-SE trending Padma fault, Karotoa-Banar Fault, Tista-Old Brahmaputra Fault, N-S fault trending Dubri-Jamuna-Madhupur fault, and E-W trending Dauki fault, makes Bangladesh more prone to earthquakes (Ansary and Sharfuddin 2002; Rahman et al. 2020). The more occurrences of an earthquake propagate, the more energy is released, proportional to the magnitude. Thus, the earthquake follows the frequency-magnitude linear association in Bangladesh, which tends to go downward with more energy release.

Another reason for linearity might be the logarithmic operation of the number of events, which is proportional to the increasing number of events. The result shows that the number of magnitudes of 4.5 is higher than those of 5, 6, and 7, respectively. The probability of earthquake occurrences of different magnitudes is defined as a and b. It is noted that the earthquake magnitude data for Bangladesh was not appropriately recorded and contained missing values. Alam (2019) mentioned that a good background of historical earthquake data is essential to evaluating the seismicity in close coincidence with the geotectonic elements. According to Godano (2014), the value of logN has gradually been lowered as seismograph networks increase their sensitivity by adding station locations and improving technology. The decrease in the b value is typically reported within less than one of the completeness threshold, where the b value is 0.392. The same model as Alam et al. (2003) envisaged the b value being 0.41, which is close to what this paper used as the b value. Our research suggests that the frequency-magnitude mean is equal (Table S1) because of the linear representation of the frequency–magnitude relationship. We demonstrate whether the Bangladesh earthquake frequency-magnitude relationship follows downward linearity with a minimum b value that appears to be the maximum energy release based on the b value.

4.2. Strain accumulation outcomes

Our findings reveal that the strain accumulation in Bangladesh is getting upward with increasing magnitude. The study also shows the seismicity through strain release, which Benioff also introduced as an indication of the seismicity of an area. The findings illustrate that the strain release depends on the magnitude; thus, the upward trend of strain releases in Bangladesh indicates many earthquake occurrences (Rasel et al. 2019; Kamal et al. 2021).

Bangladesh is prone to maximum strain release because sediments are building up there, up to 20 km thick. These sediments continuously accumulate strain due to the N-E direction in the Bengal Delta at about 5 cm/yr (Alam and Dominey-Howes 2016). Bangladesh has some faults, such as the Padma fault, which goes from NW to SE, the Karotoa-Banar fault; the Tista-Old Brahmaputra fault; the Dubri-Jamuna-Madhupur fault, which goes from N to S, and the Dauki fault, which goes from E to W. These faults may also control the release of strain. Gutenberg and Richter (1954) reported that magnitude is related to the energy radiated from an earthquake source in the form of elastic waves. This proves magnitude strain proportionality, where energy is released with increasing magnitude. The physics underlying the Gutenberg-Richter relationship is the stress-induced manipulation of the geometry of the ubiquitous distributions of stress-aligned, fluid-saturated, crack-critical micro-cracks pervading almost all rocks throughout the crust’s uppermost mantle (Crampin et al. 2013). The present study demonstrates that the amount of strain energy released is the same magnitude as the observed radiated seismic energy.

4.3. Magnitude and strain release using wavelet transform

The results clearly showed the moderate level of the earthquake on the Richter scale and the same earthquake magnitude occurrence uniformly across the range of effects. According to Simons et al. (2006), a magnitude error of ±1 can be predicted based on the wavelet approach. The amplitude and frequency content of the initial seconds of the signal did not depend on epicenter distance, so frequency information could be used to estimate earthquake magnitude directly without correcting for the distance effect (Colombelli et al. 2012). Ziv (2014) tested the idea of an empirical magnitude-scaling relationship in a study where a set of observations were plotted against magnitude. Instead, Kumar and Foufoula-Georgiou (1997) report that wavelet analysis is an excellent way to examine nonstationary signals. This study shows the red patch continuity to show the country’s proneness to the risk of earthquakes with magnitudes ranging from 5.5 to 6.8. This may be the cause of severe earthquakes in the future. Kanamori and Kikuchi (1993) showed that earthquakes could have produced tsunamis much more significantly than expected from their magnitude alone, which means slow earthquakes cannot be negligible.

It is essential for features, such as earthquake sources, propagation paths, and site conditions to be evaluated (Nolasco 2014). The mechanism surprisingly shows the cone of influence (COI), where the part of the wavelet spectrum where edge effects become essential is defined as the e-folding time for the autocorrelation of wavelet power at each scale. Donoho and Johnstone (1994) conducted a sparse representation of strain analysis using wavelet bases in the early 1990s to show the orthogonal wavelet representation. The study findings are consistent with the results. Our study computes the approximate stress-strain curves from the acceleration, where the cone of effluence is shown clearly, similar to Zeghal et al. (1996). The study demonstrated the release of energy under specific strain conditions in Bangladesh from 1992 to 2017, where energy release on a particular scale was found.

4.4. Magnitude and strain release using short-time Fourier analysis

Leakage is the main problem with the short-time Fourier transform (STFT) (Proakis and Manolakis 1996), which has been lessened by hammering style 32 windowing. According to Okaya (1992), an STFT can also be implemented by choosing frequency domain windows instead of time domain windows (Okaya et al. 1992). However, Nuttli (1983) stated that instead of sampling the time axis with moving windows, the frequency axis could be tested by a set of fixed bandwidth bandpass filters whose center frequencies are distributed uniformly along the frequency axis.

Welch (1967) used Welch overlapping segments that overlapped by 50% and were windowed with a Hamm taper, though this one contains 0% overlapping following the Hamming style. The STFT shows the continuously changing spectrum of magnitude. Gibbons et al. (2008) stated that the spectrogram produced provides an image indicating the time a burst of energy occurs on a seismogram. The STFT format of magnitude may have some discontinuity in the missing data of the Bangladesh earthquake record. The study examines how Bangladesh’s magnitude changes over time, showing that more giant earthquakes could happen. Kausel and Assimaki (2002) indicated that the low-amplitude, high-frequency components of strain induce narrower loops than the fundamental loop. However, the present study provides loops as the mode of patches, which are broader and show more energy release in the case of the Bangladesh earthquake. Some experimental studies indicate that some energy is dissipated even at shallow strain levels (Lanzo and Vucetic 1999; Haque et al. 2020), though the study shows a high strain release of a particular magnitude. Vetterli and Kovačević (1995) reported that the window function’s width determines the time localization’s accuracy using a Hamming style of 32 windows. We found a release of the strain of different magnitudes in Bangladesh, which may have presented its higher release through the red patches. Cohen (1995) stated that the STFT would not give the exact values for averages of functions of frequency or time. The findings may cause disorientation and give less clear pasteurization than the wavelet spectrum. Our study shows a large-scale release of higher-strain energy in different amounts from 1992 to 2017.

4.5. Mechanism of magnitude and strain release using multitaper model analysis

Even though there are some missing values in the study, the multitaper method helps to find a smooth figure of magnitude. Other researchers, including Thomson (1982), Park et al. (1987), and Riedel and Sidorenko (1995), found that the resistance to spectral leakage and the variance of a spectral estimate will be in tradeoff as long as just one taper is used. The smoothed spectrum estimate with a boxcar taper gives a good approximation of the parts of the spectrum with large amplitudes and low spectral leakage. Thomson (1982) compared the multitaper estimate with the smoothed direct estimates by making a synthetic seismic wave train with a certain spectrum.

The multitaper taper method makes spectra more stable at high and low frequencies, so more spectral ratios were used in this study. It is worth mentioning that earthquakes occur at considerable distances from plate boundaries in a region characterized by compressive stresses (Zoback and Zoback 1981), low strain and seismic rates, and high crustal seismic velocities (Rahman et al. 2020). As determined by multitaper analysis, the amount of stress released during an earthquake is thought to be greater than in plate boundary environments (Scholz et al. 1986; Choy and Boatwright 1995), where the smoothing analysis prevents leakage due to some missing raw data. A dependency of stress energy with seismic moment was found throughout the studied magnitude range (e.g. Boatwright 1994; Li et al. 1995), whereas the magnitude of the consequences increases the energy release. Within each data set, there is a trend towards increasing energy with increasing earthquake size. The release of energy must correlate with the change in principal stress direction and the magnitude of an earthquake (Scholz 1990). Our study delineates the continuous oscillation of strain release of different magnitudes in Bangladesh. It illustrates a particular range of magnitude, providing different releases of strain continuously following the oscillation pattern. It shows that the same event can repeatedly happen after a certain amount of time.

5. Conclusion

The present study assessed the seismicity in Bangladesh for the first time using the proposed wavelet, short-time Fourier, and multitaper models. We showed the constant values of a and b in the GR relationship, the total amount of strain released during the study period, and the spectrograms of magnitudes and strain releases. The results showed that the same 4.5–5.8 magnitude earthquakes happened from 1992 to 2017. The results of the frequency-magnitude analysis showed a linear line that illustrates the same magnitude range and the higher magnitude range that is possible in the near future. In contrast, strain release also clearly demonstrated the upward trend of the linear line. We found that the year with a larger magnitude propagates higher strain, which explains the proportional relationship between magnitude and strain release. The obtained spectrograms reveal a moderate magnitude ranging from 4.5 to 5.6 with red patches within continuity. The multitaper model analysis showed the F value and the degree of freedom as an apparent anomaly, indicating a data deficiency. The seismic linear line in the frequency-magnitude relationship suggests a possible future major earthquake occurrence. Moreover, the spectrograms also revealed that the magnitude ranges of 4.5–5.8 may trigger a major earthquake event in the future. We also observed that inactivity might be the sign of a large-magnitude earthquake or a precursor to major earthquakes. The historical earthquake catalogue shows a significant earthquake in every century. According to our findings, there is a chance of a significant earthquake in the twenty-first century. Right now, the northeastern part of Bangladesh is more likely to have an earthquake than other parts.

More research needs to be done on the wavelet, short-time Fourier, and multitaper models to determine how earthquakes, crustal stress, fault rupture, and likely strain release are related to specific faults. This could lead to identifying barriers to earthquake recurrence rates and characterizing seismic hazards. Besides, a user-friendly operational software package can broaden the applications of these models for solving seismicity problems. More case studies are needed to demonstrate its applicability in the future.

Acknowledgements

We are also greatly acknowledged Begum Rokeya University, Rangpur, for all sorts of support. We greatly acknowledge the Bangladesh Meteorological Department (BMD) for providing the necessary datasets for this research

Data availability statement

Data are available based on the reasonable request of the corresponding author.

Additional information

Funding

Authors thankfully acknowledge the Deanship of Scientific Research for proving administrative and financial supports. Funding for this research was given under award numbers RGP2/363/44 by the Deanship of Scientific Research; King Khalid University, Ministry of Education, Kingdom of Saudi Arabia.