Assessing the Impact of Brackish Water on Soil Salinization with Time-Lapse Inversion of Electromagnetic Induction Data

Abstract

Over the last decade, electromagnetic induction (EMI) measurements have been increasingly used for investigating soil salinization caused by the use of brackish or saline water as an irrigation source. EMI measurements proved to be a powerful tool for providing spatial information on the investigated soil because of the correlation between the output geophysical parameter, i.e., the electrical conductivity, to soil moisture and salinity. In addition, their non-invasive nature and their capability to collect a high amount of data over broad areas and in a relatively short time makes these measurements attractive for monitoring flow and transport dynamics, which are otherwise undetectable with conventional measurements. In an experimental field, EMI measurements were collected during the growth season of tomatoes and irrigated with three different irrigation strategies. Time-lapse data were collected over three months in order to visualize changes in electrical conductivity associated with soil salinity. A rigorous time-lapse inversion procedure was set for modeling the soil salinization induced by brackish irrigation water. A clear soil response in terms of an increase in electrical conductivity (EC) in the upper soil layer confirmed the reliability of the geophysical tool to predict soil salinization trends.

1. Introduction

Soil salinization has become one of the major environmental and socioeconomic issues globally, and this is expected to be exacerbated further by projected climatic change. Salinity is one of the main soil threats that reduces soil fertility and affect crop production because it can decrease plant growth and water quality, resulting in lower crop yields and degraded stock water supplies [1,2,3,4]. Soil salinization increases when the overexploitation of groundwater in coastal areas leads to pumping from a brackish or saline irrigation source due to saltwater intrusion. In addition, treated wastewater, increasingly used in water-scarce environments to tackle climate change, can have a significant salt concentration, depending on the treatment strategy. In such conditions, it is crucial to develop soil monitoring systems that are able to capture spatial and temporal dynamics with a high degree of accuracy.

Sensors typically used for agricultural purposes are installed in a few sparse points at no more than two or three previously defined depths. Likewise, it is unrealistic to collect over time a high number of soil samples for deriving the electrical conductivity of a saturated soil extract (ECe), which is the most useful and reliable measure of soil salinity. In the last decade, apparent electrical conductivity, defined as ECa, has been increasingly used for investigating soil properties. ECa is a measure of the bulk electrical conductivity of the soil and is influenced by various factors, such as soil porosity, the concentration of dissolved electrolytes, texture, the quantity and composition of colloids, organic matter, and water content in the soil [5]. ECa is an apparent value because it has the equivalent electrical conductivity of a homogeneous half-space whose measured response depends not only on the mentioned soil properties but also on instrumental configurations.

Geospatial ECa measurements are very fast and easy to collect because they do not require ground contact and, for this reason, they allow a range of spatial coverage from a few meters to several hectares at different depths, depending on the target to be investigated and the electromagnetic sensor to be used. In addition, collecting ECa data over time along the same transects allows temporal variations in the ECa to be correlated with soil properties that change over time. Over time, the ECa parameter has been widely used as a proxy for investigating soil salinity, as recent studies confirm. Corwin et al. [6] improved guidelines to broaden the scope of the application of ECa-directed soil sampling to map field-scale salinity on orchards under drip irrigation. Vanderlinden et al. [7] highlighted the potential of ECa to perform salinity monitoring at the field or farm scale. Emmanuel et al. [8] validated the potential of extending the electromagnetic induction (EMI) to characterize wetland soil properties, such as salinity, improving sampling plans, and interpolating soil property estimates to unsampled regions. Scudiero et al. [9] performed the ECa survey to identify locations that can be repeatedly sampled to infer the frequency distribution of soil salinity. Rodrıguez-Perez et al. [10] evaluated the usefulness of apparent electrical conductivity (ECa) data to identify variations in soil chemical and physical properties and moisture content. Lesch et al. [11] documented the spatial salinity mapping using EM survey data. Herrero and Hudnall [12] proved the value of electromagnetic induction (EMI) to map the salinity of the rootable layer.

Although ECa is widely used in agronomic applications, its precise indication remains somewhat unclear. What does ECa mean? Is it really representative of the soil properties? Can ECa be effectively related to the saturated past extract of ECe?

