Comparative Performance of Gaussian Plume and Backward Lagrangian Stochastic Models for Near-Field Methane Emission Estimation Using a Single Controlled Release Experiment

1. Introduction

Methane (CH₄) is a major component of natural gas (NG) and a potent greenhouse gas (GHG), with a global warming potential (GWP₁₀₀) 27.9 times larger than CO₂ over 100-year time period [1]. Anthropogenic activities like oil and gas (O&G) extraction, agriculture, biomass burning, and waste disposal account for between 60 and 65% of the total CH₄ emission globally, while natural emitters like termites, freshwater ecosystems, and wetlands contribute to the remaining fraction [2]. Atmospheric CH₄ concentration has increased by an average of 13 ppb yr^-1 since 2020 [3,4], contributing to changes in the global climate over the past decade. As a controllable emission source, mitigating CH₄ emissions from O&G production sites has been targeted for reducing overall greenhouse gas emissions [5,6].

Monitoring the efficacy of mitigation efforts depends on accurate quantification of facility-level emissions. Yet quantifying facility-level emission estimate remains highly uncertain due to the effects of atmospheric conditions and aerodynamically complex equipment present at O&G extracting sites [7,8,9,10,11,12]. O&G production facilities comprise thousands of components (connectors, valves, and flanges) on multiple pieces of equipment wellheads, separators, and water/condensate storage tanks, often with several to dozens of potential leak points per facility [13,14]. The potential for emissions from different heights and position on multiple components add additional complexity to quantification, therefore further challenge quantification using measurements [15,16].

Emissions on O&G production facilities typically can be emitted from vented, fugitive and combustion-related sources [14]. Gas is typically vented as the result of a process action (e.g. from pneumatic controllers on the wellhead or tank flashing) or from maintenance activities (such as the blowdown of pressurized equipment or liquid uploading). Fugitive emissions are unintended leaks and are typically emitted from flanges, valves or seals as equipment on the facility ages or becomes damaged [14,17]. Combustion-related emissions arise from burning of natural gas in turbines, engines, heaters, and flares [14,18]. Emissions from some vented and most fugitive sources (including vents, pneumatic controllers, valves, flanges, and tank hatches) are typically treated as point sources, whereas other vented and combustion emissions may behave as spatially distributed (or effective area sources) depending on measurement scale [19,20,21]. Point sources are usually easier to detect as they have higher methane concentrations at the emission point but are not necessarily larger than area sources as they emit over a smaller area.

Methane emission quantification using large-scale, remote sensing approaches, like aerial surveys, satellite observations, and tall tower inversion methods are suitable for constraining basin-wide CH₄ fluxes but their applicability in making facility-scale emission estimations are limited by both spatial resolution and temporal coverage [22,23]. These methods typically rely on spatial integration of methane enhancements over footprints spanning meters to several kilometers which prohibits detection and attribution at the equipment or facility level [22,23]. Temporally, their measurements are either very long (tall-tower multiple day averages; [24]) or very short (satellite and aerial surveys <1 minute observation), and infrequent sampling limits sensitivity to intermittent or short-term emission events [25].

To address this, near-field techniques that use near-field measurement and dispersion models are used to generate facility-level emission estimates. However, the accuracy of the emissions calculated using these dispersion approaches is limited due to meteorological variability and modeling assumptions, which are rarely valid at complex O&G sites [10,26,27].

Among commonly used dispersion modeling methods applied to make CH₄ emission estimates across a range of source types are the Gaussian plume (GP) and backward Lagrangian stochastic (bLS) approaches [27,28,29,30]. The GP method models the spread of the pollutant plume from a point source following a Gaussian profile under steady-state conditions with constant meteorological parameters [31,32,33]. Under its standard formulation, the GP approach models atmospheric transport as the dispersion and advection of a continuous point source by the mean wind field, with turbulent dispersion in the vertical and crosswind directions, based on analytical advection-diffusion solution [34,35]. This approach’s main advantages are mathematical simplicity and computational efficiency. In the O&G sector, the GP approach has been widely applied to quantify the CH₄ emissions from individual pieces of production equipment to facility-scale using mobile and fixed survey and continuous monitoring strategies [29,32,35]. The GP-based approach has been applied to both onshore and offshore production facilities, where near-field conditions and steady-state assumptions are presumed to be reasonable [26,29,31,32]. Typical estimates quantified using GP methods in O&G applications range from small (<1 kg CH₄ h^-1) equipment-level sources to larger facility-level emissions (order of 10 kg CH₄ h^-1) [29,31,32,33]. The uncertainty associated with emissions calculated using the GP approach is mostly driven by the steady state assumptions made and simplified source representations, particularly when there are intermittent emissions, multiple sources and complex aerodynamic conditions caused by structures at the facility [26,31,33]. The insensitivity of the model to changes in wind speed and direction adds additional uncertainty to the GP based estimates [11,32]. Evaluations of the GP method in both field studies and controlled experiments indicate an uncertainty of the order ±50% to ±100%, with larger errors observed within highly variable meteorological conditions or when plume intermittency/geometry is not adequately represented [29,31,32,33].

