the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
Real-time Monitoring of Petroleum Hydrocarbons in Groundwater using Hybrid Machine Learning Architectures
Abstract. Monitoring petroleum hydrocarbon (PHC) plumes in groundwater is essential for managing oil contamination but is often hindered by high costs. We evaluated machine learning (ML) frameworks that estimate concentrations of benzene, ethylbenzene, and xylenes (BEX), using affordable, in situ water quality parameters (iWQPs) as inputs: pH, dissolved oxygen, electrical conductivity, and oxidation-reduction potential. Due to a scarcity of field data, we trained and tested models on high-resolution virtual data generated by a reactive transport model. We compared a long short-term memory (LSTM) network against classical algorithms (multiple linear regression, random forest, support vector regression, XGBoost) and an LSTM-XGBoost hybrid. Model performance depended on the underlying geochemical relationship between iWQPs and BEX. Accurate predictions (R² ≥ 0.80, MAPE < 2.3 %) were achieved when iWQPs were strongly correlated with BEX degradation (e.g., as a primary electron donor); the LSTM model yielded predictions within a 5 % error margin for 70 % of the test cases. Performance declined sharply (R² < 0) during periods where iWQPs were correlated with non-volatile dissolved organic carbon, another component of dissolved PHC. Incorporating hydraulic head data improved accuracy by informing the model of groundwater flow dynamics. While the LSTM model struggled to extrapolate beyond its training data (e.g., during extreme flow events), it reliably detected the direction of concentration trends, providing a valuable trigger for adaptive monitoring. We also demonstrated how a hybrid Kalman filter could successfully capture concentration trends after source removal through recursive updating. Our proposed ML framework provides BEX level estimation for improved groundwater monitoring.
- Preprint
(1839 KB) - Metadata XML
-
Supplement
(2320 KB) - BibTeX
- EndNote
Status: final response (author comments only)
-
CC1: 'Comment on egusphere-2025-5842', Giacomo Medici, 13 Jan 2026
-
AC3: 'Reply on CC1', Chen Lester Wu, 21 Jun 2026
Thank you for these comments, and we responded to the comments below:
General comments
Good and robust research on contaminant transport in groundwater. The authors need to provide more detail before publication. See my specific comments to fix the issues.
Response:
Thank you for the comment, we provided more details for clarification and answered the questions.
Specific comments
Lines 26-27. “Groundwater contamination by petroleum hydrocarbons (PHCs) remains an environmental challenge, particularly in areas affected by historical spills or leaks”. General statement not backed up by references. Please, insert general literature on the topic.
- Agbotui, P. Y., Firouzbehi, F., Medici, G. 2025. Review of effective porosity in sandstone aquifers: insights for representation of contaminant transport. Sustainability, 17(14), 6469.
- Li, G., Huang, W., Lerner, D. N., Zhang, X. 2000. Enrichment of degrading microbes and bioremediation of petrochemical contaminants in polluted soil. Water Research, 34(15), 3845-3853.
Response:
We will add references and cite these for this general statement, thank you!
Line 88. The aim of the research is clear. But what about the 3 to 4 specific objectives? Please, describe them by using numbers (e.g., i, ii, and iii).
Response:
Thank you for the comment, we will add the specific objectives in the introduction.
Line 90. If you use MODFLOW-2005 you need much more detail on the boundary conditions.
Line 93. Provide more detail on the boundary conditions also for MT3DMS.
Response:
These were included in the manuscript as previously commented by other referees. Thank you!
Line 220. Why only MAE? What about just Mean Error and Root Mean Squared Error?
Response:
In this study, MAE, R², and MAPE were selected as baseline metrics to provide a general assessment of model performance. While additional metrics such as mean error or RMSE could certainly be included, simply adding more statistics does not necessarily lead to clearer interpretation or better decision‑making. Evaluating the practical relevance of prediction errors, particularly in the context of groundwater quality, requires expert judgment and examination of the actual predicted and observed concentration trends, which we have also presented and discussed.
Lines 490-510. I can see 5 bulletin points in your conclusions. Therefore, the specific objectives (see comment above) must be the same number to match.
Response: to be discussed
The specific objectives will be added to the introduction and will be explicitly listed. The conclusions have been revised to reflect and address these objectives, though not necessarily in a one-to-one manner, as some conclusions refer to multiple objectives. Also, the conclusions were further revised based on recommendations by reviewers.
Figures and tables
Would you like to add the flow field output of MODFLOW?
Response:
The flow field output of MODFLOW was available in the supporting information of our previous paper (Wu et al., 2024). Hence, we did not include it in this manuscript.
What about horizontal slices for the contaminant transport? I can see only a vertical one.
Response:
The horizontal slices for the contaminant transport might only make sense for 3-dimensional model. Since we considered only a 2D cross-section, simplifications were made.
Figure 1. Important figure, make it larger.
Response:
We increase the size of Figure 1, thank you!
Figure 1. What about a vertical scale in meters above the sea level?
Response:
The RTM was based on a previous modelling study of the Bemidji site (Ng et al., 2015) with the original y-axis from 418m to 424m as elevation. However, since we wanted to simplify the representation of saturated zone, we decided to use depth instead of elevation. Thus, since we did not include the unsaturated zone in our model, there is no value in the vertical scale above the sea level in the y-axis.
Figure 3. Un-clear. Please, rise the graphic resolution.
Figure 5. Make the legend larger. You can use 3 lines.
Response:
We will make sure these figures have the right resolution and legends. Currently, all figures were generated using Python with DPI of 600. As this is already a good resolution, there might be issues with uploading the figures.
Citation: https://doi.org/10.5194/egusphere-2025-5842-AC3
-
AC3: 'Reply on CC1', Chen Lester Wu, 21 Jun 2026
-
RC1: 'Comment on egusphere-2025-5842', Massimiliano Schiavo, 18 Feb 2026
-
AC1: 'Reply on RC1', Chen Lester Wu, 21 Jun 2026
Thank you for the comments. We would like to respond to each comments:
- Greetings. I have revised the paper entitled ‘Real-time Monitoring of Petroleum Hydrocarbons in Groundwater using Hybrid Machine Learning Architectures’. The paper deals with Machine Learning-based approaches to assess the contaminant plume of benzene pollutants and their fate in groundwater. The work is scientifically sound, and its methodology is robust. However, I’m recommending some adjustments before further proceeding down the publication path. Best regards.
- I cannot state the boundary conditions, extensions, and their type in your work. Hence, I cannot fully understand how you tuned your model.
Response:More information about the RTM were added for clarification:“The model also incorporated transient water table fluctuations and spatial heterogeneity in hydraulic conductivity. The synthetic aquifer domain represents a saturated porous medium inspired by the Bemidji crude oil spill site (Ng et al., 2015), with a length of 300 m and a thickness of 10 m, discretized into 150 columns (Δx = 2 m) and 50 layers (Δz = 0.2 m). The bottom boundary was assigned a no-flow condition, representing an underlying clay aquitard, and the upper boundary received a uniform recharge flux of 178 mm/yr. For the initial steady-state flow field, constant-head boundaries were imposed on the left and right sides of the model to reproduce the observed hydraulic gradient of 0.0035 m/m at the Bemidji site (Ng et al., 2015).”Stated in the earlier paper (Wu, 2024), the right boundary condition was derived from a separate 5.2‑km regional model, in which the same sinusoidal fluctuation was applied at the upstream end. The resulting damped and phase‑shifted water‑table response at 300 m was extracted and used as the time‑varying head boundary on the right side of the final model. The bottom boundary remained no‑flow, and transmissivity was assumed constant due to the small magnitude of water‑table fluctuations relative to the saturated thickness.- It is clear that geological (hence facies conductivity) heterogeneity is a key driver in your transport model. However, the methodology lacks supporting equations and a sufficient explanation. Indeed, at lines 91-93, you write “The RTM used in this study is based on the model developed in Wu et al. (2024), implemented in Python using FloPy within a Jupyter Notebook environment. Groundwater flow was simulated using MODFLOW 2005, and contaminant transport was modeled with MT3DMS, incorporating advection and dispersion processes”. Therefore, I suggest carving a brief section, maybe before the present section 2.3, to explain governing equations (see also Bedekar et al., 2016). Please incorporate a short section. I would just write a couple of lines for MODFLOW 2005 and then focus on transport:
“The MT3D model simulates dissolved solute transport in groundwater using the advection–dispersion–reaction equation. In terms of concentration per unit volume of water, the governing mass balance is∂C/(∂t ) = ∇ ⋅ (𝐃∇𝐶) − 𝐯 ⋅ ∇𝐶 + q_s/(θ ) (𝐶𝑠 − 𝐶) − 𝜆𝐶where 𝐶 is solute concentration, 𝜃 is porosity, 𝐯 = 𝐪/𝜃 is the seepage velocity derived from the Darcy flux 𝐪 (computed by MODFLOW), 𝐃 is the hydrodynamic dispersion tensor, 𝑞𝑠 and 𝐶𝑠 represent fluid sources/sinks and their concentrations, and 𝜆 is a first-order decay coefficient.When linear equilibrium sorption is included, transport is retarded and the equation becomeswhere the retardation factor iswith 𝜌𝑏 the bulk density and 𝐾𝑑 the distribution coefficient.”Response:Thank you for the suggestion, we will include this in the manuscript, along with the model description:“The MT3D model simulates dissolved solute transport in groundwater using the advection–dispersion–reaction equation (Bedekar et al., 2016). In terms of concentration per unit volume of water, the governing mass balance is∂C/(∂t ) = ∇ ⋅ (𝐃∇𝐶) − 𝐯 ⋅ ∇𝐶 + q_s/(θ ) (𝐶𝑠 − 𝐶) − 𝜆𝐶where 𝐶 is solute concentration, 𝜃 is porosity, 𝐯 = 𝐪/𝜃 is the seepage velocity derived from the Darcy flux 𝐪 (computed by MODFLOW), 𝐃 is the hydrodynamic dispersion tensor, 𝑞𝑠 and 𝐶𝑠 represent fluid sources/sinks and their concentrations, and 𝜆 is a first-order decay coefficient. “ As we did not simulate sorption of hydrocarbons in our model, we will not add the linear sorption equations.- Which numerical schemes were adopted in the transport simulations for the advection, dispersion, and reaction terms? Specifically, which advection solver was used (Upstream finite difference, TVD, MOC/MMOC, or HMOC), and how were dispersion (finite-difference formulation) and reactions handled (explicit or implicit solution)?
Response:These were also included in the manuscript for more information:“The RTM used in this study is based on the model developed in Wu et al. (2024), implemented in Python using FloPy within a Jupyter Notebook environment. Groundwater flow was solved using MODFLOW‑2005 with a finite‑difference formulation, while solute transport was simulated with MT3DMS. For the advection term, we used the Total Variation Diminishing implicit scheme in MT3DMS which also minimizes the numerical dispersion. The dispersion term was computed using the standard central finite‑difference formulation of MT3DMS, with longitudinal and transverse dispersivities and a single molecular diffusion coefficient specified in the DSP package. All geochemical reactions, including multicomponent LNAPL dissolution, aerobic and anaerobic biodegradation, mineral precipitation/dissolution, cation exchange, and gas outgassing, were solved implicitly using PHREEQC‑2 within PHT3D. PHT3D couples MT3DMS transport with PHREEQC via operator splitting, where PHREEQC applies an implicit Newton–Raphson solver and stiff ordinary differential equation integration for equilibrium and kinetic reactions.”The Jupyter Notebook containing the full reactive transport model code with this information is uploaded in the open access repository https://doi.org/10.4121/f7742f02-ee3a-4a84-adf1-625b4a9fd703 (Wu, 2024).- Can you frame your work in the wider Machine-Learning literature? There is a wide class of Genetic Algorithms that perform extremely well and could be adapted to your framework (Rajwar et al., 2023; Schiavo & Pedretti, 2026). I think the paper should deal with this branch of research, explaining (i) why your methodology is better (or not) and (ii) possible employments of metaheuristics-driven approaches.
Response:Thank you for this suggestion. Metaheuristic approaches can indeed be adapted to groundwater monitoring frameworks, particularly for tasks such as feature selection or designing optimal sensor placement strategies.However, the primary objective of this study is not to compare models or methodologies, but to evaluate the feasibility of using machine learning architectures to estimate petroleum hydrocarbon concentrations from in-situ water quality data. In this context, the focus is on assessing whether sequence learning models (e.g., LSTM based architectures) provide practical benefits for real time groundwater monitoring. Metaheuristic algorithms could certainly complement this framework, for example, by optimizing model parameters or selecting the most informative sensor combinations. Such studies still fall outside the scope of the present proof-of-concept study.To contextualize the applied methodology in our research, we clarified the primary goal and improved the analysis. This also addressed another comment we received on the comparison of different models applied:“For comparison, traditional regression models were also implemented for their interpretability and established use in water quality monitoring (Singh et al., 2021): multiple linear regression (MLR) to provide a benchmark for linear relationships (Sani Gaya et al., 2020); support vector regression (SVR) to capture nonlinear patterns (Banadkooki et al., 2020); and RF and XGBoost to handle complex feature interactions (Szomolányi and Clement, 2023). While these models can represent temporal dynamics when supplied with engineered features (e.g., lagged variables), in this study they were applied using daily iWQP measurements as independent observations. This choice allows for a focused evaluation of the extent to which explicitly modeled temporal dependencies, such as LSTM-based approaches, improve predictive performance relative to models without embedded sequence learning. However, the primary goal is to assess the feasibility of ML approaches for BEX monitoring, rather than strict model comparison.”This was also added in the final discussion of the paper:“Selecting the most suitable regression model for field application remains challenging due to deviations between real-world conditions and RTM simulations. Future work may explore other machine learning models or the integration of metaheuristic driven optimization (i.e., Schiavo & Pedretti, 2026) to further enhance model performance and operational deployment.”- Can you highlight the driving parameters for each of your three scenarios, and how they are involved in XGBoost regression? I suggest offering a framework where these parameters are well-highlighted. See e.g. Schiavo & Pedretti, 2026, Figure 1.
Response:Is our understanding correct that the term “driving parameters” refers to the sensor parameters that influences the XGBoost predictions the most?Identifying the most influential input variables can indeed be valuable when the objective of a study is to interpret model behavior or to optimize a predictive framework. However, our focus is to evaluate whether LSTM‑based sequence‑learning architectures can effectively estimate BEX concentrations from sensor data under different hydrologic stress conditions. The additional ML models (including XGBoost) were included only as comparative baselines for the LSTM performance, not as targets for detailed interpretability analysis.For this reason, we chose not to expand the scope of the study to include scenario‑specific feature‑importance analyses or a full interpretability framework for XGBoost. We will clarify this point in the manuscript to avoid the impression that the comparative models are intended to be analyzed at the same depth as the LSTM architecture. We note, however, that once this framework is validated with field data, identifying the dominant sensor drivers and conducting a full interpretability analysis will be further explored as one of the next steps for operational deployment and model transparency.- Figures 3 and 4 are very hard to read properly. I suggest, at least for the latter one, employing boxplots to show the MAE.
Response:Thank you for the suggestion, it may be hard to read the figures at first glance since it was not a standard way to visualize the results. However, we intentionally designed Figures 3 and 4 to provide a visual impression of the spread and temporal variability in model performance across different monitoring wells and sliding windows due to the sheer number of data points, rather than for exact numerical quantification. We described the figures in detail in the manuscript. We acknowledge that boxplots are an effective tool for summarizing distributions, but we believe that they would make the temporal scattering vague which is what we want these figures to convey. Since the variability across time windows is important, we prefer to preserve that scatter to be shown in the figures. For now, the current design serves this purpose effectively, though we remain open to any specific formatting recommendations to improve readability.- Suggested References:
- Bedekar, V., Morway, E.D., Langevin, C.D., and Tonkin, M., 2016, MT3D-USGS version 1: A U.S. Geological Survey release of MT3DMS updated with new and expanded transport capabilities for use with MODFLOW: U.S. Geological Survey Techniques and Methods 6-A53, 69 p., http://dx.doi.org/10.3133/tm6A53
- Rajwar, K., Deep, K., & Das, S. (2023). An exhaustive review of the metaheuristic algorithms for search and optimization: Taxonomy, applications, and open challenges. Artificial Intelligence Review, 56(11), 13187–13257. https://doi.org/10.1007/s10462-023-10470-y
- Schiavo, M., & Pedretti, D. (2026). Genetic and Iterative Metaheuristics-Informed Algorithms for Precision Shallow Groundwater Modeling and Drought Inference. Journal of Geophysical
- Research: Machine Learning and Computation, 3(1), e2025JH000854. https://doi.org/10.1029/2025JH000854
Citation: https://doi.org/10.5194/egusphere-2025-5842-AC1
-
AC1: 'Reply on RC1', Chen Lester Wu, 21 Jun 2026
-
RC2: 'Comment on egusphere-2025-5842', Anonymous Referee #2, 31 May 2026
================
General comments
================
This manuscript is a description of a numerical test to asses the applicability of ML algorithms to estimate hydrocarbon concentrations in a thin aquifer from measurements of pH, dissolved oxygen, electrical conductivity, and oxidation-reduction potential.
The work presents some flaws, which are discussed in detail in the specific comments below. This prevents the results of this work to be considered of general validity; in other words, they are strictly related to the specific test cases.
In conclusion, I think that the manuscript cannot be considered for publication in its resent version, but requires an accurate major revision.
================
Specific comments
================
- I found it quite confusing: ML architectures are “modeling” tools, they cannot be considered as “monitoring” tools. I would suggest to change the title, possibly as “Hybrid Machine Learning Architectures to Estimate Petroleum Hydrocarbons Concentrations in Groundwater from Proxy Data”. Moreover, the whole text should be modified accordingly.
- Numerical tests are conducted with a 2D model. In other words, this assumes that the horizontal component of groundwater and contaminant flow is directed along a straight line and that all physical quantities are homogeneous along the perpendicular direction. What if a more realistic 3D flow and transport setup is considered? In particular, the scenario with a single injection well creates a 3D flow, unless it is assumed that an arrays of wells is installed, perpendicular to the flow direction. The very basic physical assumptions supporting the numerical flow and transport model should be presented in a better way, and the limitations should be discussed.
- The application of the Kalman filter refers to a single specific scenario. Therefore it is overall of minor relevance for the work.
- Line 205. Assuming a Gaussian distribution for measurement errors is common and generally acceptable. However, assuming a Gaussian distribution for process noise or modeling errors is more debatable, even if it is common practice in the scientific literature, but without a proper physical support.
- Line 222. A percentage error has a very different practical relevance for low and high concentrations. For high concentrations, a small percentage error may correspond to high values: this could impact also the overcome of thresholds to declare water drinkable. For low concentrations, higher percentage errors could be acceptable, because their practical impact could be less significant. Therefore, I am not sure that the use of MAPE is the best choice.
- Figure 3. The “90 days” sequence length seems to proved the worst results, doesn’t it? I could not find any comment about this.
- Figure 5. The behavior of the LSTM predictions for the simulation time from year 35 to year 41 is very strange: the behavior of the learning time (from year 35 to year 39) is quite regular: a clean signal with annual period and small amplitude around a constant average. The prediction shows a counterphase signal with a much higher amplitude. The explanation given at lines 284ff is quite confusing. The poor fit is for the first examined period, not the third one, and it refers to years 35 to 41, not 40 to 60. May be, I misunderstood the text. Similar comments could apply to figure 8.
=================
Technical comments
=================
- Line 67. “30-360 days” is equal to -330 days. I warmly suggest to follow the instructions of section “7.7 Clarity in writing values of quantities” of the NIST “Guide for the Use of the International System of Units (SI)” downloadable at https://physics.nist.gov/cuu/pdf/sp811.pdf.
- Lines 82, 200, 445. Is Simon (2006) the best reference? Why not recalling the seminal papers by Kalman? Or some of the textbooks and papers on the use of KF in (non-linear) hydrology?
- Figure 2. I think it would be useful to add a label, e.g., “Time (years)”, for the x axis. Why using two squares for each year? Wouldn’t it be more simple and clear to draw one cell for each year?
- Figure 3. I understand that the numerical simulation of flow and transport covered 100 years. I couldn’t find this information in the text: did I miss it?
Citation: https://doi.org/10.5194/egusphere-2025-5842-RC2 -
AC2: 'Reply on RC2', Chen Lester Wu, 21 Jun 2026
Thank you for the comments. We would like to reply to each comments:
================General comments
================
This manuscript is a description of a numerical test to assess the applicability of ML algorithms to estimate hydrocarbon concentrations in a thin aquifer from measurements of pH, dissolved oxygen, electrical conductivity, and oxidation-reduction potential.
The work presents some flaws, which are discussed in detail in the specific comments below. This prevents the results of this work to be considered of general validity; in other words, they are strictly related to the specific test cases.
In conclusion, I think that the manuscript cannot be considered for publication in its resent version, but requires an accurate major revision.
Response:
Thank you for the comment. As this was part of a PhD manuscript, we received feedback and comments from the committee members during the review process. We revised the paper to improve the conclusions and tighten the justification for the methodology choices. We also improved the manuscript based on the comments from the referees.
================
Specific comments
================
1. I found it quite confusing: ML architectures are “modeling” tools, they cannot be considered as “monitoring” tools. I would suggest to change the title, possibly as “Hybrid Machine Learning Architectures to Estimate Petroleum Hydrocarbons Concentrations in Groundwater from Proxy Data”. Moreover, the whole text should be modified accordingly.
Response:
We appreciate this observation. However, the use of ML-based approaches within a monitoring framework is well-established in the literature, with several peer-reviewed studies employing similar "ML-based monitoring" terminology (e.g., Groundwater Quality Monitoring Using In-Situ Measurements and Hybrid Machine Learning with Empirical Bayesian Kriging Interpolation Method; In Situ Monitoring of Groundwater Contamination Using the Kalman Filter). While we acknowledge that ML architectures are indeed modeling tools in a technical sense, their application here is oriented toward monitoring. Specifically, monitoring through nowcasting petroleum hydrocarbon concentrations from proxy data. We believe that adopting "modeling" in the title risks introducing ambiguity, as it may imply simulation of the full physical system, which is not our intent in this work. We would also like to give an example about the use of machine learning in weather applications, where ML models are widely regarded as monitoring tools rather than simulation models. We would therefore keep the current title, as we believe it better reflects the practical purpose of our proposed framework.
2. Numerical tests are conducted with a 2D model. In other words, this assumes that the horizontal component of groundwater and contaminant flow is directed along a straight line and that all physical quantities are homogeneous along the perpendicular direction. What if a more realistic 3D flow and transport setup is considered? In particular, the scenario with a single injection well creates a 3D flow, unless it is assumed that an arrays of wells is installed, perpendicular to the flow direction. The very basic physical assumptions supporting the numerical flow and transport model should be presented in a better way, and the limitations should be discussed.
Response:
Thank you for this comment. A 2D representation may indeed oversimplify the true 3D nature of groundwater flow. In this study, however, the objective was not to reproduce the full spatial complexity of a real aquifer, but to evaluate the feasibility of using low‑cost in‑situ water quality parameters (iWQPs) to infer BEX concentrations under different hydrologic stress conditions. For this purpose, a 2D vertical slice provides a controlled environment in which the dominant hydrogeochemical processes (i.e., mass dissolution, redox evolution, mixing, and plume attenuation) are preserved, while simplifying the effects of lateral heterogeneity that could potentially reduce the signal‑to‑noise relationship between iWQPs and BEX.
The assumption is that the chemical evolution of groundwater along the flow path is primarily governed by redox zonation and reaction kinetics, which tend to be similar along parallel streamlines even in 3D systems. Furthermore, the relative changes in iWQPs (e.g., EC, DO, ORP, pH) that serve as proxies for biodegradation and plume behavior are expected to remain qualitatively consistent regardless of whether the flow is 2D or 3D. We will clarify the limitations of a 2D model in the manuscript, such as anisotropic dispersion and cross‑sectional spreading. One thing to note is that, once this framework is validated with field data, extending the approach to 3D flow and transport can be part of a future study and to quantify how spatial heterogeneity influences the robustness of sensor‑based BEX estimation.
3. The application of the Kalman filter refers to a single specific scenario. Therefore it is overall of minor relevance for the work.
Response:
The Kalman filter was kept intentionally simple and scenario‑specific to provide potential avenues for future improvements or other models that could potentially be used.
4. Line 205. Assuming a Gaussian distribution for measurement errors is common and generally acceptable. However, assuming a Gaussian distribution for process noise or modeling errors is more debatable, even if it is common practice in the scientific literature, but without a proper physical support.
Response:
Thank you for the comment. The assumption of Gaussian process and measurement noise was adopted for simplicity and to support an initial proof‑of‑concept demonstration. In future work, particularly after field validation of the full framework, non‑Gaussian error structures and more advanced filtering approaches could be examined. However, these were not part of the scope of this paper.
5. Line 222. A percentage error has a very different practical relevance for low and high concentrations. For high concentrations, a small percentage error may correspond to high values: this could impact also the overcome of thresholds to declare water drinkable. For low concentrations, higher percentage errors could be acceptable, because their practical impact could be less significant. Therefore, I am not sure that the use of MAPE is the best choice.
Response:
This is a great point, and we agree. In this study, MAPE was used only as one component of the performance evaluation, and not as a standalone indicator, precisely because percentage‑based errors behave differently across the concentration range. To avoid relying on a single metric, we also reported additional error measures and presented the actual predicted and observed (RTM-based) concentration time series. Showing the raw estimates was intentional, as it allows readers to directly assess the practical relevance of errors at both low and high concentrations.
6. Figure 3. The “90 days” sequence length seems to proved the worst results, doesn’t it? I could not find any comment about this.
Response:
At observation well X1Z2, the 90‑day sequence produced 8 test windows with R2 > 0.80, compared with 5 windows for the 360‑day sequence. However, the 90‑day sequence also generated more windows with negative R values than the 360‑day sequence. Because these mixed outcomes make it challenging to categorically state that one sequence length performs “worse” than the other, and because this comparison does not influence the main conclusion, which is that the 30‑day sequence consistently provides the highest R2 and most stable performance, we did not include a detailed discussion of the worst sequence results in the manuscript.
7. Figure 5. The behavior of the LSTM predictions for the simulation time from year 35 to year 41 is very strange: the behavior of the learning time (from year 35 to year 39) is quite regular: a clean signal with annual period and small amplitude around a constant average. The prediction shows a counterphase signal with a much higher amplitude. The explanation given at lines 284ff is quite confusing. The poor fit is for the first examined period, not the third one, and it refers to years 35 to 41, not 40 to 60. May be, I misunderstood the text. Similar comments could apply to figure 8.
Response:
We understand the source of confusion. In Figure 5, the plotted BEX concentrations represent the training period, but the LSTM is not trained on BEX itself, only on the water quality parameters (iWQPs). Although the BEX signal during years 35–39 appears relatively stable, the corresponding iWQP signals in this interval exhibit a different temporal pattern, which leads to the mismatch observed in the predictions in the test year. The explanation in our manuscript referenced the “hydrogeochemical periods” defined in our previous study (Wu et al., 2024), where the relationship between iWQPs and BEX changes over time. The terminology used (“periods”) may have contributed to the confusion, since it refers to those previously defined hydrogeochemical periods rather than the chronological segments shown in the figure. We will revise the text to specifically mention “hydrogeochemical period” instead of just “period”.
As an example:
“The window with the poor fit (R² < 0) coincided with hydrogeochemical period 3 (around year 20 to year 60), where the BEX compounds reached a state of quasi-equilibrium. In this hydrogeochemical period, the rolling Spearman's correlation between iWQPs and BEX dropped to approximately zero from around year 25 to year 40. Instead, iWQPs became highly correlated with non-volatile dissolved organic carbon (NVDOC), a fraction of the dissolving light non-aqueous phase liquid (LNAPL). NVDOC became the primary electron donor, and thus the dominant driver of geochemical change. This results in the poor performance of the LSTM model for BEX concentration estimation.”
=================
Technical comments
=================
1. Line 67. “30-360 days” is equal to -330 days. I warmly suggest to follow the instructions of section “7.7 Clarity in writing values of quantities” of the NIST “Guide for the Use of the International System of Units (SI)” downloadable at https://physics.nist.gov/cuu/pdf/sp811.pdf.
Response:
Thank you for this, we changed it into “to” and will change the others into em-dash as applicable.
2. Lines 82, 200, 445. Is Simon (2006) the best reference? Why not recalling the seminal papers by Kalman? Or some of the textbooks and papers on the use of KF in (non-linear) hydrology?
Response:
The reference to Simon (2006) was used because it served as the primary practical resource for implementing the Kalman filter in this study. However, we agree that it is appropriate to also acknowledge the foundational work. We will therefore add a citation to Kalman’s seminal paper: (Kalman, 1960).
3. Figure 2. I think it would be useful to add a label, e.g., “Time (years)”, for the x axis. Why using two squares for each year? Wouldn’t it be more simple and clear to draw one cell for each year?
Response:
Thank you for this suggestion. Figure 2 is intended as a schematic illustration of the rolling‑window concept rather than a literal representation of the dataset. The two squares per year were included purely for visual clarity to show how the training and testing windows move forward in time. For clarity, we will revise the figure to add an x‑axis label (e.g., “Time (years)”) and adjust the layout so that each year is represented by a single cell.
4. Figure 3. I understand that the numerical simulation of flow and transport covered 100 years. I couldn’t find this information in the text: did I miss it?
Response:
Similar to the comment of referee 1, the information on the reactive transport model will be included in the paper. Thank you for pointing this out.
Citation: https://doi.org/10.5194/egusphere-2025-5842-AC2
Interactive computing environment
Hybrid Machine Learning Models for Estimating Petroleum Hydrocarbon Concentration in Groundwater Chen Lester R. Wu et al. https://doi.org/10.4121/0a23147e-ba85-4ba2-a058-ba199c65d711
Virtual Experiments with Reactive Transport Modelling using FloPy: Transport and Degradation of Dissolved Petroleum Hydrocarbons in Groundwater Chen Lester R. Wu et al. https://doi.org/10.4121/f7742f02-ee3a-4a84-adf1-625b4a9fd703
Viewed
| HTML | XML | Total | Supplement | BibTeX | EndNote | |
|---|---|---|---|---|---|---|
| 1,491 | 780 | 116 | 2,387 | 247 | 103 | 108 |
- HTML: 1,491
- PDF: 780
- XML: 116
- Total: 2,387
- Supplement: 247
- BibTeX: 103
- EndNote: 108
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1
General comments
Good and robust research on contaminant transport in groundwater. The authors need to provide more detail before publication. See my specific comments to fix the issues.
Specific comments
Lines 26-27. “Groundwater contamination by petroleum hydrocarbons (PHCs) remains an environmental challenge, particularly in areas affected by historical spills or leaks”. General statement not backed up by references. Please, insert general literature on the topic.
- Agbotui, P. Y., Firouzbehi, F., Medici, G. 2025. Review of effective porosity in sandstone aquifers: insights for representation of contaminant transport. Sustainability, 17(14), 6469.
- Li, G., Huang, W., Lerner, D. N., Zhang, X. 2000. Enrichment of degrading microbes and bioremediation of petrochemical contaminants in polluted soil. Water Research, 34(15), 3845-3853.
Line 88. The aim of the research is clear. But what about the 3 to 4 specific objectives? Please, describe them by using numbers (e.g., i, ii, and iii).
Line 90. If you use MODFLOW-2005 you need much more detail on the boundary conditions.
Line 93. Provide more detail on the boundary conditions also for MT3DMS.
Line 220. Why only MAE? What about just Mean Error and Root Mean Squared Error?
Lines 490-510. I can see 5 bulletin points in your conclusions. Therefore, the specific objectives (see comment above) must be the same number to match.
Figures and tables
Would you like to add the flow field output of MODFLOW?
What about horizontal slices for the contaminant transport? I can see only a vertical one.
Figure 1. Important figure, make it larger.
Figure 1. What about a vertical scale in meters above the sea level?
Figure 3. Un-clear. Please, rise the graphic resolution.
Figure 5. Make the legend larger. You can use 3 lines.