ECa is a depth-weighted parameter giving limited information about the variation in the conductivity with depth. In fact, ECa does not provide a rigorous correlation between the Earth’s electrical conductivity (EC) structure and measured responses, which are affected by several factors, such as coil distance and orientation, sensitivity, and data error.

In the last few years, an inversion procedure has been implemented in order to produce a reliable soil EC model. The use of inversion codes has grown rapidly as the need for an effective EC distribution in the subsoil has become crucial. Several inversion methods [13,14,15] and software tools [16,17] have been created to estimate the distribution of soil EC from measured ECa data. These codes solve the complex equations of Maxwell’s electromagnetism to generate forward models used to minimize an objective function and derive models from apparent raw data.

Jadoon et al. [18] identified a quantitative distribution pattern of soil salinity by the joint inversion of multicomponent EMI measurements. Paz et al. [19] highlighted time-lapse EMI as a reliable tool for evaluating the risk of soil salinization and supporting the evaluation and adoption of proper agricultural management strategies. Farzamian et al. [20] inverted ECa data to model the salinity of an oasis in Tunisia, which may affect agricultural productivity and the sustainability of crop production. Dakak et al. [21] produced soil salinity maps at various depths through EMI inversion. Shaukat et al. [22] used EMI inversion as a robust and effective method for the risk assessment of new shrub plantations. A benefit of the inversion technique lies in its capacity to predict the depth-wise distribution of EC. This facilitates linking laboratory soil samples obtained from different depths at specific sites with the corresponding EC values at those depths. As emphasized by [20], this method permits the incorporation of all soil samples into a single calibration/validation procedure, leading to a more accurate calibration at any chosen depth. Furthermore, employing a larger number of soil samples for calibration and validation enhances the reliability of the outcomes. Many of these references show a high correlation between ECa or inverted EC against ECe, confirming their potentiality in agricultural applications.

In the proposed case study, a time-lapse inversion approach was tested on ECa data collected on an experimental field in the South of Italy where different irrigation strategies were applied during the growing season of the tomato. The primary objective of time-lapse inversion is to accurately detect variations in conductivity at specific locations across different time intervals. While independent data inversions can be conducted individually, providing insight into changes in modeling results over time through the subtraction of pixel-by-pixel values from a reference dataset, it is important to recognize that temporal changes in conductivity values may not exclusively reflect actual changes in subsurface conductivity depending on the data noise level and inversion artifacts. This is particularly true without considering the reference model and prior information [23].

Raw ECa maps, as well as those detailed on the experimental setup and soil properties, have already been published [24]. In this paper, we focused on the inversion of ECa data, both as a single snapshot and time-lapse imaging. For the proposed case study, a 2D reference transect was extracted from three different plots subjected to different irrigation strategies (irrigation with freshwater and agro-industrial treated wastewater).

The aim of the paper was (a) to highlight the capability of the inversion tool to produce a detailed 2D EC soil distribution; (b) to image the spatio-temporal EC evolution distribution over an irrigation season; and (c) to assess the impact of brackish water on soil salinization in the short term.

2. Materials and Methods

2.1. Basics of ECa Parameter

ECa is a sensor-based indirect measurement that is strictly affected by some physical and chemical properties, such as soil salinity, soil moisture, clay content, and cation exchange capacity (CEC). Archie [25] and Rhoades et al. [26] developed a theoretical basis for the relationship between ECa and soil properties, such as the soil water content, the electrical conductivity of the soil water, soil bulk density, and the electrical conductivity of the soil particles.

According to these premises, over the last few decades, ECa has been widely used as a soil quality indicator [27,28,29].

Over time, several devices have been manufactured both in the time and frequency domain for measuring ECa according to the electromagnetic induction (EMI) theory. Since only frequency domain measurements were used in this work, a brief description of this method was reported. An EMI sensor is made of two coils, including a transmitter and a receiver.

A time-varying current (I_P) circulating in the loop coil T, named the transmitter coil, generates a time-varying magnetic field, which is in phase with the current and with the same rate of change according to Ampere’s law. This field, referred to as the primary magnetic field (H_P), induces electric currents in a conductive body, thus generating secondary EM fields. A receiver coil records a signal that is the sum of the primary and secondary fields (Figure 1).