The bLS approach estimates emissions by simulating turbulent air particles backward in time from receptor to potential source using a stochastic representation of atmospheric turbulence [27,36]. Observed concentration enhancements are used to calculate surface emission fluxes or emissions by quantifying the fraction of particle trajectories that originate from emitting surface, without assuming Gaussian concentration profile [36,37]. This method assumes homogenous surface conditions and constant turbulence over the averaging period [27,36]. The accuracy of emission estimates improves with larger numbers of modelled particle trajectories and better resolved wind statistics [9,27,36,37]. Evaluation of CH₄ quantification in the O&G sector indicates that bLS provides robust estimates when turbulence is well characterized [29,38]. The bLS approach can better represent aerodynamic conditions between the source and receptor and atmospheric stability than the GP approach. However, these advantages are accompanied by a higher computational requirement and greater sensitivity to turbulence characterization [13,32]. Studies have shown that emission estimates calculated using a bLS model typically have uncertainties between ±30 and ±50%, with uncertainties rising toward a factor of two or greater when turbulence is poorly resolved [28,29,38].

Although bLS and GP approaches are widely used for CH₄ emission quantification at O&G facilities, most prior evaluations of dispersion model performance relied on aggregated facility-level emissions, multi-source configurations and extended temporal averaging [11,29,39]. These studies are valuable for estimating facility-level but have limited success in isolating the effect of short-time-scale plume intermittency, atmospheric variability, and source-receptor geometry on the performance of the model approach in near-field conditions. Consequently, the relative characteristics and performance of bLS and GP models to characterize single-emission point sources, representative of individual equipment-level emissions, are currently not well constrained for complex O&G environment.

This study aims to directly address this knowledge gap using data from controlled, single-point source CH₄ release experiments to evaluate the GP and bLS approaches’ performance under near-field O&G conditions. Using a constant, controlled emission source, we aim to systematically identify the impact of short time-scale meteorological variability and plume-geometry on the accuracy of calculated emissions. The experimental framework enables a direct comparison of bLS and GP performance under identical measurement and atmospheric conditions. In contrast to prior studies that relied on long-term averaging and source aggregation to evaluate dispersion model performance [29,33,36,38], this study isolates key uncertainty drivers, model-specific sensitivities, and failure modes in the GP and bLS approaches. As a result, the work may provide insight into the conditions under which GP and bLS yield reliable emission estimates and clarifies their relevance for CH₄ quantification in near-field O&G settings.

2. Materials and Methods

2.1. Site Description and Experiment Design

All experiments reported in this study were performed at the Methane Emission Technology Evaluation Center (METEC) test facility, operated by Colorado State University’s Energy Institute in Fort Collins, Colorado, USA. At the time when the measurements for this study were made, METEC was an experimental site that simulated NG emissions from five O&G production pads using controlled release emission points at likely/known fugitive and vented sites on tanks, separators, and well-heads (Figure 2.1). The NG used for the controlled releases was stored in compressed gas cylinders, and the emission rates controlled using precision orifices, pressure regulators, and solenoid valves [29].

Figure 1. Aerial view of METEC’s experimental layout in Ft. Collins, Colorado. Each labeled “pad” consists of oil and gas equipment fitted with controlled release points at numerous places in locations where real-world leaks are known to occur. In this study, experiments were conducted on all pads. Equipment types for releases in these experiments were wellheads (on pad 4), separators (on pads 1–5), and tanks (on pad 4) at varying rates and release heights/locations. METEC has expanded its capabilities in the time since experiments were performed and since this photo was taken.

Figure 1. Aerial view of METEC’s experimental layout in Ft. Collins, Colorado. Each labeled “pad” consists of oil and gas equipment fitted with controlled release points at numerous places in locations where real-world leaks are known to occur. In this study, experiments were conducted on all pads. Equipment types for releases in these experiments were wellheads (on pad 4), separators (on pads 1–5), and tanks (on pad 4) at varying rates and release heights/locations. METEC has expanded its capabilities in the time since experiments were performed and since this photo was taken.

