\makeatletter
\@ifundefined{HCode}
{\documentclass[screen,CRGEOS,Unicode,biblatex,published]{cedram}
\addbibresource{crgeos20250013.bib}
\newenvironment{Table}{\begin{table}}{\end{table}}
\newenvironment{noXML}{}{}
\let\citep\parencite
\let\citet\textcite
\def\xcitealp#1#2{\citeauthor{#1}, \citelink{#1}{#2}}
\newcommand*{\citelink}[2]{\hyperlink{cite.\therefsection @#1}{#2}}
\def\defcitealias#1#2{}
\let\citetalias\textcite
\let\citepalias\parencite
\def\thead{\noalign{\relax}\hline}
\def\endthead{\noalign{\relax}\hline}
\def\ndash{\text{--}}
\def\tabnote#1{\vskip4pt\parbox{.95\linewidth}{#1}}
\def\botline{\\\hline}
\def\dollar{\$}
\def\tbody{\noalign{\relax}}
\def\tsup#1{$^{{#1}}$}
\def\tsub#1{$_{{#1}}$}
\def\og{\guillemotleft}
\def\fg{\guillemotright}
\let\permil\textperthousand
\RequirePackage{etoolbox}
\usepackage{upgreek}
\usepackage{wasysym}
\usepackage{pifont}
\def\jobid{crgeos20250013}
%\graphicspath{{/tmp/\jobid_figs/web/}}
\graphicspath{{./figures/}}
\newcounter{runlevel}
\usepackage{multirow} 
\def\nrow#1{\@tempcnta #1\relax%
\advance\@tempcnta by 1\relax%
\xdef\lenrow{\the\@tempcnta}}
\def\morerows#1#2{\nrow{#1}\multirow{\lenrow}{*}{#2}}
\let\MakeYrStrItalic\relax
\def\refinput#1{}
\def\back#1{}
\csdef{Seqnsplit}{\\}
\newcommand*{\hyperlinkcite}[1]{\hyper@link{cite}{cite.#1}}
\DOI{10.5802/crgeos.310}
\datereceived{2025-01-06}
\daterevised{2025-09-10}
\dateaccepted{2025-08-18}
\ItHasTeXPublished
}
{
\PassOptionsToPackage{authoryear}{natbib}
\documentclass[crgeos]{article}
\usepackage[T1]{fontenc} 
\usepackage{pifont}
\def\CDRdoi{10.5802/crgeos.310}
 \def\citelink#1#2{\citeyear{#1}}
 \def\xcitealp#1#2{\citealp{#1}}
\makeatletter
\def\CDRsupplementaryTwotypes#1#2{}
 \def\selectlanguage#1{} 
 \def\hyperlinkcite#1#2{\citeyear{#1}}
\let\newline\break
}
\makeatother

\dateposted{2025-10-06}
\begin{document}

\begin{noXML}

\CDRsetmeta{articletype}{review}

\TopicFR{Sismologie (sources et structures)}
\TopicEN{Seismology (sources and structures)}
\TopicFR{Zone critique et socio-\'ecosyst\`emes}
\TopicEN{Critical zone and socio-ecosystems}

\title{Ambient field seismology in critical zone hydrological sciences}

\alttitle{Sismologie de champ ambiant dans les sciences hydrologiques
de la zone critique}

\author{\firstname{Marine A.} \lastname{Denolle}\IsCorresp\CDRorcid{0000-0002-1610-2250}}
\address{Department of Earth and Space Sciences, University of Washington, Seattle, WA, USA}
\email[M. Denolle]{mdenolle@uw.edu}

\author{\firstname{Qibin} \lastname{Shi}\CDRorcid{0000-0002-4211-9187}}
\addressSameAs{1}{Department of Earth and Space Sciences, University of Washington, Seattle, WA, USA}
\email[Q. Shi]{qibins@uw.edu}

\author{\firstname{Tim} \lastname{Clements}\CDRorcid{0000-0001-6632-1796}}
\address{United States Geological Survey, 350 N. Akron Rd. Moffett Field, CA 94035, USA}
\email[T. Clements]{tclements@usgs.gov}

\author{\firstname{Lo\"{i}c} \lastname{Viens}\CDRorcid{0000-0002-7975-7799}}
\address{Los Alamos National Laboratory, P.O. Box 1663 Los Alamos, NM 87545, USA}
\email[L. Viens]{lviens@lanl.gov}

\author{\firstname{Veronica}\nobreakauthor\lastname{Rodriguez-Tribaldos}\CDRorcid{0000-0002-1325-0434}}
\address{GFZ Helmholtz Centre for Geosciences, Telegrafenberg 14473 Potsdam, Germany}
\email[V. Rodriguez-Tribaldos]{verort@gfz-potsdam.de}

\author{\firstname{Fabrice} \lastname{Cotton}\CDRorcid{0000-0002-9242-3996}}
\addressSameAs{4}{GFZ Helmholtz Centre for Geosciences, Telegrafenberg 14473 Potsdam, Germany}
\email[F. Cotton]{fcotton@gfz-potsdam.de}

\shortrunauthors

\keywords{\kwd{Critical zone}
\kwd{Environmental seismology}
\kwd{Ambient noise monitoring}}

\altkeywords{\kwd{Zone critique}
\kwd{Sismologie environmentale}
\kwd{Bruit sismique ambient}}

\thanks{David and Lucile Packard Foundation, LA-UR-25-20095 by Triad
National Security, LLC, U.S. Department of Energy (no.
89233218CNA000001), Los Alamos National Laboratory Center for  Space
and Earth Science Chick Keller Fellowship,  USGS Mendenhall fellowship,
National Science Foundation's Seismological Facility for the
Advancement of Geoscience (SAGE) Award under Cooperative Agreement
EAR-1724509}

\begin{abstract} 
Passive ambient noise monitoring is an emerging tool in environmental
seismology, leveraging the ambient seismic field to assess temporal
variations in shallow subsurface properties. This review focuses on the
potential and challenges of using scattered coda waves from noise
correlation functions to monitor critical zone dynamics. The
sensitivity of seismic velocities to various environmental factors,
including precipitation, snowmelt, atmospheric pressure, and
groundwater fluctuations, underscores the method's versatility. While
coda waves excel in detecting subtle changes due to their scattered
nature, ballistic waves provide higher spatial resolution, albeit with
challenges in source stability. Advances in seismic sensing, including
distributed acoustic sensing and low-cost geophone networks, have
enabled high-resolution monitoring of hydrological processes,
subsurface deformation, and seismic hazards. Integrating seismic data
with hydrological models provides insights into water storage, pore
pressure changes, and soil moisture dynamics. However, limitations in
spatial resolution, calibration with ground truth data, and coupled
effects between environmental factors remain key challenges. This
review emphasizes the importance of interdisciplinary approaches in
refining methodologies, enhancing sensor deployments, and addressing
data gaps. Passive seismic monitoring offers opportunities to
understand critical zone processes and their broader impacts on seismic
hazards and environmental sustainability.
\end{abstract}

\begin{altabstract} 
La surveillance passive du bruit ambiant s'est impos\'ee comme un outil
de transformation de la sismologie environnementale, en tirant parti du
champ sismique ambiant pour \'evaluer les variations temporelles des
propri\'et\'es de la zone critique. Cette \'etude est une revue de la
litt\'erature sur le potentiel et les d\'efis li\'es \`a l'utilisation
du champ sismique ambient pour surveiller la dynamique de la zone
critique. La sensibilit\'e des vitesses sismiques \`a divers facteurs
environnementaux, tels que les pr\'ecipitations, la fonte des neiges,
la pression atmosph\'erique et les fluctuations des eaux souterraines,
souligne la polyvalence de la m\'ethode. La corr\'elation du champ
sismique ambient pr\'esente des opportunitit\'es de mesurer ces
changements en continue. Les ondes coda de la corr\'elation excellent
dans la d\'etection de changements subtils en raison de leur nature
dispers\'ee. Les ondes balistiques de la corr\'elation offrent une plus
grande r\'esolution spatiale, bien qu'elles posent des probl\`emes de
stabilit\'e de la source. Les progr\`es de la d\'etection sismique,
notamment la d\'etection acoustique distribu\'ee et les r\'eseaux de
g\'eophones \`a faible co\^ut, ont permis une surveillance \`a haute
r\'esolution des processus hydrologiques, de la d\'eformation du
sous-sol et des risques sismiques. L'int\'egration des donn\'ees
sismiques aux mod\`eles hydrologiques permet de mieux comprendre le
stockage de l'eau, les variations de la pression interstitielle et la
dynamique de l'humidit\'e du sol. Cependant, les limites de la
r\'esolution spatiale, l'\'etalonnage avec des donn\'ees de terrain et
les effets coupl\'es entre les facteurs environnementaux restent des
d\'efis majeurs. Cet examen met en \'evidence la n\'ecessit\'e
d'approches interdisciplinaires pour affiner les m\'ethodologies,
am\'eliorer les d\'eploiements de capteurs et combler les lacunes en
mati\`ere de donn\'ees. La surveillance sismique passive a de profondes
implications pour la compr\'ehension des processus des zones critiques
et de leurs impacts plus larges sur les ph\'enom\`enes sismiques.
\end{altabstract} 

\maketitle

\vspace*{3pc}

\twocolumngrid

\end{noXML}

\defcitealias{Solazzi21}{ibid.}
\defcitealias{Lu21}{ibid.}
\defcitealias{Schulz18}{ibid.}

\section{Introduction}\label{sec1}

\subsection{Broad overview}\label{sec1.1}
The Critical Zone (CZ) is Earth's near-surface region where freshwater,
air, soil, rock, and living organisms interact, extending from the
vegetation canopy at the surface down to freely circulating groundwater
stored in aquifers \citep{giardino2015introduction}. Monitoring soil
moisture in the CZ and understanding its dynamics is critical to
improving agricultural management, monitoring droughts, floods, and
groundwater recharge, and monitoring the stability of slope failures,
and thus has implications for food security, climate systems, and
surface hazards (floods, landslides, and earthquake-induced
liquefaction).  The spatial heterogeneity of soils' properties (e.g.,
hydraulic conductivity) and their potential for evapotranspiration,
geological structure, meteorological forcing (e.g., precipitation,
wind, and radiation), and vegetation control the spatial heterogeneity
of soil moisture \citep[e.g.,][]{Chaney13}. Conventional soil moisture
monitoring is through  (1)~in situ soil moisture probes
\citep{su2014critical} that have local sensitivity or  (2)~remote or
airborne data that use electromagnetic wave signals (amplitudes and
phase) to estimate soil moisture and have regional or global
sensitivity \citep{babaeian2019ground}. As stated by
\citet{Flinchum22}, ``a major challenge in CZ science is to obtain
measurements that have both high spatial resolution and cover large 
areas''. 

Passive ambient noise monitoring refers to the methodology of using
repeated cross-correlation of the ambient seismic field, recorded at a
single or a pair of seismic stations, to estimate the changes in
seismic properties in the subsurface. We focus this review on the
scattered waves in the coda of coherent ambient seismic waves, which
are typically surface waves in the context of soils and sediments.
``Ambient noise monitoring'' is transforming environmental seismology
as it provides new passive methods to monitor the mechanical response
of the shallow subsurface
\citep[e.g.,][]{larose2015environmental,LeBreton2021LandslideApplications,wang20}. Temporal variations in 
near-surface seismic properties stem from a range of interacting
factors whose relative importance is frequently site-dependent,
underscoring the spatial heterogeneity of soils. This review sets out
to build a unified framework for monitoring hydromechanical behavior in
the critical zone with 
\mbox{ambient-noise} seismology, examining the method's
strengths and 
\mbox{limitations.}

Relative changes in shallow seismic velocities, $\mathrm{d}v/v$
hereafter, have been observed to vary between 0.1\% and 50\% under many
tectonic and volcanic factors
\citep[e.g.,][]{wegler2007fault,Brenguier08a,brenguier2008towards,obermann2014seismic,yukutake2016determination,feng2020detecting,cabrera2023spatio,donaldson2019crustal,okubo2024monitoring,Esfahani24} 
and solid Earth tides
\cite[e.g.,][]{sens2019probing,mao2019high,takano2023tidal,Kramer23}.
In addition to the rise in analyses monitoring temporal changes in the
subsurface induced by earthquakes or volcanic activity, there are
numerous environmental factors that change seismic properties of the
subsurface: varying snow load 
\citep{wang2017seasonal,donaldson2019crustal}, rainfall load and
diffusion
\citep{wang2017seasonal,andajani2020spatial,feng2021controls,Clements23}, 
atmospheric pressure fluctuations
\citep{gradon2021seismic,Kramer23,Kramer2025}, river and lake levels
\citep[e.g.,][]{RodriguezTribaldos2021AquiferFiber}, thawing and
freezing of soils \citep[e.g.,][]{james2019insights}. 
\mbox{Figure~\ref{fig:scales}} represents the vast range of time and spatial
scales at which the hydrological phenomena have been reported to carry
seismic signatures in the literature reviewed in this\break manuscript.

\begin{figure*}
\vspace*{-2pt}
\includegraphics{fig01}
\vspace*{-2pt}
\caption{\label{fig:scales}Range of spatial and time scales of the
hydrological processes reported in the literature. Left: lateral
spatial extent ranging from centimeters to 1000~km. Right: depth
scale from centimeters to 100~m. 
Solid boxes represent the range
of time and spatial scales where seismic signatures correlate with
specific hydrological processes. Long-term climate effects have been
reported at all spatial scales and depths.}
\vspace*{-4pt}
\end{figure*}

The core principle of ambient noise monitoring is that repeated
measurements of the seismic wavefield, measured through the
cross-correlation of the ambient seismic field between channels and
stations, capture the seismic signatures of changes in effective stress
or pore pressure. Coda waves extracted from ``noise correlation
functions'' (NCFs) are commonly used in ambient noise monitoring, as
they are less affected by variations in noise sources than direct waves
\citep[e.g.,][]{hadziioannou2009stability,colombi2014temporal}.
Recently, using the direct, or ballistic, waves of the NCF, typically
either P-waves or surface waves
\citep[e.g.,][]{Garambois19,Takano19,Mordret20,Brenguier20,Sheng2024},
or specific scattered phases such as P--S \citep[e.g.,][]{Lu21}, has
helped improve the 
\mbox{spatial}
\mbox{localization} of velocity changes. However,
the use of ballistic waves for monitoring requires some knowledge about
stability and location of repeated anthropogenic sources 
\citep[e.g., trains,][]{Higueret25}. In addition, coda waves are more sensitive than
ballistic waves to small changes in $\mathrm{d}v/v$, as they sample the
medium many more times due to multiple scattering 
\citep{snieder2002coda,colombi2014temporal}. With the assumption that
the scatterers themselves do not change positions, the perturbation of
the seismic velocity averaged along the coda wave path is linearly
proportional to the accumulated phase shift: 
$\mathrm{d}v/v = - \mathrm{d}t/t$,
where $t$ is the phase lag and $\mathrm{d}t$ is the phase shift in the coda
\citep{poupinet1984monitoring}. This review focuses on the observations
made with coda waves but will highlight the potential for ballistic
wave analysis in relevant cases. Beyond resource monitoring, many of
the most damaging natural hazards, such as landslides,
earthquake-induced liquefaction, and nonlinear site amplification, are
also governed by hydromechanical processes in the Critical Zone. These
hazards arise from the same physical mechanisms discussed throughout
this 
\mbox{review:} changes in effective stress, pore pressure, saturation,
and soil structure. In particular, the response of near-surface soils
to strong shaking is modulated by water content and stress history,
linking seismic hazard to hydrological variability. By framing these
hazards as manifestations of CZ dynamics, ambient field seismology
offers a unifying observational 
\mbox{approach:} the same seismic observables
(e.g., $\mathrm{d}v/v$, attenuation ($Q$), and the
Horizontal-to-Vertical Sprectral Ratio (HVSR)) used to track
groundwater or vadose (unsaturated) processes can also be used to
detect precursors to slope failures, assess site amplification
conditions, and evaluate liquefaction susceptibility. This review aims
to bridge these application areas and demonstrate how environmental
seismology contributes to an integrated understanding of both
hydrogeophysical processes and geohazards in the CZ.\looseness=-1

\subsection{Seismic structure of the near surface}\label{sec1.2}
The seismic structure of the shallow subsurface, especially in the
critical zone, is marked by a surface low-velocity zone with shear
wavespeeds ($V_{\mathrm{S}}$) 
\mbox{ranging} from 130--1000~m/s spanning unsaturated and
saturated soil, densities ranging from 1.7--2~g/cm$^3$, and P wavespeed
($V_{\mathrm{P}}$) of 350--2500~m/s
\citep[e.g.,][]{Befus2011,Boore2016,Nagashima21}. The shallowest layers
are often highly weathered and have a greater sensitivity to
hydrological processes at and above the water table depth than more
consolidated sediments \citep{Anderson13,Okay23}. Scattering
attenuation in sediments, generated by higher porosity and crack
density in shallow sediments, is often equal to or greater than
intrinsic attenuation, which is dispersed as heat
\citep{Jemberie05,Eulenfeld16,LIN2023}. Combined intrinsic anelastic
and scattering effects may yield strong attenuation in shallow
sediments, with total quality factors $(Q)$ as low as  2--10 for S and P
waves in shallow soils \citep{Malagnini96}.

\section{Relating seismic properties to water content}\label{sec2}
Relating seismic measurements into hydrological processes requires
understanding three topics:  (1)~\textbf{hydrological data} and, when
lacking measurements, how seismologists infer it;  (2)~\textbf{rock
physics models} to estimate elastic and anelastic properties with
varying water content;  (3)~\textbf{seismic wave spatial and depth
sensitivities} depending on the wave types.\looseness=-1

\subsection{Hydrological models}\label{sec:hydro_model} %2.1
Water content, either in the form of water saturation or pore pressure,
is a key element in the hydromechanical response of the critical zone.
Due to the logistical challenges of obtaining high-resolution in situ
measurements of water content and the resulting sparsity of data
coverage, numerical and process-based models are often used to
understand hydrological processes within the critical zone. The water
budget includes sources (precipitation, snow melt, and runoff from
uphill) and sinks (drainage, runoff, evapotranspiration).

The main approach that seismologists have used to estimate pore
pressure in the fully saturated regime is through coupling infiltration
from rainfall or snow melt with simplified outflow models based on
either baseflow \citep{Sens-Schonfelder2006,feng2021controls} or the
poroelastic response to surface water input
\citep{Roeloffs96,roeloffs1998persistent,Talwani2007,Tsai2011,andajani2020spatial,wang2017seasonal}. 
For regions with a significant soil moisture reservoir, the antecedent
precipitation index has been used to represent the capacity of the
vadose zone \citep{illien2021subsurface}. \citet{Clements23} found that
baseflow and poroelastic-based models of pore pressure, along with the
cumulative deviation from the moving mean of precipitation
\citep{Smail2019CumulativeWisconsin}, equally worked as long-term
proxies of groundwater level when compared to hydrologically-induced
$\mathrm{d}v/v$.

Hydrological models that account for the full water budget with inputs
and outputs (runoff, downward drainage, and evapotranspiration) are
more conducive to modeling vadose zone dynamics. \citet{Lecocq2017}
modeled water storage as a budget between input (precipitation) and
output (river runoff and evapotranspiration) over a time period of
three decades. \citet{shen2024fiber} modeled water storage as a mass
balance between precipitation (input), evapotranspiration modulated by
surface moisture measurements and temperature (output), and deep water
drainage (i.e.\ groundwater downward draining as an output).
Additionally, groundwater levels can be inferred from river
\citep{RodriguezTribaldos2021AquiferFiber} or lake \citep{Clements23}
levels,  as surface body water levels are often hydraulically connected
to the groundwater table.\looseness=-1

\subsection{Hydromechanical models of Earth materials}\label{sec2.2}
Theoretical, laboratory, and field studies have demonstrated that both
$V_{\mathrm{P}}$ and $V_{\mathrm{S}}$ vary with water saturation in
rocks \citep{nur1969effect, oconnell1974seismic, kuster1974velocity}.
In partially saturated rocks, seismic velocity perturbations caused by
changes in water content are controlled mainly by the density effect.
Increasing bulk density through water saturation, achieved by replacing
air with water, results in a gradual reduction in both $V_{\mathrm{P}}$
and $V_{\mathrm{S}}$. When the rock achieves nearly fully saturated
conditions (80\%--90\%), the non-compressibility of the fluid in the
pores dominates the density effect, and $V_{\mathrm{P}}$ experiences a
sharp increase. Commonly, $V_{\mathrm{S}}$ continues to decrease with
saturation. In a few laboratory cases, however, studies have shown that
a small increase in $V_{\mathrm{S}}$ can also be observed when
saturation is close to 100\%, which may be due to the effect of surface
tension in small pores causing an increase in the shear rigidity of the
rock \citep{Knight90}. This increase in $V_{\mathrm{P}}$ and
$V_{\mathrm{S}}$ at full saturation will mostly be observed near the
water table and the capillary fringe. Beyond this general behavior,
studies have shown that other hydrological processes can influence the
behavior of seismic velocities in the unsaturated zone (e.g., see below
discussion of Section~\ref{sec:vadose}). For example, strong matric
suction effects can lower the effective stress and cause a consequent
increase in $V_{\mathrm{P}}$, which competes against the
$V_{\mathrm{P}}$ reduction caused by the density effect
\citep{Linneman21}. 

Below the water table, under \textbf{fully saturated} conditions, pore
structure controls the relative change in seismic velocity. For
low-porosity rocks with crack-shaped pores, $V_{\mathrm{P}}$ increases
with full saturation, whereas $V_{\mathrm{S}}$ remains unchanged with
saturation \citep{nur1969effect}. For high porosity rocks with circular
pores, $V_{\mathrm{P}}$ increases with saturation whereas
$V_{\mathrm{S}}$  decreases with saturation
\citep{toksoz1976velocities}. At shallow depths, where confining
pressure is low, increases in pore pressure thus increase
$V_{\mathrm{P}}$ and decrease $V_{\mathrm{S}}$. At the near-surface,
Rayleigh waves are far more sensitive to variations in $V_{\mathrm{S}}$
than $V_{\mathrm{P}}$ \citep{song1989sensitivity}. Thus, numerous
observations have been made of a decrease in Rayleigh wave speed
($V_{\mathrm{S}}$) with increasing groundwater levels or pore pressure.
The constitutive relations for nonlinear elasticity show that relative
change in elastic wavespeeds, $\mathrm{d}v/v$, is proportional to the
dilatational strain $\epsilon_{kk}$
\citep{Murnaghan1937FiniteSolid,Ostrovsky2001DynamicGeomaterials}, and
observations suggest that:
{\begin{equation}
\mathrm{d}v/v = -\beta \epsilon_{kk},  
\end{equation}}\unskip
where $\beta$ is in the range 
$1\times10^3$--$1\times10^{4}$ (see compilation of values by
\citet{Clements23} and a more general formulation by \citet{Wang24}).
Because pore pressure induces dilatational strains through poroelastic
effects, hydrologists have found linear relationships between
dilatational strains and hydraulic head ($S_{sk} \approx 1\times10^{-5}$,
\citet{Burbey2001StressStrainCharacterization}),  $\epsilon_{kk} =
S_{sk} \,\mathrm{d}h$, such as 
{\begin{equation}
\mathrm{d}v/v = - \beta s_{sk} \,\mathrm{d}h.
\end{equation}}\unskip
These relationships give a sensitivity of $\mathrm{d}v/v$ to water
levels such that local wavespeed changes by 1--10\% for 1~m of
water change. \citet{Rivet15} also proposes the alternative
relationship that 
{\begin{equation}
\mathrm{d}v/v= \mathrm{cov}(\mathrm{d}v/v,p)/\mathrm{var}(p) , 
\end{equation}}\unskip
where $p$ is the poroelastic pressure response to a point source of
rainfall. This relation has been used in several studies comparing
temporal changes in $\mathrm{d}v/v$ and near-surface
\citep{Delouche23,Barajas21}.

In the \textbf{partially saturated} regime, above the water table,
partial saturation effects on seismic wavespeeds are less trivial. The
first-order effect is the change in density that occurs by replacing
the air in the pores with water:
{\begin{equation}
\rho = (1 - \phi) \rho_s + \phi 
( S_w \rho_w + (1- S_w)\rho_a ),
\end{equation}}\unskip
where $\phi$ is porosity, the densities of the host rock, water, and
air are $\rho_s$,  $\rho_w$,  $\rho_a$, and saturation in the pores is
$S_w$. Capillary forces also effectively increase the rigidity and bulk
modulus of the material, such that effective $V_{\mathrm{P}}$ and
$V_{\mathrm{S}}$ increase in the vadose zone \citep{Solazzi21}.
Additionally, \citet{Sakaki10} demonstrated that water retention
creates a dynamic sensitivity to capillary pressure, which is
proportional to the rate of saturation change. This results in
increased capillary pressure during drying cycles and decreased
capillary pressure during wetting cycles. With increased saturation
with depth, the competing factors of overburden pressure with static
capillary forces yield a relatively constant
$V_{\mathrm{P}}/V_{\mathrm{S}}$ ratio down to the water table depth
\citep{Solazzi21}, below which $V_{\mathrm{S}}$ remains low in fully
saturated media, $V_{\mathrm{P}}$ jumps to high values greater than the
sound speed in the water, yielding a high
$V_{\mathrm{P}}/V_{\mathrm{S}}$ ratio, high Poisson's ratio and the
sensitivity of Rayleigh wavespeed discussed above. \citetalias{Solazzi21}
found a particularly strong sensitivity of Rayleigh wavespeed to the
partial saturation in medium-to fine-grained soils (see 
Figure~\ref{fig:profile}), though larger pore sizes in concrete samples
exhibited increased sensitivity of $\mathrm{d}v/v$ to 
saturation \citep{diewald2024impact}. 
There is  {empirical} evidence that such sensitivity remains
relatively constant for the middle range of saturation (e.g.,
30--80\%). Still, dynamic capillary effects are much stronger at lower
and higher saturation levels \citep{Dong16,diewald2024impact,Sakaki10}.
\looseness=-1

\begin{figure*}
\includegraphics{fig02}
\vspace*{-2pt}
\caption{\label{fig:profile}Generic seismic profile of a critical zone
from \citet{Solazzi21}. Left plot: saturation increases from the
surface until the water table (left). Middle plot: elastic physical
properties (density $\rho$, P wave velocity $V_{\mathrm{P}}$, S wave
velocity $V_{\mathrm{S}}$) as a function of depth. Right plot: Surface
wave dispersion as a function of frequency with (blue lines) and
without (red dashed lines) capillary forces included in the model.
\citet{Solazzi21} showed that when the $V_{\mathrm{P}}/V_{\mathrm{S}}$
ratio, which is usually a metric for water content, is about constant
in the partially saturated zone, surface wave velocities are sensitive
to saturation due to capillary suction effects.}
\vspace*{-2pt}
\end{figure*}

It should be noted that $V_{\mathrm{S}}$'s sensitivity to water level greatly
depends on the lithology: it is strong for clays but small for sand
\citep{Dong16}. Such a difference may explain why most seismological
studies exhibit spatially variable sensitivity. 

Changes in seismic attenuation, $Q$, in the partially saturated zone
are less obvious. P-wave attenuation is greatly sensitive to partial
saturation due to the contrast in the compressibility of mixed fluids
(air vs.\ water) \citep[e.g.,][]{Chapman21,Amalokwu14}, whereas S waves
are insensitive to compressible fluids. Numerous studies have
correlated values of $Q_p/Q_s \leq 1$
with the presence of air/gas in a
partially saturated medium \citep[e.g.,][]{Hudson23}. Nevertheless,
theoretical and empirical  \citep{winkler1982seismic,diewald2024impact}
studies have shown that S-wave or S-wave-dominated seismic attenuation
varies with saturation. Seismic attenuation within partially saturated
porous materials is mainly caused by the relative movement between the
solid matrix and the pore fluid. This movement generates shear stresses
in the fluid, leading to energy loss through viscous friction
\citep{mavko1979wave}. 

\subsection{Wave types and their various spatial and depth
scales}\label{sec2.3}
Direct ballistic waves can be used to probe large changes in the
subsurface, as phase shifts in the direct waves may be precise enough
to estimate the seismic velocity changes. By using direct waves, the
localization of the velocity perturbation is achieved simply through
time-lapse tomography. Using ambient noise seismology, studies have
focused on time-lapse shear wave tomography using surface wave
dispersion curves
\citep[e.g.,][]{dou2017distributed,Mordret20,Cheng22,Sun25}. The
typical frequencies used in these analyses range from 2 to 15~Hz,
providing sensitivity to the 
\mbox{shallow} 
\mbox{seismic} structure from a few
meters to approximately 
\mbox{40--50~m} in depth. The spatial resolution
of these measurements is typically on the order of tens of meters.
Ballistic wave-based tomograpy is limited by the density of the seismic
channels and the quality of the high-frequency signals. Besides surface
waves, several studies have leveraged distinct seismic phases, such as
refracted P waves \citep{Brenguier20} or converted phases \citep{Lu21},
when possible. One key element of using ballistic waves is that their
spatial and depth sensitivity depends on the source-receiver geometry
and the Earth's structure \citepalias{Lu21}. Because time-lapse tomography,
which infers absolute seismic velocities in space (laterally and in
depth) and time, rock physics models can extract temporal and spatial
variations in hydro-thermo-mechanical properties of the critical zone
\citep[e.g.,][]{Sun25}. 

When distinct seismic phases are undistinguishable or when the changes
are subtle, seismologists use coda waves instead to image changes in
the subsurface. Because scattered waves accumulate travel-time
perturbations due to changes in seismic velocities, phase shifts are
increased in late coda waves. For homogeneous velocity 
\mbox{perturbations,} such \mbox{delays~accumulate}
linearly regardless of the wave type involved
\citep[e.g.,][]{obermann2013imaging}. 

The depth sensitivity of coda waves is most often assumed to be that of
fundamental mode Rayleigh wavespeed, given a perturbation in $V_{\mathrm{S}}$ as a
function of depth and frequency. Sensitivity kernels, which link
observed wavespeed perturbations ($\mathrm{d}v/v$) to true $V_{\mathrm{S}}$
perturbations between seismic stations often rely on solutions from the
radiative transfer model for short lag-time measurements
\citep{Pacheco2005,Planes14}.  A good rule of thumb is that the
penetration depth of the Rayleigh wave is around one horizontal
wavelength \citep{Barajas22}. When seismic velocity and density are
relatively constant with depth, scattered body waves dominate the late
coda of NCFs, though  sensitivity kernels may couple surface and body
wave sensitivity at depths below the critical zone 
\citep{obermann2013imaging,Obermann2016Lapse-time-dependentMedia,obermann2019,Margerin19,Barajas22}.

\section{Review of sensing technologies}\label{sec3}
We now focus on the technologies that have served the passive
seismological community to investigate critical zone hydrology.

\subsection{Seismological observation networks}\label{sec3.1}

\subsubsection{Permanent seismic networks}\label{sec3.1.1}
The most commonly used sensors for ambient noise monitoring of the
critical zone are broadband seismometers, typically deployed directly
at the surface. Broadband seismometers have been deployed since the
1960s, although continuous digital recording began in the late 1990s.
By 2002, most seismic networks had begun continuous recording of their
entire broadband networks, which fueled the use of the ambient seismic
field for seismological research \citep[e.g.,][]{Campillo03,shapiro05}.
The deployment of broadband sensors serves two purposes. First, the
global detection and monitoring of explosions and earthquakes drive the
deployment of these sophisticated stations worldwide, particularly on
oceanic islands. With this purpose in mind, broadband sensors typically
record at high gain to retrieve weak seismic waves from distant
earthquakes and explosions. Second, regional and national seismic
networks, with the specific goal of monitoring earthquake sources, have
favored the deployment of broadband seismometers on bedrock and near
fault zones. Because seismic waves amplify in soft sediments, network
engineers have chosen to avoid these sites and focus on the station's
sensitivity to fault zone waves. These locations are not ideal for
subsurface hydrology, as the groundwater reservoirs and soils are
predominantly located in basins. However, the advantages of broadband
stations are their long-standing, multi-decadal, continuous records.
Because their sensitivity spans a broad range of frequencies,
seismologists have repurposed broadband seismometers to probe the
near-surface in the 0.1--10~Hz frequency range, covering both powerful
microseismic ambient fields and stationary anthropogenic seismic noise
\citep[e.g.,][]{hasselmann1963statistical,Diaz16}. 

Accelerometers are also used for permanent seismic monitoring and
typically record at lower digital gain to retrieve high fidelity of
strong ground accelerations and higher frequencies (above 0.1~Hz).
Regional seismic networks (e.g., the Southern California Seismic
Network SCSN, the Northern California Earthquake Data Center NCEDC, and
the Pacific Northwest Seismic Network PNSN) comprise several types of
seismometers, ranging from broadband seismometers to accelerometers,
for both monitoring microseismicity and early warning of strong
earthquakes. Regional and global networks of seismometers (e.g., GSN,
ORFEUS, Geoscope) archive their data on publicly accessible data
servers, either on-premise servers or in cloud storage
\citep[e.g.,][]{Yu21}. Most of these installations are made at the
surface only.

Several observatories \textit{pair surface with borehole
instrumentation} enabling top and bottom analysis of the near-surface
environment. After the 1995 Hyogo-ken Nanbu Earthquake, NIED (National
Research Institute for Earth Science and Disaster Prevention)
constructed an uphole/downhole observation network, KiK-net (Kiban
Kyoshin network), with approximately 700 stations
\citep{aoi2004strong}. Each KiK-net station has a borehole of 100~m or
more at depth, and accelerometers have been installed both on the
ground surface and at the bottom of the boreholes.  Such pairs of
sensors 
\mbox{enable~direct} measurement of near-surface dynamics or physical
properties between the two sensors. The KiK-net accelerometers in Japan
have been unrivaled by other downhole networks in the volume of
available acceleration seismograms, given the quality, homogeneity of
the information describing subsurface layers, and the rate of
seismicity. In addition to national observatories, another model for
focused surface-to-borehole \mbox{observations} exists in the geotechnical
community, referred to as ``liquefaction arrays''. In the United
States, such arrays include the Wildlife Liquefaction Array
\citep{Holzer07} and the ``Seattle Liquefaction Array''
\citep{steidl2017instrumented}, which feature piezoelectric pore
pressure sensors, along with several shallow borehole seismometers, to
continuously sample water pressure and earth vibration.

Near-surface structure information of the site on which these seismic
observations are made is sometimes available. For instance,  the
velocity profiles and geological information, as well as the observed
records, are widely accessible on the NIED website. 
\citet{zhu2021open} created an open-source site database that contains
site 
\mbox{characterization} parameters directly derived from available
velocity profiles, including average wave velocities, bedrock depths,
and velocity contrast. Additionally, the site database includes
topographic and geological proxies inferred from regional models or
maps. Geotechnical arrays also typically conduct surveys of the seismic
and hydrological structure \citep[e.g.,][]{youd2004geotechnical}. 

\subsubsection{Temporary seismic arrays}\label{sec3.1.2}
Most seismic imaging of the near-surface entails a small aperture array
of a few 100~m with seismic channel spacing of a few meters. When
analyzing deeper and larger groundwater reservoirs, such as aquifers in
sedimentary basins, the analysis aperture may range from 10~km to
100~km. Given these spatial considerations, seismologists use direct or
refracted P waves \citep[e.g.,][]{Flinchum22}, direct surface waves
\cite[e.g.,][]{Sobolevskaia24,Sun25}, and coda waves dominated by
multi-scattered S and surface waves \citep[e.g.,][]{shen2024fiber} to
explore the spatial and temporal variations in the CZ seismic
properties.

For temporary seismic deployments, low-cost geophones, called ``nodes''
such as the FairField 3C or the SmartSolo nodes, have also been used
for short-term experiments to image and monitor water table,
precipitation, soil moisture, or seasonal frozen ground freeze--thaw
cycles using seismic methods 
\citep{oakley2021seismic,hua2023high,gochenour2024ambient,liu2025quantitatively,Behm19,hua2023high,Grobbe21}. 
Nodes are easy to deploy in  ``large-N'' arrays, fast to install in
rapid response to earthquakes or sudden environmental events; their
disadvantage is a short battery life of a month that limits long-term
monitoring.

A novel technology transforming the field of seismology, particularly
in near-surface seismic imaging and monitoring, is distributed acoustic
sensing (DAS). DAS technology transforms traditional fiber-optic cables
into dense arrays of sensors, enabling the measurement of distributed
strain or strain-rate changes along their length
\citep[e.g.,][]{hartog2017introduction}. Current DAS technology can
continuously record the seismic wavefield at a spatial resolution of a
few meters along cable lengths that can exceed 100~km. Several studies
have also demonstrated that DAS can record strain-rate signals over a
broad range of frequencies, ranging from 0.01 to hundreds of Hertz
\citep{lindsey2019illuminating,Paitz21}, though near-source, high
ground motion observations are often  ``clipped'' due to DAS's limited
sensitivity at high strain rates \citep{lin2025analysis}. The primary
advantage of DAS lies in its capacity for continuous, real-time
monitoring across extensive areas or structures sampling the wavefield
over many closely spaced (${\sim}$3--10~m) virtual channels using
either pre-existing telecommunication fiber-optic cables or short-term
fiber deployments. DAS measures strain rate along the axial direction
of the cable, making it especially suitable for ambient noise studies,
in which Rayleigh waves traveling along the axis of the cable can be
analyzed for subsurface imaging and monitoring.  DAS arrays produce
substantial amounts of raw strain data (${\sim} TB$s), which makes
storing, sharing, and analyzing the data computationally difficult. 

\subsection{Hydrological observations networks}\label{sec3.2}
Seismologists have attempted to explain the hydrological signature by
observing time-dependent seismic wavespeeds. Occasionally, in-situ
measurements are possible near or at seismic stations, such as
groundwater wells \citep[e.g.,][]{lin2024database} or soil moisture
probes
\citep[e.g.,][]{illien2021subsurface,oakley2021seismic,Sobolevskaia24}.
Typically, the distance between the seismometer and the in-situ
hydrological sensor is up to a few kilometers, and the hydraulic
conductivity between the two sites is either unknown or highly
uncertain. Piezoelectric probes are frequently used to infer water
levels and have been employed in several local studies
\citep{voisin2017groundwater,Gaubert-Bastide22}.

Networks of groundwater wells are increasingly made available to the
public \citep{goodall2008first}, such as the National Ground Water
Monitoring Networks maintained in individual countries; however, data
rights still pose challenges for data sharing \citep{condon2021global}.
Overall, groundwater wells are either pumped for agricultural,
domestic, or industrial uses or used to monitor groundwater levels
\citep{margat2013groundwater} passively. In the past, one challenge
seismologists have faced in accessing groundwater well data to compare
with their seismic proxies has been the decentralization of groundwater
well observation networks and the non-uniformity in the temporal
sampling of the wells. 

Determining the soil moisture (SM) content of the CZ is crucial for
effective resource management and hazard mitigation during extreme
droughts and floods. There are several methods for inferring or
measuring soil moisture. \citet{babaeian2019ground} provides an extensive
review of these approaches, which are categorized as remotely sensed
(using satellite or airborne sensors), proximally sensed (e.g., GPS
reflectometry), or terrestrially sensed (e.g., time-domain
reflectrometry). These methods are only sensitive to the upper 50~cm,
limiting analysis to the upper part of the critical zone and lacking
sensitivity for deeper processes such as the weathered bedrock 
(2--30~m). Additionally, ground-based sensors are particularly local,
while remotely sensed SM measurements typically have a footprint of 10s
of kilometers: empirical downscaling strategies are often required to
bridge the scales.

When lacking in situ or remotely sensed measurements, seismologists
have turned to inferring water content from a hydrological model (see 
Section~\ref{sec:hydro_model}). For example, \citet{Lecocq2017} used
precipitation data together with a runoff, discharge, and
evapotranspiration model to estimate groundwater 
\mbox{storage.} 

\section{Review of methodologies in seismology}\label{sec4}

\subsection{Active near-surface imaging with body waves}\label{sec4.1}
Near-surface geophysics (NSG) typically refers to temporary experiments
that involve deploying sensors, transmitting a source, and inverting
data to determine water content and structural properties. NSG is a
mature field of geophysics, with a vast set of methodologies and
dedicated curriculum \citep{butler2005near,everett2013near}.
Non-invasive methods for characterizing critical zone water content
include electrical and electromagnetic surveys
\citep{siemon2009review}, magnetic susceptibility and resonance methods
\citep{hertrich2008imaging}, ground penetrating radar
\citep{klotzsche2018measuring}, gravity survey
\citep{gehman2009estimating, syed2008analysis}, and seismic, on which
this review will expand.

Seismic methods for inferring near-surface water content typically use
active seismic surveys and short-period seismometers (geophones)
\citep{Befus2011,Flinchum22}. Usually, travel time measurements of the
direct and refracted P waves are used to generate reflection and
refraction images of the near-surface \citep{Befus2011}. The spatial
scales of these surveys vary from 10 to 100s of meters, depending on
the spacing of the geophones \citep[e.g.,][]{uecker2023mapping,
Flinchum24}. Sophisticated methods, such as full-waveform inversion,
which were initially developed to image regional-scale subsurface
stratigraphy, are now being applied to the critical zone
\citep{liu2022}.

\subsection{Time-lapse surface wave tomography}\label{sec4.2}
Surface waves are dispersive, and their dispersion relation depends on
the vertical shear wave profile. Using either active or passive (e.g.,
ambient noise cross correlations) sources, surface wave imaging in
near-surface geophysics has significantly advanced the monitoring of
water content changes in the critical zone. By analyzing surface wave
dispersion, these methods utilize the frequency-dependent velocity of
Rayleigh waves to infer subsurface elastic properties, which are
sensitive to variations in water content. Dispersion inversion
techniques, such as multichannel analysis of surface waves (MASW) or
the frequency-wavenumber method \citep{Park1999}, allow the creation
of shear wave velocity profiles, and when repeated over time, they
correlate with changes in soil moisture and porosity
\citep{Socco10,Cheng22,Sobolevskaia24,Sun25} or shallow
aquifer structure \citep{gochenour2024ambient, eppinger20242d}.  

\subsection{Changes in seismic velocity $\mathrm{d}v/v$}\label{sec4.3}
The relative change in seismic velocity, $\mathrm{d}v/v$, derived from
ambient noise cross-correlation is sensitive to near-surface changes in
saturation, pore pressure, and overall effective stress.
$\mathrm{d}v/v$ monitoring of hydrological changes in the critical zone
spans from point measurements using single station cross-channel
correlations \citep{Sens-Schonfelder2006, donaldson2019crustal, Kim19,
lindner2021seasonal, feng2021controls, Clements23} to kilometer-scale
variations using inter-station cross-correlations \citep{Lecocq2017,
wang2017seasonal, Clements18, mao2022space, mao2025depth}. We
illustrate the clear evidence that $\mathrm{d}v/v$ correlates with
water content in Figure~\ref{fig:dv_ts}. 

\begin{figure*}
\includegraphics{fig03}
\caption{\label{fig:dv_ts}Temporal evolution of $\mathrm{d}v/v$ relative to
water storage seen for (a)~soil moisture (green and red lines) by
\citet{illien2021subsurface}, (b)~groundwater (blue line) by many,
including \citet{Clements18}, and (c)~in regional-scale karstic
groundwater storage (gray line) by \citet{Almagro21}.}
\end{figure*}

The temporal resolution of the $\mathrm{d}v/v$ transient depends on the
time scale of the process under study, which spans from seconds for
earthquake shaking \citep{illien2025predictable} and liquefaction
\citep{Viens_2022a} to decades for climatic processes
\citep{Lecocq2017, Clements18, mao2022space, mao2025depth}, as shown in
Figure~\ref{fig:scales}. For long-term monitoring studies and at
regional scales, NCFs are averaged in 1--30 day moving windows to
increase signal-to-noise ratio (SNR) \citep{seats2012improved}, while
other studies with short aperture may report sub-hourly temporal
resolution in $\mathrm{d}v/v$
\citep{Xie23,RodriguezTribaldos2021AquiferFiber}. The non-stationarity
of the seismic sources and their coherence is the main limiter of
temporal resolution and robustness of $\mathrm{d}v/v$ measurements.
Typically, short aperture arrays or single-station measurements can
yield high temporal resolution (e.g., minutes to hours) at high
frequency signals \citep[e.g.,][]{Xie23} while greater array aperture
and lower frequency signals tend to become stable over longer time
periods (e.g., days to months)
\citep[e.g.,][]{mao2025depth,Clements23}. Improving the quality and
temporal resolution of cross-correlation is an area of research
\citep{Baig09,hadziioannou2011improving,Viens19,yang2023optimal,He24}.
Retrieving dynamics at time scales shorter than minutes, ambient field
seismology relies on the stationarity of the ambient field. It tends to
be limited to single-channel analysis or near-sensor\break processes. 

Coda waves multiply scatter in the subsurface and thus have increased,
though still weak, sensitivity to perturbations in the medium
\citep{colombi2014temporal}. Given that the sensitivity is weak,
observed $\mathrm{d}v/v$ values typically under-predict the
perturbation at depth
\citep{obermann2013imaging,Obermann2016Lapse-time-dependentMedia,Barajas22,Yuan2021NumericalInterferometry}. 
Furthermore, $\mathrm{d}v/v$  is measured as decreasing with lag time
in the coda \citep[e.g.,][]{Gassenmeier15,Takano19}, which is also
expected from the low sensitivity of body waves to perturbations
\citep{obermann2013imaging,Barajas22}.\looseness=-1

$\mathrm{d}v/v$ can be measured in both the frequency and time domains.
The moving window cross-spectral (MWCS) technique developed by
\citet{poupinet1984monitoring} and the stretching technique developed
by \citet{lobkis2003coda} and first used in an environmental context by
\citet{Sens-Schonfelder2006} are the most commonly used methods in the
frequency and time domains, respectively. The stretching technique
works well at low SNR levels \citep{hadziioannou2009stability}, whereas
the MWCS is computationally faster to calculate \citep{Clarke2011}.
Recently \citet{okubo2024monitoring} found that the stretching and MWCS
methods may yield different results depending on which scattered wave
dominates in the coda. Additionally, \citet{Mikesell15} developed a
method based on dynamic time warping to measure local phase shift
measurements in the coda, \citet{Mao19} developed a wavelet-based
method to measure phase shifts, and
\citet{Yuan2021NumericalInterferometry} suggested combining the wavelet
and stretching techniques to improve the frequency-dependence of the
measurements relative to all previously stated methods. Overall,
studies have attempted to measure $\mathrm{d}v/v$ over a range of
seismic frequencies, which is desirable in the attempt to estimate the
depth of the\break changes.

\subsection{Changes in seismic attenuation}\label{sec4.4}
Attenuation is highly sensitive to water content, but is always
comprised of two effects, intrinsic attenuation and scattering
attenuation, which are \mbox{often} \mbox{difficult} to disentangle in observations.
A first method for measuring attenuation as related to fluids estimates
the amplitude decay as a function of distance in the frequency domain
through direct single wave measurements \citep[e.g.,][]{Hudson23} or
using the exponential decay in the spectral autocorrelation method
\citep[e.g.,][]{aki1957space,Prieto09,Nakahara12,Lawrence13}. Another
class of methods proposed by \citet{Malagnini22} utilizes a generalized
inverse technique for regression among pairs of earthquakes and
stations. A second method for attenuation measurements focuses on the
waveform envelope shape for which the radiative transfer theory
\citep{sato2012seismic,Mayor14} is frequently used to model the
dissipation of seismic energy and incorporate explicit functional forms
for intrinsic and scattering attenuation, which may be necessary to
extract the seismic signature of fluids
\citep[e.g.,][]{Hirose19,hirose2022scattering,Hirose23,vanLaaten24}.
Overall, measurements are made with absolute values of attenuation and
evaluated over time. Recently, \citet{Pinzon-Rincon25} observed  5\%
variations of seismic attenuation due to changes in shallow water
storage in the first tens of meters of the subsurface (12--23~Hz).

\subsection{Changes in ground motion site response}\label{sec4.5}
The response of the near-surface to incoming earthquake ground motions
is referred to as ``site response'' or ``site effects'', and they often
include linear effects of wave propagation in the shallow soil layers
and the damage and nonlinear deformation of these materials when ground
motion amplitudes are 
\mbox{sufficiently} high. Site effects arise due to the
presence of impedance (the product of wave velocity and density)
contrasts and subsurface topographies. Site effects also depend on the
attenuation of waves in the near-surface environments. Modeling site
effects, therefore, requires in-depth knowledge of the site's velocity
structure and the various factors controlling attenuation. 

Site response can be precisely measured using instrument-based methods,
such as spectral ratio (SR) techniques \citep{Borcherdt_1970}, the
generalized inversion technique \citep[GIT,][]{andrews1986objective},
and residual analysis \citep{Abrahamson_1992}, within ground motion
models (GMMs). SR methods require a reference-target recording pair,
such as a borehole station or a station on subsurface reference rock.
In contrast, GIT and GMM residual analysis depend on a regional seismic
network that includes the site of interest and has recorded multiple
past earthquakes with comprehensive azimuthal and distance coverage. 
SR methods have also been utilized for measuring the environmental
impacts on slope stability \citep[e.g.,][]{Borgeat24}.

Where direct measurement or modeling is not possible, amplification can
be predicted by generic prediction models linking amplifications to
various predictor variables. Predictive variables can be site
parameters from in situ geotechnical/geophysical measurements, e.g.,
shear wave velocity in the shallowest 30~m ($V_{\mathrm{S}}30$), sediment
thickness, and site resonant frequency
\citep[e.g.,][]{Hassani17,Kwak20} or indirect data from existing
regional models or maps, such as topographic slope, terrain and geology
\citep[e.g.,][]{Wald07,Weatherill20}. These indirect methods have
recently shown good performance at predicting site response due to the
emergence of quality databases and the associated use of machine
learning methods \citep{zhu2021open,Zhu21,Zhu23}.

Another metric of site response is the horizontal-to-vertical spectral
ratio (HVSR). HVSR identifies fundamental resonance frequencies of soil
layers, aiding in seismic microzonation, site amplification studies,
and subsurface characterization \citep{nakamura1989method}. It is
particularly valuable for evaluating the thickness and stiffness of
sediments, seismic hazard 
\mbox{assessment,} and earthquake-resistant design 
\citep[see {review} by][]{molnar2018application}. The influence of
environmental variations on site responses is just beginning to be
explored \citep[e.g.,][]{Vassallo22,kula2018horizontal,Colombero18}.
\citet{handel2025detecting} have also identified seasonal variations of
high-frequency ground motion characteristics in northeastern Hokkaido
and Honshu.

\section{Seismic monitoring of hydrological\newline
resources and hazards}\label{sec5}

\subsection{Vadose-zone dynamics}\label{sec:vadose} %5.1
Seismological methods offer valuable insights into the partially
saturated vadose zone, where complex water dynamics, including
infiltration, evapotranspiration, drainage, and runoff, interact with
the mechanical properties of near-surface soils and rocks. These
processes have been modeled using a range of hydrological and
poroelastic frameworks, and increasingly constrained by seismic
observations of velocity changes ($\mathrm{d}v/v$) that respond to
water saturation and stress conditions.

Several studies have used seismological data to \textbf{model
hydrological processes}. For instance, \citet{Fores18} modeled  
${\pm}$0.2\% $\mathrm{d}v/v$ variations in the 6--8~Hz frequency band
by integrating a water budget, including rainfall, evapotranspiration,
and drainage, with gravity data and Biot-Gassmann-based petrophysical
models, effectively tracking changes in water saturation in a karstic
setting. Similarly, \citet{illien2021subsurface} developed a two-layer
hydrological model incorporating an antecedent precipitation index and
rainfall thresholds to explain  ${\pm}$2\% $\mathrm{d}v/v$ variations
in the  4--8~Hz frequency range driven by transpiration and runoff.
Their model effectively captured responses to heavy monsoonal
precipitation in Nepal. In permafrost settings, \mbox{\citet{Sobolevskaia24}}
and \citet{Sun25} analyzed $\mathrm{d}v/v$ in terms of thermoelastic,
poroelastic, and saturation-driven processes caused by thawing,
precipitation, and air temperature variations.

Vadose zone seismology has been enabled by both \textbf{long-term
sparse networks} and \textbf{short-term dense deployments}. The
Cr\'epieux--Charmy water exploitation site in France exemplifies a
long-term test bed where ambient noise seismology has been combined
with hydrological observations using permanent and temporary densified
arrays, as well as piezo-electric sensors \citep{voisin2017groundwater,
Garambois19, Gaubert-Bastide22}. The Susquehanna Shale Hills Critical
Zone Observatory in the United States has supported year-round seismic
monitoring to detect shallow responses to meteorological and
hydrological forcing \citep{oakley2021seismic}. More recently, dense
distributed acoustic sensing (DAS) deployments have enabled higher
spatial resolution. For example, several studies used DAS to detect
sub-seasonal $\mathrm{d}v/v$ decreases following precipitation and
long-term increases associated with drying, with a depth sensitivity
down to 150~m in semi-arid eastern California and in permafrost
regions of Alaska \citep{shen2024fiber,Sobolevskaia24,Sun25}.

Interpreting vadose zone dynamics requires careful attention to the
differential sensitivity of seismic waves to fluid and matrix
properties. S waves are particularly responsive to changes in effective
stress and soil stiffness. \citet{Garambois19} found an approximately
linear relationship of a  1\% drop in S-wave $\mathrm{d}v/v$ and 1 m
increase in water table height at frequencies above 5~Hz, although with
deviations from that linear relation with a hysteresis observed during
wetting and drying
\citep[e.g.,][]{Garambois19,Gaubert-Bastide22,Shen25} can be explained
with capillary effects in the partially saturated zone
\citep{Dong16,Solazzi21}. In contrast, P-wave velocities are more
sensitive to changes in fluid compressibility; in the same study,
increases in P-wave velocity were attributed to greater saturation and
reduced pore compressibility. Shallow monitoring studies, such as
\citet{oakley2021seismic} and \citet{Sun25}, further highlighted that
near-surface P-wave velocities may also respond to thermoelastic
effects caused by temperature fluctuations and moisture changes in the
top few meters. DAS technologies, which primarily sense strain rates
related to mixed wave modes, offer a promising avenue to observe these
joint sensitivities over depth and time
\citep{shen2024fiber,Sobolevskaia24,Sun25}.

\subsection{Groundwater resource management}\label{sec5.2}
Groundwater monitoring with $\mathrm{d}v/v$ is a developing technique
for water resource management at the aquifer scale. There is a strong
anti-correlation 
\mbox{between} $\mathrm{d}v/v$ and relative changes in
groundwater levels, i.e., when groundwater levels decrease,
$\mathrm{d}v/v$ increases and vice-versa. \citet{Clements18} showed
that an increase of up to 16~m of groundwater level over a ${\sim} 20
\times 40$~km area during the winter of 2004/2005 was commensurate with
a  0.15\% decrease in $\mathrm{d}v/v$ in the  0.5--2.0~Hz frequency
range (sensitivity to ${\sim} 500$~m depth). Spatial changes in
$\mathrm{d}v/v$ and groundwater level were the largest in the deepest
part of the basin.  \citet{mao2022space} compared $\mathrm{d}v/v$ and
groundwater levels over the ${\sim} 40 \times 80$~km greater Los
Angeles coastal aquifers in the  0.5--2.0~Hz frequency band
(sensitivity to ${\sim} 350$~m depth). They measured a long-term
increase of up to  0.2\% $\mathrm{d}v/v$ in the central Los Angeles
(LA) basin associated with groundwater pumping from 2000--2020. Whereas
the Santa Ana basin with active groundwater replenishment had a  0.2\%
decrease in $\mathrm{d}v/v$. \citet{Almagro21} measured $\mathrm{d}v/v$
in the  0.7--0.9~Hz frequency range (depth sensitivity to ${\sim}
200$~m) over an area of ${\sim} 60 \times 60$~km in Northern Italy.
They found that  0.05\% $\mathrm{d}v/v$ variations were correlated to
$1\times10^{8}$~m$^3$ changes in total groundwater storage in a karstic aquifer. At
a much smaller spatial scale (${\sim} 500 \times 500$~m area) in the
Cr\'epieux--Charmy field site, \citet{Gaubert-Bastide22} measured a 3\%
decrease in $\mathrm{d}v/v$ in the  2--5~Hz frequency range
corresponding to a 3~m increase in groundwater level. 
Figure~\ref{fig:gwl_map} shows the progress of increasing the spatial
resolution of groundwater monitoring with $\mathrm{d}v/v$ for water
management.

\begin{figure*}
\includegraphics{fig04}
\vspace*{-3pt}
\caption{\label{fig:gwl_map}{Using seismological proxies to
groundwater levels}. Left: \citet{Clements18} first developed a smooth
map of groundwater level change in the San Gabriel basin, CA. Spatial
extent: ${\sim} 20 \times 60$~km. 
Middle: \citet{mao2022space} mapped
$\mathrm{d}v/v$ at $2 \times 2$~km with 2D sensitivity kernels in the greater
Los Angeles, CA area. Spatial extent: ${\sim} 40 \times 80$~km. Right:
\citet{Gaubert-Bastide22} mapped $\mathrm{d}v/v$ at $20 \times 20$~m in the
Cr\'epieux--Charmy field site, Lyon, France. Spatial extent: ${\sim}
500 \times 500$~m.}
\vspace*{-3pt}
\end{figure*}

Beyond measuring relative changes in groundwater levels, seismological
studies have attempted to infer groundwater dynamics. For example,
\citet{Delouche23} found that deep karstic and 
{confined} aquifers had a
rapid response to seasonal rain (30--60~days) and a reduction in
velocities. \citet{mao2025depth} found that intense rainfall events
during the 2023 winter storms in \mbox{California} replenished shallow
aquifers but not the deep aquifers that are the main contributors to
Calfornia's water storage. \looseness=-1

\subsection{Landslide monitoring}\label{sec5.3}
Landslide dynamics are profoundly affected by subsurface hydrology
\citep{Greco23}. Rapid changes in the stress field or the water
storage, e.g., due to massive rainfall events or earthquake shaking,
are \mbox{often} causes for slope failures \citep{dai2002landslide}. 
\mbox{Because} hydrological data alone do not always correlate with landslide
movements \citep{Schulz18}, additional information about the rheology
and transient processes, such as heterogeneous soil swelling in shear
zones \citepalias{Schulz18} or rheological and hydrological heterogeneity
\citep{tacher2005modelling}, is required. The pioneering work of
\citet{mainsant2012ambient} demonstrated that pore-pressure-driven
$\mathrm{d}v/v$ decreased  7\% in the 10--12~Hz frequency range in the
four days precedding a catastrophic landslide. Since then, ambient
noise monitoring has become a popular way to monitor continuously and
forecast landslides: \citet{LeBreton2021LandslideApplications} recently
provided a comprehensive review of how ambient noise monitoring has
been applied to and benefited landslide research, and in particular,
the advantages and challenges in making the $\mathrm{d}v/v$
measurements robust at short time scales for monitoring landslides.
More recently, researchers have monitored the water content in
landslides, tailing dams, and unstable rock slopes with
$\mathrm{d}v/v$, especially in response to rainfall events
\citep{liu2024monitoring,ouellet2022advanced,borgeat2025variability}.
\citet{marc2021toward} showed that the landslide susceptibility rate in
the epicentral zones of large earthquakes followed a similar recovery
as $\mathrm{d}v/v$ in the year following strong ground motion.
\citet{collins2024using} used HVSR to estimate the thickness of the
slow-moving Barry Arm landslide. Overall, ambient noise methods provide
in-situ, high temporal resolution monitoring of the coupled mechanical
and hydrological processes that drive landslides. These ground
failures, whether triggered by seismic or hydrological loading, should
be viewed as dynamic responses of the critical zone, where stress
redistribution and fluid--solid coupling govern stability.

\subsection{Seismic hazards}\label{sec5.4}
Site effects can significantly impact seismic hazards by amplifying
seismic waves and altering soil behavior in the shallow subsurface
(i.e., the upper ${\sim}$50~m). Under weak ground motions, shallow
soil layers typically behave in a linear elastic manner, with
amplification primarily dependent on the incoming seismic wave
strength. During strong shaking, soft soils may undergo substantial
dynamic strains and exhibit non-linear responses \citep{Ishihara_1996}.
This non-linearity often results in a relative de-amplification of
high-frequency ground motions and a shift of energy toward lower
frequencies  \citep{Beresnev_1996}. Numerous studies have highlighted
the potential of seismic interferometry to monitor the non-linear
behavior of shallow sediments from single \mbox{seismometers}
\citep{Bonilla_2020,Hobiger_2014,Hobiger_2016,Viens_2018,Esfahani24},
surface-borehole arrays \citep{Bonilla_2019,Nakata_2011,Qin_2020}, and
DAS arrays \citep{Viens_2022a}. $\mathrm{d}v/v$ can decrease as much as
50--60\% in the few seconds following strong motion at frequencies
above 1~Hz \citep{Bonilla_2019}. In extreme cases, phenomena such as
liquefaction \citep{Nakata_2011}, dynamic compaction
\citep{Viens_2024}, and significant variations in soil stiffness and
strength \citep{Chandra15} have been observed. Recent studies have also
highlighted the impact of environmental changes, such as rainfall and
freeze-thaw cycles, on the subsurface response to ground motions
\citep{Roumelioti_2020,handel2025detecting}. Importantly, seismic
hazard itself, via liquefaction and ground failure, is a
hydromechanical outcome of CZ processes, reinforcing the need for an
integrated perspective.

\subsection{Permafrost}\label{sec5.5}
Permafrost, which underlies vast areas of the Northern Hemisphere and
other cold regions, is warming at an increasing rate due to climate
change. The degradation of permafrost poses significant risks,
including the release of greenhouse gases such as methane and carbon
dioxide, which further exacerbate global warming, as well as threats to
ecosystems, hydrology, and infrastructure
\citep{biskaborn2019permafrost}. Conventional monitoring systems
include temperature measurements through the Global Terrestrial Network
for Permafrost and conventional active geophysics surveys, such as
resistivity and 
{ground-penetrating} radar \citep{kneisel2008advances}.
{Environmental} seismology has demonstrated great potential in monitoring
permafrost phenomena, mainly through the $\mathrm{d}v/v$ and HVSR
methods. $\mathrm{d}v/v$ monitoring constrains the timing of seasonal
freeze and thaw in the active layer, which comprises the upper layer of
soil. Thawing of the active layer in early summer induces a  5--10\%
drop in $\mathrm{d}v/v$ at frequencies of  1--30~Hz, whereas refreezing
in the late fall increases $\mathrm{d}v/v$ a similar
amount \citep{james2017improved, james2019insights, Guillemot20,
steinmann2021effect, james2021biophysical}. This  5--10\%
$\mathrm{d}v/v$ drop is much less than ${\sim}$75\% difference in
velocity between the unfrozen active layer and the frozen permafrost
body \citep{kneisel2008advances}, which suggests that most
$\mathrm{d}v/v$ studies are sensitive to depths deeper than the active
layer.  Seasonal fluctuations in the thickness of the active layer also
lead to a sharp velocity contrast at the active layer-permafrost
boundary, which can be monitored with the HVSR on land
\citep{kula2018horizontal,weber2018ambient, kohler2019potentials} and
in submarine environments \citep{overduin2015submarine,
angelopoulos2024passive}. In addition to seasonal variations in
$\mathrm{d}v/v$, \citet{albaric2021seismic} and
\citet{lindner2021seasonal} recovered long-term decreases in
$\mathrm{d}v/v$ in permafrost over 3 and 15 years, respectively,
suggesting an expansion and warming of the active layer. DAS arrays
deployed in permafrost allow for detailed inversion of permafrost
structure using tomographic methods \citep{ajo2017time} and
$\mathrm{d}v/v$ monitoring of both the timing and spatial extent of the
freeze and thaw cycle \citep{Cheng22,Sobolevskaia24}.
Recently, \citet{steinmann2022ai} employed an artificial
intelligence-based method to identify the pattern of ground freezing
from continuous ambient seismic noise in an urban environment, and
\citet{Sun25} utilized a rock-physics model to disentangle the relative
contributions of thawing, precipitation, thermoelastic, and poreelastic
stresses, each of which may dominate at specific depth ranges and
depending on the water table depth \citep{Yu25}.\looseness=-1

\begin{figure*}
\vspace*{-3pt}
\includegraphics{fig05}
\vspace*{-3pt}
\caption{\label{fig:dvv-gwl}{Change in seismic velocities as a
function of groundwater level change}. Values are shown from individual
studies that reported water levels, or pore pressure, which we
converted to water level change using a hydrostatic pressure $p = \rho
g\,\mathrm{d}H$, where $\rho = 1000$~kg/m$^3$ and  $g = 9.81$~m/s$^2$
the gravity constant. Squares represent values obtained from
inter-station correlations from regional seismic networks with station
spacing on the order of 10~km. Triangles are measurements made with
inter-station correlations of small arrays of apertures of about 1~km.
Circle markers are measurements obtained from single-station methods.
Gray triangle is a digitally extracted data from \citet{Garambois19}.
The color bar represents the middle of the frequency band used in each
study. The dashed red line is simply showing $\mathrm{d}v/v = -\mathrm{GWL}$
without any fitting.}
\vspace*{-3pt}
\end{figure*}

\subsection{Snow load and melt}\label{sec5.6}
The effects of snow load and melt on $\mathrm{d}v/v$ can be
multi-faceted. \citet{Guillemot2021} measured a  2--3\% decrease in
$\mathrm{d}v/v$ between  10--25~Hz after dry  snowfall. They
interpreted snowfall as adding new mass onto the  surface without
increasing the rigidity, thereby decreasing seismic velocity. In
contrast, vertical loads induced by the weight of snow (and or pre-snow
precipitation) have been interpreted with a positive increase in
$\mathrm{d}v/v$ at frequencies below 10~Hz
\citep{CANNATA2017165,makus2023deciphering}. Heavy, wet snow has been
interpreted as  (1)~increasing water content in the pores as correlated
with snow-water equivalent thickness \citep{hotovec2014continuous} or 
(2)~increasing pore pressure in fluid-filled cracks
\citep{Silver07,taira2016response}. 

\section{Summary of seismological evidence for ground-water
hydrology}\label{sec6}
Our review highlights the ubiquitous correlation between groundwater
levels, $\mathrm{d}v/v$, and other seismic markers. To achieve such a
correlation, studies may have had to disentangle the hydrological
signature in $\mathrm{d}v/v$ from other contributing factors (e.g.,
tectonic, thermoelastic) using modeled and/or observed data. From the
numerous studies mentioned in this review, we collected the values of
$\mathrm{d}v/v$, their corresponding $\mathrm{d}H$ groundwater level
changes, the type of seismic network configuration, and the frequencies
(i.e., depth) involved.  Figure~\ref{fig:dvv-gwl} shows how
$\mathrm{d}v/v$ has been reported to vary against groundwater levels.
The studies are compiled based on comparisons with aquifer or
water-table levels. We categorize them based on the experimental
setting: inter-station cross-correlations in regional seismic networks,
small aperture arrays, or single-station measurements.

\begin{figure*}
\vspace*{-3pt}
\includegraphics{fig06}
\vspace*{-3pt}
\caption{\label{fig:cartoon}Critical Zone processes are probed by
geophysical measurements through ground base sensors sensitive to the
subsurface (seismometers, distributed acoustic sensing) and airborne or
satellite measurements that probe the surface or very near surface
properties. Art made by Elena Hartley.}
\vspace*{-3pt}
\end{figure*}

Overall, the simple relation $\mathrm{d}v/v(\%) = - \mathrm{GWL}(m)$ is
strikingly accurate. The deviation from this relation is explained by
the following:  (1)~a choice of inter-station correlations with large
separation, and  (2)~the choice of a low-frequency band of analysis.
Both settings yield a lower sensitivity of $\mathrm{d}v/v$ to
perturbations at depth. The remaining variability is likely due to
differences in lithology, as studies range from karstic carbonate
environments to alluvial fans and large sedimentary basins. The
sensitivity of shear waves to sediments is drastically different,
especially between sand and clay \citep{Dong16}. 

\section{Future directions and research priorities}\label{sec7}

\subsection{Bridging data gaps with enhanced sensor\newline
deployment}\label{sec7.1}
Seismology is well-positioned to fill a hydrological data gap 
between remotely sensed observations and ground-based hydrological sensors,
which have a spatial sensitivity of 10--100~m and a depth 
\mbox{sensitivity}
of 1~cm to 100~m. Despite the recent progress in Earth observation,
data gaps persist in the spatial distribution of seismic sensors at
regional scales. Future studies should prioritize the co-located
deployment of seismic and hydrological sensors to improve ground truth
data and validate models. Multi-sensor data fusion will likely be a
desirable mode of analysis (see Figure~\ref{fig:cartoon}), such as
supplementing critical zone observatories \citep[e.g., CZO and
OZCAR][]{Brantley17,Gaillardet19} with added geophysical monitoring
capabilities (e.g., seismic and geodetic networks). This approach will
enhance our 
\mbox{understanding} of subsurface processes by directly
correlating seismic observations with hydrological measurements.
Emerging technologies, such as distributed acoustic sensing (DAS),
could play a pivotal role in addressing this gap, enabling
cost-effective, high-density sensor networks. 

Additionally, it should be noted that seismological data tends to be
publicly accessible, especially through open-access repositories such
as the Earthscope Data Services and the  European Integrated Data
Archive (EIDA). There is also a growing availability of hydrological
data sets, such as those provided by the United States Geological
Survey, the Hydro-Climatic Data Network \citep{lins2012usgs}, the
National Water Information System \citep{us2020usgs}, NASA's
Earth-Observing Giovanni system, and Copernicus Climate Change Service.
Open seismic and hydrological data thus allow for both timely and
long-term critical zone\break science. 

\subsection{Developing unified physical models}\label{sec7.2}
A critical challenge is the lack of a unified theoretical framework for
integrating seismic observations with hydrological and geomorphological
processes. Bridging the disciplines of continuum mechanics, granular
physics, and soil science could lead to comprehensive models that
capture the coupled effects of environmental influences. For instance,
some studies have focused on nonlinear elasticity as a linear
superposition of dilatational strains
\citep[e.g.,][]{donaldson2019crustal} or extensional strains
\citep[][]{okubo2024monitoring}, or directly modeling $\mathrm{d}v/v$
terms \citep[e.g.,][]{Sun25}, while others have focused on induced
stresses \citep[e.g.,][]{Fokker21}, and others use petrophysical models
to estimate P and S wavespeed with water content
\citep{Solazzi21,shen2024fiber,Sobolevskaia24}. 

An additional challenge is to separate elastic and poroelastic effects:
seismic wavespeeds decrease in unconfined aquifers with increasing
saturated but increase in confined aquifers due to the increased
\mbox{hydrological} load \citep{Delouche23}. Furthermore, the coupling of
thermoelastic stresses with hydrological processes is likely to play a
role, given the large thermoelastic stresses in the critical zone
\citep{sens2019probing,oakley2021seismic,Sun25}. Collaborative efforts
among seismologists, soil scientists, and hydrologists will be key to
developing a unified physical model.

\subsection{Standardizing data processing}\label{sec7.3}
Drawing unified models from the diverse literature on ambient noise
monitoring for hydrological work is challenged by the individualized
data processing scheme. Processing strategies for seismic data are
currently study-specific, limiting the generalization and
reproducibility of findings. Each study has a given sensor
configuration, seismic frequency of interest (e.g., depth sensitivity),
and lapse time in the coda (e.g., different wave types and spatial
sensitivities for early vs.\ late coda), yielding a large variability in
the proportionality between $\mathrm{d}v/v$ and groundwater levels (see
Figure~\ref{fig:dvv-gwl}). Efforts to standardize preprocessing
workflows, combined with the adoption of automated feature extraction
and machine learning techniques (e.g., recent examples to predict river
discharge or soil saturation directly from ambient noise
\citet{abi2023machine,Cunha25}), can accelerate the identification of
hydrological signatures in seismic data. These advancements will also
facilitate the development of transferable methodologies for broader
applications in CZ science.

\section{Conclusion}\label{sec8}
This review highlights the emerging role of ambient field seismology in
bridging hydrological and geophysical sciences for monitoring the
Critical Zone (CZ). Advances in ambient noise interferometry, coda wave
analysis, and time-lapse tomography have demonstrated that near-surface
seismic velocity and attenuation variations are intimately linked to
changes in water saturation, pore pressure, and effective stress.

The range of applications covered in this review, from groundwater
monitoring, vadose zone dynamics, slope stability, to seismic hazard
assessment, may appear distinct, yet they reflect shared 
\mbox{physical}
mechanisms: water-driven changes in subsurface stress and mechanical
properties. Recognizing this, we advocate for a unified conceptual
framework in which hydrological processes and geohazards are understood
as coupled systems responding to climatatic and tectonic forcings.

Seismic velocity changes do not always behave uniformly with increasing
water content: while pore pressure increases in saturated media
typically reduce velocity, partial saturation and capillary effects in
the vadose zone may increase soil stiffness and seismic speeds. These
contrasting effects underscore the importance of site-specific rock
physics, soil structure, and frequency-dependent wave sensitivity. The
interpretation of seismic observables, whether $\mathrm{d}v/v$ or $Q$,
must account for these\break regimes.

\begin{table*}
\fontsize{9}{11}\selectfont\tabcolsep2.8pt
\caption{\label{tab1}Literature review of studies that have measured
$\mathrm{d}v/v$ and  compared with  groundwater levels\vspace*{-2pt}}
\begin{tabular}{cccccccc}
\thead
Region & Wave type & Network & 
\parbox[t]{1cm}{\centering Freq (Hz)} & 
\parbox[t]{1cm}{\centering $\Delta$GWL (m)} &
\parbox[t]{1.5cm}{\centering Hydro estimate method} &
\parbox[t]{1cm}{\centering $\mathrm{d}v/v$ (\%)} & 
Reference \vspace*{2pt}\\
\endthead
\parbox[t]{1.5cm}{\centering
Southern CA, USA} & Coda & Array & 
0.2--2 & 16 & Well & 0.10 & 
\parbox[t]{4.2cm}{\centering
\citet{mao2022space}} \vspace*{2pt}\\

CA, USA & Coda & Single &
4--8 & 6 & Rainfall & 1.15 &
\parbox[t]{4.2cm}{\centering
\citet{Clements23}}\vspace*{2pt}\\

San Gabriel, CA & Coda & Array & 
2--5 & 25 & Well & 0.15 & 
\parbox[t]{4.2cm}{\centering
\citet{Clements18}} \vspace*{2pt}\\

New Zealand & Coda & Small array & 
6--8 & 2 & 
\parbox[t]{1.5cm}{\centering
Piezo, seismometers}
& 2.00 & 
\parbox[t]{4.2cm}{\centering
\citet{Voisin16}} \vspace*{2pt}\\

Lyon, France & Coda & Small Array & 
2--5 & 2 & Piezo & 2.50 & 
\parbox[t]{4.2cm}{\centering
\citet{Gaubert-Bastide22}} \vspace*{2pt}\\

Lyon, France & Coda & Small array & 
3--20 & 3 & Piezo & 2.50 & 
\parbox[t]{4.2cm}{\centering
\citet{voisin2017groundwater}} \vspace*{2pt}\\

Ketzin, Germany & Coda & Array & 
1.4--3 & 0.6 & Well & 0.60 & 
\parbox[t]{4.2cm}{\centering
\citet{Gassenmeier15}} \vspace*{2pt}\\

Texas, USA & Coda & Single & 
0.01--8 & 25 & Well& 2.00 & 
\parbox[t]{4.2cm}{\centering
\citet{Kim19}} \vspace*{2pt}\\

Germany & Coda & Array & 
0.1--0.8 & 0.1 & 
\parbox[t]{1.5cm}{\centering
Well, water balance}
& 0.01 & 
\parbox[t]{4.2cm}{\centering
\citet{Lecocq2017} (*)}\vspace*{2pt}\\

Malta & Coda & Single & 
3--10 & 0.5 & 
\parbox[t]{1cm}{\centering
Well, rainfall}
& 1.00 & 
\parbox[t]{4.2cm}{\centering
\citet{Laudi23}}\vspace*{2pt}\\

CA, USA & Coda & Small array & 
4--15 & 1.5 & River gage & 2.50 & 
\parbox[t]{4.2cm}{\centering
\citet{RodriguezTribaldos2021AquiferFiber}}\vspace*{2pt}\\

Australia & Coda & Small array & 20 & & 
0.3 & 2.00 &  
\parbox[t]{4.2cm}{\centering
\citet{olivier2017monitoring}}\vspace*{2pt}\\

Kyushu, Japan & Coda & Array & 
0.1--0.39 & 0.01 & Rainfall  & 0.02 &
\parbox[t]{4.2cm}{\centering
\citet{wang2017seasonal} (*)}\vspace*{2pt}\\

OK, USA & Coda & Array & 
0.1--1 & 0.15 & GRACE & 0.02 &
\parbox[t]{4.2cm}{\centering
\citet{Zhang23} (***)} \vspace*{2pt}\\

Greece & Late coda  & Single & 
0.3--1 & 2--5 & Wells & 0.03 &
\citet{Delouche23} \vspace*{2pt}\\

Italy & Coda & Array & 
0.1--0.9 & 100 & 
\parbox[t]{1.5cm}{\centering
Geodetic inversion} & 
0.10 &
\citet{Almagro21}(**) \vspace*{2pt}\\

Germany & Late coda & Small array & 
1.5--3 & 0.4 & Well & 0.20 & 
\parbox[t]{4.2cm}{\centering
\citet{Gassenmeier15}}
\vspace*{2pt}\\

Lyon, FR & Surface waves & Small array & 
3--40 & 2.5 & 
\parbox[t]{1.5cm}{\centering
Piezo and wells} & 1.50 &
\parbox[t]{4.2cm}{\centering
\citet{Garambois19}} \vspace*{2pt}\\

Mojave, CA & Coda & Array & 
0.1--2 & 12 & Well & 0.10 & 
\parbox[t]{4.2cm}{\centering
\citet{Tsai2011} and \citet{Meier10}} \vspace*{2pt}\\

Taiwan & Coda & Single & 
0.1--0.9 & 7 & Well& 0.40 & 
\citet{feng2021controls}
\botline
\end{tabular}
\tabnote{Regions, waves
used, array configuration, sensor configuration, frequency band of
measurements, groundwater level change $\Delta$GWL (m), 
$\mathrm{d}v/v(\%)$ values reported, and references. The ``rainfall''
method refers to diffusion, either through poroelastic or baseflow
models, of rainfall data. ``single'' refers to single station
measurements, ``small array'' refers to arrays of geophones with less
than 1~km in aperture, and ``array'' refers to regional seismic
networks with inter-station spacing of about 10~km. (*)~Converted from
the hydraulic pressure  $\mathrm{d}H = p/g \rho$, where $g=9.81$, $\rho
= 1000$~kg/m$^3$.  (**)~GWL is computed for  1000~km$^2$ basin.
(***)~\citet{Zhang23} finds the following relation: $\mathrm{d}v/v =
-1.68\times10^{-3} \Delta \mathrm{GWL} - 9.61\times10^{-3}$.}
\vspace*{-3pt}
\end{table*}

Across applications, the magnitude and timescale of seismic responses
vary widely from sub-percent $\mathrm{d}v/v$ changes over hundreds of
meters and tens of years in deep aquifers to rapid, multi-percent
$\mathrm{d}v/v$ changes in shallow landslides. While the range of
aforementioned spatial and temporal scales underscore the versatility
of ambient noise monitoring, there is a need for standardized
sensitivity metrics and frameworks to compare observations across
hydrogeologic settings.

Fully realizing the potential of ambient field seismology will require
transdisciplinary collaboration. Integrating seismic observations with
hydrological modeling, soil physics, geotechnical engineering, remote
sensing, and climate science is essential to build interpretable,
scalable, and actionable monitoring systems. Beyond collaboration, this
integration demands co-development of research questions, shared
instrumentation strategies, and open data standards across
traditionally siloed disciplines.

To translate these insights into practice, we outlined in 
Section~\ref{sec7} a set of specific future directions: from
multi-sensor data fusion and physically informed modeling to
standardization and reproducibility efforts. Together, these strategies
can help realize the promise of passive seismic methods for
understanding and forecasting critical zone dynamics in a changing
world.

\section*{Acknowledgments}
This research was partially funded by the 
David and Lucile Packard Foundation (MD). 
This manuscript has been co-authored by LV with number 
LA-UR-25-20095 by Triad National Security, LLC, under contract
with the U.S. Department of Energy (no. 89233218CNA000001). 
LV was partially supported by the 
Los Alamos National Laboratory Center for 
Space and Earth Science Chick Keller Fellowship. 
TC was partially supported by a 
USGS Mendenhall fellowship. 
The authors thank the broader community that has contributed to the
discussion with MD about this review, namely Kuan-Fu Feng, Shujuan
Mao, St\'ephane Garambois, and Christophe Voisin. The authors thank two
anonymous reviewers and Stephanie James for reviews that improved the
manuscript. Much of the open-access data presented in the publications
are hosted by the Earthscope Consortium Data Services through the
National Science Foundation's Seismological Facility for the
Advancement of Geoscience (SAGE) Award under Cooperative Agreement
EAR-1724509. In European literature, EIDA is the European Integrated
Data Archive infrastructure within ORFEUS to provide access to seismic
waveform data in European archives.

\CDRGrant[Triad]{LA-UR-25-20095}
\CDRGrant[DOE]{89233218CNA000001}
\CDRGrant[SAGE]{EAR-1724509}

\section*{Declaration of interests}
The authors do not work for, advise, own shares in, or receive funds
from any organization that could benefit from this article, and have
declared no affiliations other than their research organizations.

\back{}

\def\appendixlabel{Appendix}

\section*{Appendix}
We summarize the values of $\mathrm{d}v/v$ and groundwater changes
published in the literature in  
\mbox{Table~\ref{tab1}}. The readers will find
Python scripts to reproduce  Figure~\ref{fig:dvv-gwl} on
\url{https://github.com/Denolle-Lab/dvv-hydro-CompGeo-review-2024}.

\printbibliography
\refinput{crgeos20250013-reference.tex}

\end{document}