The measured resulting field has an imaginary part of the signal, also called out-of-phase or quadrature (Q), and a real part of the signal, the in-phase (Ph) component. Under simplified conditions, typically defined as a low induction number (LIN), the EMI sensors directly provide the subsurface apparent electrical conductivity (ECa) through Equation (1) [30]:

$${\mathrm{EC}}_{\mathrm{a} }=\frac{2}{\pi {\mathrm{fs}}^2{\mu }_0}{\left(\frac{{\mathrm{H}}_{\mathrm{s}}}{H_p}\right)}_{Qu}$$

(1)

where f is the frequency (Hz), s is the coil separation (m), μ₀ is the magnetic permeability of free space (4π × 10⁻⁷ H/m), and (H_S/H_P) Q_u is the Q_u component of the secondary H_S to the primary H_P magnetic field coupling ratio.

Conversely, the real part or the in-phase component of the measured signal is mainly affected by the magnetic permeability of the subsoil.

The ground response depends not only on the soil’s electrical conductivity but also on instrumental factors, such as coil distance, coil orientation, and frequency. In particular, the coil distance and orientation affect the depth of penetration of the electromagnetic signal. Placing the coil in a vertical position (VCP coil configuration), the topsoil layers are investigated while rotating 90° of the coils along the main axis (HCP coil configuration), and the signal investigates deeper layers.

The instrumental configuration affects the soil response, as the analysis of the cumulative sensitivity (CS) function shows. This function, defined as the ratio between the variation in the output and the variation in the input [31], quantifies how much the complex electromagnetic response recorded by the device is affected by a variation in the conductivity and/or permeability of a particular point (area or section) of the subsurface.

Figure 2a,b plots the CS distribution against depth as a function of coil and orientation distances for the CMD-Mini-Explorer (GF Instrument s.r.o, Brno), which is the EMI sensor used for collecting field data, made of a cylindrical tube that is 1.3 m long, with a 30 kHz transmitter coil and three receiver coils with 0.32 m, 0.71 m, and 1.18 m offsets, respectively. As clearly observed, the sensitivity changes significantly for three different coil distances and coil orientations.

According to [30], given a specific configuration and under LIN conditions, the effective penetration depth of an EMI sensor corresponds to CS = 0.3. Therefore, ECa-derived depths are purely indicative because they depend not only on soil properties but also on instrumental factors.

On the basis of the aforementioned considerations, the assumption that ECa cannot be used to provide quantitative information about the soil properties confirms the need to invert raw data in order to obtain a reliable EC distribution in the subsoil.

2.2. Time-Lapse ECa Dataset

The dataset used for time-lapse inversion was collected during the growing season of a tomato crop belonging to a farm located in the South of Italy. The experimental field was randomized with three different irrigation treatments: (a) plots A and B were irrigated with agro-industrial-treated wastewater with different levels of fertigation; (b) plot C was irrigated with fresh water and conventional fertigation. Pre-transplant, fertilizers were applied to the soil in all investigated plots by distributing 30 kg ha⁻¹ N and 35 kg ha⁻¹ P. Throughout the crop cycle, for plots B and C, 75 kg ha⁻¹ N, 40 kg ha⁻¹ P, and 72.5 kg ha⁻¹ K were added through fertigation. Instead, 75 kg ha⁻¹ N and 40 kg ha⁻¹ P were added through fertigation in plot A. The EC of the irrigation water was about 2000–2500 µS/cm for plots A and B and 500 µS/cm for plot C. Tomato plants were covered with an anti-hail net and grown in a net house structure. Shading nets protect the plants from sunlight, thus allowing no significant soil temperature changes throughout the irrigation season. Irrigation water was supplied every two days through an underground irrigation system.

ECa data were collected at five time points approximately in the same conditions, i.e., soon before each irrigation event, thus ensuring a homogeneous and equalized moisture distribution across the nine parcels. The initial measurement was taken before the start of the irrigation season (26th June), which served as a reference time for recording the conductivity changes. The other four datasets were collected with time intervals of about 2 weeks (10th July, 24th July, 6th August, and 31st August) during the irrigation season. The data were collected in continuous measurement mode by selecting a time of measurement equal to 1 s, meaning that the conductivity and in-phase values were measured as the average of the values measured during the selected measuring period. The measurements were collected by hand, keeping the device as close as possible to the ground, almost trailing it on the ground, in order to minimize the air layer below the sensor.