As part of the larger ADED field campaign conducted at METEC [40], an experimental schedule was generated by METEC staff to vary the location, duration and emission rate of the NG releases across the METEC site between the 19th of March and 20th of September. The location and duration were shared with us before each measurement. For each measurement, a portable methane analyzer was positioned downwind of the known release point and left in place for at least 5 minutes. The location of each measurement point was determined using a real-time kinematic (RTK) GPS unit (SparkFun RTK Facet L-Band, Boulder, CO, USA). After each 5-minute measurement, the methane analyzer was moved approximately 10 meters downwind from the emission point; this was repeated between four and seven times or until the wind changed direction

2.2. Meteorological Measurement

Meteorological data were measured at 1 Hz using the R.M. Young 81000 3D sonic anemometer (R.M. Young, Traverse City, MI, USA) permanently mounted at METEC, (Figure S1), a LI-200R-SMV (LI-COR Biosciences, NE, USA), and a Met One 597 (Met One Instruments, Grant Pass, OR, USA). The sonic anemometer was positioned 6.7 meters above ground with a maximum horizontal distance between the sonic anemometer and controlled release of 147 meters. The Met One 597 measured at 6.2 m above the ground and was located at the sonic anemometer tower. The LI-200R-SMV was located 270 m from nearest source.

The sonic anemometer provided wind speed (accuracy of ± 1%) and wind direction data (accuracy of ± 2°). The Met One 597 provided ambient air temperature (accuracy of ±2 °C), RH (accuracy of ±2%), and barometric pressure (P; hPA). Solar radiation was measured using the LI-200R-SMV model pyranometer located 250 m away from METEC over similar surface characteristics (i.e., roughness length, displacement height, topography). The solar radiation delta-T (SRDT) method was used for the PGSC classification (Table S1 and S2) [41].

2.3. Methane Mixing Ratio Measurement

Downwind methane mixing ratio data were measured using an ABB Micro-portable Greenhouse Gas Analyzer (MGGA; ABB, Zurich, Switzerland) [42]. The analyzer was connected to an inlet 1 m above the ground using 3.175 mm ID tubing (Figure S2). A sensor inlet height of 1 m above the ground is consistent with the U.S. Environmental Protection Agency’s Other Test Method (OTM)-33A near-field surface sampling guidelines [43]. The sample cell flow rate of the MGGA was 2.6 ± 0.2 min^-1, measurement accuracy of the order of 1% of the reading, and a resolution of 0.9 ppb (at 1 second) and 0.3 ppb at 3s [42].

Before each downwind measurement, the background methane mixing ratio was measured for 15 minutes and averaged. Due to in-built sampling rates, the MGGA recorded mixing ratio data at two sampling rates: 5 Hz (4.967 Hz to 6.597 Hz) for experiment 2 and 1 Hz (1.001 Hz to 1.004 Hz) for all others.

2.4. Controlled Release System

The controlled release experiment was performed between March 19^th and September 20^th, 2024. The gas composition of the controlled release NG was determined using a gas chromatograph (7890 GC, Agilent Technologies, Santa Clara, CA, USA) with nominally 85 mol% CH₄. The NG release rate ranged between 0.53 to 6.03 kg h^-1 (0.45 to 5.13 kg CH₄ h^-1) which was comparable to moderate to large fugitive leaks in upstream O&G operations [13,43]. Field studies at NG production and transmission facilities have shown that the most equipment-level methane leaks are generally under 1 kg CH₄ h^-1 [14], while component-level leaks in the storage and transmission sector typically range between 0.05 and 5 kg CH₄ h^-1 [44]. The selected release rates of NG fall in the range that represent emissions comparable to moderate fugitive leak events in upstream O&G operations. The controlled releases during our experiment were from separators, tanks, and wellheads at heights between 0.4 and 4.7 m above ground.

After each experiment, methane mixing ratio, meteorological and controlled release data were collated, time-synched and merged into a single data file. The number of measurements in each experiment is presented in Table S3. The 5-minute averages of each parameter were then determined for each measurement. Even though previous studies used >10-minute averaging periods in the past [29,36,38], a 5-minute averaging period was chosen due to practical constraints associated with conducting single controlled-release experiments during the ongoing ADED campaign. In addition, the frequent shift of wind direction at the METEC site limits the steady plume condition where longer averaging periods would have required multiple relocations of the gas analyzer, leading to additional variability in source-sensor geometry.

2.5. Modeling

2.5.1. Gaussian Plume (GP) Approach

