Hydrogeological Properties Estimation from Groundwater Level Natural Fluctuations Analysis as a Low-Cost Tool for the Mexicali Valley Aquifer

: Data on the hydrologic properties of aquifers are frequently not available or are spatially limited; additionally, their determination through aquifer tests is often logistically complicated and economically expensive. This study aimed to estimate aquifer properties by analyzing the water level response for the effects of barometric pressure and earth tide. Harmonic analysis of the time series of water level and barometric pressure recorded in three boreholes in the Mexicali aquifer provided reasonable values of porosity, speciﬁc storage, transmissivity, and compressibility of the rock materials that constitute the alluvial aquifer. The representative values of porosity (14–20%), speciﬁc storage (1.74–6.23 × 10 − 6 m –1 ), transmissivity (8.57–8.66 × 10 − 7 m 2 · s –1 ), and compressibility (3.90–8.21 × 10 − 10 Pa –1 ) were obtained. These values were consistent with the sediment types identiﬁed in the proximity of the wells analyzed. The results of this study show that the analysis of water level response to natural phenomena is a low-cost tool that provides reasonable estimates of aquifer properties. This advantage is particularly relevant in the study of aquifers where the available hydrological information is insufﬁcient.


Introduction
The detailed recording and analysis of the water-table level fluctuations in monitoring wells in response to natural and anthropogenic phenomena is an essential tool for hydrogeologists and geophysicists since their response mainly depends on the properties of the rock materials that constitute the aquifer. The water level in a well is highly sensitive to different processes such as precipitation seepage, evapotranspiration and induced recharge, regional flow and pumping, barometric pressure variations, earth tides, and seismic activity [1][2][3][4][5][6][7]. Some of these phenomena are responsible for inducing enough stress to deform the structure of the aquifer. Their effects directly modify the pore pressure and are commonly manifested as variation in the water level. In particular, periodic changes in the groundwater level are a result of two independent but related natural phenomena, barometric pressure and earth tides.
Under conditions of perfect confinement and high lateral permeability of the aquifer, fluctuations in the water level related to earth tides and barometric pressure variations are expressed through two constants, known as tidal sensitivity and barometric efficiency, respectively [8]. These parameters can be used to estimate geohydrological and poroelastic properties of aquifer materials [1,3,4,6]. However, Water 2018, 10, 586 2 of 15 just a constant is not enough to describe the frequency-dependent response of the water level to barometric pressure and earth tides for most aquifers [9][10][11].
Complementary parameters such as porosity, specific storage, compressibility, and Skempton's coefficients have been used for studying the aquifer material deformation caused by natural and anthropogenic phenomena, assuming linear and elastic deformation of the aquifer. Additionally, the response to solid-Earth tide can be used to identify the formation confinement level [12].
The conventional techniques for estimating these properties are the field aquifer test (pumping and slug tests) or laboratory tests of samples of rock material from the aquifer. However, such tests are typically expensive, time-consuming, and somewhat limited in areal extent [13]. Therefore, the assessment of aquifer properties from analysis of water table level variations as a result of barometric pressure fluctuations and earth tides provides a feasible and economical tool for improving aquifer parameter knowledge, particularly for aquifers with insufficient hydrological data.
Since 2000, the Colorado River Basin has experienced the driest 16-year period in over 100 years [14]. This period also ranks as the fifth driest 16-year period in the last 1200 years [15]. As a result, the most overallocated river system in the world is affected by a severe drought [16,17]. Therefore, the hydrologic components that constitute the Colorado River Drainage Basin require exhaustive multidisciplinary studies that provide information to evaluate environmental and economic drought effects. The Mexicali aquifer in the Delta of the Colorado River Basin has been studied and assessed in previous hydrogeological studies [18][19][20]. Mexicali Valley aquifer properties and poroelastic parameters, such as hydraulic conductivity, storage coefficient, hydraulic diffusivity, and static volume strain efficiency have been obtained, calculated, and compiled [18,[21][22][23]. However, there are many data gaps, especially concerning porosity, specific storage, transmissivity, compressibility, and Skempton's coefficients because of the high costs in terms of time and effort.
Estimates of the aquifer properties are essential in understanding the local and regional groundwater flow. Additionally, the hydrological properties derived before and during the drought are of great importance for assessing the water resources management strategies for resilience to severe drought in the Colorado River Basin. The purpose of this work was to estimate the porosity, specific storage, compressibility, storage coefficient, transmissivity, and Skempton's parameters from the analysis of groundwater response to barometric pressure and solid-Earth tide. The estimation was done to demonstrate that this is a feasible tool for filling gaps in aquifer parameter knowledge by applying this technique to existing water-level and barometric pressure data from a set of boreholes in the alluvial aquifer of Mexicali Valley.

Barometric Pressure Analysis
Water table variations show an inverse and proportional correlation to barometric pressure fluctuations. Water table level variations (WL) are related to barometric pressure changes (BP) through the barometric efficiency (BE), which can be obtained as follows [24]: Atmospheric pressure variations generate a uniformly distributed stress on the ground surface. Part of the stress is absorbed by deformation of the materials that constitute the aquifer rock and the rest is transmitted to the fluid in the porous medium [3]. If the materials exhibit a high transmissivity or specific yield, a drained condition is favored (i.e., mass transfer through flow), so that a response to atmospheric pressure fluctuations may not be observed [25]. However, it is common practice to consider the lateral flow as negligible because of the vast extension of rock formation and the almost uniform effect of atmospheric loading on the ground surface. Water level variations in the borehole can be conceptualized as changes in the aquifer pore pressure, except with a well open to the atmosphere when the barometric pressure also exerts a uniform pressure on the water surface [26]. In the case of a well penetrating an unconfined aquifer, barometric pressure variations will affect the water level within the well instantaneously. Meanwhile, if this effect is observed in the aquifer, it will show a delayed time because of the presence of air and other gases contained in the unsaturated zone. This lag time creates pressure differences between the well and the aquifer which generate outflow and inflows to the borehole, and therefore water-level fluctuations [12].
Many methods have been developed to estimate barometric efficiency. Some of them consider it independent of frequency and calculated it by using linear regression techniques [24,[27][28][29], while other methods consider barometric efficiency as frequency dependence, through transfer functions computed in just some segments or along the entire spectrum of frequencies. Additionally, barometric pressure can be evaluated alone or coupled with solid-Earth tide [8][9][10]30].

Solid-Earth Tide Analysis
Solid-Earth tides are small, periodic variations in the earth's shape due to forces of expansion and compression caused by the gravitational attraction of the celestial bodies, mainly the moon and the sun. These gravitational forces are balanced by changes in pore pressure in aquifer materials. The induced stress modifies the groundwater level in the wells located in both confined and unconfined aquifers [3,25]. Pore pressure (p) is related to applied stress (σ 33 ) using the tidal efficiency (γ e ), which is obtained as follows [1]: The fluid contained in the porous medium responds to a three-dimensional stress field, σ v . However, for the imposed deformation by solid-Earth tide or tectonic events, the three-dimensional stress field is not well-known as a priori knowledge. The vertical deformation observed on the terrain surface is approximately equal to one component of the horizontal tidal strain but with an opposite sign. The volumetric deformation of the surface rocks expressed as the sum of normal components is approximately equal to one of the horizontal tidal strain components, with an approximate amplitude of 1 × 10 −8 [26]. Thereby, it is better to examine the water level response to the areal strain (ε a = ε 11 + ε 22 ). Additionally, because σ 33 presents an opposite sign, ε v ε a , which implies that the response of water level to solid-Earth tide will be substantially lower.
Strain field due to solid-Earth tide can be estimated from the theoretical gravitational potential [31]. This theoretical strain field may differ from the real field one, mainly because of the local effects of geological faults and topography [32,33]. Geological and topographic influence is difficult to correct as a priori knowledge, and in the absence of measured strain, it is recommended that the theoretical value of strain field be used [30].

Aquifer Properties Estimation
In 1950, Jacob [34] derived expressions that relate barometric and tide efficiencies to the elasticity of the medium and can be written as: where β k is the rock-matrix compressibility (L·T 2 /M) (reciprocal of the rock bulk modulus K k ( M/L·T 2 ); ϕ is the porosity (%/100); and β f is the water compressibility ( L·T 2 /M), where units are expressed in fundamental units. If both expressions are added, the result is unity, therefore, in calculating any of these parameters, it is possible to define the other (BE = 1 − γ).
Poroelastic properties of the rock formation can be estimated because the response of the water level to atmospheric loading and solid-Earth tide stress is a result of rock aquifer deformation in undrained conditions, assuming a static-confined response of the rock formation [3,4,6,35].
Rojstaczer and Agnew [6] defined the areal strain sensitivity A s as the response of the water level W L to the areal strain ε a as: where β u is the rock-matrix compressibility in undrained conditions ( L·T 2 /M); and v u is the Poisson's ratio of the formation in undrained conditions and can be expressed as follows [36]: where B is the Skempton's coefficient, dimensionless; ρ is the fluid density (M/L 3 ); g is the gravity acceleration (L/T 2 ); v is the drained Poisson's ratio; and α is the A Skempton's coefficient, dimensionless, also known as the Biot-Willis coefficient. The Skempton's coefficients can be calculated using the following expressions: Once values for barometric efficiency and areal strain sensitive have been calculated, the compressibility, Skempton's coefficient, porosity, and specific storage (S s (1/L)) are computed as follows [6,37]: The estimation of β k , B, ϕ, and S s is through an iterative process and requires a priori values of On the other hand, if the compressibility of the rock materials that constitute the aquifer is not considered (β u = 1/K u = 0), it is possible to determine the porosity and specific storage based on the water level fluctuations due to solid-Earth tide and barometric pressure. Water level variation produced by aquifer dilatation as a result of solid-Earth tide is a function of the specific storage of the aquifer and can be calculated by measuring water level fluctuations (dh), and assuming the drained Poisson's ratio [3], the specific storage can be obtained from the following equation [4]: where, h = 0.6032 and = l = 0.0839 are the Love numbers [38]; and a is the Earth's radius (L).

of 15
The ratio of the theoretical gravitational potential dW 2 to the water level variation dh, is proportional to the ratio of the amplitude of the harmonic component of the theoretical gravitational potential A 2 (τ, θ) to the amplitude of water level variation A w (τ) at the same period (τ). This is considering that the derivative can be approximated by a finite differential, small finite change ∆W 2 and ∆h in a small period of time ∆t, the specific storage can be expressed as follows [7]: where A 2 (τ, θ) is given by: The general lunar coefficient K m , the amplitude factor b for each harmonic component with a period τ, and the latitude function f (θ) values were obtained by [7]. The classical method of representing the theoretical gravitational potential is to represent it as composed of a finite set of harmonic functions [39]. Each k-tidal harmonic component has a distinct amplitude A k , frequency f k , and phase Φ k . Melchior [40] concluded that only five of them are associated with fluctuations in groundwater level, and these harmonic components are responsible for 95% of the variation of gravitational potential [41]. They are the M 2 and N 2 semidiurnal lunar tides, the S 2 semidiurnal solar tide, the O 1 diurnal lunar tide, and the K 1 diurnal lunar-solar tide.
Estimates of amplitude A wk and phase ϕ wk for water level variations at frequencies that correspond with the components of the theoretical tidal potential are calculated from the regression coefficients a wk and b wk obtained from implementing regression techniques that minimize the mean square error, where the amplitude and phase are obtained as follows: Therefore, the areal strain sensitivity for a particular component A sk is given by: For incompressible aquifer rock materials, the change in the volume of the aquifer in response to a variation of strain could be approximated as a change in porosity [1]. This assumption is suitable for many of the types of aquifers studied in hydrogeology, except those aquifers present in low porosity rocks [3]. Thus, the porosity can be obtained as: Analyses of water level fluctuations have generally focused on the amplitude response. Cooper et al. [2] showed that the harmonic amplitude response depends on: the transmissivity (T), storage coefficient (S C ), radius of the well casing (r c ), radius of the screened or open portion of the well (r w ), periodicity of the pressure head disturbance (τ), and the inertial effects of water in the well. From these relations, a set of dimensionless parameters was derived by Hsieh et al. [42], as follows: Sr w 2 r c 2 (21) However, the amplitude of the response is generally different from that of the disturbance, and there is also a shift in phase. Hsieh et al. [42] derived expressions to estimate T from the time lag between the earth tide dilatation of an aquifer and water level response in a well. The lag in time, referred to as the phase shift (η k ) of the kth tidal constitute, is determined by [42]: The least square fitting procedure to estimate the regression coefficients, Equations (16) and (17), was applied to determine Φ tk for the O 1 and M 2 lunar components in the calculated areal strain tide. Graphs of the amplitude ratio and phase shift as a function of Equation (20) for selected values of the parameter in Equation (21) were obtained by Hsieh et al. [42]. Values of dimensionless transmissivity can be estimated and converted into standard units (L 2 /T) if the phase shift and an order of magnitude estimate of S C are calculated.

Study Area
The Mexicali Valley is located northeast of the peninsula of Baja California. Its main economic activities are agriculture and geothermal electric power generation. The primary sources of water for irrigation are surface water from the Colorado River and groundwater from the Mexicali aquifer. The geothermal steam was obtained from the geothermal reservoir (confined and deeper aquifer) at the Cerro Prieto Geothermal Field located in the central west of the Mexicali Valley.
Tectonically, both the shallow aquifer and reservoir are located in the Salton Trough, which is part of the San Andreas-Gulf of California fault system that corresponds to the boundary between the Pacific and North American tectonic plates [43]. On a local scale, the Salada-Cucapá, Cerro Prieto, and Imperial faults created a basin up to 5000 m deep filled with sediments supplied by the Colorado River, as well as transported eroded debris from the Colorado Plateau basin margins [44,45]. The fault systems are currently active and have a significant amount of associated seismic activity [46,47]. The transmissive sediments in the basin have been categorized into two main units-consolidated and unconsolidated sediments-separated by strata of a very low permeability because of hydrothermal alteration [48,49]. The consolidated deepest unit consists of sandstones that comprise the geothermal reservoir. The shallow unconsolidated unit consists of fine to coarse sands with intercalations of gravel, clays, and silts. The shallow unit has a variable thickness of 400-2500 m and contains the aquifer of the Mexicali Valley [50]. Groundwater pumping, irrigation seepage, return flows, and recharge from episodic Colorado River flows are the primary hydrologic processes affecting the potentiometric surface of the Mexicali Valley aquifer [51]. The Mexicali Valley aquifer is an unconfined and non-homogeneous aquifer of variable thickness, formed by a sequence of unconsolidated granular sediments mostly of deltaic origin. Discontinuous layers of materials with different permeabilities cause local confinement, but together behave as a single hydrogeological unit [21]. The general map of the study area including tectonic, geology, and hydrology context is shown in Figure 1.
affecting the potentiometric surface of the Mexicali Valley aquifer [51]. The Mexicali Valley aquifer is an unconfined and non-homogeneous aquifer of variable thickness, formed by a sequence of unconsolidated granular sediments mostly of deltaic origin. Discontinuous layers of materials with different permeabilities cause local confinement, but together behave as a single hydrogeological unit [21]. The general map of the study area including tectonic, geology, and hydrology context is shown in Figure 1.

Data
A set of three monitoring wells (C-03, PZ-01, and PZ-03) located in the central-west area of the Mexicali Valley were equipped with pressure transducer data loggers for recording the water level every 5 min for 83 days, from 23 May-14 August 2007, and their location is shown in Figure 1. The boreholes were drilled into Quaternary alluvial deposits, which mainly consisted of water-saturated sand, gravels, and clays. The main characteristics of the well design are presented in Table 1. The water table level fluctuations were recorded using pressure transducers (Solinst Levelogger ® , Solinst Canada Ltd, Georgetown, ON, Canada) with a measuring range of 5.00 m and accuracy of 0.05% of the full range. Barometric pressure fluctuations were also recorded every 5 min using a Solinst Barologger ® (Solinst Canada Ltd, Georgetown, ON, Canada) barometer installed in PZ-03 with a measuring range of 1.50 m and an accuracy in the order of 0.001 m. Water level and atmospheric pressure time series are shown in Figure 2a.
The sampling frequency election was based on the water level response to the seismic activity on Mexicali Valley, analyzed and reported in a previous study [23]. Detailed inspection of the water level did not show a significant response to local seismic activity (local magnitude ≤ 4.0) recorded by IRIS [52]. Therefore, in this work, we assumed that groundwater fluctuations are only caused by barometric pressure and solid-Earth tide effects.
The areal strain associated with the theoretical gravitational potential for the vicinity of each well, expressed in nanostrain (1 nstr = 1 ppb), was calculated using Some Programs for Ocean-Tide Loading (SPOTL, version 3.3.0.2) [53,54]. This computer code does not consider geological and topographic discontinuity effects for the calculation; also, the ocean tide was ignored. The areal strain is shown in Figure 2b. level did not show a significant response to local seismic activity (local magnitude ≤ 4.0) recorded by IRIS [52]. Therefore, in this work, we assumed that groundwater fluctuations are only caused by barometric pressure and solid-Earth tide effects.
The areal strain associated with the theoretical gravitational potential for the vicinity of each well, expressed in nanostrain (1 nstr = 1 ppb), was calculated using Some Programs for Ocean-Tide Loading (SPOTL, version 3.3.0.2) [53,54]. This computer code does not consider geological and topographic discontinuity effects for the calculation; also, the ocean tide was ignored. The areal strain is shown in Figure 2b.  Table 1. Summary of some basic information about the three monitoring wells. Explanation: Well-Identifier, Well-ID; Land Surface Elevation, LSE (meters above sea level, MASL); Borehole Depth, BD (meters); Water Level Elevation, WLE (MASL); radius of the well casing r c (meters); radius of the well screened r w (meters); and Length of well screened, Lws (meters).

Results
The time series of recorded groundwater level and barometric pressure, as well as calculated theoretical solid-Earth tide, were processed and analyzed by specifically written computational codes for this study. The low-frequency trend was removed from the water level and barometric pressure data to ensure that time series were stationary. Subsequently, periodic fluctuations in time series were identified using the Fourier discrete transform technique.
Amplitude spectra obtained from the water level of the well PZ-01, barometric pressure, and theoretical areal strain as a function of the frequency in cycles per day (cpd) are shown Figure 3. Groundwater level fluctuations in PZ-01 are of diurnal and semidiurnal frequency, according to the highest amplitudes on the frequency spectrum at 1.0 and 2.0 cpd, respectively. The dominant components determined in the barometric pressure spectrum are typically associated with heating and cooling of the air column near the ground surface caused by solar radiation. The five components identified in the amplitude spectrums of water level, barometric pressure, and areal strain correspond to the five major harmonics that have been reported to be responsible for 95% of the theoretical gravitational potential variation and play an important role in groundwater studies [3,41].
identified using the Fourier discrete transform technique.
Amplitude spectra obtained from the water level of the well PZ-01, barometric pressure, and theoretical areal strain as a function of the frequency in cycles per day (cpd) are shown Figure 3. Groundwater level fluctuations in PZ-01 are of diurnal and semidiurnal frequency, according to the highest amplitudes on the frequency spectrum at 1.0 and 2.0 cpd, respectively. The dominant components determined in the barometric pressure spectrum are typically associated with heating and cooling of the air column near the ground surface caused by solar radiation. The five components identified in the amplitude spectrums of water level, barometric pressure, and areal strain correspond to the five major harmonics that have been reported to be responsible for 95% of the theoretical gravitational potential variation and play an important role in groundwater studies [3,41]. The components labeled as K 1 and S 2 in the areal strain spectrum are also present and dominate the barometric pressure spectrum. Therefore, water level analysis is complicated because the groundwater response at these frequencies is caused by periodic fluctuations in atmospheric pressure, as well as gravitational potential. As a result, K 1 and S 2 are frequently ignored in the study of groundwater [26]. Moreover, water level amplitude in the N 2 frequency component is regularly much lower compared with other principal frequencies. Therefore, the signal-to-noise ratio is low and is often discarded because of being a significant source of error in the analysis. According to the previous criteria, the component frequencies O 1 and M 2 are often the only ones used for groundwater level response analysis.
The presence of components with frequencies associated with earth tide O 1 and M 2 in the water level spectrum of our data can be interpreted as a delay in the atmospheric amplitudes with frequencies K 1 and S 2 when these pass-through the vadose zone, which is like the behavior of a semiconfined aquifer. The identification of components associated with solid-Earth tide frequencies on the groundwater level spectra allows us to interpret semi-confined behavior in the proximity of the analyzed wells. The barometric efficiency was computed according to Rahi's method and listed in Table 2; the Rahi's scheme filters the water level, as well as the barometric pressure data, to remove the effects of the earth's tides from water-level fluctuations. Hence, only barometric pressure induced water-level variations and the barometric fluctuations that produce them are employed for the calculation of the barometric efficiency. Therefore, we proceed to calculate porosity and specific storage for incompressibility media, hereafter, to complement the calculation considering compressibility of the rock materials that constitute the aquifer. Table 2. Barometric efficiency (BE), water level amplitudes (A wk ), areal strain tide amplitudes (A tk ), areal strain sensitivities (A sk ), and phase shift (η k ) for each k − th tidal harmonic component (O 1 and M 2 ) for each well. Water level variation amplitudes at solid-Earth tide frequencies O 1 and M 2 were calculated. The water-level data were processed for solid-Earth tide analysis; unwanted high and low frequencies were removed using digital filters. High-frequency noise, greater than 10 cycles per day, was removed using a Chebyshev low pass filter with a cutoff frequency equal to 8.0 cpd. After that, a Chebyshev high pass filter with a cutoff equal to 0.5 cpd was used to remove low-frequency water level variations. Amplitudes A wk and phases Φ wk of the water-level filtered were computed at their exact tidal frequencies from the regression coefficients obtained by fitting the water level data to a sum of sines and cosines functions using the "t_tide" Matlab code [55]. Amplitudes A tk and phases Φ tk of the calculated areal strain tide were also computed using the regression coefficients obtained by the "t_tide" code and using Equations (16) and (17). Barometric efficiency values, amplitude, and phase for the O 1 and M 2 tidal constituents in the water level record and areal strain time-series are listed in Table 2. Areal strain sensitivity was calculated using Equation (18) and is also shown in Table 2. The tphase shift will be utilized below in the calculation of transmissivity.

Well ID
Once the barometric efficiency and areal strain sensitivity were determined, the specific storage assuming incompressible rock materials (β u = 1/K u = 0) and typical Poisson's ratio for geological materials v = 0.25 could be obtained from Equation (14). Table 3 shows the range of values of estimated S S . Table 3. Summary of aquifer properties estimated assuming incompressible rock materials (β u = 1/K u = 0): specific storage (S S ), porosity (Φ), storage coefficient (S C ), and transmissivity (T) for each O 1 and M 2 tidal harmonic component for each well. Porosity values as a function of barometric efficiency and specific storage were computed using Equation (19). In this analysis, water compressibility β f = 4.40 × 10 −10 Pa −1 ; gravitational acceleration g = 9.79 m·s −2 ; and water density ρ = 995.71 − 997.86 kg·m −3 , depending on the average temperature of the fluid recorded on each well, were considered. Values of porosity assuming incompressible rock materials are listed in Table 3.

Well ID
Using the specific storage computed and considering the saturated thickness as the length of the screened section of the borehole shown in Table 1, the storage coefficient (S C ) was computed from the relation (S C = S S ·b), and estimated values are shown in Table 3.
Determining rock-matrix compressibility (β K ) from Equation (9) involves an iterative procedure where an initial guess at rock-matrix compressibility in undrained conditions (β u ) is made. A value of β u = 2.00 × 10 −10 Pa −1 for a typical rock formation was assumed according to Rojstaczer and Agnew [6]. Equations (8) and (10) are used to compute Skempton's coefficients α and B, and then Equation (9) is used to estimate an improved value of β K . This procedure is repeated until the equations converge to a stable value of β K , α, and B.
Specific storage and porosity assuming compressibility can be obtained from Equations (11) and (12), respectively. Values for g, ρ, β f , and v described previously were assumed. Estimations of β k and B are a function of barometric efficiency and areal strain sensitivity determined.
The hydrogeological and poroelastic properties of the rock materials in each monitoring well such as rock-matrix compressibility (β k ), specific storage (S S ), porosity (Φ), storage coefficient (S C ), transmissivity (T), and Skempton's coefficients B and α, for O 1 and M 2 earth tide components, were calculated and are shown in Table 4.

Discussion
The areal strain sensitivity for each solid-Earth tide component of our data set was similar, and therefore, the specific storage value obtained is consistent. The specific storage mean value estimated for the study area of 1.74 × 10 −6 m −1 is in accordance with the type of Mexicali Valley aquifer lithology (sand-gravel-clay) and also with theoretical values reported S S = 10 −5 − 10 −6 m −1 [56]. Likewise, the obtained porosity mean value of 20% is in accordance with that of aquifers with the same lithology.
The representative compressibility of the rock materials that constitute this aquifer was calculated as 3.90 − 8.21 × 10 −10 Pa −1 . The obtained values suggest a predominance of gravel over sand-clay sediments, according to Freeze and Cherry [57]. An approximately 6% lower than average porosity was obtained. The significant reduction that occurred in proximity to well PZ-01 was probably related to a higher rigidity of rock materials in the formation surrounding the borehole. Under compressibility conditions, the average specific storage calculated was around four times higher than the value under non-compressibility conditions; however, the order of magnitude was the same. Conversely, the storage coefficient increased by one order of magnitude; therefore, the estimation of transmissivity was recalculated. The average value obtained for the Skempton's coefficients was B = 0.56 and α = 0.66. These values mean that the main induced deformation is absorbed by the fluid; and also what, the rock compressibility is insignificant, such as in the case of unconsolidated sediments. Notably, the Skempton's B coefficient calculated in this study showed better correspondence with the alluvial sediments surrounding the boreholes and is comparable with the value assumed for the same area by Sarychikhina et al. [22].
Despite not having in situ measurements of the real areal strain and considering the dependence of the barometric efficiency estimation from the method of computation used, the hydraulic aquifer properties values obtained-specific storage, porosity, and compressibility-are comparable with the typical reported values for the rock material that constitutes the aquifer [6,9,10,26,30]. Our results showed that the upper limit value of porosity was obtained considering the incompressible rock structure and the lower limit value corresponding to compressible rock structure. In contrast, the upper bound value for specific storage was obtained by considering compressible rock structure, and the lower bound is associated with the incompressibility of the rock structure. The transmissivity values obtained suggest that the method for estimating it is insensitive to the assumption of rock-matrix compressibility.

Conclusions
The analysis of the groundwater response to the deformation provided essential information on aquifer properties. Water table level and barometric pressure data were used from a set of wells drilled in the shallow aquifer of Mexicali. The methodology presented in this paper offers an accessible and inexpensive tool using just water level and barometric data recorded to obtain a reasonable approximation of hydrogeological and poroelastic parameters. In particular, aquifer rock properties, such as specific storage, porosity, and compressibility, were obtained by applying time series analysis techniques to the water table response using barometric pressure and computation of theoretical areal strain tide. These aquifer properties are needed in design strategies for sustainable exploitation of the groundwater resources of alluvial aquifers in low-income countries, where there are often high costs for hydrology tests on a vast areal extension.
The hydrogeological characterization of the Mexicali aquifer rock formation surrounding the three wells analyzed was comparable to semi-confined aquifer properties and, for the similar depositional environment, was as expected. Representative values for the porosity 14 − 20%, specific storage 1.74 − 6.23 × 10 −6 m −1 , transmissivity 8.57 − 8.66 × 10 −7 m 2 ·s −1 , and compressibility 3.90 − 8.21 × 10 −10 Pa −1 were obtained and these can be used as a benchmark because there are a lack of hydrological parameters for this aquifer area.
Future research should analyze data from a larger set of monitoring wells to extend the geohydrological and poroelastic characterization of the Mexicali Valley aquifer. Also, the research should explore geohydrological-poroelastic relationships to determine geomechanical properties associated with the rock materials that constitute the aquifer formation.
Author Contributions: Mario A. Fuentes-Arreazola carried out the analysis presented in this study and drafted the text of the manuscript. Jorge Ramírez-Hernández was fundamental in the writing of the manuscript and also formulated valuable scientific ideas. Rogelio Vázquez-González contributed the data used in this work and collaborated during the entire process with ideas, corrections, and advisory times.