ML22363A001
| ML22363A001 | |
| Person / Time | |
|---|---|
| Issue date: | 11/30/2022 |
| From: | Fagan D, Huckett J Office of Nuclear Regulatory Research, Pacific Northwest National Laboratory |
| To: | |
| References | |
| DE-AC05-76RL01830, 31310019N0001, 31310021F0022, 1b TO 31310021F0022 PNNL-33647 | |
| Download: ML22363A001 (1) | |
Text
PNNL-33647 Subsurface Radiological Survey Design and Geospatial Analysis Tool Recommendations Task 1b TO 31310021F0022 November 2022 Jen Huckett Deb Fagan Zachary Weller Moses Obiri Lisa Newburn Fred Day-Lewis David Peeler Prepared for the U.S. Nuclear Regulatory Commission Office of Nuclear Regulatory Research Under Contract DE-AC05-76RL01830 Interagency Agreement: 31310019N0001 Task Order Number: 31310021F0022 Choose an item.
DISCLAIMER This report was prepared as an account of work sponsored by an agency of the United States Government. Neither the United States Government nor any agency thereof, nor Battelle Memorial Institute, nor any of their employees, makes any warranty, express or implied, or assumes any legal liability or responsibility for the accuracy, completeness, or usefulness of any information, apparatus, product, or process disclosed, or represents that its use would not infringe privately owned rights. Reference herein to any specific commercial product, process, or service by trade name, trademark, manufacturer, or otherwise does not necessarily constitute or imply its endorsement, recommendation, or favoring by the United States Government or any agency thereof, or Battelle Memorial Institute. The views and opinions of authors expressed herein do not necessarily state or reflect those of the United States Government or any agency thereof.
PACIFIC NORTHWEST NATIONAL LABORATORY operated by BATTELLE for the UNITED STATES DEPARTMENT OF ENERGY under Contract DE-AC05-76RL01830 Printed in the United States of America Available to DOE and DOE contractors from the Office of Scientific and Technical Information, P.O. Box 62, Oak Ridge, TN 37831-0062; ph: (865) 576-8401 fax: (865) 576-5728 email: reports@adonis.osti.gov Available to the public from the National Technical Information Service 5301 Shawnee Rd., Alexandria, VA 22312 ph: (800) 553-NTIS (6847) email: orders@ntis.gov <https://www.ntis.gov/about>
Online ordering: http://www.ntis.gov Choose an item.
PNNL-33647 Subsurface Radiological Survey Design and Geospatial Analysis Tool Recommendations Task 1b TO 31310021F0022 November 2022 Jen Huckett Deb Fagan Zachary Weller Moses Obiri Lisa Newburn Fred Day-Lewis David Peeler Prepared for the U.S. Nuclear Regulatory Commission Office of Nuclear Regulatory Research Under Contract DE-AC05-76RL01830 Interagency Agreement: 31310019N0001 Task Order Number: 31310021F0022 Pacific Northwest National Laboratory Richland, Washington 99354
PNNL-33647 Acronyms and Abbreviations 2-D two-dimensional 3-D three-dimensional AI artificial intelligence AIM Aquifer Injection Modeling ArcGIS Mapping and analysis software CCM contamination concern map CORE Common Operating and Response Environment CSM conceptual site model DCGL derived concentration guideline level DOE U.S. Department of Energy DOE-EM DOE Office of Environmental Management DQA data quality assessment EPA US Environmental Protection Agency EPRI Electric Power Research Institute ERT electrical resistivity tomography FRK fixed rank kriging GEM Geospatial Extension to MARSSIM GLS generalized least square GIS geographic information system GPR ground penetrating radar GSLIB Geostatistical Software Library HH Cluster of locations with low values (HH)
HL Single high value location surrounded by low value locations (HL)
HSA historical site assessment LISA local indicator of spatial association LH Single low value location surrounded by high value locations (LH)
LL Cluster of locations with low values (LL)
MARSSIM Multi-Agency Radiation Survey and Site Investigation Manual MGK multi-Gaussian kriging ML machine learning MrDM Multi-scale Remedial Design Model MrsDM Multi-scale Remedial Sample Design Model NRC U.S. Nuclear Regulatory Commission NUREG U.S. Nuclear Regulatory Commission Regulation OLS ordinary least square PHOENIX PNNL-Hanford Online ENvironmental Information eXchange Acronyms and Abbreviations iii
PNNL-33647 PNNL Pacific Northwest National Laboratory RSSI radiological survey and site investigation RT3D Reactive Transport in Three-Dimensions SADA Spatial Analysis and Decision Assistance software SIP spectral induced polarization SOCRATES Suite of Comprehensive Rapid Analysis Tools for Environmental Sites STOMP Subsurface Transport Over Multiple Phases SVEET2 Soil Vapor Extraction End-state Tool, version 2 TRAC Tracking Restoration and Closure VSP Visual Sample Plan Acronyms and Abbreviations iv
PNNL-33647 Contents Acronyms and Abbreviations...................................................................................................... iii 1.0 Introduction .....................................................................................................................1 2.0 Compliance Survey .........................................................................................................3 3.0 Subsurface Data..............................................................................................................5 3.1 Data Types ..........................................................................................................5 3.2 Data Quality Assessment .....................................................................................6 3.3 VSP Capabilities and Advancements ...................................................................8 3.3.1 3-D Sample Surfaces and Volumes .......................................................9 3.3.2 Data Visualization ................................................................................11 3.3.3 Analysis of Multiple Contaminants/Analytes......................................... 14 3.3.4 Statistical Diagnostics..........................................................................14 3.3.5 VSP in the Cloud .................................................................................17 3.4 Additional Existing Subsurface Tools and Capabilities at PNNL ......................... 18 3.4.1 Data Visualization and Analysis ...........................................................18 3.4.2 Data Management ...............................................................................23 3.4.3 Flow and Transport Modeling ..............................................................24 3.4.4 Computational and Applied Geophysical and Geomechanics .............. 28 3.4.4.3 SIP .....................................................................................................................31 3.4.4.4 Seismic Subsurface Technologies .....................................................................32 3.4.4.5 Electromagnetic Subsurface Techniques ...........................................................33 3.4.4.6 Ground Penetrating Radar .................................................................................34 4.0 Statistical Methods for Subsurface Data Analysis ..........................................................37 4.1 Geostatistics Overview.......................................................................................37 4.1.1 Structural Analysis ...............................................................................37 4.1.2 Kriging .................................................................................................38 4.1.3 Conditional Simulation .........................................................................39 4.2 Geostatistical Methods for the Subsurface .........................................................40 4.2.1 Generalized Least Squares Regression ..............................................40 4.2.2 Local Indicator of Spatial Association ..................................................45 4.2.3 Kriging Methods ..................................................................................46 4.2.4 Geospatial Extension to MARSSIM, Multi-scale Remedial Design Model (MrDM), and Multi-scale Remedial Sample Design Model (MrsDM)....................................................................................49 4.2.5 Bayesian Methods ...............................................................................50 4.2.6 Variogram Tomography .......................................................................50 4.2.7 Artificial Intelligence/Machine Learning ................................................51 4.3 Dimensionality of Approach ...............................................................................52 Contents v
PNNL-33647 4.3.1 Layered Approach ...............................................................................52 4.3.2 Model a Complex 3-D Volume .............................................................54 4.4 VSP Advancements ...........................................................................................55 5.0 Compliance Survey Design ...........................................................................................58 5.1 Number and Location of Samples ......................................................................58 5.1.1 Model Based Survey Design ...............................................................59 5.1.2 Markov Bayes and Bayesian Elipgrid................................................... 60 5.1.3 Check and Cover .................................................................................61 5.2 VSP Advancements ...........................................................................................62 6.0 Recommendations.........................................................................................................63 7.0 References ....................................................................................................................66 Contents vi
PNNL-33647 Figures Figure 1. Subsurface Flow Diagram (Figure 3.3 from Stewart and Powers [2009]).................. 3 Figure 2. MARSSIM Roadmap (NRC 2000). ...........................................................................4 Figure 3. Quality Assurance Assessment (figure from EPA [2000] with additional recommendations for subsurface compliance survey planning in green) .................. 7 Figure 4. Basic 3-D Volumes in VSP. A polygon and a circle are shown from above in the plane view (left) and from the side in the 3-D volume view (right). ......................9 Figure 5. Current 3-D Sample Placement Methods in VSP ...................................................10 Figure 6. Example Volume Created with Rectangular Contours and Associated Depths. Contours are shown from above (left) and within the 3-D volume after applying the depths of each (right). ................................................................10 Figure 7. Example 3-D Volume Surface Created by Importing a Digital Elevation Model File ..............................................................................................................11 Figure 8. Subsurface Soil Point Data Visualization in VSP. Locations with an iodine-129 value between 0 and 4.5 are color-coded according to the legend, while locations without a valid data value (i.e., less than 0) are displayed as dark blue plus symbols. Black lines are drawn from the depth of the points to the surface to illustrate depth. ......................................................................................12 Figure 9. Currently Only 2-D Visualization of Interpolated Surface is Supported in VSP........ 13 Figure 10. Variogram Surfaces Indicating Isotropic Variation (left) and Anisotropic Variation (right).......................................................................................................16 Figure 11. Example of Three H-Scattergrams at Three H-Distances (EPRI 2016)................... 17 Figure 12. SOCRATES Component Applications ....................................................................19 Figure 13. PHOENIX Overview ...............................................................................................21 Figure 14. TRAC Web-Interface Screenshot ...........................................................................22 Figure 15. Few-Shot Learning Schematic ...............................................................................23 Figure 16. CORE Data Management Workflow .......................................................................24 Figure 17. AIM Toolbox Interface Screenshot .........................................................................25 Figure 18. Visualization of Simulation Results Using RT3D.....................................................26 Figure 19. SVEET2 Conceptual Model Schematic ..................................................................27 Figure 20. Visualizations of Simulation Results for Aqueous Phase (left) and Nonaqueous Phase (right) Transport Using the STOMP Suite of Codes ................ 27 Figure 21. Visualization of PFLOTRAN Simulation Results for Uranium Transport at the Hanford Site, Adjacent to the Columbia River in Eastern Washington .................... 28 Figure 22. Visualization of Time-Lapse ERT Images Revealing Changes in Electrical Conductivity Associated with Unanticipated Leaks During Overflow Events in Subsurface Infrastructure .......................................................................................30 Figure 23. Visualization of E4D Inversion Results Showing Time-Lapse Changes in Subsurface Electrical Conductivity Associated with Amendment Release at the Land Surface ....................................................................................................31 Figures vii
PNNL-33647 Figure 24. Visualizations of SIP Imaging Mesh (left), Imaging Results (middle), and Experimental Apparatus (right) to Understand Electrical Signatures of Biological Activity ...................................................................................................32 Figure 25. Seismic Images from the Gable Gap Area of the Hanford Site. Images include a compressional-wave velocity tomogram overlaid on migrated reflection data and images of shear-wave velocity derived from multichannel analysis of surface waves (St. Clair et al., In Review)............................................. 33 Figure 26. Electromagnetic Inversion Results and Interpreted Hydrogeologic Framework from Airborne Surveys Collected Over the Hanford Site (Jaysaval et al. 2021) .............................................................................................34 Figure 27. Inversion Results for Cross-Hole GPR Data Collected Before (left) and After (right) the Start of a Soil Desiccation Experiment at the Hanford Site (Truex et al. 2018). ............................................................................................................34 Figure 28. Methodology Suggestions as Scale of Investigation and Resolution Increase ........ 35 Figure 29. Commonly Used Variograms in GLS Estimation ....................................................38 Figure 30. Concentration of U-235 at Four Selected Wells at the Hanford Site ....................... 42 Figure 31. Comparing OLS and GLS Type 1 Error Rates ........................................................43 Figure 32. Example Well Monitoring Data (simulated) in Seven Wells that Span Three Subsurface Layers. (a) Constant mean concentration. (b) Varying mean concentration. Spatial variation is not shown. .........................................................44 Figure 33. Example LISA Application to Well Monitoring Data. a) Data and b) hotspot analysis. .................................................................................................................45 Figures viii
PNNL-33647 Tables Table 1. Data Sources Available Prior to Compliance Survey ................................................5 Table 2. Example OLS and GLS Model Output ....................................................................42 Table 3. Summary of GLS Models Used to Estimate Subsurface Means with Three Layers ....................................................................................................................44 Table 4. Recommended Updates to VSP Software for Subsurface Compliance Surveys and Geostatistical Analysis .......................................................................64 Tables ix
PNNL-33647 1.0 Introduction The U.S. Nuclear Regulatory Commission (NRC) is planning to provide guidance for the decontamination and decommissioning of subsurface soils. NRC publication NUREG-1575 (NRC 2000), Multi-Agency Site Survey Investigation Manual (MARSSIM), provides information about planning, conducting, evaluating, and documenting the final status of building surface and surface soil radiological surveys conducted to demonstrate compliance with dose or risk-based regulations or standards. It details a data quality objectives (DQOs) process to produce characterization survey results of sufficient quality and quantity (EPA 2006). Included are survey design processes and guidance on selecting appropriate measurement methods (e.g., scan surveys, direct measurements, fixed laboratory samples) and measurement systems (e.g.,
detectors, instruments, analytical methods). Data quality assessment (DQA) is another important process described therein for assessing data collected through surveys and subsequent analysis to determine whether the data quality satisfies survey objectives so that survey results can be interpreted correctly, as they apply to the decision being made (EPA 2000).
Visual Sample Plan (VSP) software is currently used when designing compliance phase surveys, but it does not include modules that explicitly consider the subsurface. Subsurface conditions lead to increased complexity relative to surface-only applications. The added third dimension, as well as multilayered geophysical phenomena, may affect potential radiological contaminant fate and transport, accessibility, and survey methods. Additionally, dose modeling that incorporates exposure pathways from possible subsurface contamination increases the complexity for establishing compliance. Subsurface survey samples are much more expensive to collect than surface samples because of the associated physical challenges, equipment, and detection capabilities required for accessing the subsurface. Such conditions demand new survey design and data analysis methods that incorporate efficiency gains to potentially reduce the total number of samples required while still providing adequate confidence in survey results and subsequent decisions. Potential gains in efficiency can be obtained by capitalizing on nonanalytical data collected during previous phases of the radiological survey and site investigation (RSSI) process, including expert judgment, prior knowledge, geographic information system data, and geostatistical tools (Stewart and Powers 2009).
This report provides a review of analytical methods applicable to subsurface soil decontamination and decommissioning. It identifies recommended updates to VSP software in support of subsurface compliance phase survey design and geostatistical analysis. This review considers subsurface complexities and practical constraints on survey sampling that will require DQA, data visualization, geostatistical modeling, and methods to determine the required number of samples and the location at which they are taken. It prioritizes VSP updates based on current capabilities, the ease of expanding them from two dimensions to three dimensions, requirements for new algorithm development, and the applicability of each method to compliance phase activities.
Making decisions about free/unrestricted release or restricted release, determining the Lower Bound of the Grey Region, and selecting which hypotheses need to be tested are not within the scope of this report. This report assumes such decisions have been made prior to the compliance survey and focuses instead on enabling VSP users to develop surveys that will implement those hypothesis tests and draw conclusions with confidence based on the collected data. Additionally, developing a methodology for dose modeling (including both surface and Introduction 1
PNNL-33647 subsurface elements for all dose pathways) is not within the scope of this report. Both the remedial action support survey and compliance survey require these activities to have been completed during previous outside efforts. Guidance on dose modeling to support development of derived concentration guideline levels (DCGL) is provided in NUREG-1757, Volume 2, Rev. 2 (NRC 2022). This report assumes that such work has been completed based on outside analyses and that appropriate decision limits and authorized limits have been determined. This report instead focuses on identifying necessary VSP improvements to enable VSP users to provide input, based on the previously completed analyses, to the compliance survey design.
Many concepts covered in this report were presented at the second NRC subsurface soils workshop on May 11, 2022 (Huckett and Fagan 2022), and additional details about methods presented at the workshop are included in this report.
The report is organized into the following sections:
- Section 2 discusses the scope of this report as it pertains to subsurface survey and geostatistical methods for the compliance phase of the RSSI process.
- Section 3 discusses types of subsurface data that could be available for the compliance phase and final status survey; the importance of DQA; current VSP capabilities for data import, manipulation, and DQA; and other existing subsurface tools and capabilities at Pacific Northwest National Laboratory (PNNL).
- Section 4 discusses geostatistical methods for evaluating subsurface data in support of compliance phase decision-making; describes approaches to such analyses, including a layered approach and a complex volume approach; and outlines required VSP advancements to support each approach. Some functionality required to implement the methods is available through libraries already in use by existing VSP algorithms (such as the Geostatistical Software Library), and the VSP code simply needs to be updated to activate the vertical dimension. Other functions will require new methods and code to be developed and programmed in VSP.
- Section 5 discusses methods for compliance survey design to determine the number and location of samples, as well as required VSP updates to apply these in the subsurface.
- Section 6 contains a summary of VSP updates. Throughout this report, recommendations are provided for methods and capabilities that should be added to VSP software to improve the ability of users to navigate subsurface investigations and decision-making. In many cases, existing resources such as those available in R statistical software (R Core Team 2021), Spatial Analysis and Decision Assistance software (SADA) (Stewart et al. 2006), and Geostatistical Software Library (GSLIB) (Deutsch and Journel 1998) could be leveraged to make these improvements in VSP. While the details are available, proposed updates are summarized in Section 6.
- Section 7 provides references.
Introduction 2
PNNL-33647 2.0 Compliance Survey The subsurface flow diagram in Figure 1 is intended to guide the RSSI process for the three-dimensional (3-D) subsurface. A detailed process, such as that outlined in the MARSSIM Roadmap, is needed for demonstrating subsurface compliance within the compliance phase.
Figure 1. Subsurface Flow Diagram (Figure 3.3 from Stewart and Powers [2009]).
For reference, MARSSIM provides information about planning, conducting, evaluating, and documenting radiological building surface and surface soil final status surveys to demonstrate compliance with dose- or risk-based regulations or standards (NRC 2000). All RSSIs complete at least the historical site assessment (HSA) and characterization survey, although the phases between these may not all be completed. The MARSSIM Roadmap (NRC 2000) is intended to guide the characterization survey process for two-dimensional (2D) surfaces, but it can be used for previous phases, if applicable. The MARSSIM Roadmap has multiple phases, as shown in Figure 2 (NRC 2000). A similar detailed roadmap corresponding to the subsurface compliance phase should be developed.
Compliance Survey 3
PNNL-33647 Figure 2. MARSSIM Roadmap (NRC 2000).
Compliance Survey 4
PNNL-33647 3.0 Subsurface Data The increased complexity of subsurface conditions results in more costly collection of subsurface survey samples, demanding new survey design and data analysis methods that can potentially reduce the total number of samples required and still provide adequate confidence in survey results and subsequent decisions. Such efficiency gains can be provided using data collected during previous RSSI phases, and based on expert judgment, prior knowledge, geographic information system data, and geostatistical tools (Stewart and Powers 2009). This section describes the types of data (and sources, such as other existing subsurface tools) that may be available for compliance survey planning and execution and DQA activities required to process and evaluate the suitability of such data. It discusses current VSP capabilities and additional capabilities recommended for future VSP development to provide better tools to VSP users.
3.1 Data Types The purpose of this section is to outline types of data expected from RSSI phases prior to the compliance survey that are available for subsurface compliance survey planning and analysis.
The compliance survey dataset can be combined with geophysical data (e.g., groundwater, soil types) for survey development purposes. Table 1 provides details about the types of geophysical data that can be used and potential uses for them. Note that this is a list of possible data types that the NRC might receive from the licensee, but not all document types will be available or required for every site and such availability will be determined site-by-site.
Table 1. Data Sources Available Prior to Compliance Survey Data Source Dataset Description(s)
Prior to RSSI Engineering drawings (facilities, structures, etc.), operations logs, geographic information system (GIS) maps, background geophysical data (surface/subsurface), water resource characterization and climate data for conceptual site model development.
HSA Risk assessment, hazard assessment, RCRA/CERCLA documentation (including soil/rock core sample data as appropriate), NEPA documentation (as appropriate), source term quantification modeling/estimates for relevant sites, contaminant fate and transport modeling.
Scoping Remedial investigation/feasibility study reports, updated conceptual site models, specification of sampling types/design/media/location, proposed statistical methods, identification, and characterization of potential contaminant plumes.
Characterization Geologic maps, soil maps, drillers logs, maps of site infrastructure, collection of groundwater levels, hydraulic tests, soil or rock cores, and development of GIS, visualizations, and maps for the site, surveillance monitoring data from previous remediation activities (if applicable),
geophysical and hydrogeological modeling results.
Remediation Characterization of plume structure and composition, conceptual site model, possibly computer models of flow and transport for the site, feasibility studies or prior relevant work demonstrating the feasibility of amendments, ongoing monitoring data to assess performance of the remedy, including routine sampling of contaminant concentration and signatures of the remedy and its effects.
Subsurface Data 5
PNNL-33647 Data Source Dataset Description(s)
Geophysical Data Borehole, cross-hole, surface, or remote sensing collection of data through electrical techniques (e.g., electrical resistivity tomography, induced polarization), electromagnetic methods (e.g., frequency and time domain electromagnetic induction, magnetotellurics, ground penetrating radar),
seismic methods (e.g., reflection seismology, seismic refraction, seismic tomography), gravity techniques (e.g., gravimetry and gravity gradiometry),
magnetic techniques (e.g., magnetometers), thermal methods (e.g.,
infrared, fiber-optic distributed temperature sensing) or multispectral/
hyperspectral methods.
Groundwater model Deterministic or stochastic subsurface numerical models of flow and transport in the vadose zone, saturated zone, or a combination, including input files, model calibration results, and predictive results. Geo-framework model describing the hydrogeology and forming the basis for a conceptual site model.
Authorized limit Authorized limit(s) based on DOE Order 458.1 (DOE 2011, 2017) or data data required to translate regulatory limits to authorized limit(s),
including hydrologic parameters (e.g., soil density, precipitation, irrigation) human health based on pre-described risk approach, and other default parameters in the RESRAD computer code.
DOE = U.S. Department of Energy; HSA = historical site assessment; NEPA = National Environmental Policy Act of 1969; RCRA/CERCLA = Resource Conservation and Recovery Act/Comprehensive Environmental Response, Compensation, and Liability Act; RESRAD = Suite of RESidual RADioactivity software codes; RSSI = radiological survey and site investigation.
The exercise of combining 3-D data from disparate sources and with varying formats, quality, and resolutions is expected to be a major undertaking and formidable challenge for any user.
However, much of this process would have been performed as part of the HSA and other RSSI phases prior to the compliance survey (see Figure 2). Additional DQA activities will be required that determine the suitability of data from prior phases for the compliance survey. These are discussed in the next section.
3.2 Data Quality Assessment The DQA process is typically applied after planning and implementation have concluded to determine whether the quality of collected data is satisfactory for subsequent analysis and decision-making purposes (EPA 2000). Given the potentially heavy reliance on previously collected data, DQA should be applied to such data during the compliance survey planning stage. The current guidance is shown in Figure 3, in which additional recommendations for the subsurface compliance phase are shown in green.
Information about previous survey and modeling efforts should be reviewed to determine its suitability for the compliance survey sampling plan, and to answer this overarching question:
Should previously collected data be used to inform compliance survey planning, sampling design, and subsequent analysis?
Answering the following questions will help address this question. Details should be included in the compliance survey data quality objective plan.
- For what purposes were data collected in previous RSSI phases (e.g., geophysical modeling, statistical modeling/kriging, hypothesis testing)?
Subsurface Data 6
PNNL-33647
- Did the key assumptions required for such purposes hold (if an assumption of independent errors was required for a hypothesis test, was such independence demonstrated)?
- What findings were made based on previously collected data? What decisions were made based on these findings?
- Was a sampling design used to collect the previously collected data? If so, what was the design?
- Were the previously collected data representative of the state of the site at that time? Are they representative of the state of the site for compliance survey purposes?
- What is the quality of the previously collected data (e.g., are there explainable/unexplainable outliers, influential points, missing data)?
Figure 3. Quality Assurance Assessment (figure from EPA [2000] with additional recommendations for subsurface compliance survey planning in green)
VSP data quality capabilities for quality assurance include outlier detection and tests of distributional assumptions, as well as retrospective power curves for the Sign and Wilcoxon Rank Sum tests. The following tools would be useful for DQA of subsurface data collected during previous RSSI phases in VSP:
- 3-D subsurface data
- Activating the z coordinate for borehole and groundwater wells to datasets, data management user interface, and visualization (maps, diagrams), activated as of July 1, 2022.
- Adding the capability to track subsurface data within a single borehole/well using (x, y, z, t) coordinates, that is tracking observations at distinct depths (z-locations) and timestamps (t) within one borehole/well located at a single (x, y) coordinate.
Subsurface Data 7
PNNL-33647
- Improving the capability to generate, filter, and track borehole/well labels.
- Having the ability to identify or import additional vertical dimension data, including soil and rock strata boundaries (e.g., hydrostratigraphic layers, lithologic layers), aquifer locations and dimensions, and layers associated with dose models for subsurface volumes.
- Data processing prior to analysis generally occurs outside of VSP via multiple tools such as spreadsheet software and databases that are flexible, available, and widely used. VSP users would benefit from additional limited data processing capabilities including:
- Unit conversion
- Ability to combine data and read metadata from disparate sources (including data files and separate metadata files) and sensor platforms and to provide details (e.g.,
instrument label, field of view, sample matrix) into a single dataset
- Ability to generate metadata for files it creates or edits including, for example, name of original dataset, coordinate system information, raster file resolution, etc.
- Ability to export 3-D data in a format that can be used by other software
- Ability to capture uncertainty and/or DQA corresponding to datasets or individual data points from disparate sources.
- Data visualization capabilities that will make it easier for practitioners to identify potential outliers, understand how data from various sources align (or not), and get a comprehensive picture of various sources of uncertainty, including:
- Plots that make potential elevated volumes or 3-D hot spots easier to identify on maps
- Capability to identify different instruments, fields of view, etc., using different symbols or color scales on site maps
- Capability to identify different analytes and/or sample matrices (e.g., groundwater, surface water, and soil types) using different symbols or color scales on site maps
- Capability to denote uncertainty and/or DQA metrics for each dataset or data point using different symbols or color scales on site maps
- Capability to automatically determine optimal symbol/color combinations given the various characteristics in a combined dataset.
These tools will also enhance the practitioners capability to perform DQA after compliance survey data collection and support subsequent subsurface analysis and decision-making processes, as well as visualization of compliance survey findings.
3.3 VSP Capabilities and Advancements VSP currently has the capability to support visualization of 3-D site volumes and data points within them. Existing capabilities and areas where further development would be necessary to fully support subsurface analysis are described in this section.
Subsurface Data 8
PNNL-33647 VSP has a wide range of capabilities for visualizing 3-D interior spaces of buildings and structures. It allows users to import schematics or floorplans, which are then converted to 3-D rotatable images. Users can map irregularly shaped rooms and room contents using object libraries or create their own. Often surface samples are collected from interior structures using swipes that cover some area (e.g., 10 cm by 10 cm square recommended for structures in MARSSIM). Sample locations and areas represented can be visualized on interior maps.
Much of this mapping capability could be extended to the subsurface for cases when such information is available. User-defined 3-D coordinates could be used to define a subsurface volume. Layer boundaries, location and depth of groundwater wells, and subsurface sample locations could be visualized in rotatable subsurface volumes. Subsurface structures, such as building footprints, piping systems, and tanks, could be mapped along with data. The ability to subset data interactively for geospatial analysis would also be useful, for example, by selecting data only at specific depths or within particular time periods. Enhanced subsurface mapping capabilities would enable different analysis methods described in this section.
3.3.1 3-D Sample Surfaces and Volumes Basic 3-D sample areas or sample volumes are currently represented in VSP using 2-D polygon contours and associated depths, where 2-D sample areas can be extended downward into the subsurface and specified depth(s) add a third dimension. Figure 4 shows an example of this with two volumes from above (left panel) and from the side (right panel).
Figure 4. Basic 3-D Volumes in VSP. A polygon and a circle are shown from above in the plane view (left) and from the side in the 3-D volume view (right).
After 3-D volumes are created in VSP, sample points can be randomly generated according to different techniques within the volumes, as shown in Figure 5: placing points randomly on different planes within the 3-D volume (left) or randomly throughout the volume (right). Historical data or coordinates of previously sampled locations can also be imported into VSP.
Subsurface Data 9
PNNL-33647 Figure 5. Current 3-D Sample Placement Methods in VSP More complex 3-D volumes can be created using elevation contours, allowing users to create surfaces and volumes with irregular sides or complex shapes. Currently, elevation contours can be drawn in VSP or imported from a map file. Figure 6 provides an example where a subsurface volume was created using two 2-D rectangular contours at depths 5 m and 10 m. The volume is shown from above (left) and in 3-D (right). Similar to VSPs room sampling capabilities, VSP can place samples on accessible surfaces within the 3-D volume (e.g., when surfaces of an excavated pit can be accessed) or within the 3-D volume (providing points to sample within the subsurface, despite the lack of accessible surface[s]).
Figure 6. Example Volume Created with Rectangular Contours and Associated Depths.
Contours are shown from above (left) and within the 3-D volume after applying the depths of each (right).
Another method of creating irregular 3-D objects in VSP is to import digital elevation model (DEM) file(s) containing subsurface depths (Dong et al. 2015). Samples can currently be placed or visualized along the irregular surface(s) of imported 3-D objects, although cannot be placed at random throughout the volume. Figure 7 shows a VSP sample area created using this method.
Subsurface Data 10
PNNL-33647 Figure 7. Example 3-D Volume Surface Created by Importing a Digital Elevation Model File Updating VSP to allow for placing sample points throughout the volume will provide value for subsurface decommissioning and the compliance phase of the process. Additionally, multiple DEM files which could represent multiple survey areas/volumes can currently be imported and visualized in VSP, either simultaneously or individually by clicking each on/off in VSPs layer control window. This might be important to allow visualization of DEM surfaces of surveyed areas during soil excavation including identification of any areas of overlap between surveys as remediation proceeds.
3.3.2 Data Visualization After data have been collected at sampled locations, contamination metric values can be plotted in VSP. Currently, the color-by-value capability allows a user to color data points based on contamination metric values (or other user-defined parameters). The depth of each point in the subsurface is represented on the map by a vertical line extending from the data location up to the surface. Figure 8 shows an example of this visualization.
In addition to visualizing individual data points, subsurface analysis could leverage geostatistical analysis (e.g., kriging, generalized linear models) to produce volumes containing interpolated data between locations where data were collected. VSP currently supports only 2-D kriging and visualizing the resulting 2-D surfaces, visualized via 2-D rasters, as shown in Figure 9.
However, we expect 3-D kriging to be a priority enhancement for subsurface surveys in the future. Thus, the current visualization capability should be extended to support 3-D visualization of interpolated data within separate layers or the whole 3-D volume. Further details about the interpolation methods and a layered versus volume approach are presented in Sections 4.2 and 4.3.
Only minor modifications to VSP would be required to allow visualization of multiple individual 2-D layers within a 3-D volume. Current capabilities would support zooming in and out, as well as rotation of individual layers. Enabling similar 3-D visualization of the whole volume would require more substantial modifications of VSP once the volume visualization has been implemented, although currently VSP does support zooming in and out as well as rotation.
Subsurface Data 11
PNNL-33647 Figure 8. Subsurface Soil Point Data Visualization in VSP. Locations with an iodine-129 value between 0 and 4.5 are color-coded according to the legend, while locations without a valid data value (i.e., less than 0) are displayed as dark blue plus symbols. Black lines are drawn from the depth of the points to the surface to illustrate depth.
Subsurface Data 12
PNNL-33647 Figure 9. Currently Only 2-D Visualization of Interpolated Surface is Supported in VSP Additional visualization and interactivity features should be added to VSP to enhance user capabilities to create and manipulate data either by manipulating existing survey area/volume data or by drawing new survey areas/volumes on the 3-D site map. Subsurface visualization should include the following;
- Ability to drape 2-D layers on a 3-D surface
- Ability to display 3-D vector layers in 3-D space (e.g., show well screen interval color-coded by value of parameter in 3-D space)
Subsurface Data 13
PNNL-33647
- Ability to apply vertical exaggeration to emphasize vertical features that are too small to see relative to the horizontal scale using a factor to specify how much greater the vertical scale is to appear than the horizontal scale (e.g., 5x, 10x, 20x)
- Ability to slice 3-D volumes vertically and horizontally, either by drawing a slice or providing 2-D coordinate inputs) and visualize data on the sliced plane including display of iso-contours
- Ability to show iso-volumes according to selected parameter values (e.g., volume greater than 5 pCi/g)
- Ability to add subsurface infrastructure and boundaries (e.g., underground building infrastructure or geologic characteristics)
- Ability to set different transparency levels for each subsurface layer or volume to aid in seeing other features in the model
- Updating VSP to enable placement of scan lanes along a DEM surface to facilitate scan survey for subsurface surfaces during remediation 3.3.3 Analysis of Multiple Contaminants/Analytes Current VSP analysis datasets include the coordinates (in local or commonly used coordinate systems), values, and units of one or more variables of interest (e.g., pCi/g or counts per minute
[cpm]). Multiple analytes are accepted, but the statistical tests and the number of sample calculations are done univariately (i.e., one analyte at a time). Adding co-kriging and multivariate analysis functionality in VSP would allow multiple analytes to be analyzed simultaneously, such that the analyses could capitalize on correlations between analytes. This would be useful for subsurface analysis, particularly in cases where data for one analyte are more readily available than another (due to time, cost, or detectability), but their behavior, presence, and/or concentrations in the subsurface are correlated. Co-kriging and multivariate analysis could make predictions/inference for the harder-to-measure analyte.
3.3.4 Statistical Diagnostics Spatial heterogeneity generally refers to the uneven distribution of a characteristic in space, such as contaminant distribution or the uniformity of a rock layer in a subsurface volume. Many geospatial and geostatistical methods require diagnostics to determine spatial heterogeneity, if correlations between metrics exist, and whether spatial variability is isotropic or anisotropic.
Isotropy is the condition that, for a single location in space, the correlation between it and neighboring locations is dependent only on the distance between the two. Anisotropy is the more complex condition where the correlation between points is a function of both the distance and direction between them.
An additional component of spatial heterogeneity is nonstationarity, of which there are two types. First-order nonstationarity occurs when the mean concentration is non-constant throughout the subsurface volume. As an example of first-order nonstationarity, consider a soluble contaminant deposited via a subsurface leak into an aquifer that results in a different mean contamination inside the aquifer compared to outside the aquifer. Second-order nonstationarity occurs when the covariance structure changes within the volume, as in the case of anisotropy. Continuing the previous example, second-order nonstationarity could occur when the aquifer crosses two distinct subsurface layers that govern contaminant fate and transport at two distinct rates or orientations in one layer compared to the other.
Subsurface Data 14
PNNL-33647 The variogram (or semivariogram) is a statistic that estimates the change in correlation between observations as a function of distance. In the case of isotropy, the variogram has the same shape for locations in every direction (angle). In the case of anisotropy, the variogram shape differs as a function of the angle within the volume.
Numerous capabilities described by the Electric Power Research Institute (EPRI) (EPRI 2016) to examine spatial heterogeneity would be valuable to incorporate in VSP for the purpose of performing such statistical diagnostics. These capabilities are summarized as follows:
- Subsurface heterogeneity assessment
- Subsurface heterogeneity introduces challenges to structural analysis and variogram estimation, as does the measurement uncertainty of subsurface data collection technologies. Subsurface heterogeneity may be induced by the presence of discontinuities (e.g., geological stratification, aquifers, engineered barriers [interfaces between natural sediment and backfill around structures during construction, building footings and zone of compression beneath structures], and fracture planes), which lead to the spatial structure changing across a site (EPRI 2016).
- Measurement techniques with different fidelities and extents of representativeness will affect compliance surveys to the extent that each is included in the planning and analysis process. Measurement methods used during early phases of decommissioning may differ from those used during the compliance survey (EPRI 2016).
- Including capabilities described in Section 3.2 for DQA will help practitioners spot disparities in variability and uncertainty between layers and measurement techniques.
- Structural (variogram) analysis along planes in several different directions will help assess whether anisotropy is present in the subsurface and should be considered. Currently, the VSP kriging module includes a variogram plotting capability. It is limited to assessing variograms in any direction but only including one variogram in a single direction in the analysis (which can also be included in the automatically generated report). The importance and implications of numerous distinct variograms are discussed in EPRI (2016),
MacCormack et al. (2017), and other work cited in Section 4.0. Recommendations for including multidirectional variograms in statistical (kriging) models to address anisotropic conditions are also included in Section 4.0.
- Updates should be considered to view variograms in multiple directions simultaneously and display all those selected for the final kriging model.
- Variogram surfaces (Figure 10) also can be used to help determine if anisotropy is present in 2-D and extended to visualize several variogram surfaces on a series of planes within the subsurface volume (EPRI 2016). Adding these capabilities to VSP would help users determine whether anisotropy is present in a subsurface volume and in which (and how many) directions.
- Capabilities to view variograms in multiple directions and view variogram surfaces are available in other statistical software (e.g., gstat and ggplot2 packages in R). While linking the variogram fits directly to the kriging models in VSP may provide a convenience, consideration should be given to how and if the costs of adding such capabilities to VSP outweigh the benefits.
- If such capabilities are added to VSP, a corresponding capability should be added to the automatically generated report to capture the diagnostics that led to subsequent sampling designs and geostatistical model selection.
Subsurface Data 15
PNNL-33647 Figure 10. Variogram Surfaces Indicating Isotropic Variation (left) and Anisotropic Variation (right)
- In addition to the methods described above, scatterplots between multiple variables and/or bivariate statistics (e.g., correlations) can help determine whether methods such as cokriging or multivariate analysis should be considered. For example, examining the relationship between gamma and alpha radiation within a subsurface volume and/or correlations between these measurements and geophysical properties can lead to selecting a co-kriging model or developing a sampling plan that capitalizes on the ease of collecting gamma measurements (EPRI 2016). Scatterplot matrices can provide powerful visual representations of such relationships by displaying scatterplots of multiply bivariate relationships in one figure. However, these capabilities are available in other statistical software (e.g., ggpairs function of the GGally package in R), thus may not need to be added to VSP (i.e., would be a lower priority item).
- H-scattergrams are an additional tool helpful in assessing spatial scale. H-scattergrams are scatterplots with the value of a variable at one location plotted against the value of the same variable at a different location, where each pair of observations is separated by distance h.
Multiple scattergrams at different h-distances can be visualized in a matrix to help determine the spatial scale of the data, where the magnitude of the correlations indicates the magnitude of the spatial scale. Figure 11 shows three h-scattergrams at different h-distances for the same metric for h-distances 1 m, 4 m, and 30 m. The larger correlations when h-distances are 1 m and 4 m and the smaller correlation when h-distance is 30 m indicate the spatial scale is larger than 1 to 4 m but smaller than 30 m (EPRI 2016).
However, these capabilities are available in other statistical software (e.g., the hscat function in the gstat package in R) and thus may not need to be added to VSP (i.e., would be a lower priority item).
Subsurface Data 16
PNNL-33647 Figure 11. Example of Three H-Scattergrams at Three H-Distances (EPRI 2016).
3.3.5 VSP in the Cloud PNNL was sponsored by the U.S. Department of Homeland Securitys Federal Emergency Management Agencys Nuclear Incident Response Team to investigate optimal use of VSP sample prioritization and geostatistical capabilities, deployed to the DOE secure Microsoft Azure cloud environment to improve access and potential future automation of some VSP functionality, as well as the interaction of VSP with other tools in that environment. PNNL identified several priorities to guide that effort. The following are relevant to subsurface compliance survey planning and geostatistical analysis:
- Provide links between VSP and other web applications
- Provide access to relevant databases also hosted in the cloud environment (e.g., those detailed below in Section 3.4), including import capabilities to load data into VSP from other web applications
- Export sample design maps into the cloud environment for use in other web applications
- Minimize/remove the requirement to install the VSP desktop application
- Provide a simplified subsurface VSP web application that requires limited user training
- Enhance computing capabilities for some functions with heavy computation costs using batch or parallel processing capabilities available in the cloud.
Subsurface Data 17
PNNL-33647 Developing a web application based on the existing VSP desktop application would require substantial web architecture development. Rather than recommending such development for VSP in its entirety, PNNL recommends a phased approach during which it is developed for select modules, identified based on their relevance to subsurface modules. PNNL recommends implementing selected modules in the desktop application using modern code and then docking them to the web application, so they remain accessible to users on both platforms.
3.4 Additional Existing Subsurface Tools and Capabilities at PNNL PNNL presented additional existing tools and capabilities at a meeting with the NRC on June 14, 2022. Subsequently, the NRC requested additional details be included in this report for consideration in future subsurface activities. Not all tools and capabilities will be directly relevant to compliance surveys but could be used in previous phases and therefore provide context to the compliance phase. This section provides such details and notes where each could be applicable in the compliance phase to show no contamination is left behind and/or where each tool or capability has the potential to integrate with VSP to provide enhanced data acquisition, geospatial analytics, or data visualization required to show compliance.
3.4.1 Data Visualization and Analysis Access to up-to-date data, data visualization tools, and applications of data analysis methods will be critical to support planning, monitoring, remediation decisions, and communication during the RSSI compliance phase. This section provides an overview of such capabilities and notes where future efforts could link these tools and capabilities with VSP to enhance subsurface sample placement and statistical analysis.
Data from the existing tools described in this section could be incorporated into geospatial modeling and survey sample decisionseither explicitly through independent model variables or implicitly by informing the conceptual site model (CSM) with subsurface layer delineation and/or anisotropy or covariance structure development. Interconnection could allow data from these tools to be imported into VSP, or vice versa and incorporating collected data or sample locations into the visualization interface. This could be accomplished by enabling VSP and these tools each with import/export mechanisms to align expected data formats and/or by creating connections between the software platforms.
3.4.1.1 SOCRATES The Suite of Comprehensive Rapid Analysis Tools for Environmental Sites (SOCRATES) 1 is a web application (see Figure 12) that provides access to environmental data, data visualization, and rapid data-driven analytics to make sense of site environmental cleanup data and support cleanup remedy decisions, optimization, and exit strategies. SOCRATES includes modules for analyzing site geo-framework/conceptual site models, groundwater levels and flow direction, contaminant concentrations, water quality parameters, plume dynamics, remediation systems, waste site risk, and remote-sensing data interpretation. SOCRATES offers an improved approach to working with disparate data types from multiple authoritative sources, with advantages including wide availability via a web browser, access to current data, availability of standard analysis methods for consistent application regardless of the analyst doing the work, and built-in functionality and logic to make analyses quick and easy. It also facilitates communication with technical staff, site managers, regulators, and stakeholders.
1 https://www.pnnl.gov/projects/socrates Subsurface Data 18
PNNL-33647 Figure 12. SOCRATES Component Applications SOCRATES currently includes data for the DOE Hanford Site and is being used to understand cleanup risks and priorities, operation and optimization of remediation systems, monitoring and sampling data and data gaps, cleanup endpoints, alerts for subsidence issues, and groundwater surface water interactions. Access is managed, allowing certain information to be publicly available and other information to be available only to site staff for use in preparing documents for review and release. SOCRATES could be expanded to provide tools to any site for decommissioning purposes in the future by establishing connections with project-specific databases. In the future, it also may be advantageous to enable interconnections between VSP and SOCRATES. For example, rather than recreating data visualization capabilities in VSP that already exist in SOCRATES, interconnection could allow VSP data to be visualized using SOCRATES modules. Further, VSP sample plan development functionality dovetails nicely with the information and desired use of SOCRATES to manage sampling.
The SOCRATES component applications provide the following capabilities:
- HYPATIA (HYdraulic Pump-And-Treat Information Analytics) provides access to and tools for analyzing pump-and-treat data, including extraction well and treatment plant chemistry and control system sensors that are often stored in separate databases.
HYPATIA also provides methods for statistical analysis of time-series data and includes calculated metrics, such as mass flow rate and injectivity for assessing pump-and-treat system performance.
Subsurface Data 19
PNNL-33647
- ARIUS (Advanced Remote-sensing Image USer interface) performs full end-to-end automated acquisition of remote-sensing datasets using a matrix of cloud-based technologies. It streamlines complex workflows and processes satellite data pertaining to variation in both ground surface elevation and temperatures. Light detection and ranging and SENTINEL-1 data are used to monitor surface displacement, whereas Landsat 8 and 9 satellites are used to monitor surface temperatures as an indirect measure of groundwater and surface water interactions.
- ARISTOTLE (Adaptive Risk-Informed System to Obtain The Likely End-state) provides geospatial visualization of waste site inventory data from multiple data sources to help assess risk and prioritize remedy actions across spatial scales. This includes estimates of chemical and radiological inventories from historical discharges important for identifying risks of residual inventories of contaminants still present in the subsurface and qualitatively and quantitatively assessing site cleanup priorities.
- ORIGEN (Online Retrieval Interface for GEologic INformation) provides access to site geology and features that can be used to create cross sections and access information about stratigraphic thicknesses, water table elevation, and well construction information.
- CRATES (Charting, Reporting And TEmporal visualizationS) is a plotting tool primarily used for visualizing groundwater concentration data and enabling users to identify well locations for plotting user-specified contaminants of concern. It includes graphics and tabulated data export features.
- GALEN (Groundwater AnaLytics for the ENvironment) provides access to multiple sources of water-level data through a single access portal and includes tools for visualization and analysis of groundwater level and flow direction over time to support site characterization and enhance remediation system design and performance monitoring.
- PLATO (PLume Analysis TOol) analyzes groundwater data to assess contaminant plume behavior, implementing data-driven, quantitative analyses based on standard statistical methods and published guidance from the U.S. Environmental Protection Agency (EPA) and the U.S. Geological Survey.
3.4.1.2 PHOENIX PNNL-Hanford Online ENvironmental Information eXchange (PHOENIX) 2 is a family of spatially enabled web applications that provide quick access to decades of valuable scientific Hanford Site data and insight through intuitive query, visualization, and analysis tools. PHOENIX provides a single access point to multiple datasets (see Figure 13) via standard web browsers, data visualization tools, and explanations of key datasets to aid understanding. By integrating previously isolated datasets and developing relevant visualization and analysis tools, PHOENIX applications enable DOE to discover new correlations hidden in legacy data to address complex issues at Hanford more effectively.
Subsurface Data 20
PNNL-33647 Figure 13. PHOENIX Overview PHOENIX could supplement VSP as a source of geospatially identified analytical chemistry data that would combine with the VSP sampling calculator to estimate the number and location of samples needed to satisfy a DQO or reduce decision uncertainty. Advanced interactions would involve augmented reality/virtual reality devices for immersive visualization/interaction with site geology (e.g., from SOCRATES), sampling locations, and analytical chemistry results.
3.4.1.3 TRAC Tracking Restoration and Closure (TRAC) 3 is a web-based application for communicating the status of and progress on groundwater cleanup at sites within the DOE Office of Environmental Management (DOE-EM) complex through a story map and associated site information and metrics. TRAC (screenshot shown in Figure 14) is intended to provide a concise narrative about groundwater cleanup information, remaining groundwater contaminant plumes, remediation technologies implemented, and waste site regulatory status. TRAC presents summary-level data across DOE-EM sites, while providing more detailed, consistent information for individual sites, down to the operable unit level. This information facilitates meaningful context for each individual site story and the collective story associated with the DOE-EM groundwater cleanup mission, thereby serving as a resource to respond to data requests and give DOE-EM leaders a comprehensive view of groundwater cleanup status across the complex. TRAC is slated for release at the end of fiscal year 2022.
Subsurface Data 21
PNNL-33647 Figure 14. TRAC Web-Interface Screenshot A version of TRAC could be implemented to provide summary-level information relevant to the NRC across multiple sites and linked to from VSP. This would provide a higher-level story map view for communicating the status of decommissioning and restoration issues.
3.4.1.4 Few-Shot Machine Learning Traditional, point-source based, destructive characterization methodologies for subsurface systems are expensive, present the potential risk for human exposure, and have long been challenged to collect sufficient data to accurately describe the complexities of the system. Data collected via borehole samples represents the state of the system at the specific place and time of data collection but may or may not be representative of the system over time or in total when sampled boreholes are not representative of the whole system. It also affords large uncertainty in forecasting the evolution of the systems. These limitations present significant uncertainties in assessing (for example) radiation dose from subsurface transport pathways.
To improve characterization while minimizing costs and exposure risk, few-shot machine learning (see Figure 15) is being advanced in conjunction with remote subsurface sensing techniques and high-performance forward prediction. This approach is being developed to reliably estimate the subsurface property distributions, including (but not limited to) permeability, porosity, and hydraulic conductivity, that control fate and transport of radioactive material, thereby addressing the scarcity of characterization data and complexities of heterogeneous subsurface systems. Further, applying machine learning approaches to a combination of discrete well or borehole datasets with modeled or surrogate data (i.e., output from flow and transport models, indirect measurements collected through electrical resistivity tomography
[ERT], spectral induced polarization [SIP], etc.), may lead to predictive models and predictions that capture relationships between multimodal measurements and variables of interest that can be used for the purpose of compliance survey design.
Subsurface Data 22
PNNL-33647 Figure 15. Few-Shot Learning Schematic By applying deep learning to integrated multimodal sensing data, the performance of machine learning approaches is being defined, through training with large and small datasets, to estimate governing system-scale subsurface parameters and their spatiotemporal evolution. These advancements will reduce the uncertainty of system-scale characterization and radiation dose assessments, minimize costs, and increase worker safety and protection of human health and the environment.
3.4.2 Data Management Environmental information and data are inherently complex because the data come in different forms and formats, and multiple entities may control the collection and management of data, resulting in variations in data quality and access. Yet access to data is vital to supporting environmental restoration decision-making. Strong expertise in data management underpins CORE (Common Operating and Response Environment) and the development of the Hanford Environmental Information and Data Index (HEIDI) for the DOE Richland Operations Office as part of their Hanford Environmental Data Management program (HEDM), which is a formal program for managing environmental data and the associated records, materials, and systems at the Hanford Site. Based on an initial evaluation of approaches (Ham and Crockett 2021),
HEIDI is being developed as a catalog of data sources that will enable long-term access and retrievability for the multiple independent sources of data that might otherwise be difficult to discover.
3.4.2.1 CORE CORE 4 software is a highly configurable system composed of broadly applicable modules for rapid setup and customization to meet a wide range of data management needs. The CORE design is derived from extensive experience in information integration, data access applications, field data collections, and operational analysis to support decision-making. The CORE data collection system uses a web-based platform to provide globally accessible data forms, fuse diverse data streams, enable context-based situation awareness, and provide programmatic and technical oversight. Using role-based access, the CORE data collection system allows diverse sets of information to be collected and controlled. The system also provides an interactive interface that streamlines communication.
4 https://www.pnnl.gov/available-technologies/common-operating-response-environment-core Subsurface Data 23
PNNL-33647 While CORE is particularly useful for test and evaluation applications, it has the flexibility for collection and management of data in a range of applications. CORE could tie into VSP work across the spectrum of planning analysis, data collection, and data management/reporting shown in Figure 16.
Figure 16. CORE Data Management Workflow 3.4.3 Flow and Transport Modeling PNNL has been developing and applying predictive models for groundwater flow and contaminant transport in hydrologic systems since the 1960s. Developments include analytic models, standard finite difference and finite element codes for fluid flow and solute transport, particle transport codes, and high-performance flow and reactive transport codes designed to run on parallelized machines. These tools are used to optimize groundwater characterization activities, predict groundwater flow velocity and direction, predict groundwater contaminant concentration and transport rates, and predict the long-term impacts and estimate risk to the public. Critical leadership has been delivered in defining the necessary level of complexity in subsurface modeling to support environmental analyses with consideration of cost-benefit tradeoffs (Freedman et al. 2017), thereby providing guidance and case studies exemplifying three primary elements: 1) modeling approach; 2) description of process; and 3) description of heterogeneity.
Although the groundwater contribution to dose assumes a single average concentration in compliance demonstrations, particularly for reactor decommissioning, data from the existing tools described in this section could be incorporated into geospatial modeling and survey sample decisionseither explicitly through independent model variables (e.g., generalized least square [GLS] regression models described in Section 4.2.1) or implicitly by informing the CSM with subsurface layer delineation and/or anisotropy or covariance structure development.
Interconnection could allow data from these tools to be imported into VSP, or vice versa and incorporating collected data or sample locations into the visualization interface. This could be accomplished by enabling VSP and these tools each with import/export mechanisms to align expected data formats and/or by creating connections between the software platforms.
Currently, users could export such data from a separate software application and merge it with Subsurface Data 24
PNNL-33647 the site data that includes contamination measurements. Users could use diagnostics like those described in Section 3.3.4 to determine if flow and transport data are correlated with the variable(s) of interest and, if so, include them in regression or kriging models (as such capabilities become available in VSP). In the future, VSP in the cloud could be developed to enable data connections more directly between VSP and other software applications.
3.4.3.1 AIM Toolbox The AIM (Aquifer Injection Modeling) Toolbox 5 software, developed by PNNL for the EPA, provides a collection of analytical solutions suitable for evaluating the potential extent of an area affected by subsurface injection operations (see screenshot in Figure 17). These operations are regulated under EPAs Underground Injection Control program and are typically related to oil/gas development, waste disposal, or subsurface mining or storage. Each of the analytical algorithms provided in the AIM Toolbox has different approaches, assumptions, and/or focuses with respect to the nature and processes in the subsurface and the nature of the injection operations. Collectively, the set of analysis algorithms provides a broader evaluation of the area that can potentially be affected by an injection operation.
Figure 17. AIM Toolbox Interface Screenshot By providing estimates of the injectate plume extent, the software supports technical aspects of planning, evaluation, and overseeing injection activities. That is, the results help assess when further regulatory controls (e.g., monitoring, reporting) may be required and, in the case of proposed disposal activities into an underground source of drinking water, the extent of the affected area that would require exemption from protection under the Safe Drinking Water Act (Safe Drinking Water Act of 1974). The AIM Toolbox software is available as a single-page web application, providing an interface for necessary inputs and visualization of results in both chart and map panes.
5 https://www.pnnl.gov/projects/aim-toolbox Subsurface Data 25
PNNL-33647 3.4.3.2 RT3D Reactive Transport in Three Dimensions (RT3D) 6 is a software package for simulating the 3-D, multispecies, reactive transport of chemical compounds (solutes) in groundwater (see Figure 18). The tool provides an easy-to-use and flexible framework that is applicable to natural attenuation, accelerated bioremediation, or other reactive transport modeling scenarios.
Predefined modules are available for common bioremediation scenarios, but users also have the flexibility to add any reaction kinetics desired/suitable to represent multiple chemical species in aqueous and adsorbed phases. RT3D focuses on reactive transport in groundwater (saturated systems) and uses the groundwater flow solution obtained from the U.S. Geological Survey MODFLOW code.
Figure 18. Visualization of Simulation Results Using RT3D 3.4.3.3 SVEET PNNL partnered with the Department of Defense and DOE-EM to develop the Soil Vapor Extraction End-state Tool, version 2 (SVEET2) 7 shown schematically in Figure 19. SVEET2 is a spreadsheet tool that allows users to easily enter site-specific information to obtain groundwater and soil gas concentration estimates resulting from a vadose zone contaminant source(s).
SVEET2 is used to quantify the impacts of a vadose zone source as part of a structured process to evaluate whether a soil vapor extraction system should continue operating, be optimized, transition to another remedy, or be terminated for site closure. SVEET2 is the result of a recent update to expand the range of permissible input parameter values, thereby enhancing the applicability of the tool for Department of Defense, DOE, and other sites. SVEET2 is not itself a numerical simulator, but rather it uses pre-modeled simulations for a generalized conceptual framework to estimate results for site-specific conditions. SVEET2 is based on 5,760 Subsurface Transport Over Multiple Phases Simulator (STOMP) simulations for combinations of parameters with a nonlinear impact on the results and scaling for parameters with a linear impact on results.
https://www.pnnl.gov/projects/remediation-performance-assessment/soil-vapor-extraction Subsurface Data 26
PNNL-33647 Figure 19. SVEET2 Conceptual Model Schematic 3.4.3.4 STOMP and eSTOMP PNNLs STOMP 8 is a suite of numerical simulators for solving problems involving coupled flow and transport processes in the subsurface. As shown in Figure 20, STOMP is based on mathematical equations that describe our understanding of hydrologic, thermal, thermodynamic, geochemical, and geomechanical processes.
Figure 20. Visualizations of Simulation Results for Aqueous Phase (left) and Nonaqueous Phase (right) Transport Using the STOMP Suite of Codes 8
https://www.pnnl.gov/projects/stomp Subsurface Data 27
PNNL-33647 STOMP has been applied to domains including environmental remediation and stewardship, geothermal resources, production of natural gas hydrates, subsurface permanent storage of carbon dioxide, and oil and gas recovery using conventional and unconventional technologies.
Some of the unique applications and features of the simulator for environmental work are its ability to model vegetated surface barriers and model flow and transport in deep vadose zones under thermally altered states. Exascale STOMP (eSTOMP) 9 is a highly scalable (parallel) version of STOMP.
3.4.3.5 PFLOTRAN PFLOTRAN 10 leverages massively parallel, high-performance computing to simulate large-scale nonisothermal multiphase flow, multicomponent reactive transport, and ERT problems in the subsurface environment (see Figure 21). The code is designed to predict the future state of environmental systems and better inform stakeholders in the regulatory decision-making process. Major contributions to PFLOTRAN development continue to be funded by DOE. The code is the primary physics simulator within the DOE Office of Nuclear Energy Spent Fuel and Waste Science Programs Geological Disposal Safety Assessment Framework.
Figure 21. Visualization of PFLOTRAN Simulation Results for Uranium Transport at the Hanford Site, Adjacent to the Columbia River in Eastern Washington 3.4.4 Computational and Applied Geophysical and Geomechanics High-resolution, real-time computational and engineering geophysics underpin the development and implementation of state-of-the-art geophysical tools and methods to improve the characterization of subsurface properties and processes. 11 Geophysical data can be collected from many different platforms (such as at the ground surface, between wellbores, and within wellbores) to interrogate subsurface variability over a variety of spatial scales and resolutions.
9 https://www.pnnl.gov/estomp 10 https://www.pflotran.org/
11 https://www.pnnl.gov/computational-and-applied-geophysical-and-geomechanics-laboratory Subsurface Data 28
PNNL-33647 The use of geophysical methods combined with conventional measurements 12,13for example, seismic source zone characterization via computational methods like moment tensor inversion and magnitude estimationprovides spatially extensive information about the subsurface in a minimally invasive manner at a comparatively high resolution. These capabilities include high-resolution geologic field mapping to develop improved local and regional fault maps that better inform geologic site characteristics. Such applied geophysical and geomechanics science capabilities can support subsurface characterization and monitoring for siting and decommissioning.
Data from existing tools described in this section could be incorporated into geospatial modeling and survey sample decisionseither explicitly through independent model variables (e.g., GLS regression models described in Section 4.2.1) or implicitly by informing the CSM with subsurface layer delineation and/or anisotropy or covariance structure development.
Interconnection could allow data from these tools to be imported into VSP, or vice versa and incorporating collected data or sample locations into the visualization interface. This could be accomplished by enabling VSP and these tools each with import/export mechanisms to align expected data formats and/or by creating connections between the software platforms.
Currently, users could export such data from a separate software application and merge it with the site data that includes contamination measurements. Users could use diagnostics like those described in Section 3.3.4 to determine if flow and transport data are correlated with the variable(s) of interest and, if so, include them in regression or kriging models (as such capabilities become available in VSP). In the future, VSP in the cloud could be developed to enable data connections more directly between VSP and other software applications.
3.4.4.1 ERT ERT provides 2-D, 3-D, and 3-D + time images of subsurface electrical structure, representative of geologic and fluid properties. Applications of ERT have demonstrated its ability to resolve subsurface geologic frameworks and detect time-lapse changes associated with natural and engineered processes including solute transport and water movement (Johnson et al., in press, 2021, 2012). ERT has been used for high-resolution site characterization, plume monitoring, and detection of leaks from subsurface infrastructure including tanks and pipes, as illustrated in Figure 22, where ERT revealed leakage associated with an unanticipated vacuum pipe overflow in the subsurface. It could also be used to discover and map subsurface infrastructure.
Advanced capabilities for near-real-time delivery of ERT results to operators in the field provide actionable information to site operators (Johnson et al. 2015).
- Collected ERT data could be used to inform sampling plans and statistical models by confirming previously identified leak areas, pointing to potential leak areas, and highlighting where underground infrastructure should be considered when sampling (e.g., structures or piping prohibit bore hole sampling). These data could be used to put constraints on the extent of contamination predicted by incorporating the coordinates into a statistical model with an indicator variable (e.g., forcing predicted concentration levels to be zero where underground structures are detected). If relationships between the measured phenomena and contamination are understood prior to statistical model fitting, the VSP user could specify contamination likelihoods as part of the Bayesian Elipgrid procedure (if included in future VSP capabilities). ERT data could be included in statistical models to improve the 12 https://www.pnnl.gov/applied-subsurface-science-and-characterization-laboratory 13 https://www.youtube.com/watch?v=7Zr94UkpAUI Subsurface Data 29
PNNL-33647 accuracy of model predictions (kriging) and uncertainty estimates. Since sampling plans will consider expected contamination levels and uncertainty, potentially aiming to reduce uncertainty through additional data collection, such information could prove valuable if it reduces uncertainty associated with some areas of the subsurface volume.
- ERT data could be entered directly (through future VSP capability) to delineate subsurface infrastructure and boundaries (e.g., underground building infrastructure or geologic characteristics)either through drawing/specifying the coordinates or through direct data import/entry and coding as impermeable or solid for the purposes of informing sampling plans (e.g., label an underground pipe as impermeable such that any plan to drill a borehole in/near it would need to consider and negotiate its presence). See Section 3.3.2 for discussion of VSP capability recommendations.
Figure 22. Visualization of Time-Lapse ERT Images Revealing Changes in Electrical Conductivity Associated with Unanticipated Leaks During Overflow Events in Subsurface Infrastructure 3.4.4.2 E4D E4D 14 is a scalable 3-D geophysical modeling and inversion software code designed for subsurface imaging and monitoring. E4D uses geophysical measurements to generate images that serve as indicators of subsurface hydrostratigraphic boundaries, contaminant distributions, injection delivery extent, and dynamic biogeochemical processes. The E4D software15 provides customization for a wide range of field and laboratory applications, using static and time-lapse ERT, SIP, or seismic or radar travel-time tomography. E4D can be coupled with flow and transport models (e.g., PFLOTRAN, eSTOMP) for cost-effective, pre-field planning and 14 https://www.pnnl.gov/projects/e4d 15 https://www.youtube.com/watch?v=vrcb2N_TZ8Q Subsurface Data 30
PNNL-33647 physically based field-scale feasibility assessments. E4D leverages distributed-memory high-performance computing, providing solutions for analyzing large geophysical datasets that involve high computational demands. E4D is highly scalable and can be used on systems ranging from laptop computers to state-of-the-art supercomputing facilities. Web services can be integrated with E4D for real-time processing of electrical geophysical data to deliver near-real-time results in the form of interactive, online visualizations. The example in Figure 23 shows a visualization where E4D was used to monitor amendment release and its subsequent infiltration.
Figure 23. Visualization of E4D Inversion Results Showing Time-Lapse Changes in Subsurface Electrical Conductivity Associated with Amendment Release at the Land Surface
- Collected E4D data could be used to inform sampling plans and statistical models though quantified subsurface contaminant distributions and/or hydrostratigraphic boundaries. Such measures should be considered in survey planning to reduce uncertainty associated with contaminant locations in the case of the former and to put constraints on where contamination is predicted in the latter (e.g., forcing predicted concentration levels to be zero at points past hydrostratigraphic boundaries). Contaminant distributions from E4D could be used to inform VSP user specifications of contamination likelihoods in a Bayesian Elipgrid procedure (possible future VSP capability). Such E4D data could also be included in statistical models as an independent variable to predict directly measured contaminant and quantify the uncertainty in those predictions. Reduced uncertainty could result in fewer new sampling locations being required to reach acceptable levels of uncertainty. E4D data could be used as input to future VSP capability to delineate subsurface contours of likely contamination or boundary coordinateseither through drawing/specifying the coordinates or through direct data import/entry and coding.
3.4.4.3 SIP Like ERT, SIP is an electrical method, but whereas ERT represents a single (zero) frequency SIP represents subsurface electrical properties across a spectrum of frequencies. Compared to ERT, which is primarily sensitive to pore fluids, SIP provides additional information about soil and rock materials and important phenomena that occur at the mineral/fluid interface. SIP has been used to monitor subsurface biological activity and its effects on soils (as shown in Figure 24), mineralization and precipitation processes, and sequestration of contaminants (Johnson et Subsurface Data 31
PNNL-33647 al. 2022; Johnson et al. 2015). SIP has also been used to estimate permeability and reactive surface area (Robinson et al. 2018; Revil et al. 2015). Advanced capabilities for SIP measurement, simulation, and data analysis have been developed and are being used to quantitatively monitor groundwater remediation and understand the role of immobile porosity in prolonging cleanup times.
- Similar to the data types above, SIP quantification of contaminant sequestration, permeability, and reactive surface areas may be correlated with subsurface contaminant distributions and therefore has the potential to improve statistical prediction (kriging) models (e.g., through independent variables in a regression model). To the extent correlations between these metrics and subsurface contamination are understood prior to statistical model fitting, they could be used by VSP users to specify contamination likelihoods through the Bayesian Elipgrid procedure (if included in future VSP capabilities). Such data could also be incorporated as input to future VSP capability to delineate subsurface contours of likely contamination or boundary coordinateseither through drawing/specifying the coordinates or through direct data import/entry and coding, as recommended in Section 3.3.2.
Figure 24. Visualizations of SIP Imaging Mesh (left), Imaging Results (middle), and Experimental Apparatus (right) to Understand Electrical Signatures of Biological Activity 3.4.4.4 Seismic Subsurface Technologies Seismic methods are used to characterize subsurface geologic structure (as shown in Figure 25, which shows seismic images of the Hanford Site). Seismic methods are valuable in developing stratigraphic models and hydrogeologic framework models, and they are extensively used in estimation of geotechnical properties that are important for assessment of earthquake risk and the design of buildings in areas prone to seismic hazards. Ongoing development of advanced capabilities in distributed acoustic sensing and automated seismic sources are expected to enable high-resolution seismic monitoring of subsurface processes including water table fluctuations and saturation changes in moisture conditions (Fokker et al. 2021; Linneman et al. 2021). Compared to other geophysical methods such as electromagnetic and ERT, seismic methods are well suited to investigations over a broad range of scales, including deeper subsurface investigations (10-100 m depth).
Subsurface Data 32
PNNL-33647 Figure 25. Seismic Images from the Gable Gap Area of the Hanford Site. Images include a compressional-wave velocity tomogram overlaid on migrated reflection data and images of shear-wave velocity derived from multichannel analysis of surface waves (St. Clair et al., In Review).
Similar to the previous data types, seismic and acoustic data that quantify stratigraphical and hydrogeological characteristics in the subsurface volume may be correlated with subsurface contaminant distributions. Including such metrics in a statistical prediction (kriging) model (i.e.,
as independent variables) would improve the accuracy and uncertainty of predictions to the extent they are indeed correlated. If such correlations are understood prior to statistical model fitting, VSP users could use that information to specify contamination likelihoods in a Bayesian Elipgrid procedure (if included in future VSP capabilities). Such data could be imported or entered to delineate subsurface contours of likely contamination or boundary coordinates, either through drawing/specifying the coordinates or through direct data import/entry and coding.
3.4.4.5 Electromagnetic Subsurface Techniques Electromagnetic techniques are similar to ERT in that they image subsurface electrical structure, but where ERT relies on conduction (Ohms Law), electromagnetic techniques rely on induction (Faradays Law); consequently, electromagnetic techniques do not require direct coupling of instrumentation to the ground surface and are therefore amenable to more rapid data collection from vehicles, boats, drones, and helicopters/planes. Mobile data collection allows for cost-effective coverage of long transects or large areas. Airborne electromagnetic techniques have been used to image subsurface structure at the Hanford Site (Figure 26; Jaysaval et al. 2021) and as have ground and waterborne electromagnetic techniques for diverse characterization problems (CHPRC, 2010; Mangel et al. 2022). These techniques are also used for location of underground targets including pipes and tanks. PNNL has expertise and has developed advanced capabilities for electromagnetic modeling and inversion (Jaysaval et al. 2022, 2021, 2016, 2015, 2014). Measurements are problematic in the presence of power lines or unknown buried metal infrastructure.
Subsurface Data 33
PNNL-33647
- Similar to ERT, electromagnetic data that quantifies subsurface electrical structures could be used in statistical prediction (kriging) models and Bayesian Elipgrid procedures to inform sampling and compliance survey designs and planning.
Figure 26. Electromagnetic Inversion Results and Interpreted Hydrogeologic Framework from Airborne Surveys Collected Over the Hanford Site (Jaysaval et al. 2021) 3.4.4.6 Ground Penetrating Radar Ground penetrating radar (GPR) is a high-frequency electromagnetic method capable of submeter resolution of soil structure. GPR is used extensively for nondestructive geotechnical, concrete, and pavement assessments and to locate buried targets and infrastructure. In most earth materials, GPR signals penetrate only a few meters to a few tens of meters; hence GPR is used for shallow characterization and shallow targets. Time-lapse GPR has been used to monitor infiltration, solute transport, and saturation changes, and to detect leaks from subsurface infrastructure. PNNL has used ground-based and cross-well GPR for diverse problems including soil desiccation monitoring, as in the cross-well example shown in Figure 27 (Truex et al. 2018).
- Similar to the technologies above, GPR data that quantifies subsurface infrastructure and leak detection could be used in statistical prediction (kriging) models and Bayesian Elipgrid procedures to inform sampling and compliance survey designs and planning.
Figure 27. Inversion Results for Cross-Hole GPR Data Collected Before (left) and After (right) the Start of a Soil Desiccation Experiment at the Hanford Site (Truex et al. 2018).
Subsurface Data 34
PNNL-33647 3.4.4.7 Hydrogeophysics Toolbox The concept of a Hydrogeophysics Toolbox was first advanced in near-surface geophysics over 20 years ago and is well established (e.g., Day-Lewis et al. [2017]). The Toolbox comprises borehole, cross-hole, and airborne geophysical methods (Figure 28). It was built on the premise that geophysical method selection should be based on the following:
- 1. Project objectives: For a given method to support project objectives, it must provide information relevant to the problem.
- 2. Site conditions: For a method to work at a given site, conditions there cannot preclude data acquisition or mask targets (e.g., presence of underground infrastructure, lithology, or soil type).
- 3. Resolution, spatial coverage, and depth of investigation requirements: What is achievable with a given method must also align with those required by the project objectives.
- 4. Potential synergies between complementary methods: Some project objectives will call for use of multiple methods in concert, capitalizing on the synergy between methods.
Figure 28. Methodology Suggestions as Scale of Investigation and Resolution Increase Although these concepts are widely accepted, their implementation is limited by disconnects between practitioners and those deciding to contract, fund, or require geophysical investigations. The following articles about method selection tools have been published to address these disconnects, improve method-selection practices, and maximize return on geophysical investments;
- The Fractured Rock Geophysical Toolbox Method Selection Tool (Day-Lewis et al. 2016) and the Groundwater/Surface-Water Method Selection Tool (Hammett et al. 2022) guide users to enter information about project goals and site characteristics and then provide guidance on methods likely to address stated goals under site-specific conditions.
Subsurface Data 35
PNNL-33647
- The Scenario Evaluator for Electrical Resistivity tool (Terry et al. 2017a, b) is a simplified image-appraisal tool that allows users to assess the ability of ERT to resolve hypothetical targets in the subsurface, given user-specified characteristics of the targets (shape, depth, contrast) and survey design.
Additional method selection and image-appraisal tools could be designed to address the needs of the NRC related to topics including siting, leak detection, remediation monitoring, and closure or decommissioning.
Subsurface Data 36
PNNL-33647 4.0 Statistical Methods for Subsurface Data Analysis The Task 1a report describing statistical methods for continuous data collection (Fagan et al.
2022) discusses the need to account for spatial autocorrelation when making statistical decisions using geospatial data. If spatial autocorrelation is not accounted for when present, statistical hypothesis tests and statistical interval calculations may not meet acceptable decision error tolerances. Additionally, subsequent identification and delineation of potentially elevated areas of residual contamination may not reflect the true state of a site. These considerations apply to subsurface volumes as well and have additional complexity with a third dimension.
Section 4.1 provides an overview of geostatistical methods, applicable in three dimensions for subsurface applications.
As the complexity of the subsurface structure(s), geophysical characteristics, fate and transport mechanisms, and contamination increases, the complexity of geostatistical methods required to assess compliance will increase as well. We attempt to cover the range from simple to complex methods in this section.
4.1 Geostatistics Overview EPRI (2016) provides a helpful guide to understanding geostatistics (also referred to as geostatistical analysis) and its role in RSSIs at nuclear sites. Key concepts presented by EPRI are summarized here to provide context for more specific methodological approaches considered in the subsequent sections. This report assumes, however, that the reader has a basic understanding of geostatistics, including knowledge of concepts such as variogram models and their parameters (i.e., range, sill, nugget), kriging models (i.e., models used to predict or interpolate observations between locations where data were collected), and simulation (i.e., using parametric or nonparametric distributions or models to generate plausible observations, with variability). Further information on these topics can be found in numerous resources, including Cressie (1993), Deutsch and Journel (1998), Goovaerts (1997), Griffith (2017), and other references used throughout this report.
Geostatistics can be considered a class of methods that enable inference and simulation of spatially (and perhaps temporally) referenced variables based on observed measurements from a set of sampled locations. It is different from traditional statistical methods in that it explicitly models spatial correlation. EPRI (2016) considers three major elements of geostatistics:
structural analysis, kriging, and conditional simulation.
4.1.1 Structural Analysis Structural analysis includes estimation and inference of spatial structure from available data, usually using a (semi)variogram function to describe dependence and relative correlation between two observations as a function of the distance between them. A sample variogram is a function of the distances separating sampled locations, which are calculated based on x-y location pairs in 2-D or x-y-z location triplets in 3-D. We expect the reader to be familiar with the sill, nugget, and range components of the variogram (details can be found in EPRI [2016] and other references cited in this section of the report). Figure 29 depicts some of the six commonly used variogram models.
Statistical Methods for Subsurface Data Analysis 37
PNNL-33647 Figure 29. Commonly Used Variograms in GLS Estimation Variogram distances can also incorporate directional information between locations. Spatial structure is called isotropic when the variogram is equivalent regardless of the orientation (direction) between points, but anisotropic when it differs depending on its orientation. An example of anisotropy on the surface is when directional contaminant deposition results due to wind, where concentrations are more correlated in the same direction as the wind than they are in the direction perpendicular to the wind. Anisotropy also may be present in the subsurface within an individual layer or across layerswhen conditions result in directional contaminant deposition. Anisotropy may occur multi-directionally, for example, when contamination fate and transport properties differ between distinct soil layers. This phenomenon can be incorporated into geospatial models using multiple variograms (of the form and function shown in Figure 29).
It is critical to assess whether anisotropy is present and account for it in modeling and survey designs, if so. Anisotropy can increase or decrease the number of sample locations, depending on the number of distinct variograms. It is also important to account for in statistical modeling and analysis because it will affect uncertainty estimates, statistical hypothesis test conclusions, and subsequent decisions. In structural analyses, anisotropy can be diagnosed by examining multiple directional variograms and/or employing statistical hypothesis tests (Weller and Hoeting 2016). Anisotropy can be modeled by allowing variogram ranges to vary in different directions or by estimating coordinate transformations for which isotropy holds (Schabenberger and Gotway 2017). Horizontal dependence within a layer is likely to be stronger than vertical dependence across layers (EPRI 2016).
4.1.2 Kriging Kriging is a tool used to predict (interpolate) quantities of the variable of interest over an entire site (area in 2-D or volume in 3-D) based on a set of sample observations at locations within the site. While some kriging methods ignore measurement uncertainty, others are based on spatial models with small- and/or large-scale spatial variability terms. Kriging uncertainty, or prediction variance for predicted values at each location, is governed by the number of observations, distances between sample locations, variability in observed values, and variogram differences (unidirectional or multidirectional). However, kriging uncertainty only represents observed and estimable geophysical properties of a site, where the variogram provides the estimation Statistical Methods for Subsurface Data Analysis 38
PNNL-33647 mechanism. Estimating the kriging uncertainty provides estimates of the probability of exceeding a threshold, which can in turn be used to optimally place new sample points and calculate the number of sample locations required for decision-making (EPRI 2016).
Different kriging methods are used to address different challenges. Universal kriging can be used when there is nonstationarity in the mean (e.g., when mean concentration is a function of location). Indicator kriging is suggested in the literature to address non-normality by transforming a continuous variable of interest into a categorical variable (e.g., above or below a threshold) and corresponding nonparametric estimation (EPRI 2016). However, PNNL does not recommend transforming continuous measurements into discrete categories for this purpose alone, because several alternative methods can handle the continuous but non-normal data.
Indicator kriging should be used when presence/absence data or categorical data are available.
Cokriging can be used to study multiple interrelated variables and may be useful when measurements of one variable are easier/cheaper to obtain than others (e.g., when gamma radiation measurements are easier to collect than alpha radiation measurements or modeled data are more feasible to collect than field measurements), but the two are correlated so that inference about the more challenging data can be made based on the more accessible data.
Spatial regression models, also referred to as kriging with external drift, can be used for kriging when explanatory variables (covariates) are available at locations for which the response variable is not (EPRI 2016).
4.1.3 Conditional Simulation Conditional simulation is used to generate individual realizations of the variable of interest, where each realization includes the same statistical and geostatistical properties determined during structural analysis. The statistical properties are those represented by the variogram structure of the observed data, and the geostatistical properties are those represented by the physical characteristics of the site. Simulations model the impact of statistical and geophysical properties at each location by conditioning on observed data at nearby locations. Conditional simulation makes use of the variogram structure and accounts for prediction error (variability) associated with the interpolated (kriged) estimates, providing a set of equally probable predicted observations at each location. Conditional simulation may be useful (or required) in the following applications (EPRI 2016):
- Nonlinear mapping functions: When the variable of interest is modeled through a nonlinear relationship with the interpolated (kriged) values (e.g., via an inverse transformation to achieve normality), analytical uncertainty estimates may not be practical or possible to calculate.
- Spatial continuity: While interpolated (kriged) estimates provide the best estimates of the variable of interest in many cases, it does not do so when local variability is present but unobserved. Conditional simulation may provide insight into unobserved local variability.
However, some prior knowledge (from data or subject matter expertise) should be incorporated to ensure it is represented accurately.
- Probabilistic assessment of global quantities: Conditional simulation may be used to develop probability distributions for interpolated (kriged) values based on derived values (e.g., total source term over an area or volume).
- Complex dynamical simulation: When site survey results will be used to inform fate and transport models, conditional simulation may be useful for generating probabilistic estimates of model inputs using Monte Carlo framework.
Statistical Methods for Subsurface Data Analysis 39
PNNL-33647 Conditional simulation generates realizations of the spatial process while honoring the observations (Schabenberger and Gotway 2017). It is designed to preserve the characteristics of the data (e.g., mean, distribution, spatial correlation) in the realizations. Sequential (Gaussian or indicator) simulation, conditioning by kriging, simulated annealing, and simulation from convolutions are all methods that can be used to produce realizations of a random field that can then be used to assess uncertainty. Gaussian and multi-Gaussian simulation use the simple kriging mean and variance to define a point-by-point cumulative distribution function of the spatial process. Rather than kriging values at pixels where no data were observed, values from the cumulative distribution function are simulated via a Monte Carlo sampling algorithm. Multiple simulations then result in a pixelwise empirical distribution of the spatial process. Pixelwise empirical distributions estimated from conditional simulations lead to data combinations that can be used to estimate statistical quantities (e.g., mean, upper percentile) and uncertainties.
Historically, geostatistical simulation approaches are used to simulate geophysical phenomena, such as properties of flow and transport, rather than to predict concentrations of contaminants.
Combining the output of simulation models with observed concentrations is usually done by way of cokriging, multi-Gaussian kriging, or stratification (Goovaerts 1997, Goovaerts 2000; Stewart et al. 2006). Recent research includes exploring numerical variograms (Pannecoucke et al.
2020) estimated from simulation data and alternative methods for estimating kriging variance, including simulation uncertainty (Silva 2021).
Many conditional simulation algorithms exist, including sequential (multi)Gaussian simulation, sequential direct simulation, simulated annealing, and filter simulation, and the area continues to be an area of active research (Ortiz 2020; Metahni et al. 2019). A number of these algorithms are implemented in existing software products. For example, several R packages, including geoR (Ribiero and Diggle 2001), varycoef (Dambon et al. 2021), and RandomFields (Schlather et al. 2015), perform multi-Gaussian simulation. Additionally, the Geostatistical Software Library (GSLIB) contains code for sequential Gaussian simulation (see Section V.2.3 in Deutsch and Journel 1998). Several VSP geostatistical methods are based on GSLIB implementation.
4.2 Geostatistical Methods for the Subsurface In this section, we review numerous geostatistical methods for analysis of collected subsurface data and subsequent decision making in the context of the compliance phase of nuclear site decommissioning. Modeling and remediation of subsurface contamination can vary widely in their complexity. Factors such as the depth of contamination, composition of the subsurface (e.g., how many layers, what are layers made of), contaminant transport for both current and future scenarios, and exposure pathways all contribute to the complexity of subsurface applications. Keeping this complexity in mind is important when considering appropriate methods for both modeling and decision-making. Low-complexity sites may require only small extensions of the methods used for surface decommissioning while highly complex sites may require extensive modeling, the use of new modeling techniques, and/or complex decision frameworks.
4.2.1 Generalized Least Squares Regression GLS is a method used to estimate the unknown parameters in a linear regression model in the presence of spatial or temporal autocorrelation. GLS includes the commonly used ordinary least square (OLS) as a special case, as its name suggests. It is employed in circumstances where the OLS estimator is not the best linear unbiased estimator due to heteroskedasticity and dependence (in contrast to fundamental assumptions of homoskedasticity and independence, or Statistical Methods for Subsurface Data Analysis 40
PNNL-33647 absence of correlation). An advantage of GLS is that it does not make distributional assumptions about the data, but only requires assumptions about the first and second moments (Schabenberger and Gotway 2017), although this can introduce challenges to developing confidence intervals representing the uncertainty of predicted values at unobserved locations.
GLS accounts for spatial structure in data by incorporating it into the estimation of the covariance matrix based on a fitted variogram model of the data, and it has the flexibility to capture different spatial structures across layers or a grouping variable, making it a suitable method for subsurface investigation. In practice, the variance-covariance matrix is computed by first fitting the appropriate variogram model to the empirical variogram obtained from the data.
The spatial form is converted to covariances using the corresponding covariance function after selecting the variogram model based on available data or prior information. If the spatial structure is exponential, for example, an exponential variogram can be fit and the spatial structure captured using an exponential covariance function.
Contaminant fate and transport and expected impacts of geophysical characteristics on contamination in the subsurface could be incorporated into GLS methods through explanatory variables in the regression model.
VSP would require algorithm development and coding to add GLS capabilities for model fitting, estimation, prediction, and hypothesis testing.
4.2.1.1 Hypothesis Testing As described in the Task 1a report (Fagan et al. 2022), GLS regression can be used within the MARSSIM framework to accurately estimate average activity over a site or to compare the average activity to a threshold value using a statistical hypothesis test in surface applications.
Similarly, GLS can be used for hypothesis testing in subsurface applications. GLS is recommended in favor of OLS, which cannot incorporate the non-Gaussian spatial structure of the data. For example, consider the data in Figure 30 showing U-235 concentrations observed in four wells near the 241-BX Tank Farm at Hanford.
Assume we want to compare mean U-235 concentration at the site (estimated using well observation average ) to a historical average of 5 pCi/g and apply the following statistical hypothesis testing framework 16:
- : 5 (mean concentration of U-235 does not exceed historical average)
- : > 5 (mean concentration of U-235 is greater than historical average).
16 The null hypothesis that the mean concentration of U-235 does not exceed the historical average is consistent with alternative Scenario B in NUREG-1575 (i.e., that the survey unit meets the release criteria). The default scenario in MARSSIM, is Scenario A, which assumes the survey unit does not meet the release criteria.
Statistical Methods for Subsurface Data Analysis 41
PNNL-33647 Figure 30. Concentration of U-235 at Four Selected Wells at the Hanford Site In Table 2, both OLS and GLS models are used to estimate the means, standard errors, t-values, and p-values corresponding to this test. The OLS model overestimates average U-235 concentration and results in a p-value less than 0.05, leading to rejection of the null hypothesis and a conclusion that the site U-235 levels exceed the historical average. This represents a Type I error that erroneously leads to the conclusion that the site is dirty when the site is actually clean. The GLS model, on the other hand, leads to a smaller mean and larger p-value, and fails to reject the null hypothesis, and thus leads to the conclusion that the site mean does not exceed the historical average.
Table 2. Example OLS and GLS Model Output Model Standard Error t-value p-value OLS 8.777 0.325 11.625 <0.0001 GLS 6.423 1.840 0.773 0.2198 These are results for one realization of this example (based on n = 150 observations). When multiple simulations are used, expected Type I error rates can be calculated. Figure 31 presents average Type I error rates for the same example but varying numbers of observations. When using the OLS model, the Type I error rates are always higher than the GLS model and increase as the number of observations increases. The reason for this is that pseudo replicates (observations that are not statistically independent but treated as if they are) produce additional observations that lead to "more significant" results, but these results are incorrect due to the fact that the spatial dependence was not taken into consideration. It is counterintuitive and highlights where practitioners can go wrong when spatial dependence is present but not accounted for.
Error rates for the GLS model, on the other hand, remain roughly constant and near the nominal level ( = 0.05) regardless of the number of samples. This shows that the likelihood of making an erroneous conclusion is greater using an OLS model and increases as sample sizes increase, where using the GLS model, which addresses spatial dependence and heteroskedasticity results in the expected nominal error rates.
Statistical Methods for Subsurface Data Analysis 42
PNNL-33647 Figure 31. Comparing OLS and GLS Type 1 Error Rates 4.2.1.2 Mean Estimation In the subsurface, the following four scenarios were considered for mean estimation:
- 1. Mean and spatial structure are the same in all layers (but spatial parameters vary).
- 2. Mean is the same but spatial structure varies (and spatial parameters vary).
- 3. Mean varies but and spatial structure is the same in all layers (spatial parameters vary).
- 4. Mean and spatial structure vary (and spatial parameters vary).
Figure 32 shows the observed data (simulations) for seven wells in a three-layer subsurface.
Observations are stacked vertically within each well, with white space representing boundaries between each layer occurring at depths around 95 and 120 units. Figure 32a depicts the first two scenarios with observations when the mean is the same in all layers (50 units). Figure 32b depicts the second two scenarios with observations when the mean differs between layers (from 30 to 60 to 100 units as depth increases).
Statistical Methods for Subsurface Data Analysis 43
PNNL-33647 Figure 32. Example Well Monitoring Data (simulated) in Seven Wells that Span Three Subsurface Layers. (a) Constant mean concentration. (b) Varying mean concentration. Spatial variation is not shown.
GLS regression was applied to the simulated data from each scenario to examine its capabilities for estimating mean concentrations and standard errors. It was applied to each layer separately and to the combined volume. Table 3 shows the results. Regardless of whether layers are modeled separately or as a combined volume the mean values are estimated accurately.
However, when the mean varies across layers the model of the complex volume is not precise, as reflected by its high standard error in those scenarios.
Table 3. Summary of GLS Models Used to Estimate Subsurface Means with Three Layers Mean is Equal across Layers Mean Varies between layers True True Spatial (simulated) Estimated Mean (simulated) Estimated Mean Layer Structure Mean (standard error) Mean (standard error) 1 47.1 (1.33) 100.0 97.1 (1.33) 2 48.1 (1.40) 60.0 58.4 (1.40)
Same in 3 50.0 51.4 (0.77) 30.0 31.2 (0.77) all layers Complex 49.8 (0.72) 63.3.0 64.3 (18.51) volume 1 48.1 (1.33) 100.0 97.1 (1.33) 2 Varies 48.4 (1.70) 60.0 58.1 (1.70) 3 between 50.0 51.2 (0.75) 30.0 31.5 (0.74)
Complex layers 48.9 (0.82) 63.3 64.5 (18.64) volume In practice, mean values, spatial structure, and spatial parameter values are likely to be unknown prior to data collection. As a result, histograms of contamination metrics and variogram plots should be used to help determine which strategy is optimal (i.e., modeling the layers independently or as a complex volume), given observed subsurface complexity. If prior information is available from previous RSSI phases, then models can be used to estimate uncertainty and inform survey sample design for the compliance phase. Additional considerations to determine the dimensionality of the approach are discussed in Section 4.3.
Statistical Methods for Subsurface Data Analysis 44
PNNL-33647 4.2.2 Local Indicator of Spatial Association A local indicator of spatial association (LISA), also known as the Local Morans I statistic, is a method presented by Anselin (1995) to identify local clusters and local spatial outliers. The method can be applied in 2-D surface and 3-D subsurface applications to discover hot spot areas or volumes, respectively. In the subsurface, it can be used to detect statistically significant volumes of high contamination within otherwise low contamination volumes and, conversely, it can detect cool spots, or volumes of significantly low contamination compared to surrounding high contamination volumes.
A LISA value is calculated for each observed location. Positive LISA values indicate neighboring locations have similarly high or low contamination. Negative LISA values indicate neighboring locations are higher/lower concentration. Each LISA value has a corresponding p-value, where lower p-values indicate statistically significant clusters of locations. There are four types of statistically significant clusters:
- 1. Several locations with high values (HH)
- 2. Several locations with low values (LL)
- 3. Single high value location surrounded by low value locations (HL)
- 4. Single low value location surrounded by high value locations (LH).
Clusters of high values (HH and HL) indicate prospective areas of relatively heightened contamination (i.e., hot spots), while low values (LL and LH) indicate potential areas of comparatively low activity (i.e., cool spots).
To illustrate the application of LISA to subsurface compliance survey analysis, it was applied to a simulated dataset shown in Figure 33 (data were simulated based on concentrations observed in eight wells near Hanford's 241-BX Tank Farm). The data are shown in Figure 33a with high values represented by red points. Locations with no statistically significant spatial grouping, cold spot clusters, and spatial low spot outliers are illustrated by blue points. The hot spot analysis is shown in Figure 33b, in which significant hot spot clusters are indicated in red and cold spot clusters in blue. Based on this analysis, conclusions could be made for each well.
Figure 33. Example LISA Application to Well Monitoring Data. a) Data and b) hotspot analysis.
Statistical Methods for Subsurface Data Analysis 45
PNNL-33647 In this example, LISA was applied to the observations at sampled locations. However, additional research is recommended to determine how the LISA method could be expanded to apply to kriged or modeled volumes so that it can be used to identify hot spot volume boundaries that include locations between wells. Additional research would be required to determine if and how contaminant fate and transport or the impacts of geophysical characteristics on contamination in the subsurface could be incorporated into LISA.
VSP would require algorithm development and coding to add LISA capabilities for hot spot detection. Additional research would be required to determine if and how data from previous parts of the RSSI process and/or outside data (e.g., contaminant fate and transport datasets, modeled geophysical subsurface characteristics) could be incorporated into LISA. If they could, then it seems promising for developing an initial sampling design based on such data.
4.2.3 Kriging Methods Kriging methods are a valuable tool for geostatistical modeling (Gelfand and Schliep 2016), with standard methods including simple, ordinary, block, and indicator kriging. These assume an underlying spatial Gaussian process for which the conditional expectation of a new observation is assumed to be linear with respect to the observed data. Kriging or prediction models must address change-of-support challenges when the size, dimension, shape, and/or orientation of observations or measurements on different variables differ from the support required for the kriged (predicted) surface or volume. For example, when observed values are collected at different resolutions or with different fields of view. EPRI (2016) recommends block kriging and/or downscaling to address direct field measurement support typically being smaller than scanning technology support, and the support associated with a survey unit being much larger than either of those (EPRI 2016).
Gaussian models provide convenient frameworks for creating confidence and prediction intervals for values at unobserved locations (e.g., testing for hot spots by estimating exceedance at a location), estimation and testing of the site mean, and Gaussian simulation.
One of the criticisms of kriging is that it results in over-smoothed local variability, or predictions that are less variable than is realistic (or observed) in the environment (Diggle and Ribeiro 2010). Research has shown, however, that more complex methods are not always the best optioneven for complex sitesand that careful consideration should be given to selecting the right method depending on the desired statistics (Metahni et al. 2019).
Additional research will be required to determine if and how contaminant fate and transport and geophysical characteristics in the subsurface could be incorporated into each kriging method.
VSP currently has 2-D kriging capability for simple, ordinary, and indicator kriging (Fagan et al 2022; Matze et al. 2014). Fixed rank kriging (FRK) were included in a recent VSP release (September 2022). Existing 2-D kriging methods in VSP will be applicable for sites that model the subsurface as a collection of 2-D layers, although a set of assumptions must hold in those cases.
- Subsurface volume is modeled as a collection of 2-D layers with isotropic contamination where previously collected contamination data include one data point per x-y location in the z-direction (i.e., data were only collected at one depth and that depth is the same for every x-y location).
Statistical Methods for Subsurface Data Analysis 46
PNNL-33647
- Subsurface volume is modeled as a collection of 2-D layers with isotropic contamination.
Previously collected contamination data are observed at multiple z-locations per x-y location (e.g., multiple data points along the z-axis within a borehole or well at each x-y location) but observations are isotropic and stationary and can be averaged at each x-y location.
The existing 2-D FRK capability will be applicable in the following cases:
- Subsurface volume is modeled as a collection of 2-D layers with anisotropic contamination in the x-y direction where previously collected contamination data include one data point per x-y location in the z-direction (i.e., data were only collected at one depth and that depth is the same for every x-y location).
- Subsurface volume is modeled as a collection of 2-D layers with anisotropic contamination in the x-y direction. Previously collected contamination data are observed at multiple z-locations per x-y location (e.g., multiple data points along the z-axis within a borehole or well at each x-y location) but observations are isotropic and stationary in the z-direction and can be averaged at each x-y location.
Existing 2-D kriging methods will need to be updated with the capability to model anisotropy or account for it by including independent variables (FRK intrinsically models the anisotropy but cannot currently incorporate independent variables) in the model in the following cases:
- Subsurface volume is modeled as a collection of 2-D layers with anisotropic contamination in the x-y direction. Data were collected at one or multiple z-locations per x-y location and can be averaged within each x-y location. Supplementary data are available within each layer (e.g., data from models or measurements such as those described in Section 3.4) and will be used as independent variables in the kriging (prediction) model to account for the anisotropy.
All VSP kriging methods (including FRK) will need to be updated to accommodate 3-D data when the subsurface volume is modeled as a complex 3-D volume or as a collection of 3-D layers with anisotropy in the z-direction. For example, one layer may represent the vadose zone and another an aquifer; layers may represent geologic formations or soil horizons; or distinct layers informed by different dose models may each be modeled as a separate 3-D layer. Each physical layer could potentially be further discretized into multiple layers for statistical sampling and analysis based on mean concentrations, variograms, and/or other statistical properties.
4.2.3.1 Multi-Gaussian Kriging VSP implements multi-Gaussian kriging (MGK) as described by Deutsch and Journel (1998).
MGK uses the Gaussian conditional cumulative distribution function at each location (pixel) where predictions are desired using the simple kriging mean and variance structure. Data are converted to z-scores, kriged, and then back-transformed, resulting in predictions that are influenced by local data values rather than simply the distance between data. MGK is an intrinsic model, meaning predictions tend toward a local rather than global mean. This improves predictions from a geoscience point of view (Goovaerts 1997). Although there is no assumption of mean stationarity, VSP recommends removing large-scale trends prior to kriging using MGK.
The MGK module in VSP contains several options for maps of concentration and uncertainty.
Concentration maps can be produced for mean, median, or user-specified percentiles and uncertainty maps for conditional variance, interquartile range, and reference uncertainty index.
Additionally, VSP provides an option for generating maps showing the probability that a potential contaminant will exceed a given threshold.
Statistical Methods for Subsurface Data Analysis 47
PNNL-33647 VSP would require updates to the current 2-D MGK module so that it could accommodate 3-D data.
4.2.3.2 Kriging Non-Gaussian Data In many applications a Gaussian assumption may not be appropriate; for example, when concentration measurements are zero or positive, have a skewed distribution, or take the form of presence/absence measurements. A common way to address such situations is to transform the data so that the Gaussian assumption is met with the transformed data, and then perform estimation and kriging on the transformed data. Cressie (1993) suggests using trans-Gaussian kriging, a general approach to developing optimal kriging predictors in these scenarios. A common example is lognormal kriging in which a logarithmic transformation is applied to the data and the Gaussian assumption holds in the transformed data. Trans-Gaussian kriging explores different ways to appropriately back-transform kriged estimates and uncertainty to the original unit of measure. Multi-Gaussian kriging is a related approach that considers the change-of-support problem, where non-Gaussian point measurements are used to estimate block values. If the data represent concentration measurements, these methods can be used to estimate site average concentrations and identify hot spots.
Spatially explicit generalized linear models and indicator kriging are other ways to model non-Gaussian data. These models are commonly used for spatially referenced binary (0/1) outcomes or count data. They assume a data-generating model (e.g., Bernoulli, Poisson) and use a link function to model the mean of the data-generating process. They also typically include spatial random effects to account for unobserved covariates or spatial structure.
Indicator kriging can be used to krige with binary (0/1) outcomes and is a nonparametric approach that can be used to estimate the probability of exceedance at a location. If the data represent concentration measurements, these methods can be used to estimate site averages and identify hot spots. However, if the data only consist of 0/1 outcomes, it may not be possible to estimate site means using these methods.
VSP would require updates to the current 2-D non-Gaussian module so that it could accommodate 3-D data.
4.2.3.3 FRK FRK is one solution when large amounts of data make standard kriging calculations slow or intractable, although other solutions are available for the large data problem (Heaton et al.
2019). FRK also allows data with different fields of view (i.e., datasets with different supports) to be combined into a single model. While the volume of subsurface data may not warrant use of FRK in the context of computational speed (since we expect small datasets), it might be useful to predict values at locations with subsurface data collected with different fields of view, such as geophysical and fixed laboratory data (Cressie and Johannesson 2008; Zammit-Mangion and Cressie 2021). An additional benefit of FRK is that it employs a flexible covariance structure that does not assume stationarity or isotropy, which could be very useful for modeling complex and changing spatial structure and anisotropy in subsurface applications.
Application of this method to the subsurface should be explored, especially with respect to its flexible covariance assumptions. VSP would require updates to the current 2-D FRK module so that it can accommodate 3-D data.
Statistical Methods for Subsurface Data Analysis 48
PNNL-33647 4.2.4 Geospatial Extension to MARSSIM, Multi-scale Remedial Design Model (MrDM), and Multi-scale Remedial Sample Design Model (MrsDM)
Stewart (2011) developed a decision framework called the Geospatial Extension to MARSSIM (GEM) using geostatistical simulation to model the uncertainty (probability) about exceeding a DCGL anywhere within the study area that can integrate different types of data, including field methods, laboratory methods, etc. The method was also summarized by Gogolak (2022). The GEM is not limited to any particular geostatistical simulation methodology. The GEM case study uses sequential indicator simulation (Goovaerts 1997), citing its ability to encode values as 0 or 1 as a strength for combining metrics from different sources with varying accuracy and precision. Stewart (2011) notes however that other methods are viable to integrate field measurement data, such as sequential Gaussian simulation. GEM simply requires a model that can model the uncertainty around exceeding DCGLs and can integrate various types of data.
Although multiple sources and types of information can be used to formulate distinct CSMs, the variation between them and the uncertainty associated with each is not routinely accounted for.
GEM addressed this limitation by using a stochastic conceptual site model to delineate the probability of complying with the regulatory limit(s). GEM uses a model-based decision rule(s) called the regulatory limit rule that requires this probability to be less than a specified limit and geostatistical simulations site to calculate the probability of exceedance. When the test fails and the site is determined to be out of compliance, GEM provides a method of delineating the spatial extent of the subsurface volume(s) to be remediated using the Multi-scale Remedial Design Model (MrDM) approach. When the identified subsurface volumes prove too costly for remediation in its entirety, the Multi-scale Remedial Sample Design Model (MrsDM) approach is proposed for determining additional sample locations that could decrease modeled uncertainty and therefore decrease MrDM remedial volumes. MrDM and MrsDM can be implemented in an iterative fashion to determine if and how much additional subsurface volume requires remediation (Stewart 2011). Check and cover sample designs are also described in SC&A (2022), Stewart (2011) and NUREG/CR-7021 (Stewart and Powers 2009). This method is discussed further in Section 5.0 in this report.
Because GEMs MrDM and MrsDM are optimized under the paradigm that a site has failed compliance (using observed and inferred concentrations) and are used to seek remediation volumes, it may be better geared toward remediation phase planning and decision making than toward the compliance phase. Additional research would be required in the following areas:
- To determine if and how contaminant fate and transport or of geophysical characteristic impacts on contamination in the subsurface could be incorporated into GEM, MrDM, and MrsDM.
- To incorporate quantify uncertainty of the estimates from each step and set of results and incorporate that into subsequent decision making.
- In cases where remediation is pursued and excavation occurs, resulting in exposed surface.
Research will be required to determine if and how GEM can be used to determine vertical and lateral boundaries (i.e., how deep should the survey/remediation go, cutoff depth conditions, and/or when to continue to keep excavating to measure contamination)
VSP (Matzke et al. 2014) includes a module that identifies sampling redundancy using global kriging weights to identify the relative importance of samples in mapping contaminant plumes and to identify sample locations or wells that could be removed from the sampling schedule.
VSP would require new code development to implement GEM, MrDM, and MrsDM for the Statistical Methods for Subsurface Data Analysis 49
PNNL-33647 purposes of characterization an/or remediation surveys. In cases where these methods were used for the characterization survey and collected data indicate the site is ready for compliance phase, the findings (including estimated contamination and uncertainty) could be used to design the compliance survey. In cases where characterization survey results indicate remediation will be required, MARSSIM provides guidance on next steps.
If resources have been developed but unreleased in previous SADA versions, they could be useful for future VSP development purposes.
4.2.5 Bayesian Methods The use of Bayesian geostatistical methods has rapidly increased over the last two decades due to increases in computational power and availability. An advantage of using Bayesian geostatistical methods is that they account for uncertainty in covariance parameter estimates, and they provide a framework for incorporating prior knowledge into estimates of concentration, soil properties, etc. A disadvantage is that they can be computationally burdensome; however, several methods have been developed to reduce computational time (e.g., Banerjee et al. 2008; Heaton et. al. 2019). Empirical Bayesian kriging (Gribov and Krivoruchko 2020) has increased in popularity due to its availability in ArcGIS Pro software. One disadvantage is that it does not model anisotropy or provide cokriging capability (ArcGIS Pro 2.8).
Spatially explicit Bayesian regression models may be useful in subsurface applications because they can combine concentration measurements and covariate information (e.g., soil types, outputs from hydrologic models) to predict concentrations at unobserved locations while accounting for spatial correlation. Unlike traditional frequentist approaches, these methods quantify uncertainty in the estimates of the spatial covariance parameters and propagate this uncertainty into interpolated (kriged) values. New Bayesian methods for identifying locations that exceed a threshold, particularly those discussed in French and Hoeting (2015) are worth researching for applicability to the subsurface.
An area requiring further research is the incorporation of prior belief of contamination in a fully Bayesian paradigm. The Markov Bayes method uses conditional probability to incorporate prior belief and measurements to update the contamination concern map (CCM) but does not provide a posterior distribution of contamination (and therefore lacks uncertainty estimation as part of the analysis). A fully Bayesian approach would provide uncertainty as a function of observations based on the number of and spatial proximity of hard data points and the prior belief map.
Contaminant fate and transport and geophysical characteristics could be incorporated into Bayesian methods through explanatory model variables and/or prior distribution(s).
Research to identify the best Bayesian methods to include in VSP would be required, and then subsequent algorithm development and coding that would be needed to accommodate Bayesian analysis for subsurface applications.
4.2.6 Variogram Tomography Most geostatistical models assume second-order stationarity for the covariance function. If anisotropy is assumed to exist, the anisotropy is assumed to be globally constant also; that is, the direction and magnitude of anisotropy is the same across the domain. This assumption is inappropriate for complex geological environments in the subsurface when highly nonlinear patterns of contamination exist. Further, the random field representing a subspace can be either Statistical Methods for Subsurface Data Analysis 50
PNNL-33647 non-homogeneous or homogeneous and it is important to accurately characterize it in subsurface modeling to avoid biased parameter estimates. Here are potential approaches for addressing these.
- Incorporating locally varying anisotropy: Imposing constant anisotropy (i.e., second-order stationarity) in models with known locally varying spatial features can compromise the accuracy of model estimates. Rather than using an approach that incorporates multiple variograms, which is potentially computationally expensive depending on the data size, incorporating locally varying anisotropy can achieve decreased computational requirements as well as increased accuracy (Lillah et al. 2012). Boisvert and Deutsh (2011) provide an algorithm that uses GSLIB in combination with other executable packages for kriging with locally varying anisotropy.
- Kriging with external drift: Kriging with external drift accounts for nonstationarity in the mean by modeling it as a linear combination of covariates (Wadoux et al. 2018; Goovaerts 2000) where the linear combination is allowed to vary locally using either a frequentist (traditional) approach (Dambon et al. 2022) or a Bayesian approach (Gelfand et al. 2003). Wadoux et al.
(2018) proposes an extension of kriging with external drift which accounts for nonstationarity in both the mean and variance. Similar to the mean model, the variance is modelled as a linear combination of covariates. Such extensions of kriging with external drift could be explored further for subsurface modeling. Where GLS methods assume a global functional relation for the mean, kriging with external drift allows the effect of each explanatory variable to vary spatially as a function of location.
- Variogram matrix functions for vector random fields with second-order increments (needs further investigating): Developing novel covariance or variogram matrix functions for vector random fields with second-order moments can also yield suitable models to model subsurface (Ma 2011; Du and Ma 2012). For example, Du and Ma (2012) proposes variogram matrix functions for vector elliptically contoured random fields. Further research studies on this can be explored further for subsurface modeling.
VSP updates would require further research to identify which of these methods to include, and then subsequent algorithm development (or library identification) and coding to accommodate their use in subsurface applications.
4.2.7 Artificial Intelligence/Machine Learning The development and use of machine learning (ML) and artificial intelligence (AI) methods have dramatically increased over the last 15 years. One of the reasons for the popularity of AI/ML methods is their flexibility to discover and model complex and nonlinear relationships in massive datasets. AI/ML methods have been applied to subsurface data to perform tasks such as delineating layers (Wohlberg et al. 2005), clustering observations (Romary et al. 2015), and mapping contaminant plumes (Tao et al. 2019). As one example, recent extensions to random forest algorithms have been developed for both global (Hengl et al. 2018) and local (Georganos et al. 2021; Ancell et al. 2021; Benito 2021) spatial regression problems.
AI/ML methods could be used to predict contaminant levels at unmeasured locations by combining concentration measurements with other subsurface measurements or model outputs (e.g., groundwater transport, soil properties, layer information). While traditional AI/ML methods are particularly suited to scenarios with large amounts of data available to train and test the models, they are not appropriate when very few observations are available and they often lack interpretability. Another challenge with using these methods is that it can be difficult to quantify Statistical Methods for Subsurface Data Analysis 51
PNNL-33647 uncertainty in order to create confidence intervals or perform hypothesis tests. In Section 3.4.1.4, few-shot machine learning was discussed as an approach being advanced in conjunction with remote subsurface sensing techniques and high-performance forward prediction for the purposes of improving site characterization while minimizing costs and exposure risks. An approach is currently being developed to reliably estimate subsurface property distributions, including permeability, porosity, and hydraulic conductivity, that control fate and transport of radioactive material. Additionally, such an AI/ML approach could be used to combine information from discrete well and borehole datasets with modeled data (e.g., from flow and transport models or indirect measurements collected through ERT). The resulting predictive models and predictions may capture relationships between these disparate datasets for the purpose of compliance survey design.
Currently, PNNL is applying deep (few shot) learning to such multimodal sensing data, to help investigate subsurface parameters that govern subsurface systems. Findings from current research may provide guidance leading to recommendations about AI/ML methods suitable for inclusion in VSP for characterization and/or compliance phase survey design. Typically, AI/ML methods lack uncertainty estimates and therefore would not be appropriate for hypothesis testing or decision-making purposes. We recommend a scoping study to follow up on current research and other recent methodological developments to identify which AI/ML methods are best suited to characterization or subsurface compliance phase applications, make a recommendation, and plan the path forward to implementation in VSP.
4.3 Dimensionality of Approach Two paradigms can be considered for subsurface compliance surveys:
- Layered Approach: When contaminant distributions are homogeneous within layers of the subsurface, a layered approach may provide simplification to an otherwise complex conceptual site model. A layered approach that considers each individual layer as a decision unit (e.g., with different acceptable limits [DCGLs]and/or exposure pathways) may effectively reduce the 3-D subsurface volume problem to a set of distinct 2-D problems.
- Volume Approach: When spatial dependence exists between layers, layers are not well-defined, there is heterogeneity in the effects of geophysical properties on contaminant fate and transport within and between layers, and/or to ensure or confirm that a layered approach does not underestimate contamination. Simplification of these problems may not be feasible, so sampling and statistical analysis for compliance surveys should consider the entire complex volume as a single decision unit to account for the complex processes governing the end state. Further, data may not be available to model individual layers with rigor and hence a complex volume approach could be the best option due to data volume.
- Section 4.2.1 provides a good example of differences between layered and volume approaches.
4.3.1 Layered Approach A layered approach when layers are distinct, well-defined, and homogeneous with respect to the geophysical properties governing contaminant fate and transport, sampling and statistical analysis of compliance can consider each layer as a separate decision unit, potentially with a unique acceptable limit. Geophysical models and/or dose models can be used to identify layers and acceptable limits (e.g., DCGLs) for each layer. Further, if exposure pathways applicable to each layer are examined, it may be determined that some layers need not be considered in the Statistical Methods for Subsurface Data Analysis 52
PNNL-33647 compliance phase (i.e., when exposure risk is very low or nonexistent). Additional complexity in decision statements will be required, however, to address potential outcomes where data in one layer support decisions that no further remediation is required but support further remediation in another.
This simplified approach that considers each individual layer as a decision unit effectively reduces the 3-D problem back to a set of distinct 2-D problems. Although each layer represents a 3-D volume, the homogeneity of it implies that each observation on the layer can be treated as a point on an x-y plane. This in turn allows surface, or MARSSIM-like methods, to be implemented with each layer. In a MARSSIM-like approach, horizontal layers of the subsurface could be classified as Class I, II, or III, depending on their probability of containing contamination. The horizontal extent of potentially impacted volume(s) would need to be determined. PNNL anticipates that subsurface geostatistical tools like those reviewed in Section 3.4 will help define these boundaries, in addition to the probability of contamination classification. Layers or volumes could also be classified by their contribution to dose. A combination of probabilities and contribution to dose could also be used to define decision units for evaluation against dose-based criteria or action level and areas requiring additional sampling to verify a final status decision.
Considerations for using this approach include:
- A layered approach will ignore spatial dependence between layersit should not be used if non-negligible vertical spatial correlation is present and/or impacts contaminant fate and transport.
- Due to physical constraints of collecting samples via wells and/or boreholes, samples in deeper layers cannot be obtained without also sampling/accessing the layers above them.
Consequently, the total number of sample locations and sample placement for the volume may be governed by the layer with the highest sample size (particularly when a deeper layer requires a higher sample size than a shallower layer).
- Different DCGLs may be applicable to distinct subsurface layers. Decision statements and alternative decisions may differ from layer-to-layer if findings and results of statistical testing differ.
- Risk information must be incorporated into the data analysis when depth discrete data is needed to assess risk. For example, surface or shallow subsurface concentrations may be most limiting for radionuclides dominated by an external dose pathway. Layer concentrations should not be averaged if doing so leads to underestimated concentrations.
- Averaging across the z-coordinate results in a loss of fine-scale variation; possible that an elevated area is missed because it gets averaged out.
- A simple layered approach will therefore only be applicable for sites that model the subsurface as a collection of 2-D layers that meet the following assumptions:
o Contamination is isotropic and data have been collected at one depth per x-y location (e.g., at one depth within each borehole and the same depth in every borehole).
o Contamination is isotropic and data were previously collected at multiple z-locations per x-y location (e.g., multiple data points along the z-axis within a borehole or well at each x-y location) but observations are similar and can be averaged at each x-y location.
Statistical Methods for Subsurface Data Analysis 53
PNNL-33647 o Contamination is isotropic in the x-y direction within one or more 2-D layers and a kriging (prediction) model will be used to account for this anisotropy (a different model could be used within each layer) and previously collected contamination data include one data point per x-y location in the z-direction or multiple isotropic/similar points that can be averaged for each x-y location.
- When there is z-direction anisotropy, a complex 3-D volume approach should be taken (i.e., data should not be averaged at individual x-y locations) where the subsurface is considered as a complex volume and/or comprised of multiple layers, each of which is itself a complex 3-D volume.
As discussed in Fagan 2022, VSP uses the kriging algorithms available through the GSLIB (Deutsch and Journel 1998), which includes ordinary, simple, and indicator kriging. Using GSLIB, VSP can produce maps of kriged 2-D surfaces, conditional variance, interquartile range of predictions, and reference uncertainty index. VSP can also produce contours on kriged surfaces that represent the probability of exceeding a specified threshold or an upper confidence or tolerance limit. In cases where data are collected at multiple z-locations within a layer, the GSLIB algorithms used for 2-D layers could be extended to 3-D volumes by activating the z-coordinate, which is available in GSLIB.
Current implementation of the GSLIB in VSP uses isotropic variogram models. Anisotropic models may be required in a layered approach when spatial dependence in one or more layers is not uniform but directional in the x-y plane. GSLIB includes anisotropic extensions to the isotropic models. Adding these 2-D modeling and hypothesis testing capabilities in VSP would add value for both surface and subsurface approaches where anisotropic contamination distributions are present in the layers.
Current visualization and user interface capabilities in VSP would need to be updated to show independent layers and results simultaneously and enable toggling layers on and off. Compute time may increase as the number of layers and size of data in each increase. Parallel or cloud computing in VSP may be required if it becomes too cumbersome.
4.3.2 Model a Complex 3-D Volume When the subsurface is too complex to model independent layers, a complex approach to statistical sampling and analysis in the subsurface may be required. Such complexity could arise when spatial dependence exists between layers, layers are not well-defined, and/or there is heterogeneity in the effects of geophysical properties on contaminant fate and transport, both within and between layers. Simplification may not be feasible, so sampling and statistical analysis for compliance surveys should consider the entire complex volume as a single decision unit to account for the complex processes governing the end state. Unfortunately, such models and analyses will be more complex to implement, understand, and communicate. Therefore, excellent communication and training should be prioritized in future efforts to add this capability in VSP. Further, VSPs user interface, help files, and automatically generated reports should focus on clear communication of complex topics, for the purpose of communicating within and between practitioners, regulatory bodies, and stakeholders.
Considerations for using a complex volume approach include:
- Spatial dependence between layers cannot be ignored when non-negligible vertical spatial correlation is present and/or impacts contaminant fate and transport.
Statistical Methods for Subsurface Data Analysis 54
PNNL-33647
- The total number of sample locations in the subsurface volume may be governed by the deepest regions of the subsurface. In cases where deeper locations require sampling and data collection, data will need to be collected at the same x-y location but at shallower z-locations due to physical constraints of collecting samples vertically, via wells and/or boreholes.
- Concerns of cross-contamination between depths/layers (similar to the layered approach).
- A single DCGL may apply to the entire subsurface volume. Decision statements and alternative decisions may apply to the entire volume. Alternatively, while a single DCGL could be applicable to the entire vadose zone, individual layers could include DCGLs that cannot be diluted or averaged across the larger volume (e.g., the surface soil layer dominates the dose and surface concentrations should not be averaged with lower concentration residual radioactivity at lower depths). In such cases, both a layered approach and a complex volume approach should be used to ensure that compliance survey design incorporates all relevant information and that all required estimates can be obtained with desired margins of error for statistical decision making.
- Collecting data from all numerous z-locations will result in information about fine-scale variation, lowering the likelihood that an elevated area is missed due to averaging. For groundwater-dependent pathways and large-scale excavation scenarios (e.g., basement excavation or a large construction project), total inventory and volume-wide averages will likely suffice due to mixing of residual radioactivity in a well or due to large-scale soil disturbance. For other intruder scenarios involving disturbance of relatively small soil volumes (e.g., well drilling scenario), maximum concentrations for individual volumes may be more appropriate. Additional research is needed to determine the importance of elevated areas with respect to the likelihood of exposure and dose, and how this depends on the size of the elevated area.
Implementation of the GSLIB in VSP currently only uses the 2-D functionality but could be updated to include 3-D functionality. Similar to the layered approach, anisotropic models will likely be required in a complex volume approach when spatial dependence is not uniform but directional in the x-y plane and in the z-direction. Current visualization capabilities in VSP can show 3-D volumes; however, user interface updates discussed in Section 3.3 would be required.
4.4 VSP Advancements Many capabilities for visualization and preliminary analyses exist in other subsurface visualization and analysis tools, including those detailed in Section 3.2 and other publicly available sources. VSP development should consider integration with such tools, either through the current desktop software or in future migration of VSP to the cloud.
Several statistical methods described above would enhance VSPs benefits in the subsurface.
This report proposes prioritizing methods based on existing VSP capabilities and ease of expanding them from 2-D to 3-D, applicability to the compliance phase, and requirements for new algorithm development.
- High priority: Enhance current 2-D capabilities to 3-D and enable anisotropy
- Enhance 3-D data handling, computation, and user interface
- Enhance 3-D data visualization, including observations, multidirectional variograms, and results of statistical models and analysis Statistical Methods for Subsurface Data Analysis 55
PNNL-33647
- Add anisotropy capabilities (variogram modeling and visualization) to all 2-D kriging algorithms currently included in VSP. 2-D FRK already handles anisotropy and will not need to be updated.
- Update all 2-D kriging capabilities to 3-D (including FRK), including anisotropy.
- All of the above for both layered and complex volume approaches Discussion: VSP already uses GSLIB which supports 3-D kriging and anisotropy, so the VSP updates to enable these features would only require updates to data management and user interface components of VSP and would not require algorithm development.
Enhancements of 3-D visualization can expand on existing 3-D visualization capabilities in VSP. Overall, these improvements can be classified as upgrades and expansions to existing VSP capabilities rather than requiring significant new development of computational functions.
- Medium priority: Develop algorithms and code to support methods with promise for subsurface compliance phase applications:
- GLS
- LISA hot spot detection
- Bayesian Elipgrid
- Bayesian methods
- Expanding current Gaussian and MGK capabilities from 2-D to 3-D.
Discussion: Implementing GLS in VSP can leverage semivariogram calculation and variogram modeling already available in VSP; additional updates would be required to incorporate explanatory variables beyond location coordinates in GLS and/or Bayesian regression models. VSP statisticians and software developers would need to review computational details for the new methods and evaluate which computational methods to add to VSP. Where possible, methods from GSLIB, R, boost, SADA, or other computational libraries would be used.
- Low priority: Develop (or use existing) algorithms and code to support methods with promise for other RSSI phase applications
- Geostatistical simulation including GEM, MrDM, and MrsDM
- AI/ML
- Variogram tomography Discussion: VSP does not currently contain these capabilities. There may be existing SADA code that implements GEM, MrDM, and MrsDM that VSP could leverage to implement these methods. There are numerous powerful and mature commercial off the shelf AI/ML frameworks available, so implementation of such methods would likely be best accomplished by developing external packages or tools and then integrating them with VSP (either through libraries that can be ported into VSP or through establishing an interconnection between VSP and outside tools). Variogram tomography could likely be implemented by expanding current use of the GSLIB library, as could kriging with external drift.
Additionally, several functions in SADA (Version 5) described in the SADA user guide (Stewart et al. 2006) could be added to VSP, in collaboration with SADA developers.
- 3-D visualization capabilities that allow the user to:
Statistical Methods for Subsurface Data Analysis 56
PNNL-33647
- View sampling locations and corresponding data from different angles (i.e., rotating the subsurface volume)
- Visualize interpolated data within the subsurface volume produced by inverse distance weighting interpolation
- Computing and visualizing directional 3-D variograms, which would be especially useful for evaluating the presence of anisotropy in the data
- Search for 3-D hot spots by estimating the probability of discovering a hot spot shape specified by a user-defined ellipsoid (i.e., Bayesian Elipgrid).
To integrate these capabilities into VSP, review of SADA source code will be required to evaluate the feasibility of integration versus re-coding the algorithms. VSP is currently architected to rely on several FORTRAN-based dynamic link libraries (DLLs) including a customized version of GSLIB, so wrapping key SADA functionality into a DLL may be a viable option to provide a technical solution for accessing those capabilities. An alternative approach would be to refactor and reimplement the algorithms in C++ for direct integration into the VSP source code and then use the SADA code for validation and verification to ensure the VSP implemented version is correctly performing the analyses. Selection of which approach to use for which method would be based on code complexity and estimated levels of effort required.
While the visualization capabilities would most likely not need to rely on SADA code for VSP implementation, PNNL would consider visualization in SADA as a guide for best practices and lessons learned, seeking to provide similar or improved functionality and user experience with 3-D visualizations.
Statistical Methods for Subsurface Data Analysis 57
PNNL-33647 5.0 Compliance Survey Design Sampling plans for the compliance survey in subsurface applications should be developed with specific survey objectives in mind. These typically include demonstrating that the potential dose is below the release criterion for each survey unit (e.g., considering the DCGLw) and/or that the potential dose from small areas of elevated activity is below release criteria (e.g., DCGLemc)
(Barr et al. 2022). In surface applications, the Wilcoxon Rank Sum or Sign tests are frequently used to demonstrate DCGLw compliance, and scanning surveys are typically used to demonstrate compliance with DCGLemc. The required coverage of scanning in such surveys is determined by the class designation of the survey unit.
There are several challenges associated with determining the number and location of samples for subsurface applications. An understanding of fate and transport of contaminants from data collection and modeling efforts will be a critical part of compliance survey planning. Lack of scanning survey data and difficulty identifying the boundaries of a subsurface volume represented by a sample, however, present challenges. Another difficulty is that different subsurface layers may require different numbers of samples to achieve desired levels of statistical confidence for declaring that each is in compliance. This may not be a challenge when excavation is used for remediation, but when borehole or other costly sampling methods are required, a method for choosing an appropriate sample size will need to be identified. Finally, there may be different DCGLs for specific depths and thicknesses in the subsurface due to different exposure pathways (NRC 2022).
5.1 Number and Location of Samples As is the case in surface applications, the number of samples in subsurface compliance surveys should be a function of desired levels of statistical confidence and the statistical methodology that will be used to determine if decision-making criteria have been met. The complexity of subsurface contamination and challenges with sampling may preclude the use of the simple nonparametric Wilcoxon Rank Sum or Sign tests, therefore alternative methods will need to be identified.
- Various criteria can be used to determine the number of samples to collect in the subsurface.
- The number of samples can be chosen to minimize prediction error (e.g., mean squared prediction error).
- The number of samples can be chosen to minimize prediction uncertainty (e.g., average kriging standard error).
- The number of samples can be determined based on an assessment of the value of each additional sample, quantified by either prediction error or prediction uncertainty.
The choice of survey locations should be made so that the resulting sample data provide representative information which are sufficient to demonstrate that local areas of elevated activity are unlikely to exist, with high confidence. While previous site activity or decommissioning activities may have established boreholes or wells where samples have previously been collected and can continue to be collected, compliance survey planning (e.g., DQA activities discussed in Section 3.2) should establish whether each sampling location should be used as part of the compliance survey or not. Since nonrandom sampling locations can lead to biased assessments (e.g., if they were preferentially selected in proximity to leaks Compliance Survey Design 58
PNNL-33647 and/or using judgment sampling) (EPRI 2016), additional random samples should be placed to ensure that site estimates and resulting decisions are not biased. Further, local areas of elevated concentrations and classification like that used in the surface MARSSIM approach could vary by layer. Therefore, it will be important to consider designs that provide adequate coverage of each area and layer within the subsurface volume to ensure decisions reflect these conditions.
General guidance for sampling in the subsurface is provided in EPRI (2016). The following highlight main points therein:
- Sample spacing should consider spatial scale (distance between two points at which spatial correlation is negligible).
- If spacing is greater than the spatial scale, then inference will be difficult and uncertainty estimates will not reflect the spatial structure (EPRI 2016).
- If spacing is less than the spatial scale, redundancy can occur and there may be little added value of geostatistics (and additional sample points) (EPRI 2016).
- Edge and/or boundary sampling should be located throughout the sample domain. Low coverage at borders of the domain can lead to increased prediction uncertainty due to extrapolation (EPRI 2016).
- Prior information, although it is valuable and should be used to guide sample designs, if samples are preferentially clustered around regions of concern, declustering methods should be considered that assign weights to data points to mitigate biasing effects (EPRI 2016).
- Dimensionality sampling must be done in 3-D and particular attention should be paid to potential anisotropy (EPRI 2016).
5.1.1 Model Based Survey Design Given the complexity of subsurface investigations and the cost and difficulty of subsurface sampling, we anticipate decisions about compliance to be predicated on model-based inference rather than design-based inference (de Gruijter et al. 2006; Stewart 2011). Models will combine concentration measurements with spatial location information and/or subsurface features (e.g., soil type, transport model output) in a spatially explicit regression framework to predict concentrations and the uncertainty in those estimates at unobserved locations (Lark 2012).
Previously collected measurement data and geophysical models will be important for understanding spatial correlation, spatial variation, and ultimately for identifying areas with high prediction error and uncertainty, which will inform model-based compliance survey designs.
Section 4.0 cataloged numerous geostatistical methods to estimate contamination and uncertainty at unsampled locations. The number and location of samples required to meet survey goals will depend on which of these is the planned approach. The approach, in turn will depend on the conditions of the subsurface (i.e., isotropy, anisotropy) and whether a 2-D layer or 3-D complex volume approach is used. In sites with anisotropy.
There are several methods for optimizing the design of surveys for 2-D applications (Myers 1997). In de Gruijter et al. (2006), authors discuss model-based design approaches for estimating global (e.g., mean) and local (e.g., hot spots) quantities. Zhu and Stein (2006) use simulated annealing or a two-step algorithm for finding sampling designs that provide good point and interval predictions. Krause et al. (2008) use entropy and information gain to optimize sampling locations under the assumption of a Gaussian process. The extension of these Compliance Survey Design 59
PNNL-33647 methods to the subsurface is not immediately obvious because sampling a specific subsurface location may require excavation of, or borehole drilling through, layers above. If excavation is being used, it may be possible to use the 2-D approaches. When excavation is not being used, a possible approach is using the maximum of the sample sizes determined for each layer.
Additional samples may be necessary for verifying the absence of hot spots in certain layers.
Some approaches have been developed for optimizing the placement of additional boreholes in subsurface applications. Morshedy and Memarian (2015) examined this problem for the purpose of minimizing kriging variance, or prediction error. Zhang et al. (2022) examined the effectiveness of different sampling designs on the performance of 3-D soil mapping. Because existing sampling locations are likely to exist, frameworks for determining the number and location of additional sampling locations may be desirable. Judgement samples may also be desired if there are locations that have previously been identified as having elevated measurements.
PNNL recommends a scoping study and report to research and develop model-based survey sample design approaches corresponding to OLS for isotropic cases and to GLS for anisotropic cases, when the study goals include estimating global (e.g., mean) or local quantities (e.g., hot spots). The scoping study would review de Gruijter et al. (2006), Zhu and Stein (2006), and Krause et al. (2008), and related/more recent publications; determine the best option(s) for licensees; and outline a path forward for implementation. A fully Bayesian methodology should also be investigated as part of this scoping study. Particularly, Bayesian model-based approaches that provide a measure of uncertainty in posterior estimates of contamination and that use continuous measurements should be reviewed. Flow and transport model data and geophysical characteristic metrics like those discussed in Section 3.4 could be incorporated into such models through a Bayesian prior distribution and/or through independent variables in a Bayesian regression or kriging model. For example, the approaches described in French and Hoeting (2016), Hewitt et. al. (2019), or Cressie and Suesse (2020) could all be used to identify regions that exceed the DCGLemc.
5.1.2 Markov Bayes and Bayesian Elipgrid SC&A (2022) recommends Bayesian Ellipgrid for initial survey design, before a detailed CCM is available, and Markov Bayes for secondary survey design to indicate where additional data need to be collected. Markov Bayes and Bayesian Elipgrid are not strictly Bayesian methods but could use their conditional probability framework to incorporate soft data described by Stewart et al. (2006) as well as hard data (observations). SC&A proposed using Bayesian Elipgrid as the initial survey design, based on knowledge about hot spot sizes that may exist and then after the data are collected according to this design, using Markov Bayes for secondary survey design.
Bayesian Elipgrid is an extension of the (non-Bayesian) Elipgrid algorithm used to create sampling designs and is intended to detect regions of elevated contamination (SC&A 2022).
Unlike Elipgrid, Bayesian Elipgrid assumes elevated zones exist with user-specified probabilities, specified based on knowledge about the site and/or results of previous RSSI phases. Both Elipgrid and Bayesian Elipgrid can be used to determine the number and location of samples. While they have been applied to surface applications, they may not be appropriate for subsurface applications because both assume a uniform probability that contamination exists. This assumption may not hold in the subsurface, especially with anisotropy expected in subsurface volumes. In these cases, calculations assuming isotropy when in fact anisotropy is Compliance Survey Design 60
PNNL-33647 present, calculations will likely overestimate the number of samples required, resulting in too many samples (Stewart and Powers 2009). Assuming isotropy could also lead to sampling that does not capture, for example, contamination concentrated along transport pathways. Such considerations should be examined when using this method for survey design.
Once data are collected using the initial survey design developed using Bayesian Elipgrid, SC&A (2022) recommends using Markov Bayes to create a probability map (e.g., probability of contamination exceeding a threshold) by combining the prior belief (soft data) with observations collected through the survey (hard data). Markov Bayes is a distribution-free method that honors the observations. It can be used to create point estimates of the probability of elevated contamination in areas or volumes of the subsurface (SC&A 2022). Note that Markov Bayes is not technically a fully Bayesian methodit does not rely on a statistical prior distribution or likelihood function to derive a posterior distribution accounting for distributional assumptions. Markov Bayes, therefore, does not provide a measure of uncertainty in resulting probability maps. Further, information is lost when hard data are converted to zeros and ones via thresholding (SC&A 2022). There is no distinction between, for example, a measurement well below the threshold and a measurement slightly below the threshold, because both are coded as zero. 17 SC&A recommends using the combination of Bayesian Elipgrid and Markov Bayes for the purpose of characterization or initial scoping survey designs early in the RSSI process. We agree that these methods may be appropriate at these stages but note that they would be less appropriate in the compliance phase, primarily because the Markov Bayes results do not include a measure of uncertainty and therefore could not be used to support statistically-based decisions about a site. VSPs current 2-D module for detecting hot spots implements the Davidson (1995) Elipgrid algorithm but does not include the Bayesian Elipgrid or Markov Bayes modules. We recommend adding both in 2-D and in 3-D, which would require new code development. Resources may exist in GSLIB and [possibly] in a SADA development version that has not been released, either or both could be leveraged for VSP development purposes.
5.1.3 Check and Cover NUREG/CR-7021 and SC&A (2022) propose check and cover as a possible metric for choosing the sample size. The goal of this approach is to minimize the maximum sum of value-weighted distances, and it provides a metric to characterize the value of each additional sample (Stewart and Powers 2009). A desirable feature of this approach is that it can incorporate information from the CCM to choose the number and location of samples. A drawback of this approach is that it does not account for uncertainty in the CCM or utilize a metric (e.g., statistical power, kriging variance) related to the decision that will be made with the resulting data. While it can provide a sampling plan that covers much of the site and areas of potentially elevated concentrations, it does not optimize sampling locations with respect to the decision that will be made about the site. Thus, check and cover is most appropriate for sampling campaigns early in the investigation and remediation process (Stewart 2011) but not in the compliance phase. This method should remain a low priority to add to VSP until other methods relevant to the compliance phase have been added.
17 The Markov assumption also needs to be checked when applying this method (Goovaerts 1997; Goovaerts and Journel 1995), a critical and non-trivial assumption.
Compliance Survey Design 61
PNNL-33647 5.2 VSP Advancements New survey design methods required for the subsurface will be required in VSP. Each method will need to be developed and implemented, including new design dialogs where users will enter inputs and subsequent computational results would be shown that determine the required number of samples and sample placement. User input will be required to control how samples should be placed within the 3-D volume (e.g., at specified depths, randomly throughout the subsurface volume, or at multiple depths at sampled x-y locations). Help files and report generators would be required to document the basis for sampling plans and specific parameters that were chosen to arrive at the final survey designs.
Current 3-D placement and visualization methods in VSP (described in Section 3.3) could be used to display the resulting designs, although improvements to data filtering and visualization would be needed to streamline the user experience of working with complex 3-D data and subsequent sampling plans. An immediate need will be to update current sample definition in VSP to include not only the x-y data collection location but also to enable observations at one or more z-coordinates to be associated with the sample location so that different samples can be placed and identified along the vertical direction when well or borehole sampling will be performed.
Implementing new survey design methods in VSP will require additional research to determine the best algorithm to use and development of the code required to implement it. Where possible, methods from GSLIB, R, boost, SADA, or other computational libraries would be used.
Compliance Survey Design 62
PNNL-33647 6.0 Recommendations Table 4 contains the VSP updates recommended for subsurface compliance phase survey sample design and geostatistical analysis.
Recommendations 63
PNNL-33647 Table 4. Recommended Updates to VSP Software for Subsurface Compliance Surveys and Geostatistical Analysis Priority Task ID VSP Improvement Adding the capability to track subsurface data within a single borehole/well using (x, y, z, t) coordinates, that is tracking observations at distinct depths (z locations) and timestamps (t) within 1
one borehole/well located at a single (x, y) coordinate (i.e., back-end software development required to handle 3-D-4-D data).
2 Improve VSP's capability to generate, filter, and track borehole/well labels (i.e., back-end software development required to handle 3-D-4-D data).
1 3-D Data Management 3 Add the capability for unit conversion Ability to combine data and read metadata from disparate sources (including data files and separate metadata files) and sensor platforms and to provide details (e.g., instrument label, field of 4
view, sample matrix) into a single dataset.
5 Ability to capture metadata including uncertainty and/or DQA corresponding to individual datasets or individual data points from disparate sources in a combined dataset.
2 3-D Sample Point Placement Upgrade 6 Update VSP to allow placing sample points throughout the volume will provide value for subsurface decommissioning and the compliance phase of the process.
Scoping report to determine what is required to finalize Lag-k Report: NRC report formatting requirements, address previous technical comments, Information Release; determine if other 7
3 Lag-k Report methods should be considered for this type of data; include a summary and recommendation in Task 1a report.
8 Implementation of lag-k (or alternative recommended) method.
Extend 3-D mapping so that user-defined 3-D and vertical dimension information and data could be used to define subsurface volumes, soil and rock strata boundaries, location and depth of 9 groundwater wells, aquifer locations and dimensions, and layers associated with dose models for subsurface volumes. User-defined information could be imported, entered manually via the data interface.
10 Allow user to specify, visualize, drape multiple individual 2-D layers separately (as planes) within a 3-D volume (e.g., observed, modeled and/or interpolated data).
11 Enable display of 3-D vector layers in 3-D space (e.g., show well screen interval color-coded by value of parameter in 3-D space).
4 3-D Data Visualization & DQA 12 Enable 3-D volume visualization (not as layers but as a complex volume), including observed and/or interpolated data.
Allow users to specify and apply vertical exaggeration to emphasize vertical features that are too small to see relative to the horizontal scale using a factor to specify how much greater the 13 vertical scale is to appear than the horizontal scale (e.g., 5x, 10x, 20x).
Enable the user to control how much of the 3-D volume is displayed and its transparency level (i.e., to display and apply sampling or analysis modules) and to set different transparency levels 14 for each subsurface layer or volume to aid in seeing other features in the model 15 Plots that make potential elevated volumes or 3-D hot spots easier to identify on maps (e.g., iso-contours) 5 LISA Analysis 16 Add LISA capabilities for hot spot detection.
17 Add the capability to identify different instruments, fields of view, etc., using different symbols or color scales on site maps 18 Add the capability to identify different analytes and/or sample matrices (e.g., groundwater, surface water, and soil types) using different symbols or color scales on site maps 6 Enhanced 3-D Data Management 19 Add the capability to denote uncertainty and/or DQA metrics for each dataset or data point using different symbols or color scales on site maps 20 Add the capability to automatically determine optimal symbol/color combinations given the various characteristics in a combined dataset.
Enable the user to slice 3-D volumes vertically and horizontally, either by drawing a slice or providing 2-D coordinate inputs) and visualize data on the sliced plane including display of 21 isocontours 22 Include ability to show isovolumes according to selected parameter values (e.g., volume greater than 5 pCi/g)
Enhanced 3-D Data Visualization 7 Update visualization to include rotatable subsurface volumes and structures including building footprints, piping systems, and tanks as separate data layers that can be turned on/off in the Capabilities 23 VSP interface Provide the ability to subset data interactively (e.g., by selecting data only at specific depths or observed at specific times) for visualization and geospatial analysis. Build on data filtering 24 capability in Analyze Data module so that filtered data is displayed on the map and in the coordinate view Variogram analysis & visualization Add the capability to calculate and view variogram/semivariogram plots in multiple selected directions on same plot, each view providing the estimated semivariogram values at the sill, range, 8 25 (multiple directions, 2-D and 3-D) nugget for a specific direction (in up to four directions)
Recommendations 64
PNNL-33647 Priority Task ID VSP Improvement Add the capability to calculate and view multiple variogram/semivariogram surfaces, where a surface corresponds to the semivariogram values at all distances within a plane (slice) of the 3-D 26 volume, providing information about whether anisotropy exists or not.
Add the capability to display H-scattergrams, or scatterplots with the value of a variable at one location plotted against the value of the same variable at a different location, where each pair 27 of observations are separated by distance h. Provide the capability to view multiple scattergrams at different h-distances to help determine spatial scale.
Update existing 2-D Gaussian and MGK modules to accommodate 3-D data. Update existing 2-D kriging modules (ordinary, simple, indicator kriging) to accommodate 3-D data. Update existing 9 3-D Kriging & Anisotropy 28 2-D kriging modules and new 3-D kriging models to include anisotropy using variograms in multiple directions. (See also #19-20.) Include automatic calculation and fit of variogram models in multiple directions 10 GLS Analysis 29 Add GLS modeling capabilities, for model fitting, estimation, prediction, and hypothesis testing.
Scoping report to document deep-dive literature review, identify method(s), present method to NRC, determine details of tasks for implementing at least one new subsurface model-based 30 survey design method (e.g., based on uncertainty estimates resulting from GLS model or posterior distribution from hierarchical Bayesian model) in VSP (not including implementation).
Enhanced Subsurface Sampling 11 Document details of method, tasks required for implementation, determine if and what existing open-source statistical software code is available, and chart a path forward.
Capability Implement new survey design methods to determine the number and placement of new sample locations in VSP will require additional research to determine the best algorithm for this use 31 and then development of the code required to implement it. (Estimate could change based on scoping report.)
12 3-D Elipgrid 32 Update current Elipgrid algorithm in VSP to accommodate 3-D data Scoping report to identify method, determine details of tasks for implementing at least one Bayesian method in VSP (not including implementation). Document details of method, determine if Bayesian Geostatistical & Sampling 33 13 and what existing open-source statistical software code is available, chart a path forward, present method to NRC.
Method 34 Implement Bayesian method(s) from scoping report.
35 Enable placement of scan lanes along a DEM surface to facilitate scan survey for subsurface surfaces during remediation 36 Add the capability to analyze collection time and location information to calculate velocity and identify significant time gaps or deviations from planned survey velocity.
Enhanced Continuous Data 14 37 Add the capability to analyze survey altitude to evaluate the average altitude and flag or adjust for deviations from the expected survey altitude.
Management Improve hot spot detection and boundary identification by adding alternative formulations of upper tolerance limits (UTLs) and adding upper simultaneous limit (USL) methodology accounting 38 for spatial autocorrelation (i.e., to determine effective sample size of correlated data leading to improved estimation of standard error and confidence bounds).
39 Display scatterplot matrices between multiple variables and/or bivariate statistics (e.g., correlations) to determine whether cokriging or multivariate analysis should be considered.
15 Cokriging in 2-D and 3-D Add cokriging and multivariate analysis functionality in VSP would allow multiple analytes to be analyzed simultaneously, such that the analyses could capitalize on correlations between 40 analytes, particularly useful in cases where data for one analyte are more readily available than another but their behavior/presence/concentrations is correlated.
Scoping report to document details of tasks required to implement GEM (MrDM and MrsDM) in VSP (not including implementation). Collaborate with Rob Stewart to gather details of each 41 method, determine if and what SADA code is available, and chart a path forward.
16 GEM, MrDM, and MrsDM Implement GEM, MrDM, and MrsDM. If resources exist in development but unreleased SADA versions, they could be useful for development purposes. (Estimate could change based on 42 scoping report.)
43 Identify and implement which AI/ML method(s) are best suited to subsurface compliance phase applications.
17 AI/ML for Subsurface Identify AI/ML methods best suited to subsurface in compliance phase survey design and statistical analysis, determine details of tasks required to implement at least one AI/ML method in 44 VSP (not including implementation). Determine if and what software code has already been developed, if and what existing open-source statistical software code is available and chart a path forward. Document details in scoping report. (Estimate could change based on scoping report.)
18 Check & cover for subsurface 45 Add check and cover capability as VSP sampling module 19 2-D and 3-D Bayesian Elipgrid 46 Implement Bayesian Elipgrid in 2-D and 3-D 20 Markov Bayes 47 Implement Markov Bayes methods to create a probability map, including adding functionality for users to combine add prior belief (soft data) with collected observations (hard data).
21 3-D FRK 48 Update existing 2-D FRK module to accommodate 3-D data.
Recommendations 65
PNNL-33647 7.0 References Ancell, E. & Bean, B. (2021). Autocartspatially-aware regression trees for ecological and spatial modeling. Accessed July 15, 2022. arXiv preprint arXiv:2101.08258 Anselin, L. (1995). Local indicators of spatial associationLISA. GEOGRAPHICAL ANALYSIS, 27, 2,93-115.
ArcGIS Pro 2.8. What is empirical Bayesian kriging? ArcGIS Pro Documentation. Accessed July 15, 2022. https://pro.arcgis.com/en/pro-app/2.8/help/analysis/geostatistical-analyst/what-is-empirical-bayesian-kriging-.htm Banerjee, S., Gelfand, A.E., Finley, A.O., Sang, H. (2008). Gaussian predictive process models for large spatial data sets. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 70(4), 825-848.
Banovac, K.L., Buckley, J.T., Johnson, R.L., McCann, G.M., Parrott, J.D., Schmidt, D.W.,
Shepherd, J.C., Smith, T.B., Sobel, P.A., Watson, B.A., Widmayer D.A., Youngblood, T.H.
(2006). Consolidated Decommissioning Guidance Decommissioning Process for Materials Licensees. NUREG-1757, Volume 1, Revision 2. U.S. Nuclear Regulatory Commission, Washington, D.C.
Barr, C.S., Clark, S. Chapman, G.C., Clements, J.P., Esh, D.W. , Fedors, R.W., Huffert, A.M.,
Kauffman, L.A., LaFranzo, M.M., Lowman, D.B., McKenney, C.A., Parks, L.L., Schmidt, D.W.,
Schwartzman, A.L., Watson, B.A. (2022). Consolidated Decommissioning Guidance Characterization, Survey, and Determination of Radiological Criteria. NUREG-1757, Volume 2, Revision 2. U.S. Nuclear Regulatory Commission, Washington, D.C.
Benito, M.B. (2021). spatialRF: Easy spatial regression with random forest. R package version 1, No. 0. https://cran.r-project.org/web/packages/spatialRF/spatialRF.pdf Boisvert, J.B. & Deutsch, C.V. (2011). Programs for kriging and sequential Gaussian simulation with locally varying anisotropy using non-Euclidean distances. Computers and Geosciences, 37(4), 495510. https://doi.org/10.1016/j.cageo.2010.03.021 CH2M HILL Plateau Remediation Company. (2010). Testing Ground Based Geophysical Techniques to Refine Electromagnetic Surveys North of the 300 Area, Hanford, Washington.
Richland, Washington, USA. SGW-47996, Rev. 0.
Cressie, N. (1993). Statistics for spatial data (Revised edition.). Wiley-Interscience Publication.
Cressie, N. & Johannesson, G. (2008). Fixed rank kriging for very large spatial data sets.
Journal of the Royal Statistical Society: Series B (Statistical Methodology), 70(1), 209226.
Cressie, N. & Suesse, T. (2020). Great expectations and even greater exceedances from spatially referenced data. Spatial Statistics, 37. doi:10.1016/j.spasta.2020.100420 Dambon, J., Sigrist, & F., Furrer, R. (2021). varycoef: An R Package for Gaussian Process-based Spatially Varying Coefficient Models. https://doi.org/10.48550/arxiv.2106.02364.
References 66
PNNL-33647 Davidson, J.R. (1995). ELIPGRID-PC: Upgraded Version. ORNL/TM-13103. Oak Ridge National Laboratory, Oak Ridge, TN.
Day-Lewis, F.D., Slater, L.D., Robinson, J., Johnson, C.D., Terry, N., & Werkema, D. (2017).
An overview of geophysical technologies appropriate for characterization and monitoring of fractured rock sites. Journal of Environmental Management, 204(Part 2), 709720.
https://doi.org/10.1016/j.jenvman.2017.04.033 Day-Lewis, F.D., Johnson, C. D., Slater, L. D., Robinson, J. L., Williams, J. H., Boyden, C. L.,
Werkema, D., & Lane, J. W. (2016). A Fractured Rock Geophysical Toolbox Method Selection Tool. Ground Water, 54(3), 315-316.
Myers, D. E. (1999). GSLIB: Geostatistical Software Library and Users Guide, Second Edition.
Computers & Geosciences, 25(3), 309-312.
de Gruijter, J., Bierkens, M., Brus, D., & Knotters, M. (2006). Sampling for natural resource monitoring. Springer. https://doi.org/10.1007/3-540-33161-1 Diggle, P.J. & Ribiero, R.J., Jr. (2010). Model-based Geostatistics. Springer, New York.
DOE. (2011, 2017). Radiation Protection of the Public and the Environment. DOE Order 458.1.
Dong, M., Neukum, C., Hu, H., & Azzam, R. (2015). Real 3-D geotechnical modeling in engineering geology: a case study from the inner city of Aachen, Germany. Bulletin of Engineering Geology and the Environment, 74(2), 281300. https://doi.org/10.1007/s10064-014-0640-6 Du, J. & Ma, C. (2012). Variogram matrix functions for vector random fields with second-order increments. Mathematical Geosciences, 44(4), 411425. https://doi.org/10.1007/s11004-011-9377-y EPA. (2006). Guidance on Systematic Planning Using the Data Quality Objectives Process, EPA QA/G-4. EPA/240/B-06/001. U.S. Environmental Protection Agency, Washington, D.C.
EPA. (2000). Guidance for Data Quality Assessment, EPA QA/G-9. EPA/600/R-96/084. U.S.
Environmental Protection Agency, Washington, D.C.
EPRI. (2016.) Guidance for Using Geostatistics to Develop Site Final Status Survey Program for Plant Decommissioning. EPRI, Palo Alto, CA.
Fagan, D., Obiri, M., Weller, Z., Newburn, L., & Bunn, A. (2022). Recommendations for VSP enhancements for continuously collected survey data. Fulfillment of Task 1a TO 31310021F0022. PNNL-32664.
Freedman, V.L., Truex, M.J., Rockhold, M.L., Bacon, D.H., Freshley, M.D., &Wellman, D.M.
(2017). Elements of complexity in subsurface modeling, exemplified with three case studies.
Hydrogeology Journal, 25(6), 18531870. https://doi.org/10.1007/s10040-017-1564-6 Fokker, K., Ruigrok, E., Hawkins, R., & Trampert, J. (2021). Physics-Based Relationship for Pore Pressure and Vertical Stress Monitoring Using Seismic Velocity Variations. Remote Sensing, 13(14), 2684. https://doi.org/10.3390/rs13142684 References 67
PNNL-33647 French, J.P. & Hoeting, J.A. (2016). Credible regions for exceedance sets of geostatistical data.
Environmetrics, 27(1), 414.
Gelfand A.E., Kim, H.J., Sirmans, C.F., & Banerjee, S. (2003). Spatial Modeling with Spatially Varying Coefficient Processes. Journal of the American Statistical Association, 98(462), 387 0396. doi:10.1198/016214503000170 Gelfand, A.E. & Schliep, E.M. (2016). Spatial statistics and Gaussian processes: A beautiful marriage. Spatial Statistics, 18(Part A), 86104.
Georganos, S., Grippa, T., Gadiaga, A.N., Linard, C., Lennert, M., Vanhuysse, S., Mboga, N.,
Wolff, E., & Kalogirou, S. (2021). Geographical random forests: a spatial extension of the random forest algorithm to address spatial heterogeneity in remote sensing and population modelling. Geocarto International, 36(2), 121-136.
Gogolak, C. (2022). Memo/technical report summarizing Geospatial Extension to MARSSIM (GEM) method presented in Robert Stewart dissertation. Provided by NRC.
Goovaerts, P. (1997). Geostatistics for natural resources evaluation. United Kingdom: Oxford University Press.
Goovaerts, P. (2000). Geostatistical approaches for incorporating elevation into the spatial interpolation of rainfall. Journal of Hydrology, 228(1), 113129.
Gribov, A., & Krivoruchko, K. (2020). Empirical Bayesian kriging implementation and usage.
Science of the Total Environment, 722.
Griffith, D. (2017). Spatial Statistics and Geostatistics: Basic Concepts. In: Shekhar S., Xiong H.
& Zhou X. (eds) Encyclopedia of GIS. Springer. https://doi.org/10.1007/978-3-319-17885-1_1650 Hammett, S., Day-Lewis, F.D., Trottier, B., Barlow, P.M., Briggs, M.A., Delin, G., Harvey, J.W.,
Johnson, C.D., Lane, J.W., Jr., Rosenberry, D.O., & Werkema, D.D. (2022). GW/SW-MST: A Groundwater/Surface-Water Method Selection Tool. Groundwater, 60(6), 784.
https://doi.org/10.1111/gwat.13194 Ham, K.D. & Crockett, D.L. (2021). Evaluation of Data Catalog Software for Hanford Site Environmental Datasets. Pacific Northwest National Laboratory Richland, Washington 99354.
PNNL-31960 Rev 0 and DVZ-RPT-066 Rev 0.
Heaton, M.J., Datta, A., Finley, A.O., Furrer, R., Guinness, J., Guhaniyogi, R., Gerber, F.,
Gramacy, R.B., Hammerling, D., Katzfuss, M., Lindgren, F., Nychka, D.W., Sun, F., & Zammit-Mangion, A. (2019). A case study competition among methods for analyzing large spatial data.
Journal of Agricultural, Biological and Environmental Statistics, 24, 398-425.
https://doi.org/10.1007/s13253-018-00348-w Hewitt, J., Fix, M.J., Hoeting, J.A. & Cooley, D.S. (2019). Improved return level estimation via a weighted likelihood, latent spatial extremes model. Journal of Agricultural, Biological and Environmental Statistics, 24(3), 426443.
References 68
PNNL-33647 Hengl, T. & R.A. MacMillan. (2019). Predictive Soil Mapping with R. OpenGeoHub foundation, Wageningen, the Netherlands. www.soilmapper.org.
Hengl, T., Nussbaum, M., Wright, M.N., Heuvelink, G., & Grler, B. (2018). Random forest as a generic framework for predictive modeling of spatial and spatio-temporal variables. PeerJ, 6, e5518. https://doi.org/10.7717/peerj.5518.
Huckett, J. & D. Fagan. (2022). Statistical Methods for Subsurface Surveys to Support Decommissioning. Presentation delivered at the NRC subsurface soils workshop on May 11, 2022. PNNL-SA-173062.
Jaysaval, P., Shantsev, D., & de la Kethulle de Ryhove, S. (2014). Fast multimodel finite-difference controlled source electromagnetic simulations based on a Schur complement approach. Geophysics, 79(6), E315E327. https://doi.org/10.1190/geo2014-0043.1 Jaysaval, P., Shantsev, D., & de la Kethulle de Ryhove, S. (2015). Efficient 3-D controlled-source electromagnetic modelling using an exponential finite-difference method. Geophysical Journal International, 203(3), 15411574. https://doi.org/10.1093/gji/ggv377 Jaysaval, P., Shantsev, D., de la Kethulle de Ryhove, S., & Bratteland, T. (2016). Fully anisotropic 3-D EM modelling on a Lebedev grid with a multigrid preconditioner. Geophysical Journal International, 207(3), 15541572. https://doi.org/10.1093/gji/ggw352 Jaysaval, P., Robinson, J.L., & Johnson, T.C. (2021). Stratigraphic identification with airborne electromagnetic at the Hanford Site, Washington. Journal of Applied Geophysics, 192.
https://doi.org/10.1016/j.jappgeo.2021.104398 Jaysaval, P. & Johnson, T.C. pGEMINI: Parallel Geophysical Electromagnetic Modeling and Inversion for Natural and Induced sources D Forward modeling for active source, Computational Geosciences, Under Review.
Johnson, T., Strickland, C., Thomle, J., Day-Lewis, F., & Versteeg, R. (2022). Autonomous time-lapse electrical imaging for real-time management of subsurface systems. The Leading Edge 41 (8), 520-528. doi: https://doi.org/10.1190/tle41080520.1 Johnson, T.C., Slater, L.D., Ntarlagiannis, D., Day-Lewis, F.D., & Elwaseif, M. (2012).
Monitoring groundwater-surface water interaction using timeseries and time-frequency analysis of transient three-dimensional electrical resistivity changes. Water Resources Research, 48, 7.
DOI:10.1029/2012WR011893 Johnson T.C., Versteeg, R., Day-Lewis, F.D., Major, W., & Lane, J.W. (2015). Time-Lapse Electrical Geophysical Monitoring of Amendment-Based Biostimulation. Groundwater, 53(6),
920-932. doi:10.1111/gwat.12291 Johnson, T.C., Burghardt, J., Strickland, C., Knox, H., Vermeul, V., White, M., Schwering, P.,
Blankenship, D., Kneafsey, T., & EGS Collab Team18. (2021). 4D proxy imaging of fracture dilation and stress shadowing using electrical resistivity tomography during high pressure injections into a dense rock formation. Journal of Geophysical Research: Solid Earth, 126(11),
e2021JB022298. https://doi.org/10.1029/2021JB022298 18 https://publications.mygeoenergynow.org/grc/1033768.pdf References 69
PNNL-33647 Krause, A., Singh, A., & Guestrin, C. (2008). Near-optimal sensor placements in Gaussian processes: Theory, efficient algorithms and empirical studies. Journal of Machine Learning Research, 9(2), 235-284.
Lark, R.M. (2012). Towards soil geostatistics. Spatial Statistics, 1, 92-99.
https://doi.org/10.1016/j.spasta.2012.02.001 Linneman, D.C., Strickland, C.E., & Mangel, A.R. (2021). Compressional wave velocity and effective stress in unsaturated soil: Potential application for monitoring moisture conditions in vadose zone sediments. Vadose Zone Journal. 20(5). https://doi.org/10.1002/vzj2.20143.
Lillah, M. & Boisvert, J. (2012). Inference of 2-D and 3-D locally varying anisotropy fields for complex geological formations. Ninth International Geostatistics Congress, Oslo, Norway June 11-15.
Ma, C. (2011). A class of variogram matrices for vector random fields in space and/or time.
Mathematical Geosciences, 43, 229242, 2011.
MacCormack, K., Arnaud, E., & Parker, B.L. (2018). Using a multiple variogram approach to improve the accuracy of subsurface geological models. Canadian Journal of Earth Sciences, 55, 786-801.
Mangel, A., Linneman, D., Sprinkle, P., Jaysaval, P., Thomle, J., & Strickland, C. (2022).
Multifrequency electromagnetic geophysical tools for evaluating the hydrologic conditions and performance of evapotranspiration barriers, Journal of Environmental Management, 303.
https://doi.org/10.1016/j.jenvman.2021.114123 Matzke, B.D., Wilson, J.E., Newburn, L.L., Dowson, S.T., Hathaway, J.E., Sego, L.H., Bramer, L.M., & Pulsipher, B.A. (2014). Visual Sample Plan Version 7.0 Users Guide. PNNL-23211, Pacific Northwest National Laboratory, Richland, WA. Available at https://vsp.pnnl.gov.
Metahni, S., Coudert, L., Gloaguen, E., Guemiza, K., Mercier, G., & Blais, J-F. (2019).
Comparison of different interpolation methods and sequential Gaussian simulation to estimate volumes of soil contaminated by As, Cr, Cu, PCP and dioxins/furans. Environmental Pollution, 252(Part A), 409-419. https://doi.org/10.1016/j.envpol.2019.05.122 Morshedy, A.H., & Memarian, H. (2015). A novel algorithm for designing the layout of additional boreholes. Ore Geology Reviews, 67, 34-42.
Myers, J. (1997). Geostatistical error management: quantifying uncertainty for environmental sampling and mapping. John Wiley & Sons.
NRC. (2000). Multi-Agency Radiation Survey and Site Investigation Manual (MARSSIM).
NUREG-1575, Rev. 1; EPA 402-R-97-016, Rev. 1; DOE/EH-0624, Rev.1; U.S. Nuclear Regulatory Commission, Washington, D.C.
Ortiz, J.M. (2020). Introduction to sequential Gaussian simulation, Predictive Geometallurgy and Geostatistics Lab. Queens University, Annual Report, 7-19.
References 70
PNNL-33647 Pannecoucke, L., Le Coz, M., Freulon, X., & de Fouquet, C. (2020). Combining geostatistics and simulations of flow and transport to characterize contamination within the unsaturated zone.
Science of the total environment, 699. doi.org/10.1016/j.scitotenv.2019.134216 R Core Team. (2021). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. https://www.R-project.org/
Revil, A., Binley, A., Mejus, L., Kessouri, P. (2015). Predicting permeability from the characteristic relaxation time and intrinsic formation factor of complex conductivity spectra.
Water Resource Research, 51(8). https://doi.org/10.1002/2015WR017074.
Ribiero, P.J., Jr., & P.J. Diggle. (2001). geoR: A package for geostatistical analysis. R-NEWS 1, 2, 1609-3631.
Robinson, J., Slater, L., Keating, K., Robinson, T., Weller, A., Rose, C., & Parker, B. (2018). On Permeability Prediction From Complex Conductivity Measurements Using Polarization Magnitude and Relaxation Time. Water Resources Research, 54(5), 3436-3452.
https://doi.org/10.1002/2017WR022034 Romary, T., Ors, F., Rivoirard, J., & Deraisme, J. (2015). Unsupervised classification of multivariate geostatistical data: two algorithms. Computers & Geosciences, 85,96-103.
SC&A, Inc. 2022. Guidance on Surveys for Subsurface Radiological Contaminants. White Paper. Prepared for U.S. Nuclear Regulatory Commission Office of Nuclear Regulatory Research.
Schabenberger, O. & Gotway, C.A. (2017). Statistical methods for spatial data analysis: Texts in statistical science. Chapman and Hall/CRC.
Schlather, M., Malinowski, A., Menck, P., Oesting, M., & Strokorb, K. (2015). Analysis, simulation, and prediction of multivariate random fields with package RandomFields. Journal of Statistical Software, 63(8), 1-25.
Silva, V.M. (2021). Estimating the uncertainty of kriging estimates: a practical review and the proposal of two novel approaches. Geociencias UNESP (Universidade Estadual Paulista), 40, 3, 641-650. https://doi.org/10.5016/geociencias.v40i03.15899 St. Clair, J., A. R. Mangel, T. Johnson. In Review. New insights from Legacy Seismic Data Regarding Basalt elevations and Permeability on the Hanford Site. Environmental and Engineering Geoscience.
Stewart, R. (2011). A Geospatial Based Decision Framework for Extending MARSSIM Regulatory Principles into the Subsurface. PhD dissertation. University of Tennessee Knoxville.
Stewart, R. & Powers, G. (2009). A Subsurface Decision Model for Supporting Environmental Compliance. Report for U.S. Nuclear Regulatory Commission. NUREG/CR-7021. U.S. Nuclear Regulatory Commission, Washington, D.C.
Stewart, R.N., Welsch, C., Purucker, T., & Stewart, D. (2006). An Introduction to Spatial Analysis and Decision Assistance (SADA): Environmental Applications for Version 5.
References 71
PNNL-33647 Tao, H., Liao, X., Zhao, D., Gong, X., & Cassidy, D. (2019). Delineation of soil contaminant plumes at a co-contaminated site using BP neural networks and geostatistics. Geoderma, 354, 113878. https://doi.org/10.1016/j.geoderma.2019.07.036 Terry, N., Day-Lewis, F., Robinson, J., Slater, L., Halford, K., Binley, A., Lane, J., Jr., &
Werkema, D. (2017a). The scenario evaluator of electrical resistivity (SEER) survey design tool.
Groundwater, 55(6), 885890. https://doi.org/10.1111/gwat.12522.
Terry, N., Day-Lewis, F., Robinson, J., Slater, L., Halford, K., Binley, A., Lane, J., Jr., &
Werkema, D. (2017b). The Scenario Evaluator for Electrical Resistivity (SEER) Survey Design Tool v1.0: U.S. Geological Survey Provisional Software Release, 01 May 2017.
https://doi.org/10.5066/F7028PQ1.
Truex, M.J., Chronister, G.B., Strickland, C.E., Johnson, C.D., Tartakovsky, G.D., Oostrom, M.,
Clayton, R.E., Johnson, T.C., Freedman, V.L., Rockhold, M.L., Greenwood, W.J., Peterson, J.E., Hubbard, S.S., & Ward, A.L. (2018). Deep Vadose Zone Treatability Test of Soil Desiccation for the Hanford Central Plateau: Final Report. Pacific Northwest National Laboratory, Richland, WA. PNNL-26902.
Safe Drinking Water Act of 1974, as amended. 42 U.S.C. §300f et seq.
Wadoux, A.M.C., Brus, D.J., & Heuvelink, G.B. (2018). Accounting for non-stationary variance in geostatistical mapping of soil properties. Geoderma, 324, 138147.
https://doi.org/10.1016/j.geoderma.2018.03.010 Weller, Z.D. & Hoeting, J.A. (2016). A review of nonparametric hypothesis tests of isotropy properties in spatial data. Statistical Science, 31(3), 305-324. DOI: 10.1214/16-STS547 Wohlberg, B., Tartakovsky, D.M., & Guadagnini, A. (2005). Subsurface characterization with support vector machines. IEEE transactions on geoscience and remote sensing 44(1), 4757.
Zammit-Mangion, A. & Cressie, N. (2021). FRK: An R package for spatial and spatio-temporal prediction with large datasets. Journal of Statistical Software, 98, 148.
Zhang, Y., Saurette, D.D., Huq Easher, T., Wenjun, J. I., Adamchuk, V.I., & Biswas, A. (2022).
Comparison of sampling designs for calibrating digital soil maps at multiple depths. Pedosphere, 32(4), 588-601.
Zhu, Z., & Stein, M.L. (2006). Spatial sampling design for prediction with estimated parameters.
Journal of Agricultural, Biological, and Environmental Statistics, 11(1), 24-44.
References 72
PNNL-33647 Pacific Northwest National Laboratory 902 Battelle Boulevard P.O. Box 999 Richland, WA 99354 1-888-375-PNNL (7665) www.pnnl.gov l www.nrc.gov