The GP approach is used to model the dispersion of pollutants from a continuously emitting point source. The model assumes three-dimensional steady-state dispersion in the downwind direction following a normal (Gaussian) distribution in the vertical (z, m) and crosswind/lateral (y, m) directions from the centerline of the plume (x, m) [31,32,34,35]. The horizontal (σ_y, m) and vertical dispersion coefficients (σ_z, m) were obtained from distances downwind and PGSC (Tables S3 and S4). The mixing ratio C (μ g m⁻³) at any given point downwind in 3-dimensional space (x, y, z) can be calculated from the mean wind speed (U, m s^-1) and emission rate (Q, g s^-1) from the source of height (H, m) (Equation 1). The equation 1 was applied in its inverse form substituting mixing ratio and solving algebraically for unknown emission rate Q.

$$C\left(x,y,z\right)={\displaystyle \frac{Q}{2\pi U{\sigma }_z}}{e}^{-{\displaystyle \frac{y^2}{2{\sigma }_y^2}}}\left[\left(e^{-{\displaystyle \frac{{\left(z-H\right)}^2}{2{\sigma }_z^2}}}\right)+\left(e^{-{\displaystyle \frac{{\left(z+H\right)}^2}{2{\sigma }_z^2}}}\right)\right]$$

(1)

The Gaussian Plume approach makes several assumptions about the behavior of the emission source and plume, these include: 1. the source emits gas at constant rate; 2. the sensor is positioned more than 100 m from source; 3. the terrain between the source and the sensor is homogeneous; 4. meteorological conditions are horizontally uniform and steady over the averaging period; 5. dispersion in crosswind and vertical directions follows a Gaussian distribution; and 6. effects of complex surface roughness and flow recirculation are negligible [45]. One major shortcoming of this approach is that emissions can always be calculated even if any of the assumptions are violated [26,45].

The model is mathematically simple and efficient in terms aof computation and implementation. However, it is not clear how useful this approach is for calculating emissions in the real-world, as the model is based on a steady-state assumption of emission with constant meteorological conditions over homogenous terrain. Therefore, it is unlikely that the GP approach can be used for conditions at oil and gas sites which include flow obstructions and are often located in heterogenous surroundings [38,45].

2.5.2. Backward Lagrangian Stochastic (bLS) Approach

The backward Lagrangian Stochastic (bLS) approach tracks the trajectories of released particles in the atmosphere [27,36,37,38] inversely from the detector to the emission source using stochastic processes that simulate the turbulent nature of the dispersion [27]. The bLS model used in this study was WindTrax (version 2.0.9.7) developed by Thunder Beach Scientific (Nova Scotia, Canada) [37]. Data input to the model comprised 5-minute averages of wind speed, wind direction, air temperature, barometric pressure, PGSC, CH₄ mixing ratios, background mixing ratios, roughness length, meteorological measurement location/heights, source location/height, and measurement location/height [Table S4].

The initial number of particles used in the inversion for each 5-minute measurement was 500,000. In the advent of a quantification error, the number of particles was increased in steps of 1 million from 1 million to 20 million until the error was resolved. Despite increasing particle count, 18 samples did not resolve the errors using 20 million particles. The average particle count used for the simulations was 6.4 million.

2.6. Performance Analysis

The calculated emissions were compared against the known controlled-release values to evaluate models’ performance in varying meteorological conditions and plume geometries. Metrics used to assess model performance for this study were: 1. calculated emissions within a factor of two of the known emission rate (FAC2) [equation S1]; 2. the geometric mean (GM) bias calculated by using geometric mean of the ratio of estimated to known emission [equation S2]; 3. the median factor of error (MFoE) calculated by taking exponential of median absolute logarithmic error between estimated and known emission [equation S3] [47,48]; 4. mean absolute error (MAE) calculated as the mean of absolute difference between estimated and known emission rate [equation S4]; 5. the root mean square error (RMSE) calculated as a square root of mean of the squared differences between estimated and known emission [equation S5]; and 6. mean absolute bias (MB) calculated as the mean of the signed difference between estimated and known emission [equation S6].

The FAC2 measures the proportion of emissions calculated by the model that are within a factor of two of the known emission rates, i.e. +100%, -50%. This was used because factor of two criterion is widely adopted in model evaluation in atmospheric dispersion as a practical benchmark for acceptable agreement subject to measurement uncertainty and meteorological variability [47,48]. The GM identifies over/underestimation of model prediction, returning a value of one for an unbiased model, a value less than one for a model that underestimates the emission, and a values greater than one for a model that overestimates emissions. The MFoE calculates a value representing the multiplicative spread of error (always ≥1). MFoE values close to 1 indicate a model that calculates emissions with a small deviation from the known value. These three indicators, when taken together, are used to understand the veracity of the model’s calculated emissions in varying environmental conditions.