During the last EMI campaign, on 31st August, six soil samples were collected at three points, each belonging to a single plot at two different depths, 0–0.2 m and 0.2–0.4 m, from the ground surface in order to provide the ground truth for soil salinity. EC was measured with a multi-range Cryson-HI8734 electrical conductivity meter (Crison Instruments, S.A., Barcelona, Spain). Details of the sampling procedure and soil analysis are described in [24]. In this paper, only the ECe values were extracted in order to provide a correlation function with the inverted EC. In order to test the inversion procedure and compare the findings, nine different 2D EMI transects corresponding to the locations of soil sampling (see Figure 3) were extrapolated from the three plots.

2.3. Time-Lapse Inversion Procedure

The time-lapse inversion procedure is aimed at estimating the electrical conductivity variations over time along the reference transects. The EM4Soil v4.5 software package [16], based on Occam’s regularization [32], was used to invert the ECa data. The code uses a quasi-2D inversion, assuming that below each measured location, a 1-dimensional variation in calculated soil conductivity (EC, dS/m) is constrained by variations under neighboring locations [13].

In the time-lapse inversion procedure, the choice of some parameters is crucial for providing accurate inversion models: (a) the starting model for solving the inverse problem; (b) the inversion algorithm to be used; (c) the spatial damping factor (λ); and (d) the temporal damping factor (α) [33].

Two inversion algorithms, S1 and S2 [34], provide a different level of constraint to the model parameters. In the S1 algorithm, the corrections to the model parameters at each iteration were calculated by solving the system of equations:

$${(\mathrm{J}}^{\mathrm{T}}\mathrm{J}+\lambda {\mathrm{C}}^{\mathrm{T}} \mathrm{C}) {\delta }_{\mathrm{p}}{=\mathrm{J}}^{\mathrm{T}}\mathrm{b}$$

(2)

where δ_p is the vector comprising the corrections of the parameters (logarithm of conductivities, p_j) of an initial model; b is the vector containing the differences between the logarithm of the observed and calculated apparent conductivities. J is the Jacobian matrix, and λ is the aforementioned damping factor.

Conversely, the S2 option has one more constraint (Equation (3)) and can produce smoother results than S1.

$${(\mathrm{J}}^{\mathrm{T}}\mathrm{J}+\lambda {\mathrm{C}}^{\mathrm{T}} \mathrm{C}) \delta {\mathrm{p}=\mathrm{J}}^{\mathrm{T}}\mathrm{b} +\lambda {\mathrm{C}}^{\mathrm{T}} \mathrm{C} (\mathrm{p} -{\mathrm{p}}_0)$$

(3)

where p₀ refers to a reference model.

The spatial smoothing or damping factor λ [35] determines the amplitude of the parameter corrections in the space domain, which controls the balance between the data fit and model roughness. The selection of λ is determined either through the “L curve” method [36] or through trial and error to determine the value that most accurately represents the expert expectations based on the study site and data fit. A smaller damping factor tends to refine more detailed model parameters, especially in areas with larger expected spatial variability [37]. The temporal damping factor α is a regularization factor that gives the weight for minimizing the temporal changes in the conductivity along the time [33]. As the α value increases, the resulting reference models from the inversion become more similar. A value of zero indicates that no temporal constraints are applied, resembling a traditional independent inversion.

A linear solution, based on the ECa cumulative response (CF) and a non-linear solution (FS), was used for forward calculations in order to convert depth-profile EC to ECa [30]. In our study case, both CF and FS modeling, as well as S1 and S2 algorithms, were tested with no substantial differences observed in the model output. The outputs obtained with CF modeling and the S2 algorithm are shown in the section Results section. A uniform starting model considering the average ECa value for each plot (ECa = 1.00 dS/m for plots A and B; ECa = 0.60 dS/m for plot C) was considered to solve the inversion problem, while λ and α parameters were set to 0.07 and 0.05, respectively. These parameter sets were selected after conducting several tests and comparing the results in terms of inversion misfit.