The MAE measures the average absolute magnitude of error between the estimated and true emission rates, regardless of direction. The RMSE predicts the square root of the average of the squared errors, which is highly sensitive to outliers. Unlike MAE and RMSE, MB shows the bias direction of the error. The positive error in MB means that the model is overestimating, and the negative means that it is underestimating. These three indicators, when taken together, are used to understand magnitude, variability and directional bias of the model estimation errors.

2.7. Statistical Analysis

Numerical performance metrics alone are insufficient to explain the variation in the model’s predictive capacity and emphasize physical interpretation [47,48]. To better understand the model’s performance, statistical analysis was conducted to evaluate the effect of other factors like release rates, signal strength, wind meander, and plume geometry. In each case, the effect of each factor was evacuated using the FAC2 metric to identify optimal conditions for each environmental or geometric factor.

Signal strength: Signal strength is defined as the 5-minute averaged mixing ratio reported by the methane analyzer. Signal strength bins were grouped into four evenly populated quartiles: ≤ 2.3 ppm, 2.3–2.9 ppm, 2.9–4.0 ppm, and >4.0 ppm.

Release rates: Release rate bins were created to determine if a model performs consistently across varying source strengths. Release rate bins (kg h^-1) were defined as: low for emission rates less than 1.8 kg h^-1; medium between 1.8-3.4 kg h^-1; and high between 1.8-3.4 kg h^-1.

Atmospheric stability: Atmospheric stability was categorized using PGSC and derived from 5-minute averages of solar radiation and wind speed (Table S1 and S2). All measurements at METEC were unstable daylight conditions, specifically classes A (extremely unstable), B (moderately unstable), and C (slightly unstable).

Wind meander: Wind meander captures the directional variability of the wind direction field [18,30]. In this study, we define wind meander as the standard deviation of wind direction during the 5-minute averaging interval and indicates whether dispersion is governed by ambient airflow or turbulence (mechanical or convective). Strong variation in wind direction can disrupt the steady state assumption of GP, thereby reducing the reliability of model estimates [36,48]. Wind meander data categorized calculated emissions into seven 5° bins (±0°-5°, ±5°-10°, ±10°-15°, ±15°-20°, ±20°-25°, ±25°-30° and > ±30°.

Plume centeredness: The plume centeredness score is the ratio of crosswind distance (y) to the horizontal dispersion coefficient (σ_y), a measure of lateral plume spread [32]. This score indicates the relative distance of a measurement point to the centerline of the downwind plume. Larger values indicate the measurement position was near the plume edge with smaller values corresponding to the plume axis [49]. Centeredness score was categorized into four bins: 0-0.5; 0.5-1; 1-2; and >2.

Plume height: The height match score is determined by using the ratio of the absolute difference between source and sensor height to the plume’s vertical dispersion coefficient, σ_z, and used to indicate the extent of vertical overlap between the sensor and plume heights. In this analysis, four height match score bins (0-0.5, 0.5-1, and 1-2) are used to indicate the relative vertical distance between the analyzer inlet and the plume centerline and are expressed as multiples of plume height standard deviation (σ_z).

Plume centerline distance: The plume centerline distance represents the downwind distance between the source and the analyzer inlet measured along the plume’s trajectory. Dispersion and dilution increase with downwind distance. At smaller downwind distances the plume is likely to be more concentrated, while farther downwind, mixing is enhanced and signal strength is smaller [34,49]. To evaluate this, the plume centerline distances were binned into three categories: 0-40, 40-80, and >80 m.

3. Results

3.1. Measurement Data

3.1.1. Methane Mixing Ratio

Between March and September 2024, 121 5-minute methane mixing ratio measurements were made downwind of single point source releases. These data comprise 5-minute averaged CH₄ mixing ratios measured at 1 m above ground level. The observed average CH₄ mixing ratios ranged from a minimum of 1.8 ppm to a maximum of 21.1 ppm. The median of mixing ratio across the measurement set was 2.8 ppm, and the interquartile range spanning 2.3 to 4 ppm.

3.1.2. Meteorology

The 5-minute wind speeds during the measurements ranged between 0.9 and 7.9 m s^-1, with a median of 3.1 m s^-1 and an average of 3.3 ± 1.4 m s^-1. The data was filtered for wind speeds <1 m s^-1, as low speeds do not provide adequate transport conditions for a well-defined plume and are not typically included in downwind dispersion analysis [36,38,50]. The 5-minute averaged wind direction data included sampling winds from nearly all directions (between 18 to 342°). The measurements were conducted in air temperatures between 3.7 and 33.4 °C. The average 5-minute relative humidity during measurements varied between 10% to 59% (Table S5) across all measurements, reflective of the dry atmospheric conditions of Northern Colorado.

3.2. Emissions

3.2.1. Calculated Emissions

Emissions for each downwind measurements were calculated using both the GP and bLS approaches. Emissions calculated using the GP approach range from 4.2×10^-27 to 14.7 kg CH₄ h^-1, with median value of 0.92 kg CH₄ h^-1 and mean of 2.03 kg CH₄ h^-1. Corresponding bLS estimates range from 0.043 to 22.52 kg CH₄ h^-1, with median value of 1.42 kg CH₄ h^-1 and mean of 2.7 kg CH₄ h^-1.

3.2.2. GP and bLS Performance Analysis

Overall, 25% and 52% of emissions calculated using the GP approach and the bLS model were within the FAC2 bounds, respectively (Table 1). This percentage was calculated from the total emissions within the FAC2 bounds and the number of samples, with a total number of samples for GP and bLS of 121 and 103, respectively. The geometric means (GM) of the GP approach are 0.08 and bLS is 0.63. The median factor of error (MFoE) of the GP and bLS models are 3.5 and 1.9, respectively. The MAE and RMSE were lower for bLS calculated emissions (1.5 and 2.6 kg h^-1, respectively) compared to GP predictions (2.4 and 3 kg h^-1, respectively). The MB of the bLS model has a smaller residual bias (0.04 kg h^-1) than the GP approach (-0.8 kg h^-1).

When the emission rate calculated using the GP approach is plotted against known emission rates (Figure 2A), the low R² value (0.05) reflects the weak linear relationship between estimated and known emissions, while the gradient (m = 0.45) shows the calculated emissions are largely underestimated. The line of best fit between emissions calculated using the bLS approach to the known emission rates has a gradient of 1.26 and R² of 0.27 (Figure 2B).

3.3. Statistical Analyses

For both modeling approaches, the highest number of calculated emissions between the FAC2 bounds was for a signal strength greater than 4 ppm, although the bLS modeled data was in better agreement (30% for GP and 69.2% for bLS). For the GP approach, the FAC2 percentage suggested a better agreement when the signal strength was greater than 2.9 ppm (Figure 3; Table S6). For the bLS modeled data, the FAC2 percentage is higher in all signal strength bins compared to GP, ranging from 27% at the lowest signal strength to 69% in the highest bin.

4. Discussion

4.1. Comparative Performance of GP and bLS

In this study, emissions calculated from measurements downwind of controlled-single source release using a bLS approach was more representative than emissions calculated using a GP approach. Overall, the percentage of emissions between the FAC2 bounds (± factor of 2) using the GP approach was 25% compared to 52% for emissions calculated using the bLS model (Table 1). Similarly, the multiplicative predictive error was higher for the bLs data (3.5) compared to the GP data (1.9) as was the geometric mean bias (bLS data underpredicted by 37% while the GP data were an underprediction of 92%). For this dataset, both modeling approaches underpredict emissions, but emissions calculated using the bLS model were notably closer to known controlled release rates. Error-based metrics, e.g. MAE, RMSE, and MB, also indicate the better agreement between the bLS emission data and known release rates. The MAE and RMSE were lower for bLS emission estimates (1.5 and 2.6 kg h^-1, respectively) compared to GP predictions (2.4 and 3 kg h-¹, respectively). The MB of bLS data has a small residual bias (0.04 kg h^-1), while the underprediction of the GP data is much larger (-0.8 kg h^-1).

We suggest the difference in the performance by the two modeling approaches is due to underlying fundamental physical assumptions. The GP approach is based on the steady-state transport of a plume under uniform wind conditions and homogenous terrain; these assumptions are frequently violated in O&G sites. As a result, GP modeled emissions are strongly influenced by deviation from idealized plume geometry and meteorological stability. The bLS simulates individual random particles driven by turbulent atmospheric variability, allowing it to capture the varying nature of particle trajectories within a plume. This stochastic nature of the bLS model adapts to variation in wind conditions that leads to better near-field estimations. Unlike bLS, the GP approach depends on parameterized dispersion coefficients to represent average plume dispersion, restricting the ability to capture and simulate turbulence variations directly. As a result, bLS achieves better agreement with true emission rates as it can simulate the behavior of the plume more realistically than the GP approach (Figure 2 and Figure 4).

Overall, the bLS model had greater predictive reliability compared to the GP model, with approximately twice as many calculated emissions falling between the FAC2 bounds. While the GP approach was faster and more convenient for calculating estimations, emissions were comprehensively underestimated due to wind variability and complexity in plume geometry. Conversely, the bLS model, though computationally more intensive, more accurately estimated the emissions under unstable near-field O&G conditions.

Practically, the performance advantages of the bLS approach must be interpreted in the context of its operational constraints. While bLS yields more accurate near-field estimates, it requires high particle count, stochastic trajectory simulations to internally represent turbulent transport, and final inversions, making it computationally demanding and impractical for real-time monitoring applications. In contrast to bLS, the GP method, although prone to underpredicting emission estimates, remains practical for real time application due to its ease of implementation and computational simplicity.

4.2. Statistical Analyses

The statistical analysis was used to identify the geometric and environmental conditions associated with improved agreement with known release rates and to provide physical interpretation for the performance differences of GP and bLS approaches. Rather than relying solely on aggregate performance metrics, model agreement with known rates was assessed using FAC2 criterion across binned subsets of signal strength, release rate, plume centeredness, height match score, wind meander and atmospheric stability.

The bLS model exhibits higher FAC2 percentage in all signal strength bins compared to GP, with percentages ranging from 26.9% at the lowest signal strength to 69.2% in the highest bin. This suggests bLS is more effective with higher signal strength. We would expect a higher FAC2 percentage with increasing signal strength, but the pattern of the results is not entirely consistent. Weaker signals typically reflect that the plume overlaps the sensor partially, limiting the detection signal, leading to inaccuracy in emission estimation [36,38]. Collectively, the result indicates that sufficient signal strength is essential for a reliable estimation, emphasizing the importance of sensor placement strategies that maximize mixing ratio enhancement above background levels.

The release rate binning describes the contrasting responses of model to source strength. Results reinforce the assertion that the performance of the GP model does not systematically improved with increasing emission magnitude. By comparison, bLS shows clear dependence of release rate magnitude on its performance, with FAC2 percentage increasing from 45% for low emissions to 80% for high emissions. These findings imply that higher release rate alone does not necessarily improve GP-based quantification, whereas bLS performance improves when sufficient source strength is present. From the measurement perspective, sensor placement determines whether increased release rates translate into improved model performance, especially for bLS, since only configurations that measure the plume allow higher emission magnitudes to convert into higher signal strength.

Atmospheric stability strongly influences both models’ output, as it is incorporated into the formulation [34,36,49]. The bLS model comprehensively outperforms the GP model in each PGSC category based on FAC2 percentage. But, in the outcomes from both the models there is limited variation in FAC2 percentage between stability classes. Therefore, variation in stability has not systematically translated into improved performance.

Wind meander score indicates whether dispersion is governed by steadier airflow or by turbulence (possibly from rapidly shifting weather or from nearby equipment obstructions, or a combination of both). Strong variation in wind direction can disrupt the steady state assumption of GP, thereby reducing the reliability of model estimates [48,49]. The bLS model showed greater resilience as wind variability increases. The FAC2 percentage for bLS exceeded 50% in four meander bins in the bLS analysis. This increased robustness reflects the particle-based formulation of bLS, in which the wind direction variability is explicitly represented through stochastic transport. The analysis reflects that wind meander is critical factor in dispersion modeling, highlighting the importance of consistent plume interception through optimized sensor placement and selection of sampling durations that characterize moderate directional variability.

When the plume centeredness score is small (0-0.5), a small deviation in wind direction may cause large swings in GP estimations [49], especially for shorter averaging intervals that may miss short-term deviations in the signal, much like the 5-minute measurements used in this study. The bLS approach reflected greater robustness to centeredness but is less sensitive towards varying bins. Overall, the result underscores the importance of lateral source-sensor alignment for better estimations using GP and bLS.

The FAC2 percentage showed a decreasing pattern with increase in the height match score (Table 3.7). But the pattern is completely opposite for bLS model as the FAC2 percentage increased with increasing height match score. But this description is deceiving because the sample distribution across the height match score bin is highly inconsistent. There are only 4 and 1 samples in the score range of 1-2 for GP and bLS among 121 and 103 samples, respectively. This may not represent the analytical conclusion of models performing better at this height score. However, in 0-0.5 match score range, where maximum samples are available (98 and 90 for GP and bLS, respectively), the bLS had a better FAC2 percentage.

Dispersion and dilution increase with plume centerline distance, which can influence the model’s predictability. At a shorter distance, the plume is more concentrated, while with increasing centerline distance, mixing is enhanced, and signal strength is also smaller [34,49]. Disregarding the bins with smaller sample sizes (>80 m; only 7 out of 121 and 103 for GP and bLS, respectively), the highest FAC2 percent for both modeling approaches occur at distance of between 40 and 80 m downwind of the source. This indicates the favorable range of centerline distance for sensor placement likely between 40 and 80 m downwind of an O&G facility. We suggest 40 to 80 m downwind is the distance at which the plume sufficiently disperses and concentration stabilizes, leading to better detection, whereas at shorter distances (<40 m), the plume is still developing (with small lateral size) and sensitive to small turbulence (especially in unstable conditions).

Overall, the statistical analyses show that near-field methane quantification performance by GP and bLS approach is driven by combined effects of signal strength, sampling geometry, and plume interception rather than by individual meteorological factors alone. For both modeling approaches, reduced agreement was associated with poor plume alignment and low signal strength, suggesting estimates based on weak plume interception may needed to be filtered from quantitative analysis. The result also suggests that prioritizing downwind sensor placement directly downwind of likely emission source and at distances greater than 40 m may result in more accurate emissions calculated using GP and bLS-based approaches for O&G facilities.

4.3. Limitation and Scope

This study was constrained by several limitations that pre-defined the scope of the experiments and gives bounds to which the results should be interpreted. First, the study is based on only one experimental site (METEC) representative of mock upstream O&G in dry windy conditions in Northern Colorado and featured a facility with single controlled releases. As such, the study is unlikely to fully represent the dispersion complexities of different aerodynamic and emissive environments. Second, all the measurements for this study were performed under daytime unstable atmospheric conditions (PGSC classes A-C) which limits the representation of stable or neutral conditions which are common on nighttime operations or during synoptic transitions like frontal passage.

Third, emission estimates calculated by the models were derived from relatively short-duration measurements (5- minutes). Although, this choice may enhance operational relevance for leak detection and quantification, it increases sensitivity to plume source-sensor geometry and wind meander. The use of a longer averaging interval may result in different performance behavior, particularly for GP in which the formulation relies on steady-state assumptions.

Fourth, the bLS model is influenced by stochastic sampling limitations and computational constraints. Under more turbulent atmospheric conditions, bLS generally offers more reliable estimates compared to GP, provided that sufficient particle counts and accurate meteorological inputs are used. However, calculating emissions using the bLS model takes longer than the GP approach due to its computational demands. Although particle counts were progressively increased to improve numerical stability, some samples failed to provide the estimates even at maximum particle number employed (defined as 20 million in this study).

Finally, this study does not attempt optimize model formulations or develop new dispersion algorithms. Instead, both the models were implemented using established and commonly accepted configurations to evaluate model performance in near-field conditions using the same underlying dataset. Future studies could expand the analysis presented here to include complex terrain, longer averaging durations, multiple simultaneous sources, and stable/neutral atmospheric conditions. While it is likely that the results presented here are suitable for use in leak detection technology in the Denver-Julesburg Basin, it is currently unknown if it will be useful in more vertiginous or heavily vegetated oil and gas producing regions of the Marcellus or Uinta Basins.

5. Conclusions

This study evaluated the performance of backward Lagrangian Stochastic (bLS) and Gaussian Plume (GP) dispersion approaches by comparing calculated CH₄ emissions to controlled release single-point emissions in near-field conditions at a simulated O&G production facility. The emission estimates were calculated using 5-minute averaged downwind mixing ratio measurements taken at the METEC facility and were compared against the true emission rates across a range of emission strengths, meteorological and stability conditions, and plume geometries. The bLS approach showed better predictive performance compared to those calculated using the GP approach. The data calculated by the bLS model also had a lower multiplicative error and reduced bias relative to GP. Other error-based metrics further confirmed the bLS model performed better, as it yielded lower RMSE and MAE than GP.

Statistical analyses of the emission data show that the lateral and vertical alignment of source and sensor plays a critical role in emission estimations as measurements made closer to the plume centerline and at a distance between 40 and 80 m downwind yielded the best FAC2 agreement. High wind meander degraded ability of both approaches to generate representative emissions particularly with the GP approach as it potentially violates the modelling approach’s assumption of steady-state emissions. Data suggest emissions calculated by the bLS model are comprehensively in better agreement but the computational demands of the modeling approach and integration into fenceline systems limit real-time applicability.

Overall, the study clarifies the conditions under which GP and bLS modelling approaches produce reliable near-field CH₄ emission estimates which may provide practical guidance for measurement design and interpretation in O&G applications. Realistically the bLs model can be used to calculate emissions within a factor of two 52% of the time, which may also improve as calculations are made using high concentration data that are measured directly downwind at distances greater than 40 m of the O&G production facility. Data suggests emissions calculated by the GP approach are likely an underestimate and that filtering for signal strength or instrument positioning will improve the estimates by a small extent.