In the pre-processing stage, a single dataset containing all the ECa readings obtained over space and time was defined as input for the EM4SOIL code.

3. Results

The inversion was performed for the nine selected EMI transects. As the findings of the inversion results were consistent with the plots subjected to the same irrigation strategies, for simplicity, we show a single time-lapse EMI inversion for plots A, B, and C.

The inversion results were visualized as static 2D EC images and normalized time-lapse EC differences. For a clear comparison of the soil response to different irrigation strategies, both results are scaled with the same color range.

Figure 4 shows the inverted EC for the transect corresponding to plot A, irrigated with brackish water, at five different time points.

At the reference time t1 (26th June), a two-layer model was identified as follows: (a) a low conductive top layer (EC < 0.60 dS/m at z < 0.5 m from ground surface); (b) a high conductive bottom layer (0.20 < EC < 0.60 dS/m). During the irrigation season, the top resistive layer was conductive, while the bottom conductive layer experienced a decrease in conductivity.

Variations in EC over time are clearly observed in Figure 5. After two weeks from the start of irrigation (Figure 5a), slight changes could be observed, with positive changes in the yellow area and negative changes in the green area. Over time, EC increasing in the topsoil layer and decreasing in the bottom layer became more prevalent, as observed in Figure 5b. This trend intensifies in Figure 5c, which corresponds to the EC distribution after 41 days from the start of irrigation and persists until the end of the irrigation season (Figure 5d).

To evaluate the accuracy of the inverted model, the observed vs. calculated fit was visualized for each coil configuration at every observation point. (Figure 6). The RMSE between calculated and observed data provides the error model. Generally, the datasets align along the fit line, although some misfits are evident. This trend is not surprising given the level of accuracy of the ECa data. In fact, the estimated RMSE of the geophysical model is 0.15 dS/m, which is the misfit between the observed data (field data) and the theoretical one (calculated data). In terms of percentage value, it is roughly equivalent to about 16%, expressed as the ratio between the RMSE in dS/m and the average value of the observations.

Similar trends were recorded in the transect belonging to plot B. Figure 7 and Figure 8 correspond to the static 2D images and normalized time-lapse EC differences, respectively. Small differences are visualized in the bottom layer, where a markedly lower decrease is recorded. In addition, the main increase in EC in the top layer occurs in Figure 7d, corresponding to the EC distribution after 66 days and, hence, later compared to plot A.

According to the same procedure used for plot A, the RMSE was estimated through the misfit between observed and calculated data, corresponding roughly to about 16%.

On the other hand, a different soil response was observed in plot C, irrigated with freshwater. Figure 9 and Figure 10 show less marked variations in EC, both positive and negative. Given the poor salinity in the freshwater used for irrigation, the surficial increase can be attributed solely to the increase in moisture content in the soil during the irrigation season, while the decrease in EC observed in the cross-section can be attributed to the root water uptake.

The estimated RMSE for plot C is 9 dS/m, corresponding to about 17%.

Figure 11 shows the calibration function of inverted EC vs. EC, whose values are derived from [24]. To plot this graph, we extracted the inverted EMI values for all nine EMI transects. Only one EC point was removed from the dataset due to an unexplained drift from the general trend, probably as an error in soil sampling analysis.

Some differences between the two datasets can be observed. For low salinity values (less than 1.2 dS/m), EMI overestimates Ece, while for high salinity values, EMI slightly underestimates ECe. This is not surprising because, in soils with low salinity, other factors can increase the EMI signal (clay content, soil moisture, etc). On the other hand, overestimation could be due to data noise or inversion artifacts that could mainly affect the EC distribution for very shallow soil.

Apart from these differences, the correlation is clear, confirming the assumption that the inversion process improves the soil modeling.

4. Discussion

This study was conducted on tomato plants moderately tolerant to salinity [38,39]. Nevertheless, different levels of salt in the soil or in the irrigation water can induce changes in plant morphology and physiology and address severe consequences on crop yield [40,41,42]. In addition, the trend of accumulating salts during the irrigation season may have negative implications not only for soil health but also for groundwater quality. In fact, excess salt can be removed by winter rainfalls and pushed into the vadose zone until the aquifer is reached, causing the oversalinization of the water as a resource [43]. With these premises, the implementation of effective tools capable of monitoring the salinization dynamics in the soil plays a crucial role.

Two weeks after the start of the irrigation, slight changes in EC were observed (Figure 5a, Figure 8a and Figure 10a). The positive changes identified the soil-wetting effects caused by the increase in soil moisture. Conversely, decreases in EC observed below the top layer could be correlated to the root water uptake, which causes negative conductivity changes, as observed in previous works [44,45,46,47].

The changes in EC emphasize the soil–plant–water interaction during the irrigation season, as Figure 5, Figure 8, and Figure 10 highlight. In particular, at the end of the irrigation season, the inversion of the ECa data clearly distinguished the salt accumulation in the topsoil layer in plots irrigated with brackish water (plots A and B) and the root water uptake in the bottom layer of the three plots (Figure 12). For the higher EC differences observed in plots A and B compared with the observations in plot C, those irrigated with fresh water (plot C) had reduced root water uptake activity due to the high osmotic pressure that inhibits the water flow from the soil to the plant, as observed in [48,49,50,51,52]. In the context of saline or brackish groundwater management, this evidence should be taken into account in order to balance water requirement and consumption and, hence, to increase water saving and protect the soil. The lower the salt water required, the lower the accumulation of salts while preserving the crop yield.

The model error recorded in the three plots reflects a certain level of noise in the raw ECa data. Such noise cannot be estimated with the CMD Mini-Explorer when the most commonly “continuous mode” of the measurement is used, i.e., the ECa data are collected when the instrument moves in the field. This is a critical aspect when the visualization of the results is based on ECa maps. Although ECa produces a rapid and fast visualization of the soil properties, it includes undetermined systematic or random errors. Conversely, inverting data allowed us to estimate a model error, remove spikes or noisy channels, and invert again the filtered datasets. In addition, as ECa is a depth-weighted parameter, it could not be correlated with the ECe or EC obtained through other sensors due to the different resolutions, depth of investigation, and sensitivity. In fact, according to the manufacturer’s indications, the upper soil layer (0–0.50 m) is potentially investigated with three ECa measurements (VCP_0.32m, VCP_0.71m, and HCP_0.32m), leading to an ambiguous correlation function.

On the other hand, using the inverted EC data allowed us to plot a single EC vs. ECe calibration function because the inverted model provides an estimation of EC distribution at different depths. This is essential when a limited number of soil samples are available as “ground truth” because all inverted EC data can be used in a single calibration function, such as the proposed case study. It is worth mentioning that no temperature correction was applied either on the modeling results after the inversion procedure or over the raw data. For a more precise soil salinity assessment, such an option should be considered in future studies.

5. Conclusions

This study made use of ECa data collected and processed using time-lapse methods to assess the impact of brackish water used as an irrigation source during the tomato growth season. In three plots where different irrigation strategies were used, repeated ECa data were collected at five time points. This approach provided promising answers to the research questions highlighted in the aims of the study. These findings have highlighted how inverted EC allows us to accurately identify soil salinization when different irrigation strategies are used.

In fact, a clear soil response to different irrigations of water, brackish and freshwater, was observed.

Compared with the traditional raw ECa visualization, the processing of the data through time-lapse inversion allowed for a higher level of detail in the soil properties to be visualized. Although this study was conducted in the short term of a single irrigation season, a significant increase in EC in the upper layer could have strong implications in terms of the accumulation of salt. At the same time, the activity of water uptake from the roots was imaged, confirming the versatility of the geophysical tool in the agronomic investigation. The inversion of the ECa allowed for a single correlation function, EC vs. ECe, to be defined, although based on a few points, by comparing the EC data extracted in correspondence to the sampling points. Depending on the extension of the area to be investigated, a significant number of data points can strengthen the calibration function in order to accurately convert the geophysical outcome into hydrological properties of interest.

The capability of producing an accurate correlation function through the inverted model represents an added value with respect to the use of the ECa, which is a depth-weighted parameter and could address meaningless correlations with point scale measurements.

As electromagnetic data are increasingly widespread in the scientific landscape, it is strongly recommended that the inversion procedure be routinely used in the comprehension of soil properties and dynamics.