Learning earthquake ground motions via conditional generative modeling

Domain: biology
Published: 2026-03-16
Editorial outcome: major_revision

Abstract
Predicting high-fidelity ground motions for future earthquakes is crucial for seismic hazard assessment and infrastructure resilience. Conventional empirical simulations suffer from sparse sensor distribution and geographically localized earthquake locations, while physics-based methods are computationally intensive and require accurate representations of Earth structures and earthquake sources. We propose an artificial intelligence (AI) spectrogram generator, Conditional Generative Modeling for Ground Motion (CGM-GM). CGM-GM leverages earthquake magnitudes and geographic coordinates of earthquakes and sensors as inputs, when postprocessed with phase information, capturing spatially continuous Fourier amplitude spectra (FAS) as well as properties such as P and S arrivals, and waveform durations, without explicit physics constraints. This is achieved through a probabilistic autoencoder that extracts latent distributions in the time-frequency domain and variational sequential models for prior and posterior distributions. We evaluate the performance of CGM-GM using small-magnitude earthquake records from the San Francisco Bay Area, a region with high seismic risks. Here, we report that CGM-GM demonstrates potential for complementing physics-based simulations and non-ergodic empirical ground motion models, as well as shows promise in seismology and beyond.

Introduction
The accurate prediction of ground motion waveforms and their characteristics for future earthquakes is crucial for assessing seismic hazards and ensuring the safety and resilience of critical infrastructure. However, it is challenging and resource-intensive to obtain comprehensive ground motion observation across a wide geographic area. Furthermore, predictions of earthquake rupture processes and estimates of the Earth’s elastic model remain to exhibit significant uncertainties. The development of precise and robust ground motion prediction methodologies has long been of great interest in seismology and earthquake engineering to complement the limited recorded data.

Existing ground motion simulation studies branch into two streams: stochastic and physics-based approaches. The first stream, stochastic methods, is rooted in the modulation of Gaussian white noise to reproduce the desired ground motion characteristics [1–4]. They provide a computationally efficient framework to synthesize ground motion data by calibrating stochastic process-based models to match the historical recordings. However, potential limitations exist in representing spatial continuity and physical phenomena. The second stream, physics-based methods, is based on the numerical solution of wave equations [5–9] while considering comprehensive physical characteristics, including fault ruptures [10–13], heterogeneous earth media, and site-specific effects. Although recent advancements in high-performance computing enable simulating high-frequency waveforms of large-magnitude earthquakes (e.g., up to 10 Hz) [14], physics-based methods are computationally demanding. For example, ground motion simulations of the San Francisco Area over a domain of 120 km\times90 km\times35 km require 128 NVIDIA A100 GPU nodes and take 6 hours to compute up to a frequency of 5 Hz. Simulations at higher frequencies tend to be computationally prohibitive and these data are typically complemented by stochastic simulations [15, 16]. Furthermore, physics-based simulations face challenges from significant uncertainties in wave theory, subsurface elastic models, and source characteristics.

More recently, machine learning (ML) and artificial intelligence (AI) have shed new light on this classic task, primarily through their capability of accelerating earthquake modeling processes. One representative line of work involves the use of neural operators for modeling seismic wave propagation [17, 18] and extending their application to complex 3D elastic wave phenomena [19, 20]. Although these data-driven ML methods show efficiency, by avoiding the stringent time-step constraints in traditional time-domain physics-based numerical approaches, they require a large amount of high-fidelity data. Researchers have also resorted to incorporating physical constraints into ML models, such as physics-informed neural networks (PINNs) [21]. In particular, by leveraging physical principles as a prior, PINNs have been used for predicting ground motions with a limited amount of training data [22–24]. However, PINNs are known to exhibit fundamental failure modes in network optimization [25–27]. Moreover, due to their specific design of objective functions, PINNs encounter significant limitations in generalizing to different initial conditions and subsurface elastic models. PINNs can be regarded as physics-based and still suffer from uncertainties in the problem setup. In addition, the spectral bias of fully-connected neural networks used by both neural operators and PINNs typically constrains the resolution to low frequencies [23], making broadband synthesis of waveforms challenging. Furthermore, graph neural networks (GNNs) have also been employed to earthquake ground motion synthesis [28–30], leveraging their ability to capture spatial correlations. These models are capable of predicting ground motion intensities at more distant stations within the seismic network.

Generative modeling has emerged as an alternative approach for scientific modeling to capture the complexities of natural phenomena, as demonstrated for fluid dynamics [31, 32] and molecular science [33, 34]. It is an inherently stochastic method that incorporates random noise as inputs and uses probabilistic processes to generate diverse and realistic data. A detailed discussion of commonly used generative models can be found in Supplementary Section 1.2. In the context of earthquake ground motions, generative modeling can produce various waveforms while capturing model uncertainties. This capability is critical due to the significant variability and unpredictability of real-world earthquake ground motions, influenced by factors such as seismic sources, propagation paths, and site conditions. Generative models, which are not governed by specific wave equations or subsurface models, need to learn meaningful wave physics (i.e., governing equations) and site/source conditions from sparse and irregular sensor and earthquake distributions. Generative Adversarial Networks (GANs) and their variants [35–39] have demonstrated success in simulating ground motions, which paves the way for generative modeling of ground motions. Their ability to produce broadband waveforms makes them a valuable tool in generative modeling. However, their training process can be sensitive to design choices, often requiring careful tuning to mitigate challenges such as mode collapse and instability [40]. These prior works focus on using distances, magnitudes, and near-surface velocity structures (Vs30) as inputs, allowing for capturing ergodic stochastic characteristics. Nevertheless, to date, generative AI models have not been demonstrated to capture spatially variable wave propagation effects.

In this paper, we demonstrate that generative modeling can be more robust and powerful by introducing Variational Autoencoder (VAE) and incorporating geospatial coordinate information. This approach presents spatial continuity in median spectral components of the generated ground motions. To achieve this goal, we introduce a Conditional Generative modeling framework for Ground Motion (CGM-GM). Recent work [41] indicates that dynamic VAE models are a flexible, user-friendly, and robust alternative to GAN models for time series generation, mitigating common issues such as mode collapse and training instability. More broadly, this feature enhances the accessibility of this framework for seismologists who may not have extensive experience with generative modeling methods. The primary advantage of dynamic VAE models is their use of variational sequential architectures (e.g., recurrent neural networks) in both prior and posterior distributions, which facilitates the capture of the temporal evolution (dynamics) of time series. Our main methodological contribution lies in the design of dynamic VAE models for learning time-frequency information and the conditional embedding of physical parameters, including earthquake magnitudes, depths, and geospatial coordinates. This strategy enables the CGM-GM framework to implicitly learn the underlying physics and spatial heterogeneity from observation data, even without explicitly incorporating physical principles. As a result, our generative models can produce realistic simulations that are aware of as a combination of source, path, and site effects.

We demonstrate the effectiveness of our method by focusing on the San Francisco Bay Area (SFBA). Despite the relative infrequency of large-magnitude earthquakes, the densely populated SFBA remains a region of significant public interest due to its history of such events. This is different from previous studies [35, 36, 38] that have leveraged a much larger number of large-magnitude earthquakes recorded in Japan. Instead, our work uses small-magnitude earthquakes that are crucial for characterizing earthquake ground motion, especially for the linear effects of the path and site. Studying small-magnitude earthquakes allows us to estimate the region-specific variations of the ground motion, which can help reduce seismic hazard in areas where no large earthquake records are available [42, 43]. This task presents a significant challenge since our model needs to generalize well with lower signal-to-noise ratio (S/N) data. Moreover, the seismic data are recorded at a limited number of stations per event due to the released seismic energy and attenuation, which leads to a relatively sparse seismic network compared to wavelengths.

Here, we show the performance of our proposed model in producing Fourier Amplitude Spectra (FAS) maps and their spatial heterogeneity as well as reproducing properties such as P and S arrivals, and waveform durations. These features distinguish our framework from existing generative models in the context of ground motion generation. Additionally, CGM-GM results exhibit good agreement between generated samples and ground truth data across the entire frequency range (2,15 Hz), including waveform shapes, peak ground velocity (PGV) distributions, FAS, and arrival time of earthquake ground motions in the SFBA. We recognize that specifying only magnitude and source location does not uniquely define an earthquake event, and limits the capability to generate a set of waveforms from a single event. Nevertheless, we demonstrate that the model is capable of capturing spatial continuity in the median generated FAS.

Results
Datasets

The dataset for the study of small-magnitude earthquakes in 1990-2022 in the San Francisco region is downloaded from the Northern California Earthquake Data Center (NCEDC) database. The stations of interest are chosen within a 50 km radius from the Hayward fault, and the events with magnitude M < 4 recorded within a rupture distance Rrup < 100 km from the selected stations are included. We focus on two horizontal components H1 and H2 of particle velocity considering their direct importance for earthquake hazard analysis. The H1 and H2 components correspond to the East-West (E-W) and North-South (N-S) directions, respectively. After the data selection process illustrated in the Methods section, we retain 5,108 and 5,301 recordings available for H1 and H2 components, all within the frequency range of [2,15] Hz.

Conditional generative model

We present our CGM-GM framework based on a conditional dynamic VAE framework, as shown in Figure 1(a-d) for generating realistic ground-motion time-series data. The objective is to develop a generative model for predicting ground motions D* at unobserved sources and site locations, with varying earthquake magnitudes, based on actual sparse seismic recordings D. Firstly, in our model, we use Short-Time Fourier Transform (STFT) [44, 45] to decompose the time sequence data x1:T into amplitude and phase information in each time window, where T is the length of the ground motion time series. STFT is an effective technique to extract the time and frequency information, which is also considered in the previous implementation of the GAN model for ground motion generation [36]. Our dynamic VAE model is trained on amplitude information A1:\tau, which is in the time-frequency domain with a new time sample \tau used in STFT (\tau < T). The key methodological contribution lies in the specific design of the prior and posterior distributions in VAE models, where we incorporate temporal dynamics using recurrent neural networks (RNNs) for learning time-frequency information [41, 46]. Next, we employ two strategies to obtain phase information and apply inverse STFT for waveform reconstruction. Specifically, during training, we use the true phase to recover ground motion data and construct a waveform loss to better capture waveform shapes. In the generation stage, phase information is estimated using phase retrieval methods. Furthermore, to generate the earthquakes in which we are interested, we integrate physical parameters (i.e., earthquake magnitudes and depths, geospatial coordinates of sensors and earthquakes) as conditional variables into the CGM-GM framework using Multilayer Perceptron (MLP) layers. These conditional parameters are fundamental for understanding earthquake applications since ground motions vary spatially due to structural heterogeneity and the magnitudes reflect the energy released by rupturing. Figure 2(c) illustrates that we can obtain spatial variations of ground motions at a given scenario of earthquake locations and magnitudes. Another motivation for using geospatial coordinates is to enable neural networks to implicitly learn the physical interactions, such as path, source, and site effects. The rupture distance Rrup and incident angles Aepi can be computed based on the coordinate information (i.e., latitude, longitude, and depth) of earthquake hypocenters and stations of interest. In this study, we have not incorporated variations in focal mechanisms due to our focus on the Hayward fault, where the majority of seismic events exhibit similar focal mechanisms. A more detailed discussion of this aspect is presented in the Discussion section.

Waveforms and spatial continuity

Based on earthquake magnitudes and geospatial coordinates of earthquake sources and stations, the generator part of our CGM-GM framework produces physically consistent ground motion data. Figure 1(e) shows three representative comparisons between the generated ground motion waveforms and the true recordings. The selected cases involve seismic waveforms of similar magnitudes but with varying rupture distances, epicenter-station azimuths, and earthquake depths. Furthermore, we produce two random generations (blue) with the defined conditional variables to compare them with the corresponding observed data (red) in the time domain. The generated ground motion sequences effectively capture waveform shapes, frequency contents, peak values (e.g., waveform amplitude), and the arrival time, even for events with different rupture distances, depths, and magnitudes. For instance, in the first earthquake (M = 2.51), the model successfully captures the moderate amplitude P-wave packet around 14 seconds, followed by the large amplitude S-wave and surface wave packets starting at 22 seconds. Noticeably, the peak amplitudes of wavefields are well generated. More detailed investigations of these peak amplitudes and their spectra are presented in subsequent sections. Another advantage of generative modeling for ground motion waveforms is to perform uncertainty analysis, as discussed in Supplementary Section 2.3. For various earthquake scenarios, the mean curves of generated waveforms can capture dynamic characteristics of ground motion data, and the uncertainty regions show a good coverage of the ground truth. This indicates that our CGM-GM framework is effective and robust for generating earthquake ground motions.

A key aspect of our proposed generative model lies in its ability to approximate ground motions for arbitrary (future) earthquakes and sensor locations. To demonstrate it, we compute FAS maps across a specific region in the SFBA. FAS provides valuable insights into the spatial variability of frequency-dependent ground motion and informs seismic hazard assessments, structural design, and risk mitigation strategies. Specifically, using our CGM-GM framework, we generate the ground motions within a selected sub-region of the SFBA, which spans longitudinally from -122.50\degree to -121.35\degree and latitudinally from 37.27\degree to 38.07\degree. A uniform 100\times100 spatial grid is sampled on a geographic map, yielding 10,000 station coordinates. We use the source parameters of an existing earthquake event that occurred at 3:16 AM on October 21, 2011. Our model is not intended to reproduce ground motions of this specific event; rather, the source parameter information is used as a representative case of an earthquake along the Hayward fault. Utilizing these conditions, 10,000 ground motion instances are generated and the FAS values are computed at each spatial location. The FAS map produced by our generative model at 10 Hz is shown in Figure 2(c). It corresponds to the median across 30 realizations, providing a stable and representative depiction of the model output. The generated median FAS map presents spatial continuity even between station pairs separated by large distances (e.g., in the southeast of the map). Moreover, we observe the FAS decay with respect to distances and its variation with azimuths, which validates the effectiveness of capturing spatial heterogeneity in ground motions by embedding geospatial coordinates of sources and stations into the CGM-GM framework. These characteristics are consistently seen across frequencies ranging from 2 to 15 Hz and under various earthquake scenarios, as detailed in Supplementary Section 2.6.

We conduct a comparative study between the generated results and baseline models specifically tailored for the SFBA [47–49], to assess the performance of our CGM-GM. The selected baseline models integrate methodologies from both ML and seismology. The first method is the empirical ground motion model (GMM), which is built from the observations and widely employed to predict ground motion intensities for seismic hazard analysis [47]. We use a state-of-the-art non-ergodic GMM [49] that incorporates location-specific effects to accurately represent the ground motion intensities. It is specifically built for the SFBA from the same dataset used in this study. The second one, termed CGM-baseline, uses the same CGM-GM architecture but only includes three conditional variables: earthquake magnitudes, source depths, and epicentral distances. Figures 2(a-b) illustrate the FAS results of the non-ergodic empirical GMM and CGM-baseline. The median FAS map generated by the CGM-baseline reveals a radial pattern, attributed to the constrained conditional embedding and the implicit assumption of a radially homogeneous Earth subsurface across the spatial domain. Hence, the contours of this FAS map are circular, showing geometrical spreading and average attenuation effects. In contrast, the CGM-GM and non-ergodic GMM capture more local features by incorporating location-specific factors, which provide a more nuanced representation of the spatial variability in ground motions. For instance, both models predict larger motions in the southern region near San Jose compared to the CGM-baseline predictions. This is reasonable considering that the soft bay mud in the area amplifies the ground motions. However, certain discrepancies are observed in the spatial distributions of FAS maps derived from our generative model and the non-ergodic GMM. The CGM-GM predicts slightly larger motions in the northwest region near San Francisco than the non-ergodic GMM, which is also plausible due to the presence of Bay Mud in the area. The CGM-GM predictions exhibit significant variations across the San Andreas and Hayward faults. The sharp lithological changes lead to substantial differences in ground motions on either side of the faults [50, 51]. The first-order agreement between the CGM-GM and non-ergodic GMM demonstrates the validity of our framework in learning hidden ground motions. In Supplementary Section 2.11, we further provide a detailed analysis of the spatial correlation in the FAS values generated by the CGM-GM model.

To further evaluate the accuracy of all models, we investigate their performance against true ground motion recordings for all earthquake events. Figure 2(d-e) illustrate the averaged FAS differences (Residuals) in the form of a natural logarithm between the true recordings and simulated samples from the non-ergodic GMM, the CGM-baseline, and the CGM-GM model across the entire frequency band of [2,15] Hz. This comparison includes the mean values and their associated uncertainties, represented by one standard deviation (std). Note that our model does not explicitly differentiate between epistemological uncertainty (due to data scarcity) and aleatory variability (arising from randomness of earthquakes) [52]. Overall, our generative model and the non-ergodic GMM show similar residuals. At frequencies between 3 and 11 Hz, the CGM-GM slightly outperforms the non-ergodic GMM, though its performance declines at frequencies above 11 Hz. The uncertainty ranges for both models overlap significantly. These findings indicate that the CGM-GM performs comparably to the state-of-the-art non-ergodic GMM. In comparison to the CGM-baseline, our CGM-GM demonstrates reduced misfits and uncertainty over the entire frequency, confirming the importance of incorporating spatial heterogeneities into generative models.

Amplitude spectra

We further assess the performance of our CGM-GM in the frequency domain, specifically analyzing FAS values and their frequency distributions versus rupture distances. To present a comprehensive evaluation, we select diverse ranges of earthquake magnitudes and geospatial coordinates, generating 100 random samples for each set of conditional variables. For instance, as shown in the top-left part of Figure 3(a), we select the range of earthquake magnitudes in [1.75, 2.25], the rupture distances in [10, 20] km, and the earthquake depths in [5, 10] km, where 63 ground motion sequences in our field dataset. We obtain the statistical values (e.g., 15th percentile, mean, and 85th percentile) for ground truth and generation results using 63 and 6,300 samples, respectively. Figure 3(a) shows the variations of FAS values across diverse earthquake magnitudes and rupture distances. Generally, the FAS results simulated by our CGM-GM model align well with the true FAS curves, especially for mean curves. However, we also observe that the CGM-GM model tends to slightly under-predict at the low-frequency part ([2, 3] Hz) and over-estimate at the high-frequency region ([13, 15] Hz). This is due to the inductive bias of VAE models, where the learned amplitude information is over-smoothing in the logarithmic space. Similarly, this phenomenon is also seen in Figure 3(b), which illustrates a heatmap of frequency distributions versus rupture distances. The error is computed with the division of the natural logarithm of generated FAS f and ground truth f*, which is given by ln(f/f*). The heatmap of the generation results exhibits a generally good alignment with that of the ground truth data. The discrepancies primarily manifest in low-frequency and high-frequency parts. However, the magnitude of logarithmic division errors remains within an acceptable range (approximately [-1.1, 1.2]). Additionally, we provide the comparisons between the generated and true FAS values across different earthquake depths in Supplementary Section 2.4, and discuss the model performance for the H2 component in Supplementary Section 2.5.

P and S arrival times and durations

To further evaluate the quality of synthesized waveforms, we analyze the statistical properties of amplitudes and arrival times of the generated ground motions in the time domain. Statistical analysis of ground motion data is essential for developing predictive models. In this study, we focus on investigating the relationship between rupture distances and three evaluation metrics: the PGV distribution, the peak arrival time (i.e., arrival times of the wavelet of the PGV), and the arrival times of P- and S- waves. Specifically, we implement 100 realizations of ground motion generation with consistent conditional variables from the SFBA dataset. All the statistical analyses are based on these multiple-run generations. Firstly, the PGV distribution across the rupture distances is visualized in Figure 4(a). The solid lines and dashed lines denote the mean curves and the boundaries of mean \pm standard deviation (std), and the dark points present the corresponding PGV values of data samples. We observe that the generated samples effectively capture the distribution of PGV values, including magnitudes and associated uncertainty regions, across diverse rupture distances. It implies that the generations can match the ground truth data well. Moreover, we analyze the arrival of seismic waves. As shown in Figure 4(b), the generated results present a high correlation with the ground truth data in terms of the arrival time versus rupture distances. To gain a deeper insight into the performance of our generative model, we further investigate the performance of capturing the arrival time of P and S waves in the EW direction. Specifically, we use the PhaseNet [53] to pick the arrival time of two types of waves with a probability larger than 50%. Figure 4(c) and (d) show the cross-validations between the generated and true arrivals, including 3,138 and 909 samples for P and S waves, respectively. The result demonstrates a good agreement and validates the effectiveness of our proposed generative modeling method. Additionally, we investigate the durations of the simulated ground motion waveforms. We compute a significant duration D5-95 based on area intensities [54], which is another important ground motion intensity metric in engineering and seismic hazard analysis. Figure 4(e) presents the cross-validations of the durations in the EW direction, which exhibits greater scatter than arrival time results illustrated in Figures 4(c) and (d). It indicates the difficulties of accurately generating and measuring seismic durations. Moreover, we observe a bias primarily for durations within the [0, 3] second range, corresponding to smaller-magnitude earthquakes. Due to lower S/N ratios in such cases, it is challenging to precisely measure the durations from actual seismic data since the relatively high noise level can artificially reduce the measured duration compared to the actual signal duration. Note that the generative model does not explicitly simulate noise, resulting in inherently higher S/N ratios. This potentially amplifies the discrepancy. Nonetheless, despite these limitations, the overall duration trends between observed and generated data demonstrate a good general correlation.

Generalization to larger-magnitude scenarios

Large-magnitude earthquakes produce intense ground shaking across wide areas, posing significant hazards to infrastructure and communities. In this part, we show the magnitude scaling performance of our CGM-GM model compared to empirical GMMs, including a non-ergodic GMM for small-magnitude earthquakes in Section Empirical GMMs and the LA21 model [55] for large-magnitude earthquakes in California. Since our model is only trained on small-magnitude earthquakes (M < 4), directly generalizing to large-magnitude scenarios will be challenging. Therefore, we fine-tune our pre-trained model on two larger-magnitude earthquake data: (i) the 2014 Napa earthquake with a magnitude of 6.02 and depth of 11.0 km, and (ii) the 2018 Berkeley earthquake with a magnitude of 4.38 and depth of 12.3 km. We present the scaling performance up to a magnitude of 8, occurring at a hypocenter depth of 7.94 km, with geographic coordinates of latitude 37\degree51.61N and longitude 122\degree15.61W (the same location as shown in Figure 2). The magnitude is varied from 1 to 8, and waveforms are generated at 73 real-world stations.

The CGM-GM predictions exhibit a reasonable increase in FAS with increasing magnitude, but deviations from the LA21 model remain challenging. First, the CGM-GM predictions remain lower than those predicted by the LA21 model for magnitudes between 4 and 7. This discrepancy may arise because our model is trained mostly on small-magnitude events, where the earthquake source rupture processes can be approximated as point sources. In contrast, large-magnitude events require consideration of finite rupture planes, which introduce deviations in the predictions. Furthermore, the CGM-GM does not account for nonlinear effects, which are typically weak in small-magnitude events but become significant at higher shaking intensities. As a result, the model produces a near-linear increase in FAS with magnitude, unlike the LA21 model, which captures the saturation behavior observed in large-magnitude events. As shown in Supplementary Section 2.10, the fine-tuning strategy improves the model performance by enhancing FAS values, which further supports this observation. To enhance model performance, we plan to explore fine-tuning using larger datasets and incorporate out-of-distribution generalization methods [56]. More details can be found in Supplementary Section 2.10.

Methods
In this section, we present the technical details of our generative model in the context of ground motion simulation. The entire framework is shown in Figure 1(a-d), including the forward and inverse STFT, the network design of dynamic VAE, the sequential model prior, and the embedding of conditional variables. Here, we use the dataset from the H1 component as an illustrative example, and the results on the H2 component are provided in Supplementary Sections 2.2 and 2.5.

Forward and inverse STFT

STFT has been widely applied in audio processing [65, 66] and seismic data analysis [67–69]. It defines a valuable category of time-frequency distributions [70] that describe the amplitude and phase relationships with respect to time and frequency for a signal. This is achieved by repeatedly applying the Fourier transform within specific time windows. Typically, we use sliding time windows with overlaps to capture signals throughout the entire time domain. Therefore, leveraging STFT in generative modeling [71, 72] facilitates the extraction of time and frequency information A from the ground motion sample x. Although the time-frequency resolution of STFT is fixed in the entire time and frequency domains with the chosen time window length, the inverse STFT is relatively stable compared to other 2D spectral decomposition methods such as continuous wavelet transform. The inverse STFT is employed to reconstruct waveforms xp from amplitude spectrograms Ap. During the generation of artificial ground motions, phase retrieval methods are used to estimate the phase information. To be more concrete, these approaches estimate the missing phase information from available amplitude measurements and then recover the timing and shape of seismic waveforms in earthquake ground motion analysis. Various mathematical techniques are leveraged for phase retrievals, such as iterative algorithms and optimization frameworks. For instance, the Griffin-Lim algorithm [73] is widely used thanks to its simplicity and effectiveness. It iteratively refines estimates of seismic wave phases to minimize the discrepancy between the original and reconstructed signals. On the other hand, the Alternating Direction Method of Multipliers (ADMM) [74] is a versatile optimization technique that has gained traction in various scientific domains, including ground motion synthesis [36]. ADMM leverages a convex relaxation framework to decompose complex optimization problems into simpler subproblems and solve the augmented Lagrangian function iteratively. In this paper, we leverage the dynamic VAE model to learn amplitude information and employ the Griffin-Lim algorithm for phase retrieval due to its simplicity. This strategy is chosen over generating phase information with the VAE model, as phase signals are complex and amplitude evolution is relatively smooth over time.

Dynamic VAEs

We consider a dynamic VAE architecture since it is specifically designed to model sequential data with temporal correlations by extending from standard VAEs [41, 46, 75]. Although VAE models have achieved great success in image processing tasks, the absence of explicit temporal modeling in standard VAEs hinders their effectiveness in tackling time series and audio data [76].

The dynamic VAE models focus on a sequence-to-sequence mode for encoding and decoding. Namely, the latent variable is constructed in the form of a temporal sequence instead of a static vector. Let us consider a time sequence x1:T = {xt \in R^N}_{t=1}^T, where T is the sequence length. Dynamic VAE models typically yield a sequence of latent variables z1:T = {zt \in R^l}_{t=1}^T. Therefore, the joint distribution of latent and observed sequences can be reformulated as

p_\theta(x1:T, z1:T) = \prod_{t=1}^T p_\theta(xt | x1:t-1, z1:t) p_\theta(zt | x1:t-1, z1:t-1). (1)

Eq.(1) is a generalized version that describes the generative process in dynamic VAEs. Researchers usually resort to state-space models (SSMs) to simplify the dependencies in conditional distributions of Eq.(1) [76]. One of the most commonly used SSM families is RNNs, which are specifically designed network architectures for handling sequential data and capturing temporal dependencies among data points. However, a significant challenge in training vanilla RNNs is the vanishing and exploding gradient problem. The reason behind this phenomenon is that the gradients can either shrink exponentially (vanishing gradients) or grow exponentially (exploding gradients) when RNNs propagate information through time. This instability hinders effective training and prevents the network from learning long-term dependencies. To alleviate such issues, gating-based RNNs, such as LSTM [77] and GRU [78], are proposed to control the flow of information through the network.

By incorporating the auto-regressive recurrence into dynamic VAE models, the generative process can be simplified as

p_\theta(x1:T, z1:T) = \prod_{t=1}^T p_\theta(xt | zt) p_\theta(zt | z1:t-1). (2)

Moreover, the approximate posterior can be reformulated as

q_\phi(z1:T | x1:T) = \prod_{t=1}^T q_\phi(zt | z1:t-1, x1:t), (3)

where q_\phi(z1:T | x1:T) works as an inference model for the latent sequence z1:T from the observed sequential data.

The training of a dynamic VAE model involves optimizing the evidence lower bound (ELBO) on the marginal likelihood of the observed data x1:T. For a given data sample x1:T, the marginal likelihood is defined as

log p_\theta(x1:T) = D_KL(q_\phi(z1:T | x1:T) || p_\theta(z1:T | x1:T)) + L(\theta, \phi; x1:T) \ge L(\theta, \phi; x1:T), (4)

where D_KL(.) denotes the Kullback-Leibler (KL) divergence between the approximate and true posterior distributions. D_KL(.) is a non-negative term. Therefore, L(\theta, \phi; x1:T) represents the ELBO, which is given by

L(\theta, \phi; x1:T) = E_{q_\phi(z1:T | x1:T)} [log p_\theta(x1:T | z1:T)] - D_KL(q_\phi(z1:T | x1:T) || p_\theta(z1:T)). (5)

The first and second terms on the right-hand side (RHS) are reconstruction loss and a KL-divergence term, respectively. The KL-divergence works as a regularizer for \phi that promotes the approximate posterior q_\phi(z1:T | x1:T) to closely resemble the prior p_\theta(z1:T).

Network design

As shown in Figure 1(b), our dynamic VAE model is trained on amplitude information A1:\tau (\tau < T) to facilitate the learning of time-frequency features of ground motion data. In the CGM-GM framework, we incorporate a GRU layer [78] into the Encoder to learn the dynamics. Subsequently, two MLP layers yield the posterior mean \mu_t and variance \sigma_t at each time stamp t. Hence, the reparameterization trick to produce the posterior sequence z1:\tau can be written as

zt \sim N(\mu_\phi(A_t), \sigma^2_\phi(A_t)). (6)

The sampled latent sequence is fed into the Decoder part to obtain the reconstructed amplitude information Ap1:\tau. In the training stage, we directly leverage the true phase information to reconstruct the time series xp1:T. In the generation stage, we conduct the inverse STFT procedure to generate the artificial ground motion sequence xp1:T by utilizing the phase retrieval method, i.e., the Griffin–Lim algorithm [73]. Furthermore, we design a sequence of Gaussian distributions to serve as a dynamic model prior [79], as shown in Figure 1(c). The mathematical formulation is given by

p_\theta*(zt | z1:t-1) = N(\mu_\theta*(z1:t-1), \sigma^2_\theta*(z1:t-1)). (7)

The model prior distribution (i.e., mean and variance) is learned through a sub-network architecture \theta*, where another GRU layer and two MLP layers are employed for capturing the underlying dynamics. Moreover, to better capture the waveform shapes, we incorporate another waveform loss to the total loss function apart from the reconstruction loss of amplitude information and KL-divergence in Eq.(5). The waveform loss is constructed by calculating the difference between the true waveforms and the recovered counterparts using the true phase information in the training. The final loss function is defined as

L(\theta, \phi; x1:T) = E_{q_\phi(z1:\tau | A1:\tau)} [log p_\theta(A1:\tau | z1:\tau)] - \beta \cdot D_KL(q_\phi(z1:\tau | A1:\tau) || p_\theta(z1:\tau)) + \alpha \cdot E_{q_\phi(z1:\tau | x1:T)} [log p_\theta(x1:T | z1:\tau)]. (8)

Here, \alpha and \beta are two weighting coefficients for waveform loss in the time domain and KL-divergence term. The incorporation of waveform loss facilitates the capturing of peak values and the generation of realistic earthquake ground motion shapes.

Generally, we clarify the training (reconstruction) and inference (generation) processes of our CGM-GM as follows. In the training stage, the dynamic VAE model is optimized to reconstruct input waveforms from their encoded latent representations z1:\tau. This process enables the latent space to learn a meaningful representation of the underlying seismic data. The prior distribution acts as a regularization term through the KL-divergence, enforcing the latent space distribution to remain close to the prior. The reconstruction loss and the KL divergence regularization guide the model in learning a smooth and structured latent space and facilitate effective sampling during the generative process. Unlike the training stage for reconstruction, the inference process (generation) does not rely on existing input waveforms to synthesize ground motion data. During inference, ground motion waveforms are generated by directly sampling from the prior distribution and passing it through the decoder to produce new waveforms. By incorporating conditional variables (e.g., geospatial coordinates, earthquake magnitudes, and depths) into the model prior and the decoder, the entire framework ensures the generation of consistent and realistic ground motion waveforms for arbitrary earthquake scenarios.

Embedding conditional variables

We design an embedding module, consisting of a stack of MLP layers, to integrate the conditional variables into the framework. The illustration of embedding conditional variables is presented in Figure 1(d). To fully encode physical knowledge into the generative process, we apply this shareable embedding module to the encoder, the decoder, and the model prior. For the embedded physical variables, we use earthquake magnitudes and the geospatial coordinates of earthquake sources and stations. The primary rationale for incorporating geospatial coordinates is that many physical properties, such as source and site effects, are inherently coordinate-based and crucial for ground motion simulation. Neural networks can implicitly learn representations of the underlying physics and kinematics with the coordinates as inputs [21, 80]. Note that, in previous papers [35, 36], the researchers also use the information Vs30. Given that Vs30 is an empirical parameter based on geological coordinates, the embedding network should be capable of implicitly capturing these velocity properties with coordinate information as inputs. Vs30 is a proxy correlated with site-specific ground-motion properties, and it is often used in the development of ground-motion models. However, estimations of Vs30 are not always available at all sites and can present significant uncertainties. Thus, we choose not to integrate Vs30 into our model as a parameter and instead estimate local site conditions directly from ground motion data with conditional generative modeling.

For comparison, the CGM-baseline employs earthquake magnitudes, source depths, and epicentral distances as conditional variables. The details of conditional variables for CGM-GM and CGM-baseline are listed in Table 1.

Data selection

A total of 626,423 recordings are collected that spans from 10 seconds prior to the event time to 60 seconds after the event time, leading to each recording of 70 seconds with 7,000 time steps per component. We perform a selection procedure to ensure only recordings with an acceptable S/N ratio are used. Specifically, the first 10 seconds of each recording are considered as noise, while the subsequent 60 seconds are analyzed as earthquake signals. The Fourier transforms of both the noise and the signal are calculated up to the Nyquist frequency of 50 Hz, and the Fourier amplitudes are smoothed using the Konno-Ohmachi window procedure [81] with a smoothing parameter bexp = 20 [81]. The S/N ratio is computed for each recording by comparing the amplitude spectra of the signal and noise. We keep those recordings with an S/N ratio exceeding 3 across the frequency range of 2 to 15 Hz. Hence, approximately 15,000 recordings per component are retained. Additionally, we conduct a rigorous visual inspection for each time series to exclude those with equivalent noise levels in the first 10 seconds and the last 60 seconds. This process results in 5,108 and 5,301 recordings available for horizontal components H1 and H2, respectively. After selecting the final recordings, a bandpass filter is applied over the frequency range of interest in [2, 15] Hz.

Figure 6(a) and (b) present the distribution of earthquake depth and the scatter plot of the corresponding magnitudes and rupture distances, respectively. The locations of selected stations and events for the H1 component are shown in Figure 6(c) and (d). Note that our methods are applicable to diverse earthquake ground motion datasets (i.e., different magnitude ranges) and geographical regions. More complete information about the dataset can be found in [82].

Implementations

This part includes the data preprocessing and the training implementations. The last 60 seconds of the SFBA dataset are selected for ground motion generation. The dataset is then split into [80%, 20%] for training and testing, respectively. Namely, the number of training samples is 4086 for the experiments. Firstly, we use STFT to obtain the time-frequency information from the ground motion data since the CGM-GM model focuses on learning the amplitude spectrograms. The window length and hop length are set as 160 and 46, respectively. Therefore, the amplitude spectrogram has a size of 81\times131, where 81 is the frequency range and 131 denotes the sequence length. The window length of 160 is equivalent to 1.6 seconds. We use signals between [2, 15] Hz for our analysis, and this window length contains three wavelets for the lowest frequency to resolve signals at this frequency. The hop length is more arbitrary, and Welch [83] recommends using a length shorter than half of the window length. Due to the computational cost, we use 46 in this study. We add a minimum threshold of 10^-10 for amplitude spectrograms and then convert them into the logarithmic space. Furthermore, the logarithmic time-frequency coefficients are normalized between 0 and 1. For conditional variables, we consider the earthquake magnitudes and geospatial coordinates of sources and stations. Each variable is scaled within [0, 1] separately.

For the model parameters in CGM-GM, we consider a 3-layer MLP with feature sizes of [32, 32, 16] for embedding conditional variables. The label embedding module, as shown in Figure 1(d), is kept fixed for different parts of the network. For the Encoder component, we use one GRU layer with a hidden dimension of 144 and two independent 2-layer MLPs with dimensions of [64, 32] to obtain the sequential latent variables. Note that the outputs from the encoder part, \mu and \sigma, have the same sequence length \tau as the input amplitude spectrogram A1:\tau. For the decoder, a stack of MLP layers with feature sizes of [128, 128, 81] is employed for reconstructing the amplitude information. For the sequential prior model, we use one GRU layer with a hidden dimension of 32 and two independent linear layers to get the sequential prior variables. This model architecture contains 0.17 million parameters.

We use the Adam [84] optimizer to train the proposed method for 5,000 epochs. The learning rate is set as 8\times10^-4 initially and decays every 100 steps with a ratio of 0.99. The weight decay parameter is set as 5\times10^-5 and the batch size is 128. Furthermore, we leverage the grid search to select the hyper-parameters in the loss function. The optimal parameters are defined as \alpha = 0.5 and \beta = 0.003.

Empirical GMMs

Empirical ergodic and non-ergodic GMMs for velocity FAS values are developed to evaluate the ground motions generated by our CGM-GM. Specifically, ergodic GMMs focus on the average scaling of ground motions, and non-ergodic GMMs account for the spatial distribution of ground motions due to path effects related to the 3-D velocity structure. For the ergodic GMM, the functional formulation of the natural log of FAS values, i.e., ln(Y), is given by

ln(Y; M, Rrup, Dhyp) = \alpha0 + \alpha1 M + \alpha2 ln(Rrup) + \alpha3 Dhyp + \deltaBe + \deltaS2Se + \deltaPe + \deltaWSe, (9)

where \alpha1 M, \alpha2 ln(Rrup), \alpha3 Dhyp denote magnitude scaling, geometrical spreading, and depth scaling terms, respectively. \deltaBe, \deltaS2Se, \deltaPe represent the source, site, and path effects with zero-mean normal distributions. \deltaWSe is the within-site residual, which is also assumed to be normally distributed. The coefficients are derived through linear regression, ensuring a smooth spectrum and imposing physical constraints on the coefficients [85]. Moreover, the total standard deviation of ergodic GMMs is defined as,

\sigma_fas = \sqrt{\tau^2 + \phi^2_S2S + \phi^2_SS}, (10)

where \tau, \phi_S2S, \phi_SS denote the standard deviation of the mixed-effects coefficients \deltaBe, \deltaS2Se, \deltaWSe, respectively. The parameters of the ergodic GMM in the SFBA at different frequency bands (2, 5, 10, and 15 Hz) are detailed in Table 2.

For the non-ergodic GMMs, the source, site, and path terms are considered spatially dependent, which are the functions of the coordinates of sources and sites. Therefore, the non-ergodic GMM is re-written as [49],

ln(Y; M, Rrup, Dhyp, ..., te, ts) = LANNP26_Adj-erg(M, Rrup, Dhyp) + \deltaL2L(te) + \deltaS2S(ts) + \deltaP2P(te, ts) + \deltaB0 + \deltaWSP_es, (11)

where LANNP26_Adj-erg is the median from the ergodic model as shown in Eq. 9, \deltaL2L(te), \deltaS2S(ts), \deltaP2P(te, ts) are the median shifts of the source, site, and path terms. te and ts represent the earthquake and site locations. Additionally, the term \deltaB0 + \deltaWSP_es denotes the aleatory variability apart from the systematic source, site, and path effects. The Gaussian Process (GP) regression is leveraged to fit the available ground motion data within the SFBA dataset by providing the medians and epistemic uncertainty [49]. The non-ergodic GMMs simulate FAS maps thanks to the capability of spatial interpolation of GP. Furthermore, let the \tau0 and \phi_SP represent the standard deviations of \deltaB0 and \deltaWSP_es, respectively. The total aleatory variance of non-ergodic GMMs is formulated as,

\sigma^2_NE = \tau^2_0 + \phi^2_SP. (12)

We use the same earthquake scenario and 100\times100 spatial coordinates of stations to synthesize the non-ergodic FAS maps.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "Review of manuscript of “Learning Physics for Unveiling Hidden Earthquake Ground Motions via Conditional Generative Modeling”\nThe manuscript presents the results of conditional generative models for ground motion simulation. The approach is based on a dynamic variational autoencoder, which does not have the drawbacks and difficulties of CGANs. The authors demonstrate promising results of this approach for scenarios with small magnitude in the San Francisco Bay area. The author investigated the potential of the proposed approach for capturing different characteristics of ground shaking.\nI think the proposed approach and the results are interesting. In particular, the proposed method has a significant application in the improvement of shaking maps in early warning systems and seismic hazard analysis for regions with low magnitude. I believe this type of generative model is part of the next generation of ground motion modeling and simulations. I suggest publishing the manuscript after revision.\nHere are my comments:\n1) Page 4, Line 8: In signal processing, it is not common to use the Short-Term Fourier Transform. The STFT stands for the short-time Fourier transform.\n2) I enjoyed Figure 04 and how the proposed method captures the arrival time. However, it is important to investigate the duration of ground shaking (e.g., significant duration) of ground-shaking signals between real and generated data. This investigation is very important for investigating shaking in geotechnical and seismic hazard analysis.\n3) In the studied range of magnitude, the source is assumed to be a point source, and in this case, the hypocenter and rupture distances are the same. You mentioned Rhyp but consider it a rupture distance. It is better to be consistent and use one of them.\n4) You considered Rhyp and depth together as conditions. I completely understand why depth is important, but this parameter is part of Rhyp, too. Why did you not use Repi instead of Rhyp?\n5) Figure 02(d,e): The Y-axis limits for these two plots are not consistent, making it difficult to compare them.\n6) Figure 03: The label of the colorbar for Error is missing (ln(FAS/hat{FAS}).\n7) The author did not explain why they consider epicenter-station azimuths as a condition parameter in the generative model. Does this parameter capture the anisotropy of media around the source and the path effect?\n8) The author should mention that the proposed method and incorporated parameters can not capture the site effects.\n9) Page 11: The author mentioned that using the time domain instead of the time-frequency domain can improve the generative model since it no longer requires phase retrieval. I agree with the intrinsic bias of phase retrieval because of the ill-posedness of this problem. However, I think using the time domain directly is more challenging since the waveform is nonstationary and requires decomposing the signal into different frequencies. It is just my opinion.\n10) I think the introduction and literature on ground motion simulations are incomplete. There are approaches that simulate ground shaking in the time-frequency domain based on a stochastic approach (Sabetta et al. 2021).\n11) Figure 1, 6, 7, and 8: The time axes in simulated records are different, which is misleading. I suggest modifying it to the same maximum time for all seismic records or mentioning in the capture that the records are limited to a specific duration.\nReference:\nSabetta, F., Pugliese, A., Fiorentino, G. et al. Simulation of non-stationary stochastic ground motions based on recent Italian earthquakes. Bull Earthquake Eng 19, 3287–3315 (2021). https://doi.org/10.1007/s10518-021-01077-1\n(Remarks on code availability)\nI was curious about the code and the approach. The code does work and the results are reproducable.", "reviewer_2": "The authors present a Variational Autoencoder (VAE) referred to as conditional generative modeling framework for ground motion (CGM-GM) that learns sequences of short-term Fourier Transforms of earthquake time series, conditioned on the coordinates of the earthquake source and receivers, and other physical parameters such as the earthquake magnitude and hypocentral depth.\nThe authors argue, correctly, that the stochastic approach of generating ground motion time series is computationally efficient but lacks the underlying physics framework and leads to waveforms that do not represent realistic arrivals of the different seismic phases. On the other hand, physics-based simulations are physically sound but are computationally taxing and carry significant uncertainties in the representation of the source, path and nonlinear, heterogeneous site effects alike.\nMachine learning approaches such as the model presented by the authors have emerged in the past 5 years, taking advantage of the advances in artificial intelligence algorithms and the proliferation of simulated and recorded ground motion data. The authors, however, use erroneous statements to demonstrate the need for their model compared to alternative machine learning methodologies for ground motion generation. Furthermore, their results do not support their claims on being able to simulate path effects and high frequency ground motions. I will elaborate below:\n1.The authors claim that the model can generate [spatially continuous earthquake ground motion waveforms]; however, the methodology described in the text does not support this claim adequately. From the framework outlined, it is unclear how the model learns the spatial correlations between different stations for a given event. The inference process described suggests that the VAE generates 1D time histories independently for each station, and there is no detailed description of how a conditional 1D model can capture the complex 3D spatial-temporal functions inherent in earthquake data.\n2. If the model is indeed learning the path effects, this should be demonstrated by evaluating multiple scenarios similar to the one in Figure 2c, and evaluating the uncertainty along paths where there are no data. In the Gaussian Process framework that non-ergodic models are build upon, lack of data translate into having the full epistemic uncertainty in the adjustment factors. If the uncertainty that the authors get is constant in areas with no data, then what they are learning is not path effects, but some averaged effects for a range of sources (and not the true wave propagation physics along the underlying medium).\n3. The paper doesn’t clearly explain how the VAE model generates synthetic ground motion samples. As far as I understand it, there are two ways of using VAE to generate time series:\n3.1 The first method involves directly sampling from the latent space, which is typically assumed to be a standard Gaussian distribution, followed by passing these samples through the decoder. This approach aligns with the sampling processes used by other generative models such as GANs and Diffusion Models. If this is the case, then the model is at its core like the GAN models that the authors criticize in the introduction.\n3.2 The second method involves reconstructing data from latent samples, where an input x (waveform sample) is first encoded to a latent representation z, followed by decoding z back to x. This method, however, does not allow for the generation of ground motion for arbitrary (future) earthquakes and sensor locations, as it relies on having existing data to reconstruct from. If this is the case, then the model is not able to generate time series at any arbitraty spatial location, let alone at spatially continuous wavefields.\n4. Following up on 3.1, the paper dismisses the use of GAN citing that \"GAN models cannot capture 3D wave propagation effects due to limited input conditional variables\" is not a limitation of GANs. Rather, the ability of GANs to model complex phenomena like 3D wave propagation largely depends on the design of the input space and the conditioning of the model. The inclusion of more or appropriate conditional variables could potentially mitigate this limitation, suggesting that the model's design, rather than the model type, influences its capability. This is one of many erroneous statement examples\n5. VAEs are known to generate samples of lower quality compared to GAN and Diffusion models. I cannot see the advantage of using the conditional VAE model compared to other-state of-the-art generative models. While GAN and Diffusion models are shown to have universal approximation ability under mild conditions, the universal approximation ability of VAE is constrained by the reliance on the Gaussian assumption in the latent space. The training objective of VAE relies on KL divergence, which encourages the latent space to match a predefined distribution (typically Gaussian). This regularization can sometimes be too strong, leading to posterior collapse and not covering all modes of the data distribution. To circumvent issues with training stability of GANs, Diffusion models such as DDPM have significantly improved their training and inference speeds and can be viewed as an extension of VAEs that also mitigates earlier concerns about training costs. Lastly, the conditional VAE model relies on a phase retrieval method that has been shown to perform worse that learning waveforms directly in the time domain, similarly to the approach explored by GAN models.\n6. The spike-like issues observed in Figures 6 and 7 appear to introduce high-frequency artifacts, which are not observed in previous machine-learning-based ground motion synthesis models. The authors should conduct an ablation study to determine whether these artifacts are attributable to the chosen model itself or are a result of the conditional variables used. This study could help clarify the source of the artifacts and guide improvements in the model design or variable selection. Knowing that results deviate from the median of GMMs and NGMMs in the high frequency regime, as evidenced also by the statistics in Figure 2 d and e, makes the generation of FAS maps at 10Hz irrelevant; the maps that would be useful to see are the <1Hz maps and their uncertainty, that show the capability of the model to capture path effects, if it does, as the authors claim.\n7. The selection of baseline models for validation in this study appears to be inadequate. Fine-tuned GANO models such as the recently published version of [32] should be included in the comparison.\n8. The model that is trained on small magnitude data, should also show the capability for reasonable scaling, with increased uncertainty, for larger magnitude events. If it doesn’t, then there is no value in having a model like this for California.\n(Remarks on code availability)", "reviewer_3": "The authors of the manuscript entitled \"Learning Physics for Unveiling Hidden Earthquake Ground Motions via Conditional Generative Modeling\" proposed a dynamic VAE, featured by variational sequential architectures (e.g., RNNs) in both prior and posterior distributions, which facilitates the capture of the temporal evolution (dynamics) of time series. Thanks to this Conditional Generative Modeling for Ground Motion (CGM-GM) strategy, the authors synthesize high-frequency and spatially continuous earthquake ground motion scenarios in the San Francisco Bay Area. Interestingly, CGM-GM embeds earthquake magnitude and hypocenter coordinates as well as sensors' geographic coordinates, so to reproduce the natural ground motion variability. CGM-GM was trained with small-magnitude data, which is remarkable because they challenges small Signal-to-Noise Ratios.\nThe manuscript is very well written. The presented results convincing. CGM-GM framework is effective and robust for generating earthquake ground motions. Its lightweight architecture highly appreciated. The most interesting aspect of our proposed generative model lies in its claimed ability to approximate ground motions for arbitrary (future) earthquakes and sensor locations.\nSome comments:\n- In the Introduction, the authors overlooked some recent literature on earthquake ground motion generation using neural operators and graph neural networks. These represent a credible alternative to the CGM-GM presented in the manuscript and it would interesting to highlight what are the advantages and limitations of adopting CGM-GM instead of one of these other approaches.\n- What steps and data would the authors deem necessary to apply CGM-GM to another arbitrary region of the world? The title entices the readership but it seems that CGM-GM would be hardly adoptable in regions with poor seismic record and/or varying focal mechanisms. It is in this reviewer's opinion that the title should reflect the claim the authors make that their primary interest is in generating waveforms along the Hayward Fault, with no domain shift for the earthquake rupturing processes.\n- Did the authors consider using ground motion numerical simulation to expand their training database? If so, how would they mix synthetics with records within the learning framework?\n- As per this reviewer's understanding, each 3-component single-station time-history is encoded separately, along with its geographical coordinate. How is spatial correlation encoded into the latent space, if so? Peak values contour plots seem rather realistic, but what about the whole time series? This question arises from the fact that the phase spectra and the ground motion coherence are of paramount importance in earthquake engineering and seismic hazard modelling.\n- The readership would benefit from a Figure showing extensive goodness-of-fit metrics comparisons. This reviewer suggests to adopt the modified Anderson's criteria, found in Olsen et al., 2010. Other GoF metrics would probably mislead the readership, because CGM-GM is conceived for stochastic generator of ground motions, not for historical ground motion reconstruction.\nMinor comments and questions are found in the attached pdf.\n(Remarks on code availability)\nData report is not available, nor data themselves. Testing failed because of the lack of test dataset.\nI suggest to add some demonstrative notebook/tutorial."}, "author_response": "Dear Editor and Reviewers,\nWe sincerely appreciate the insightful and constructive comments from the reviewers, which have helped improve our manuscript. We are pleased that all reviewers acknowledge the significance of our research and see the merit of our manuscript. In particular, we thank the reviewers for recognizing the effectiveness (Reviewers 1 and 3), robustness (Reviewer 3), and computational efficiency (Reviewers 2 and 3) of our proposed method in synthesizing arbitrary small-magnitude earthquakes. We extend our gratitude to Reviewer 3 for commending the overall quality of the manuscript.\nOur point-by-point responses are provided below marked in blue and the corresponding manuscript revisions are highlighted in red. Additionally, we have updated all figures in the revised manuscript, as we made minor adjustments to the training hyperparameters to ensure optimal network performance. These modifications did not change the interpretation of the results. We hope our manuscript will be found worthy of publication in this revised form.\nI. RESPONSE TO REVIEWER 1\nGeneral comment: The manuscript presents the results of conditional generative models for ground motion simulation. The approach is based on a dynamic variational autoencoder, which does not have the drawbacks and difficulties of CGANs. The authors demonstrate promising results of this approach for scenarios with small magnitude in the San Francisco Bay area. The author investigated the potential of the proposed approach for capturing different characteristics of ground shaking.\nI think the proposed approach and the results are interesting. In particular, the proposed method has a significant application in the improvement of shaking maps in early warning systems and seismic hazard analysis for regions with low magnitude. I believe this type of generative model is part of the next generation of ground motion modeling and simulations. I suggest publishing the manuscript after revision.\nReply: Thank you very much for your positive feedback and recommendation for publishing after revision. We will address your concerns point-by-point as follows.\nComment 1: Page 4, Line 8: In signal processing, it is not common to use the Short-Term Fourier Transform. The STFT stands for the short-time Fourier transform.\nReply: Thank you for pointing this out. We have updated the terminology to the commonly-used “short-timeFouriertransform”, assuggested. ThemodificationscanbefoundonPage4oftherevisedmanuscript.\nComment 2: I enjoyed Figure 04 and how the proposed method captures the arrival time. However, it is important to investigate the duration of ground shaking (e.g., significant duration) of ground-shaking signals between real and generated data. This investigation is very important for investigating shaking in geotechnical and seismic hazard analysis.\nReply: We agree with the reviewer that the shaking duration is also a critical factor in evaluating the simulated ground motions. We have included additional evaluations of the shaking durations for all simulated ground motion waveforms, as presented in Figure 4(e) of the revised manuscript (Page 11). This figure provides cross-validations of the shaking durations in the EW direction, demonstrating a good agreement between the generated samples and the ground truth. The duration in the NS direction shows a similar trend. The corresponding discussions have been added in Section P and S arrival times and durations (see Page 10).\nFIG. 1. The cross-validation of the shaking durations in the EW direction. The density is normalized between 0 and 1.\nComment 3: In the studied range of magnitude, the source is assumed to be a point source, and in this case, the hypocenter and rupture distances are the same. You mentioned Rhyp but consider it a rupture distance. It is better to be consistent and use one of them.\nReply: Thank you for this suggestion. You are correct that the source is treated as a point source within the studied magnitude range. We have revised R to R to represent the rupture distance hyp rup throughout the paper. Meanwhile, please also see the updates in Figures 1e, 3a, and 7–12.\nComment 4: You considered Rhyp and depth together as conditions. I completely understand why depth is important, but this parameter is part of Rhyp, too. Why did you not use Repi instead of Rhyp?\nReply: We agree that depth is inherently included in R . To ensure clarity and consistency, we hyp have replaced “rupture distance” with “epicentral distance” as one of the conditional variables in the CGM-baseline model (see Page 6). We have also updated the results with the latest setup in the revised manuscript (see Figure 2).\nComment 5: Figure 02(d,e): The Y-axis limits for these two plots are not consistent, making it difficult to compare them.\nReply: Thankyouforbringingthistoourattention. WehaveadjustedtheY-axislimitsinFigure2(d) and (e) to ensure consistency. The modifications can be found on Page 7 of the revised manuscript.\nComment 6: Figure 03: The label of the colorbar for Error is missing (ln(FAS/hatFAS).\nReply: Thanks for pointing it out. We have included the label of the colorbar for Figure 3(b). Please see the modifications on Page 9.\nComment 7: The author did not explain why they consider epicenter-station azimuths as a condition parameter in the generative model. Does this parameter capture the anisotropy of media around the source and the path effect?\nReply: Thankyouforthisquestion. Firstofall,wewouldliketoclarifythatepicenter-stationazimuths were not included as a conditional parameter in our generative model. In the CGM-GM framework, we utilize earthquake magnitude, earthquake depth, and the geospatial coordinates of earthquake sources and seismic stations as conditional variables. For comparison, the CGM-baseline employs earthquake magnitudes, sourcedepths, andepicentraldistancesasconditionalvariables. Thedetailsofconditional variables for CGM-GM and CGM-baseline are listed in TABLE (I). Our coordinate-based parameters allow us to effectively capture underlying heterogeneity. While we do not have evidence of strong anisotropy in the SFBA (for example, the USGS San Francisco Bay region 3D seismic velocity model (SFVM) does not include anisotropy, but reasonably explains observed waveforms), our conditional variable setup allows us to capture anisotropy.\nSecondly, we would like to elaborate on the rationale behind using two different sets of conditional variables. The variables selected for the CGM-baseline model represent primarily 1D phenomena, whereas those for the CGM-GM model reflect 3D wave physics and subsurface properties, enabling a more comprehensive representation of seismic events. We have included an additional introduction in Section Embedding conditional variables (see Page 17).\nTABLE I. Conditional variables used in CGM models.\nModel Conditional Variables\nCGM-baseline earthquake magnitudes, source depths, and epicentral distances\nCGM-GM earthquake magnitudes, geospatial coordinates of earthquake sources and stations\nComment 8: Theauthorshouldmentionthattheproposedmethodandincorporatedparameters can not capture the site effects.\nReply: Thanks for this suggestion. We would like to clarify that the primary objective of our proposed method is to generate spatially continuous waveforms. Our model captures wave physics by implicitly learningacombinationofsource, site, andpatheffectswithintheconditionalgenerativeAIframework. We didn’t explicitly learn these effects individually. The site effect is regarded as a part of this learned representation. Whileourcurrentframeworkdoesnotseparatesiteeffectsasindependentcomponents, addressingthislimitationanddevelopingmethodsfortheirexplicitcharacterizationwillbeanimportant directionforfutureresearch. Wehaverevisedtherelevantsectionsintheupdatedmanuscripttoclarify this point. Please refer to Section DISCUSSION on Pages 12-13.\nComment 9: Page 11: The author mentioned that using the time domain instead of the time- frequency domain can improve the generative model since it no longer requires phase retrieval. I agree with the intrinsic bias of phase retrieval because of the ill-posedness of this problem. However, I think using the time domain directly is more challenging since the waveform is non- stationary and requires decomposing the signal into different frequencies. It is just my opinion.\nReply: This is an insightful comment. We agree with the reviewer that directly modeling waveforms in the time domain presents significant challenges, particularly in capturing long-range dependencies in time series data. For instance, recurrent neural networks often encounter vanishing and exploding gradientproblems, whichhindertheirabilitytolearncomplexpatternsandlong-termdependencies[1]. However, there is still some prior work in machine learning that successfully handles audio waveforms [2]andlong-sequencetimeseriesdata[3]inthetimedomain. Thus,whilemodelingdirectlyinthetime domain is indeed challenging, we believe it represents a worthwhile direction for further exploration.\nReferences:\n[1] Naiman, Ilan, et al. “Utilizing Image Transforms and Diffusion Models for Generative Modeling of Short and Long Time Series.” in The Thirty-Eighth Annual Conference on Neural Information Processing Systems, 2024.\n[2] Van Den Oord, Aaron, et al. “Wavenet: A generative model for raw audio.” in The 9th ISCA Speech Synthesis Workshop, 2016.\n[3] Zhou, Linqi, et al. “Deep latent state space models for time-series generation.” in International Conference on Machine Learning, 2023.\nComment 10: I think the introduction and literature on ground motion simulations are incom- plete. There are approaches that simulate ground shaking in the time-frequency domain based on a stochastic approach (Sabetta et al. 2021).\nReference: Sabetta, F., Pugliese, A., Fiorentino, G. et al. Simulation of non-stationary stochastic ground motions based on recent Italian earthquakes. Bull Earthquake Eng 19, 3287–3315 (2021). https://doi.org/10.1007/s10518-021-01077-1\nReply: Thank you for recommending this paper. We agree that it provides valuable insights and will enrich the current literature review on stochastic methods for simulating ground motions. Accordingly, wehaveincludedthispaperintherevisedmanuscript. PleaseseePage2(reference[4])fortheupdated citation.\nComment 11: Figure 1, 6, 7, and 8: The time axes in simulated records are different, which is misleading. Isuggestmodifyingittothesamemaximumtimeforallseismicrecordsormentioning in the capture that the records are limited to a specific duration.\nReply: Thank you for this comment. In the figures mentioned, the time axes reflect the motion duration for each seismic record, which may vary depending on the specific event. We aim to highlight the most relevant portions of the waveforms for analysis. However, we acknowledge that this could potentiallycausemisunderstanding. Toaddressthis, wehaveaddedclarificationsinthefigurecaptions (see Pages 5 and 33), specifying that the displayed waveforms are limited to their respective motion durations.\nComment 12: Remarks on code availability: I was curious about the code and the approach. The code does work and the results are reproducible.\nReply: Thank you for reviewing our code. We have added additional tutorial examples, including guidance on training models and performing evaluations from various perspectives. Please see the updated materials in our GitHub repository: https://github.com/paulpuren/cgm-gm.\nII. RESPONSE TO REVIEWER 2\nGeneral comment: The authors present a Variational Autoencoder (VAE) referred to as con- ditional generative modeling framework for ground motion (CGM-GM) that learns sequences of short-term Fourier Transforms of earthquake time series, conditioned on the coordinates of the earthquake source and receivers, and other physical parameters such as the earthquake magnitude and hypocentral depth.\nThe authors argue, correctly, that the stochastic approach of generating ground motion time series is computationally efficient but lacks the underlying physics framework and leads to wave- forms that do not represent realistic arrivals of the different seismic phases. On the other hand, physics-based simulations are physically sound but are computationally taxing and carry signif- icant uncertainties in the representation of the source, path and nonlinear, heterogeneous site effects alike.\nMachinelearningapproachessuchasthemodelpresentedbytheauthorshaveemergedinthepast 5 years, taking advantage of the advances in artificial intelligence algorithms and the proliferation of simulated and recorded ground motion data. The authors, however, use erroneous statements to demonstrate the need for their model compared to alternative machine learning methodologies for ground motion generation. Furthermore, their results do not support their claims on being able to simulate path effects and high frequency ground motions.\nReply: Thank you for your thoughtful feedback. We appreciate the reviewer’s recognition of our general arguments and the context of our proposed approach. We understand that certain aspects of theoriginalmanuscriptmayhaveledtosomemisunderstandings. Weprovideadetailedpoint-by-point response below to address your concerns.\nComment 1: The authors claim that the model can generate [spatially continuous earthquake ground motion waveforms]; however, the methodology described in the text does not support this claim adequately. From the framework outlined, it is unclear how the model learns the spatial correlations between different stations for a given event. The inference process described suggests that the VAE generates 1D time histories independently for each station, and there is no detailed description of how a conditional 1D model can capture the complex 3D spatial-temporal functions inherent in earthquake data.\nReply: Thank you for raising this important point. To address your concerns, we would like to clarify how our CGM-GM framework captures spatial correlations and generates spatially continuous earthquake ground motion waveforms from two perspectives.\nCapturing 3D spatiotemporal patterns:Whilethegeneratedwaveformsare1Dtimeseries, thelearning process is inherently 3D due to the conditional variables embedded in the model. The embedding module incorporates the 3D setup of geospatial coordinates, including longitude, latitude, and depth of the source and site locations. This allows our CGM-GM model to learn and capture the complex spatiotemporal patterns of earthquake data. To be more concrete, let us consider a general case of predicting seismic waveforms at a station x generated from a source x , which is given by r s\nF(t,x ,x ,m) = S(t,x )∗R(t,x )∗G(t,x ,x ), (1)\ns r s r s r\nwhere F(·) denotes the prediction function that generates waveforms. S(·), R(·), and G(·) represent the source effects, site effects, and Green’s function (or path effects), respectively. ∗ is the convo- lution. In our conditional generative modeling framework, this function F(·) is approximated by the neural networks. Specifically, we incorporate physical variables (i.e., earthquake magnitudes, depths, source and receiver locations) as conditional variables v. In our implementation, x represents source s parameters, including magnitudes, source locations, and depths. The CGM-GM framework incorpo- rates source effects, such as source mechanisms and their coupling with near-source geology. Note that this leads to implicitly assuming a specific source mechanism at each source location, which is reasonable. The conditional variable of a receiver location corresponds to x and we incorporate the r site effect R(x ,t). We do not separately model S(·), R(·) and G(·), and hence these effects are r learned as a combination through the neural network. We interpret that using many earthquake events allows us to learn R(·) at a specific location, and using many stations lets us learn S(·) for a specific event. Similarly, by using many earthquake-station pairs, we learn Green’s function G(·). We have in- cluded additional clarifications in the revised manuscript. Please see Section A.1 Details of conditional variables on Pages 30-31.\nSpatial continuity: As shown in Figure 4 of the revised manuscript (Page 11), our generated 1D time series exhibit waveform characteristics, including arrival times, peak ground velocities (PGVs), and durations. They closely match the observations across numerous realizations. Furthermore, we demonstrate the reproduction of key features, such as P and S arrival times through cross-plots of their values. These results indicate that 3D wave propagation, source, and site effects are properly captured in a non-ergodic manner. The spatial continuity of ground motions generated by our model is demonstrated in the FAS map at 10 Hz, as shown in Figure 2(c). The generated FAS map exhibits spatially continuous patterns, even between station pairs separated by significant distances (e.g., in the southeast of the map). This continuity reflects the capacity of our CGM-GM framework to capture underlying spatial heterogeneities and wave physics inherent in earthquake ground motions. However, we acknowledge that the continuity of details of waveforms (reflections and scatterings) is still under investigationandsubjectofthefuturestudy. AdetaileddiscussiononthisaspectisprovidedinSection Waveforms and spatial continuity. Please see Pages 6-8 of the revised manuscript.\nComment 2: If the model is indeed learning the path effects, this should be demonstrated by evaluating multiple scenarios similar to the one in Figure 2c, and evaluating the uncertainty along paths where there are no data. In the Gaussian Process framework that non-ergodic models are build upon, lack of data translate into having the full epistemic uncertainty in the adjustment factors. If the uncertainty that the authors get is constant in areas with no data, then what they are learning is not path effects, but some averaged effects for a range of sources (and not the true wave propagation physics along the underlying medium).\nReply: First of all, we would like to clarify that our proposed model captures a combination of site, source, and path effects within a conditional generative framework. The main goal of our approach is togeneratespatiallycontinuousgroundmotionsratherthanexplicitlylearningtheseeffectsasseparate components. Tofurtherillustratethis, wepresentapseudo-patheffect, definedastheresidualbetween the predicted FAS map and the summation of the ergodic effect, site effect, and source effect obtained from the non-ergodic Ground Motion Model (GMM). The comparison of path effect from non-ergodic GMM and pseudo path effect from CGM-GM is presented in Fig. 2. The results indicate that the pseudo-patheffectderivedfromourmodelgenerallyalignswithpredictionsfromthenon-ergodicGMM. We recognize the significance of further disentangling site, source, and path effects and will explore this in future work. We have included more clarifications on this aspect in the revised manuscript. Please see Section DISCUSSION on Page 13.\nFurthermore,weemphasizethatourgoalisnottoreplacetheGaussianProcessframeworkusedinnon- ergodic GMMs with generative modeling. Instead, we view generative modeling as a complementary approach that can provide additional insights and serve as a promising tool alongside existing GMM methodologies.\nTheCGM-GMmodeldoesnotseparateepistemicandaleatoryuncertainties. Oneofourrecentstudies [1] showed that the ground motion variability depends on the hyperparameter settings. We consider that exploring uncertainties is a necessary future direction and warrants a thorough investigation, considering the importance of epistemic and aleatory uncertainties in seismic hazard analysis. Our additional tests conducted for this revision (magnitude scaling for Comment 8 from Reviewer 2 and limited data for Comment 2 from Reviewer 3) provide important initial insights into uncertainties of the CGM-GM. First, the comparisons with the GMMs at the available stations show that the aleatory variability of our model matches that of the ergodic GMM. Second, the increase in the range of predictions with decreasing training data shows that epistemic uncertainty is captured in the CGM- GM. These results are an initial path to future works on uncertainty quantification that will make the CGM-GM more practical and important for future seismic hazard analysis. We have added more clarifications in Section Waveforms and spatial continuity of the revised manuscript (Pages 7-8).\nNon-ergodic GMM CGM-GM\nFIG. 2. A comparison of path effect from non-ergodic GMM and pseudo-path effect from CGM-GM.\nReferences:\n[1] Nakata, Rie, et al. ”Simulating seismic wavefields using generative artificial intelligence.” The Leading Edge 44.2 (2025): 123-132.\nComment 3: The paper doesn’t clearly explain how the VAE model generates synthetic ground motion samples. As far as I understand it, there are two ways of using VAE to generate time series:\n• Thefirstmethodinvolvesdirectlysamplingfromthelatentspace, whichistypicallyassumed to be a standard Gaussian distribution, followed by passing these samples through the de- coder. This approach aligns with the sampling processes used by other generative models such as GANs and Diffusion Models. If this is the case, then the model is at its core like the GAN models that the authors criticize in the introduction.\n• The second method involves reconstructing data from latent samples, where an input x (waveform sample) is first encoded to a latent representation z, followed by decoding z back tox. Thismethod, however, doesnotallowforthegenerationofgroundmotionforarbitrary (future) earthquakes and sensor locations, as it relies on having existing data to reconstruct from. If this is the case, then the model is not able to generate time series at any arbitraty spatial location, let alone at spatially continuous wavefields.\nReply: First of all, we would like to clarify that our intention was not to criticize the application of GANs to earthquake ground motion generation but rather to mention the inherent challenges in their trainingprocesses. WeacknowledgethatGANshaveachievedsignificantsuccessinseismicapplications buttheirtrainingoftenrequirestrickstoaddressinstabilityissues[1]. VAEsofferamorestabletraining framework and are relatively more accessible to scientists without extensive deep learning expertise, and they complement previous generative approaches. VAEs provide a straightforward implementation while retaining the ability to capture complex patterns through their probabilistic framework. This aligns with one of our goals of making generative modeling methods more approachable to the broader scientific community. To avoid confusion, we have rephrased this part and added more clarifications in Section INTRODUCTION (see Page 3).\nSecondly,weprovideclarificationonthetrainingandinferenceprocessesofthegeneralVAEframework.\nTraining Stage (Reconstruction): This corresponds to the “second method” mentioned in this re- viewer’s comment. During training, the VAE is optimized to reconstruct input waveforms from their encoded latent representations. Specifically, a given input x (e.g., a seismic waveform sample) is encoded into a latent variable z, and the decoder reconstructs x from z. This process enables the latent space to learn a meaningful representation of the underlying seismic data. Additionally, we define a prior distribution specifically for seismic data. While a standard Gaussian prior is commonly used, we employ a sequential Gaussian prior, which better aligns with the temporal dynamics inherent in seismic waveforms. The prior distribution acts as a regularization term through the KL divergence, enforcing the latent space distribution to remain close to the prior. The reconstruction loss and the KL divergence regularization guide the model in learning a smooth and structured latent space and facilitate effective sampling during the generative process. Meanwhile, by incorporating conditional variables (e.g., geospatial coordinates, earthquake magnitudes, and depths) into the model prior and the decoder, the entire framework ensures the generation of consistent and realistic ground motion waveforms for arbitrary earthquake scenarios.\nInference Stage (Generation): Unlike the training stage (reconstruction), the inference process (gener- ation)doesnotrelyonexistinginputwaveformsxtosynthesizegroundmotiondata. Duringinference, ground motion waveforms are generated by directly sampling from the prior distribution and passing it through the decoder to produce new waveforms. The conditional variables are integrated simulta- neously to ensure that the generated waveforms are physically meaningful. This part aligns with the “first method” mentioned in this reviewer’s comment.\nThe training and inference processes of VAEs are illustrated in Figure 1(b-c) and introduced in the Dynamic VAEs and Network design sections. We have also included more clarifications in the revised manuscript. Please see Page 16.\nReferences:\n[1] Arjovsky, Martin, Soumith Chintala, and L´eon Bottou. “Wasserstein generative adversarial net- works.” International conference on machine learning. PMLR, 2017.\nComment 4: Following up on 3.1, the paper dismisses the use of GAN citing that “GAN models cannot capture 3D wave propagation effects due to limited input conditional variables” is not a limitationofGANs. Rather, theabilityofGANstomodelcomplexphenomenalike3Dwaveprop- agation largely depends on the design of the input space and the conditioning of the model. The inclusion of more or appropriate conditional variables could potentially mitigate this limitation, suggesting that the model’s design, rather than the model type, influences its capability. This is one of many erroneous statement examples\nReply: We would like to reiterate that it was not our intention to dismiss the use of GANs for ground motion generation (as noted in our response to Comment 3). In the original manuscript, we aimed to mention that previous applications of GANs in this context have often employed a different set of con- ditional variables, which may limit their ability to model 3D wave propagation effects comprehensively. We agree with the reviewer that the ability of GANs, or any generative model, to capture complex phenomena could be improved by including appropriate and sufficient conditional variables.\nTo address any potential misunderstanding, we have revised the relevant sections in the manuscript to clarify our position. We now emphasize that the limitations mentioned are tied to the specific design choices in prior work, rather than the GAN framework itself. Please see Page 3.\nComment 5: VAEs are known to generate samples of lower quality compared to GAN and Diffusion models. I cannot see the advantage of using the conditional VAE model compared to other-state of-the-art generative models. While GAN and Diffusion models are shown to have universal approximation ability under mild conditions, the universal approximation ability of VAE is constrained by the reliance on the Gaussian assumption in the latent space. The training objectiveofVAEreliesonKLdivergence, whichencouragesthelatentspacetomatchapredefined distribution (typically Gaussian). This regularization can sometimes be too strong, leading to posterior collapse and not covering all modes of the data distribution. To circumvent issues with training stability of GANs, Diffusion models such as DDPM have significantly improved their training and inference speeds and can be viewed as an extension of VAEs that also mitigates earlierconcernsabouttrainingcosts. Lastly, theconditionalVAEmodelreliesonaphaseretrieval methodthathasbeenshowntoperformworsethatlearningwaveformsdirectlyinthetimedomain, similarly to the approach explored by GAN models.\nReply: Thank you for your comment. We would like to provide clarifications from three perspectives: the choice of dynamic VAEs, the research focus of our paper, and the use of phase retrieval.\nThe choice of dynamic VAEs: We agree that vanilla VAE models can sometimes produce lower-quality samples compared to GANs or diffusion models. However, in our paper, we employ a dynamic VAE framework,whichhasbeendemonstratedtoperformwellfortimeseriesgeneration(seethisrecentwork [1]). In this context, the dynamic VAE framework has shown competitive results compared to models such as TimeGAN [2]. Moreover, VAEs are lightweight, computationally efficient, and straightforward to implement, making them especially accessible for a wider range of domain scientists who may not have extensive experience with generative models. Due to these benefits, VAEs remain an active area of research in scientific machine learning [3,4,5].\nEach generative model has distinct strengths and limitations. For example, GANs are often challenged by training instability, while Diffusion models generally require significant computational resources. To provide future clarity, we have added a detailed comparison of commonly used generative models in Section A.2 Comparisons of generative models. Please see Page 31. For your convenience, we also present the summary in Table II.\nTABLE II. Summary of GenAI models in terms of generation fidelity, sampling efficiency, training stability, and compression (✓: Good, ✗: Suboptimal).\nModel High fidelity Fast sampling Training stability Compression\nVAE ✗ ✓ ✓ ✓\nGAN ✓ ✓ ✗ ✗\nDiffusion model ✓ ✗ ✓ ✗\nThe research focus of this paper: We would like to emphasize that the primary focus of this work is not on the generative model itself but more on the ability to learn spatially continuous ground motions. VAEs provide a suitable and accessible framework for this purpose. While GANs and Diffusion models offer their own strengths, our choice of a dynamic VAE was motivated by its stability, ease of use, and effectiveness for time series learning [1]. Additionally, we are actively exploring Diffusion models in follow-up work (see Bi. et al. [6]) to address some of the limitations mentioned.\nThe use of phase retrieval: We recognize that phase retrieval methods can introduce numerical errors. However,directlylearningwaveformsinthetimedomainalsoposessignificantchallenges,asalsonoted by Reviewer 1. There is always a tradeoff between computational difficulty and model performance. Additionally, we have included the discussion of phase retrieval methods in Section DISCUSSION and pointed out that we will explore learning in the time domain as a future work. It is worth noting that many recent studies [5,7] have successfully used phase retrieval methods in generative modeling, validating its practicality in the context of ground motion generation.\nReferences:\n[1] Naiman, Ilan, et al. “Generative modeling of regular and irregular time series data via Koopman VAEs.” in The Twelfth International Conference on Learning Representations, 2024.\n[2] Yoon, Jinsung, Daniel Jarrett, and Mihaela Van der Schaar. “Time-series generative adversarial networks.” Advances in neural information processing systems 32 (2019).\n[3] Eckmann, Peter, et al. “Limo: Latent inceptionism for targeted molecule generation.” Proceedings of machine learning research 162 (2022): 5777.\n[4] Solera-Rico, Alberto, et al. “β-Variational autoencoders and transformers for reduced-order mod- elling of fluid flows.” Nature Communications 15.1 (2024): 1361.\n[5] Ning, Chunxiao, and Yazhou Xie. “Convolutional variational autoencoder for ground motion clas- sification and generation toward efficient seismic fragility assessment.” Computer-Aided Civil and Infrastructure Engineering 39.2 (2024): 165-185.\n[6] Bi, Zhengfa, et al. “Advancing data-driven broadband seismic wavefield simulation with multi- conditional diffusion model”, arXiv.\n[7]Esfahani,RezaDD,etal. “TFCGAN:NonstationaryGround-MotionSimulationintheTime–Frequency Domain Using Conditional Generative Adversarial Network (CGAN) and Phase Retrieval Methods.” Bulletin of the Seismological Society of America 113.1 (2023): 453-467.\nComment 6: Thespike-likeissuesobservedinFigures6and7appeartointroducehigh-frequency artifacts, which are not observed in previous machine-learning-based ground motion synthesis models. The authors should conduct an ablation study to determine whether these artifacts are attributable to the chosen model itself or are a result of the conditional variables used. This study could help clarify the source of the artifacts and guide improvements in the model design or variable selection. Knowing that results deviate from the median of GMMs and NGMMs in the high frequency regime, as evidenced also by the statistics in Figure 2 d and e, makes the generation of FAS maps at 10Hz irrelevant; the maps that would be useful to see are the < 1Hz maps and their uncertainty, that show the capability of the model to capture path effects, if it does, as the authors claim.\nReply: We have examined the waveforms generated by both the CGM-baseline and CGM-GM models andfoundthatthespike-likeartifactsappearinbothcases. Oneoftherepresentativeexamplesisshown in Fig. 3. Our analysis suggests that these artifacts are primarily due to the phase retrieval methods usedduringwaveformreconstruction. Sincephaseretrievalestimatesphaseinformationindirectlyfrom amplitude information, it can introduce inaccuracies that may distort waveform shapes, particularly in the high-frequency range. However, we emphasize that these artifacts do not significantly impact the model’s ability to generate realistic Fourier Amplitude Spectra (FAS) maps, as the overall spectral characteristics remain well-preserved. Hence, comparing the FAS at 10 Hz, where the CGM-GM demonstrates the best fit to the observations, is reasonable.\nFIG. 3. A representative example of comparing waveforms from CGM-GM and CGM-baseline.\nWe agree that evaluating lower-frequency FAS maps and their associated uncertainties would provide additional insights for path effects. However, the frequency range used in this study is limited to [2,15] Hz, and hence we cannot compute the FAS map at 1 Hz without building a model using wider frequency data. We selected [2,15] Hz to satisfy the S/N ratio over an entire frequency range at all used stations. Small-magnitude earthquakes typically exhibit very low S/N ratios at lower frequencies, whichwouldsignificantlyreducethesizeandreliabilityofthetrainingdataset. Nevertheless, extending the CGM-GMframeworkto incorporate lower-frequency data remains an important directionforfuture work. This could be achieved by fine-tuning the model using stations with high-quality low-frequency waveforms or by developing a framework and training strategies that effectively handle varying data distributions across different frequency ranges. For the discussion of uncertainty, please refer to our response to Comment 2.\nAdditionally, we present the FAS map at 5 Hz in Section B.6 FAS Maps (Pages 34 and 38) and provide the FAS map at 2 Hz below in Fig. 4. The selection of 10 Hz in the main text was based on its effectiveness across our CGM-GM model and all baseline models, ensuring a consistent and meaningful comparison. We have included additional FAS maps at 2 Hz in Figure 14 of the revised manuscript (Page 38).\na CGM-GM b Ergodic GMM c Non-ergodic GMM\nFIG. 4. Comparative results of FAS maps at a frequency of 2 Hz between our generations and the empirical GMMs. This seismic event is defined with a magnitude of 3.84 and a depth of 7.94 km. The epicenter (red star) is located at a geographic position with a latitude of 37◦51.6′N and a longitude of 122◦15.6′W. a, b, and c show the FAS maps of our generative model, ergodic GMM, and non-ergodic GMM..\nComment 7: The selection of baseline models for validation in this study appears to be inade- quate. Fine-tuned GANO models such as the recently published version of [32] should be included in the comparison.\nReply: Thank you for this suggestion. We have carefully considered the fine-tuned GANO model and conducted experiments to compare its performance with our proposed CGM-GM framework. We generallyfollowedthedefaultnetworksetupintheGitHubrepository(https://github.com/yzshi5/GM- GANO) and modified the number of conditional variables from 4 (original setup in GANO) to 6 (ours). We trained the network for 500 epochs with a learning rate of 0.0001. Several other important parameters are listed here: n critic is 10, λ grad is 10, and the batch size is set as 64.\nHowever, in our case, where the dataset consists of small-magnitude earthquakes with a relatively low signal-to-noise (S/N) ratio, we found that the GANO model did not perform well. Specifically, the generated waveforms exhibited noticeable discrepancies from the ground truth, as shown in both visual inspections of sample waveforms and quantitative statistical analyses of peak ground velocities and arrival times. The detailed results are shown in Fig. 5. These findings suggest that the model may struggle to generalize effectively under conditions of low S/N ratio. While we acknowledge the advancements of fine-tuned GANO models, our results indicate that their effectiveness may depend on the characteristics of the dataset, particularly for smaller-magnitude events. We appreciate the reviewer’s suggestion and welcome further discussion on this matter.\na\nb c\nFIG. 5. Model performance of GANO on our small-magnitude earthquake datasets. a shows the two randomly generated samples (blue) compared to the corresponding ground truth waveforms (red). b and c present the comparisons between ground truth and the generations of peak ground velocities and arrival times, respectively.\nComment 8: Themodelthatistrainedonsmallmagnitudedata, shouldalsoshowthecapability for reasonable scaling, with increased uncertainty, for larger magnitude events. If it doesn’t, then there is no value in having a model like this for California.\nReply: To address your concern, we show the magnitude scaling performance of our CGM-GM model compared to empirical GMMs. Since our model is only trained on small-magnitude earthquakes (M < 4), directly generalizing to large-magnitude scenarios will be challenging. Therefore, to mitigate the challenge, we fine-tune our pre-trained model on two larger-magnitude earthquake data: (i) the 2014 Napa earthquake with a magnitude of 6.02 and depth of 11.0 km, and (ii) the 2018 Berkeley earthquake with a magnitude of 4.38 and depth of 12.3 km.\nWe present the scaling performance up to a magnitude of 8, occurring at a hypocenter depth of 7.94 km, with geographic coordinates of latitude 37◦51.6′N and longitude 122◦15.6′W (the same location as shown in Figure 2 of the manuscript). The magnitude is varied from 1 to 8, and waveforms are generated at 73 real-world stations. In Fig. 6, we present a comparison of FAS for small-magnitude earthquakes and spectral acceleration (Sa) values for large-magnitude earthquake scenarios. For small- magnitudeevents,weshowthemeanandonestandarddeviationofourCGM-GMpredictionsalongside the non-ergodic GMM. The predictions from our model generally align with those of the non-ergodic GMM, with an increased standard deviation around M1 due to the limited number of observational samples available for training, as illustrated in Figure 6(b) of the revised manuscript (Page 18). For largemagnitudescaling,weevaluatetheSaresultsofourmodelcomparedtotheergodicASK14model of the Next Generation Attenuation-West2 project [1] as there are no non-ergodic models developed specifically for the SFBA.\nThe ASK14 model requires site and earthquake scenario parameters that may not always be readily available. To compute distances between the fault plane and the stations, we assume finite faults are centered at a latitude of 37◦51.6′N and a longitude of 122◦15.6′W, which corresponds to the hypocenter location used in small earthquake comparisons. We determine fault widths and lengths based on [3]. For site parameters, we utilize V values ranging from 400 to 1000 m/s and Z values s30 1.0 between 0.1 and 2 km [2]. We then calculate Sa predictions at 10 Hz for these site conditions and take the mean values. Stations within rupture distance R of 10 to 50 km are selected for this analysis. rup\nIn Figure 6, we present comparisons for rupture distances ranging from 10 to 50 km. We plot the mean and standard deviation of the generated samples. For the GMM predictions, we calculate the logarithmic mean (µ ) of ground motion estimates, along with the mean of the standard deviation mean predictions.\nThe CGM-GM predictions exhibit an increase in Sa with increasing magnitude but remain lower than those predicted by the ASK14 model for magnitudes between 5 and 7. This discrepancy may arise because our model is trained mostly on small-magnitude events, where the earthquake source rupture processes can be approximated as point sources. In contrast, large-magnitude events require consider- ation of finite rupture planes, which introduces deviations in the predictions. As shown in Fig. 7, the fine-tuning strategy improves the model performance by enhancing Sa values, which further supports this observation. Furthermore, the range of ground motion predicted by the CGM-GM model is slightly narrower than that of the ASK14 model, as the generated samples represent a single realization of an earthquakescenario. Theseresultssuggestthatacomprehensiveinvestigationintothevariabilityofthe generated ground motions is an essential direction for future research. To enhance model performance, we plan to explore fine-tuning using larger datasets and incorporate out-of-distribution generalization methods. We have included additional scaling results in the revised manuscript. Please see Section Generalization to larger-magnitude scenarios on Pages 10-12.\nAdditionally, we show the comparison between the pre-trained and fine-tuned versions of our CGM-GM new model. We found that the fine-tuning strategy can effectively reduce the variations of Sa values under large-magnitude scenarios. Please see Section B.9 Details of fine-tuning and GMM comparisons on Pages 40-41.\na b\nFIG. 6. Comparison of the magnitude scaling performance of our CGM-GM model against empirical GMMs over rupture distances between 10 and 50 km. a illustrates the magnitude scaling of FAS (10 Hz) for small-magnitude earthquakes (M1–4), benchmarked against the non-ergodic GMM. b depicts the magnitude scaling of Sa (10 Hz) for large-magnitude earthquakes (M4–8), compared with the ASK14 model.\nSecondly, we respectfully disagree with the reviewer’s comment on the value of small-magnitude earth- quakes in California. Accurately modeling small-magnitude earthquakes remains essential for seismic research, hazard assessment, and engineering applications in California and beyond. Small-magnitude earthquakes play a crucial role in evaluating the seismic response of infrastructures, particularly for serviceability-level assessments, performance-based earthquake engineering, and understanding long- term seismic risk in urban areas. Additionally, the study of small-magnitude earthquakes allows us to\nFIG. 7. Comparison of the magnitude scaling performance between our pre-trained CGM-GM model and the fine-tuned CGM-GM model over rupture distances ranging from 10 to 50 km.\nestimate the region-specific variations of the ground motion, which can lead to reducing uncertainties in the seismic hazard in the region [4,5].\nWe believe that developing a model for small-magnitude earthquakes in California is a meaningful and importantendeavor, andwehaveaddedmoreclarificationsonthisaspectinSectionINTRODUCTION on Page 3. The significance of simulating small-magnitude earthquakes outside California is also recognized by Reviewer 1. We cite Reviewer 1’s related comment as follows:\n“In particular, the proposed method has a significant application in the improvement of shaking maps in early warning systems and seismic hazard analysis for regions with low magnitude”\nReferences:\n[1] Abrahamson, Norman A., Walter J. Silva, and Ronnie Kamai. “Summary of the ASK14 ground motion relation for active crustal regions.” Earthquake Spectra 30.3 (2014): 1025-1055.\n[2] Hirakawa, Evan, and Brad Aagaard. “Evaluation and updates for the USGS San Francisco bay region 3d seismic velocity model in the east and North Bay portions.” Bulletin of the Seismological Society of America 112.4 (2022): 2070-2096.\n[3]Leonard, Mark. “Earthquakefaultscaling: Self-consistentrelatingofrupturelength, width, average displacement, and moment release.” Bulletin of the Seismological Society of America 100.5A (2010): 1971-1988.\n[4] Lavrentiadis, Grigorios, Norman A. Abrahamson, and Nicolas M. Kuehn. “A non-ergodic effective amplitude ground-motion model for California.” Bulletin of Earthquake Engineering 21.11 (2023): 5233-5264.\n[5] Chiou, Brian, et al. “Ground-motion attenuation model for small-to-moderate shallow crustal earthquakes in California and its implications on regionalization of ground-motion prediction models.” Earthquake spectra 26.4 (2010): 907-926.\nIII. RESPONSE TO REVIEWER 3\nGeneral comment: The authors of the manuscript entitled ”Learning Physics for Unveiling Hidden Earthquake Ground Motions via Conditional Generative Modeling” proposed a dynamic VAE, featured by variational sequential architectures (e.g., RNNs) in both prior and posterior distributions, which facilitates the capture of the temporal evolution (dynamics) of time series. Thanks to this Conditional Generative Modeling for Ground Motion (CGM-GM) strategy, the authors synthesize high-frequency and spatially continuous earthquake ground motion scenarios in the San Francisco Bay Area. Interestingly, CGM-GM embeds earthquake magnitude and hypocenter coordinates as well as sensors’ geographic coordinates, so to reproduce the natural groundmotionvariability. CGM-GMwastrainedwithsmall-magnitudedata, whichisremarkable because they challenge small Signal-to-Noise Ratios.\nThemanuscriptisverywellwritten. Thepresentedresultsareconvincing. CGM-GMframeworkis effective and robust for generating earthquake ground motions. Its lightweight architecture highly appreciated. The most interesting aspect of our proposed generative model lies in its claimed ability to approximate ground motions for arbitrary (future) earthquakes and sensor locations.\nReply: Thank you for your positive feedback. We appreciate the reviewer for acknowledging the effectiveness, robustness, and lightweight design of our proposed method. We address your concerns point-by-point as follows.\nComment 1: In the Introduction, the authors overlooked some recent literature on earthquake ground motion generation using neural operators and graph neural networks. These represent a credible alternative to the CGM-GM presented in the manuscript and it would interesting to highlight what are the advantages and limitations of adopting CGM-GM instead of one of these other approaches.\nReply: In the original manuscript, our focus was on discussing the general application of generative modeling for earthquake ground motion generation. We agree with the reviewer that a more detailed introductiontoothermachinelearningapproaches, suchasneuraloperatorsandgraphneuralnetworks (GNNs), would enhance the quality of the manuscript. We have expanded our discussion to include recent works on earthquake ground motion generation using neural operators [1-3] and GNNs [4-6].\nNeural operators have been applied in two main paradigms for producing seismic waveforms: (i) mod- eling wave propagation and (ii) generative AI. While our original manuscript briefly discussed neural operators for wave propagation (see Page 2), we have incorporated additional references [2-3] high- lighting recent advances in neural operators for learning 3D wave phenomena, which show promise in simulating higher-dimensional seismic events. Additionally, neural operators have been employed in a generative framework [1] for generating broadband ground motion waveforms from low-frequency components. This method considers two separate networks to learn amplitude and phase information, avoiding numerical errors from phase retrieval methods. This approach demonstrates significant im- provements over traditional methods but its potential limitations lie in the high computational cost of training two independent network modules (ignoring relationships between amplitudes and phases) and the inherent instability of the GAN architecture. Additionally, GNNs have also been employed for ground motion synthesis, leveraging their ability to capture spatial correlations. These models are capable of predicting ground motion intensities at more distant stations within the seismic network.\nWehaveincludedtheseexpandeddiscussionsandthecorrespondingreferencesintherevisedmanuscript (see Pages 2-3). We appreciate the reviewer’s feedback and welcome further recommendations if there are additional critical works we may have missed.\nReferences:\n[1] Aquib, Tariq Anwar, and P. Martin Mai. “Broadband Ground-Motion Simulations with Machine- Learning-Based High-Frequency Waves from Fourier Neural Operators.” Bulletin of the Seismological Society of America 114.6 (2024): 2846-2868.\n[2] Zou, Caifeng, et al. “Deep neural Helmholtz operators for 3-D elastic wave propagation and inversion.” Geophysical Journal International 239.3 (2024): 1469-1484.\n[3] Lehmann, Fanny, et al. “3D elastic wave propagation with a factorized Fourier neural operator (F-FNO).” Computer Methods in Applied Mechanics and Engineering 420 (2024): 116718.\n[4] Clements, T., et al. “GRAPES: Earthquake early warning by passing seismic vectors through the grapevine.” Geophysical Research Letters 51.9 (2024): e2023GL107389.\n[5] Bloemheuvel, Stefan, et al. “Graph neural networks for multivariate time series regression with application to seismic data.” International Journal of Data Science and Analytics 16.3 (2023): 317- 332.\n[6] Murshed, Rafid Umayer, et al. “Real-time seismic intensity prediction using self-supervised con- trastive GNN for earthquake early warning.” IEEE Transactions on Geoscience and Remote Sensing (2024).\nComment 2: What steps and data would the authors deem necessary to apply CGM-GM to anotherarbitraryregionoftheworld? ThetitleenticesthereadershipbutitseemsthatCGM-GM would be hardly adoptable in regions with poor seismic record and/or varying focal mechanisms. It is in this reviewer’s opinion that the title should reflect the claim the authors make that their primary interest is in generating waveforms along the Hayward Fault, with no domain shift for the earthquake rupturing processes.\nReply: This is a thoughtful observation. The purpose of this paper is not to build a pre-trained model for Hayward fault earthquakes but to demonstrate that the conditional generative modeling framework with an appropriate parameter set can capture ground motion waveforms. We address your concerns regarding applying our method to different regions from three perspectives: requirements for implementingtheCGM-GM,theresultsofusinglimitedseismicrecords,andvaryingfocalmechanisms.\nThe details are outlined below. The corresponding revisions in the manuscript can be found in Section B.8 Limited seismic records (Pages 37-40) of Supplementary Material.\nRequirements for applying the CGM-GM: The main data components include recorded seismic wave- forms and their corresponding conditional variables, i.e., earthquake magnitudes, source depths, and geospatial coordinates of sources and stations. Seismic stations of interest can be selected based on a specific range of magnitudes and hypocentral distance. This model does not necessarily require Vs 30 or focal mechanisms, though they can be included as a part of conditional parameters. This feature makes the application of CGM-GM more practical, especially in scenarios where such information is either unavailable or subject to uncertainty regarding its accuracy. The implementation process can be summarized into two parts: (i) preprocessing the datasets and (ii) applying the CGM-GM framework. The preprocessing procedures are described in Section Data selection (Pages 15-16). We perform the Fourier transforms to both the noise and the signal and smooth the Fourier amplitudes using the Konno-Ohmachi [1] window procedure. Our goal is to retain recordings with a signal-to-noise (S/N) ratio exceeding 3 across the frequency range of [2,15] Hz. To apply the CGM-GM model, we first extract time-frequency amplitude information using the Short-Time Fourier Transform (STFT). The dynamic VAE model is then trained on amplitude data alongside the conditional variables. Once the model is well trained, we generate synthetic ground motions by sampling from the prior distribution and conditioning on specific input variables.\nLimited seismic records: To assess the robustness of our CGM-GM model under limited availability of seismic datasets, we evaluate its performance on additional scenarios with reduced station and source coverage. Specifically, we design two cases of data sparsity by independently reducing the number of stations and sources in the original training dataset. We introduce a random removal ratio λ, selected from [0.01,0.1,0.2,0.3,0.4,0.5,0.9], to systematically reduce data availability. For each removal scenario, we train three independent models using three different random seeds, ensuring variability in the selected waveform samples and facilitating a fair comparison. For instance, applying λ = 0.3tostationremovalresultsinadatasetcontaining3880waveformsfortraining. Thevariationin random seeds ensures different subsets of waveforms are removed in each instance. To quantify model performance, we compute the average absolute difference between the predicted and true Fourier Amplitude Spectra (FAS) values across all data samples. The results of this evaluation under reduced data availability are presented in Fig. 8. Moreover, Fig. 9 presents a comparison of FAS maps (10 Hz) for the non-ergodic GMM, CGM-GM trained with the full dataset, and CGM-GM trained under data-limited scenarios, where source removal ratios are defined as 0.3, 0.5, and 0.9. The seismic event used for producing the FAS maps is consistent with the event in Figure 2 of the manuscript, which is characterized by a magnitude of 3.84 and a depth of 7.94 km. The epicenter, denoted by a red star, is located at a geographic position with a latitude of 37◦51.6′N and a longitude of 122◦15.6′W. Our model demonstrates robust performance, maintaining reasonably good performance and successfully generating realistic FAS maps even with a 30% reduction in data availability. These results validate the effectiveness of our approach in handling limited seismic records.\n1.6\n1.3\n1.0\n0.7\n0.4\n10 2 10 1 100\neulav\nSAF\nlaudiseR\nReducing Sources Reducing Stations\n1.6\n1.3\n1.0\n0.7\n0.4\n10 2 10 1 100\nFIG. 8. Model performance of our CGM-GM on poor seismic records. We reduce the number of earthquake sources and stations in the original dataset. The results are based on three random runs. The purple curves represent the average absolute difference between predicted and true FAS values across all waveform samples. The shading regions denote the coverages of the mean ± one std.\nVarying focal mechanisms:IntheCGM-GMframework,weusesourcelatitudes,longitudes,anddepths as conditional variables. This approach allows us to handle various focal mechanisms while implicitly assuming that focal mechanisms remain consistent (although magnitudes may vary) for earthquakes occurring at the exact same geospatial location. In our dataset, the predominant focal mechanisms are characterized by strike-slip faulting, although detailed parameters can vary [2]. Furthermore, in one of our follow-up studies, we utilized seismic records from the Geysers geothermal field [3,4]. Earthquakes inthisregionexhibitarangeofmechanisms,includingnormal,reverse,andstrike-slipfaulting[5]. With\nNon-ergodic GMM CGM-GM (full) CGM-GM (0.3) CGM-GM (0.5) CGM-GM (0.9)\nFIG. 9. A comparison of FAS maps at 10 Hz for the non-ergodic GMM, CGM-GM trained with the full dataset, and CGM-GM trained under data-limited scenarios, where source removal ratios are defined as 0.3, 0.5, and 0.9.\nappropriate parameter settings, the CGM-GM model demonstrated the ability to effectively capture this diversity. We have included additional descriptions in Section DISCUSSION (Page 12).\nReferences:\n[1]KatsuakiKonnoandTatsuoOhmachi. “Ground-motioncharacteristicsestimatedfromspectralratio between horizontal and vertical components of microtremor.” Bulletin of the Seismological Society of America, 88(1):228–241, 02 1998.\n[2] Hirakawa, Evan, and Brad Aagaard. ”Evaluation and updates for the USGS San Francisco bay region 3d seismic velocity model in the east and North Bay portions.” Bulletin of the Seismological Society of America 112.4 (2022): 2070-2096.\n[3] Nakata, Rie, et al. “Simulating seismic wavefields using generative artificial intelligence.” The Leading Edge 44.2 (2025): 123-132.\n[4] Bi, Zhengfa, et al. “Advancing data-driven broadband seismic wavefield simulation with multi-conditional diffusion model.” arXiv preprint arXiv:2501.14348 (2025).\n[5] Lin, Guoqing, and Bateer Wu. ”Seismic velocity structure and characteristics of induced seismicity at the Geysers geothermal field, eastern California.” Geothermics 71 (2018): 225-233.\nComment 3: Did the authors consider using ground motion numerical simulation to expand their training database? If so, how would they mix synthetics with records within the learning framework?\nReply: We have indeed considered the possibility of incorporating numerical simulation data into our training dataset. However, there are certain limitations associated with using these datasets, as they are typically constrained to low-frequency components (< 1 Hz). Furthermore, the discrepancies between simulation and observation need to be carefully examined. Therefore, we chose not to include simulation datasets in this work.\nComment 4: As per this reviewer’s understanding, each 3-component single-station time-history is encoded separately, along with its geographical coordinate. How is spatial correlation encoded into the latent space, if so? Peak values contour plots seem rather realistic, but what about the whole time series? This question arises from the fact that the phase spectra and the ground motion coherence are of paramount importance in earthquake engineering and seismic hazard modelling.\nReply: Thank you for your insightful question. In our current framework, we do not explicitly incorporatespatialcorrelationinthelatentspace. Itcouldbeachieved, e.g., byincorporatingaMat´ern covariance function, as stated in Section DISCUSSION (see Page 12). However, spatial correlation is implicitly learned through the inclusion of geographical coordinates as conditional variables. We acknowledge that ensuring phase continuity is inherently more challenging than maintaining spatial continuity of amplitudes, due to the rapid variability of phase characteristics. Addressing this challenge is an important direction for future research, and we are actively exploring it in a separate study. We have added further clarifications in the revised manuscript. Please see Page 13.\nComment 5: The readership would benefit from a Figure showing extensive goodness-of-fit metrics comparisons. This reviewer suggests to adopt the modified Anderson’s criteria, found in Olsen et al., 2010. Other GoF metrics would probably mislead the readership, because CGM- GM is conceived for stochastic generator of ground motions, not for historical ground motion reconstruction.\nReply: Thank you for this suggestion. We have reviewed the literature and identified both An- derson’s criteria [1] and Olsen’s method (modified Anderson’s criteria) [2] as established approaches for computing goodness-of-fit (GoF) metrics for broadband waveforms. Additionally, we note that a recent generative modeling study [3] has applied Olsen’s method in this context. We acknowledge that presenting GoF metrics provides a more comprehensive evaluation of the generated waveforms and enhances the readers’ understanding of model performance. However, we were unable to find open-source implementations for directly computing GoF metrics using Olsen’s criteria. Given this, we proceeded with the original Anderson’s criteria, for which we had existing implementation resources. Anderson’s criteria have ten distinct characteristics, each evaluated on a scale from 0 to 10. A score of 10 indicates perfect agreement. The scores for individual parameters are averaged to determine the overall GoF. A score below 4 denotes a poor fit, 4–6 represents a fair fit, 6–8 indicates a good fit, and scores exceeding 8 are an excellent fit. The GoF analysis has been performed across different rupture distances, earthquake depths, and magnitudes, with results presented in Fig. 10. Our findings indicate that the generated waveforms exhibit a good agreement with the ground truth, with mean GoF values generally around 6. We have included further clarifications in Section B.7 Evaluations of goodness-of-fit of the revised manuscript (Pages 36-37). Additionally, if the reviewer could provide a reference to an open-source implementation of the modified Anderson’s method, we would be happy to further evaluate it and incorporate it into our analysis.\nFIG. 10. An overview of GoF values between our generations and ground truth observations across various magnitudes, distances, and earthquake depths. The dashed curves and shading regions denote the mean values and the coverages of the mean ± one std, respectively.\nReferences:\n[1] Anderson, John G. “Quantitative measure of the goodness-of-fit of synthetic seismograms.” Pro- ceedings of the 13th world conference on earthquake engineering. Vol. 243. Earthquake Engineering Research Institute, 2004.\n[2]Olsen,KimB.,andJohnE.Mayhew. ”Goodness-of-fitcriteriaforbroadbandsyntheticseismograms, with application to the 2008 Mw 5.4 Chino Hills, California, earthquake.” Seismological Research Letters 81.5 (2010): 715-723.\n[3] Aquib, Tariq Anwar, and P. Martin Mai. “Broadband Ground-Motion Simulations with Machine- Learning-Based High-Frequency Waves from Fourier Neural Operators.” Bulletin of the Seismological Society of America 114.6 (2024): 2846-2868.\nComment 6: Minor comments and questions in the attached pdf.\nOn Page 7, “Hence, the contours of this FAS map are circular, showing geometrical spreading and average attenuation effects.” Very interesting. Some similar results were found by Paolucci et al., 2018. https://doi.org/10.1785/0120170293. They found that the contours were following the fault geometry. Did the authors observe something similar when addressing finite-fault scenarios?\nReply: Thank you for referencing the paper by Paolucci et al. (2018) [1]. In their paper, Figures 5(a) and 5(b) illustrate heterogeneous patterns capturing physical features and “circular patterns”, respectively. Their findings, including the alignment between maps of peak ground motion values and fault geometry, are highly relevant to our work.\nIn our manuscript, Figure 2(b) presents a radial decay pattern and average attenuation effects due to an implicit assumption of a homogeneous Earth subsurface across the spatial domain, when using earthquake magnitudes, source depths, and rupture distances as conditional variables.\nAdditionally, our CGM-GM predictions vary across the San Andreas and Hayward fault locations, as illustrated in Figure 2(c) of our manuscript. The sharp lithological changes lead to substantial differences in ground motions on either side of the faults. This alignment with fault geometry has also been observed in physics-based earthquake simulations of the San Francisco Bay Area [2]. This indeed align with the findings of [1] and also with the physics-based simulation results of for example [2]. Therefore, our CGM-GM model exhibits the capability to capture spatial heterogeneity in ground motions. By embedding geospatial coordinates of sources and stations as conditional variables, our model can effectively capture variations in combined site, source, and path effects, providing insights into more complex spatial features. We have expanded the discussion in the revised manuscript to include this relevant study [1,2]. Please see Section Waveforms and spatial continuity on Page 8.\nReferences:\n[1]Paolucci, Roberto, etal. “Broadbandgroundmotionsfrom3Dphysics-basednumericalsimulations using artificial neural networks.” Bulletin of the Seismological Society of America 108.3A (2018): 1272-1286.\n[2] Rodgers, Arthur J., et al. “Broadband (0–5 Hz) fully deterministic 3D ground-motion simulations of a magnitude 7.0 Hayward fault earthquake: Comparison with empirical ground-motion models and 3D path and site effects from source normalized intensities.” Seismological Research Letters 90.3 (2019): 1268-1284.\nOn Page 8, “This is due to the inductive bias of VAE models, where the learned amplitude information is over-smoothing in the logarithmic space.” Have the authors tested beta-VAE in this sense?\nReply: This is a good comment. We have tried beta-VAE [1] before, and the performance was comparable to the results presented in our manuscript. While beta-VAE can help mitigate the over- smoothing issue to some extent, it is primarily designed to improve the disentanglement performance in the latent space rather than specifically address over-smoothing. The over-smoothing issue in VAEs is more directly related to the fixed variance parameter in the decoder’s output distribution [2]. This fixed variance can lead to oversmoothness and posterior collapse, where the latent space becomes uninformative. Addressingthischallengerequiresstrategiesthatexplicitlytargetthedecoder’svariance parameterization. We have acknowledged this limitation in the manuscript and discussed potential methodsforfurtherimprovingmodelperformanceinsectionDISCUSSION (Page13). Forconvenience, the relevant paragraph is attached below:\n“Another direction is exploring specific techniques, such as generalized variance parameterization [48] and heavy-tailed distributions [49], to mitigate the over-smoothing issue in VAE models.”\nReferences:\n[1] Higgins, Irina, et al. “beta-vae: Learning basic visual concepts with a constrained variational framework.” ICLR, 2017.\n[2]Takida,Yuhta,etal. “PreventingoversmoothinginVAEviageneralizedvarianceparameterization.” Neurocomputing 509 (2022): 137-156.\nOn Page 11, “A potential limitation lies in the lack of the constraints of path effects.” Highly appreciated comment from the authors! “The implicit learning of paths can be achievable since the training waveform data naturally contains the path effect.” I suggest the work of Gatti and Clouteau 2020, https://doi.org/10.1016/j.cma.2020.113421.\nReply: Thank you for suggesting this interesting paper by Gatti and Clouteau (2020) [1]. Their approach is highly related to our proposed future research direction, which combines physics-based numerical simulations with experimental databases and facilitates learning the path effects. Moreover, given the practical constraints of limited data availability, we will also explore methods for explicitly incorporating path effects in the generative framework in future studies to ensure accurate modeling. We have included additional discussions on this topic in the revised manuscript and have cited the referenced paper accordingly. Please see section DISCUSSION on Page 13.\nReferences:\n[1] Gatti, Filippo, and Didier Clouteau. “Towards blending physics-based numerical simulations and seismic databases using generative adversarial network.” Computer Methods in Applied Mechanics and Engineering 372 (2020): 113421.\nComment 7: Data report is not available, nor data themselves. Testing failed because of the lack of test dataset. I suggest to add some demonstrative notebook/tutorial.\nReply: Thank you for your suggestion. We have uploaded the data report to the DesignSafe [1] (see link: https://www.designsafe-ci.org/data/browser/public/designsafe.storage.published/PRJ- 4573). The data report is also updated in the revised manuscript (see Page 20). Additionally, we have included tutorial examples, including guidance on training models and performing evalua- tions from various perspectives. Please see the updated materials in our GitHub repository: https: //github.com/paulpuren/cgm-gm.\nReferences:\n[1] Lacour, M., R. Nakata, N. Nakata, N. Abrahamson (2024). Earthquake Dataset for the Study of SmallMagnitudeEarthquakesintheSanFranciscoRegion. DesignSafe-CI.https://doi.org/10.17603/ds2- necm-5q32."}, {"round": 2, "reviewer_comments": {"reviewer_1": "I am satisfied with the revisions to the manuscript. I recommend this article for publication.", "reviewer_2": "I would like to start by thanking the authors for putting effort to respond, to the extent that they have, to my questions, and improving the manuscript.\nI still think that the work is valuable but its value is significantly overstated, to the point of even erroneously representing the contribution. I will summarize here, again, the main issues I have with this work. Some of the improvements that the authors have made have offered clarity while other improvements have reinforced false claims of this work. I insist that the authors should clarify the true contributions and remove all the comparisons and claims that point towards any other directions.\nSpecifically:\n1. The model that the authors are presenting here is one dimensional. This means that each event-station pair is trained independently of any other station, despite the fact that the coordinates of the source and station are introduced as conditional variables. There is no global variable for the VAE to determine that two adjacent stations with the same source coordinates, from the same event, come from the same earthquake. In this sense, the model IS NOT LEARNING wave propagation.\n2. Since the model is not learning wave propagation, the title LEARNING PHYSICS is false. Its learning a mapping from one point to another, independent to each other, but no coherency; and since the medium from the source to the surface is fixed, the high frequencies (here 10 Hz) are learning some shallow crustal amplification functions which one could similarly capture by replacing the source-station coordinates with Vs30. Instead, the authors are trying to learn a 3D problem with a 1D model which, not only is impossible, but also leads to erroneously spiked strong motion records like the ones they show in Figure 7.\n3. Another way to see what I am saying is to observe the incoherent pattern that emerges in the longer period range, at 2Hz for example, where 'path' and 'site' effects become muddled. There, Figure 14 shows that adjacent stations respond independently of each other for the same 'event'. In other words, there is no spatial continuity. The authors are just querying the model at adjacent stations, whose waveforms, however, even for the same location of source and magnitude of earthquake, will not necessarily be spatially correlated like they are originating from the same event.\n4. In that sense, the best comparison the authors can make is against the ergodic GMMs in Nor Cal. These treat each GM time series separately with some conditional variables. The non-ergodic GMMs implement within event variability terms by studying the residuals relative to the ergodic terms of all the event terms at the same time across the strong motion network using GPs. Why dont they show in Figure 2 comparison against the ergodic GMMs? This is a much more fair comparison since the non-ergodic GMM model is a higher dimensional model than the cGM-GM. And definitely statements like ‘outperforming a state-of the art non-ergodic GMM’ should be deleted since they are misinforming the audience that the model is capturing effects that it is not.\n5. Despite the fact that cGM-GM uses event and station coordinates, because the stations 'dont know' that their motions come from the same event, it cannot be used to plot scenaria. In that sense, all the scenaria plots like Figure 2b and 2c should be replaced with some median amplitude of many realizations, to be compared with median amplitudes of ergodic GMM plots.\n6. Page 7, line 3 of the text: [… implicit assumption of a homogeneous Earth…] should be replaced with the [… implicit assumption of a radially homogeneous Earth…]. We should also note here that this is a limitation of the baseline model that the authors selected; a more fair comparison would be to select a baseline model conditional on magnitude, distance, and a site response proxy (the same they used in their empirical GMM maps) in which case they would recover most of the capabilities that they are able to replicate with their cGM-GM model, possibly with a better time series synthesis capability.\n7. Top of page 8: again erroneous statement: the agreement of the non-ergodic median here is purely phenomenological and if they plot the ergodic median the agreement will likely be better. Also, they should not be plotting scenaria, but median of scenaria, since their model is also ergodic (in the sense that it is learning independent pairs of source-stations and not within event variability like non-ergodic empirical models).\n8. Page 10: bottom of first paragraph: Figure 4e shows bias in the arrival times, not a good correlation as the authors claim. The heat map is clearly off center.\n9. Figure 5b and description of Figure on page 10: The SA trends at 0.1sec are completely off: the empirical model saturates at M8 due possibly to nonlinear site response and the cGM-GM increases almost exponential (due to the fictitious spikes we described above). The fact that the lines cross at some point doesn’t mean there is agreement.\n10. The 2Hz and 5Hz median predictions vs. the ergodic GMM maps should be plotted along with the 10Hz. Engineers don’t care about frequencies as high as 10Hz nearly as much as they care about low frequencies.", "reviewer_3": "The authors satisfactorily addressed almost all of the comments raised by this reviewer in the first round, except Comment 2. The sense of the question was if CGM-GM could be adopted zero-shot. As per the explanation given, the suspect is that the CGM-GM is not fit for zero-shot. A comment on this would be highly appreciated.\nThe transparency of showing the goodness of fit in Fig.10 is highly appreciated.\nFilippo Gatti"}, "author_response": "Learning Physics for Unveiling Hidden Earthquake Ground\nMotions via Conditional Generative Modeling\nNCOMMS-24-45416A\nPu Ren, Rie Nakata, Maxime Lacour, Ilan Naiman, Nori Nakata, Jialin Song, Zhengfa Bi, Osman\nAsif Malik, Dmitriy Morozov, Omri Azencot, N. Benjamin Erichson, and Michael W. Mahoney\nDear Editor and Reviewers,\nWe sincerely appreciate the reviewers for their constructive comments and valuable suggestions. We are pleased that all reviewers acknowledge the value of our research. Our point-by-point responses are provided below marked in blue, and the corresponding manuscript revisions are highlighted in red. Notably, we have made several clarifications as follows:\n• Our model is not designed to directly reproduce 3D wave propagation physics; instead, it learns\na stochastic mapping that captures wave propagation characteristics between earthquake sources\nand stations. Please see the responses to Comments 1 and 2.\n• Our model captures spatial continuity, and we have provided additional supporting evidence.\nPlease refer to our response to Comment 3 for further details.\n• We clarify that specifying only magnitudes and source locations does not uniquely define a\nsingle earthquake event. We have updated all relevant FAS maps with median values from 30\nrealizations. Please see Figures (2, 13-16, 23) of the revised manuscript and the responses to\nComments 1 and 5.\nWe hope our manuscript will be found worthy of publication in this revised form.\nI. RESPONSE TO REVIEWER 1\nGeneral comment: Iamsatisfiedwiththerevisionstothemanuscript. Irecommendthisarticle\nfor publication.\nReply: Thank you very much for your positive feedback and recommendation for publishing.\nII. RESPONSE TO REVIEWER 2\nGeneral comment: I would like to start by thanking the authors for putting effort to respond,\nto the extent that they have, to my questions, and improving the manuscript.\nI still think that the work is valuable but its value is significantly overstated, to the point of even\nerroneously representing the contribution. I will summarize here, again, the main issues I have\nwith this work. Some of the improvements that the authors have made have offered clarity while\nother improvements have reinforced false claims of this work. I insist that the authors should\nclarify the true contributions and remove all the comparisons and claims that point towards any\nother directions. Specifically:\ncomments help us gain a deeper understanding of the strengths and remaining challenges of condi-\ntional generative models for ground motion waveform synthesis, particularly in how our model learns\nand represents wave propagation. Also, it is challenging to write interdisciplinary papers that make\ncontributions both to machine learning methods and to scientific domains; and so we appreciate the\nfeedback, whichhelpsustoaddressconfusionpointsandtoclarifythecontributionforbothaudiences.\nWe address each of your concerns in detail below.\nComment 1: The model that the authors are presenting here is one dimensional. This means\nthat each event-station pair is trained independently of any other station, despite the fact that\nthe coordinates of the source and station are introduced as conditional variables. There is no\nglobal variable for the VAE to determine that two adjacent stations with the same source coor-\ndinates, from the same event, come from the same earthquake. In this sense, the model IS NOT\nLEARNING wave propagation.\nReply: As the reviewer points out, our model generates time series data based on earthquake source\nand sensor locations. We would like to clarify what physics we are trying to learn: our objective is\nnot to simulate or learn physical wave propagation through a 3D Earth volume explicitly; instead, our\nmodel learns a stochastic mapping of a net effect of physical waveforms between source and sensor\nlocations. We have clarified this, as we outline in detail below. By leveraging numerous earthquake-\nsensor pairs, our approach captures wave propagation characteristics along multiple paths through\na complex neural network-based machine learning (ML) model rather than relying on a Gaussian\nprocess framework. Analogous to a non-ergodic Ground Motion Model (GMM), when our ML model\nprovides reasonable results for a broad set of earthquake-sensor pairs, it demonstrates the capability\nto generalize ground motion characteristics across arbitrary source-receiver configurations within the\nregion of interest. Furthermore, by conditioning the model on continuous geographic coordinates,\nour approach facilitates the inference/generation of spatially continuous waveforms, as evidenced by\nthe smooth spatial patterns observed in the Fourier Amplitude Spectrum (FAS) maps across multiple\nfrequencies(e.g.,2Hz,5Hz,and10Hz). RegardingyourconcernaboutspatialcontinuityinComment\n3, wehaveprovidedadditionalclarificationsandresultsinthecorrespondingresponse. Pleaseseemore\ndetails on Pages 4-6 of the response letter.\nA primary benefit of this stochastic mapping framework, compared to explicitly modeling full 3D\nwave propagation, is significantly reduced complexity. Learning 3D propagation effects would require\nsubstantially larger models and considerably greater computational costs. Additionally, our model\navoids explicit assumptions regarding spatial correlations. Different from Gaussian Process-based non-\nergodicGMMs,wherespatialcorrelationstructuresarepredetermined,ourapproachlearnsthesespatial\ndependencies directly from observed data, providing increased flexibility.\nTo help address the reviewer’s concern, we have revised the manuscript to clarify the following points:\n• We are not learning to directly reproduce 3D wave propagation physics.\n• Our model learns a stochastic mapping that captures wave propagation between earthquake\nsources and sensors.\n• The updated median FAS maps exhibit that our model produces spatially continuous patterns.\nNamely, despite the conditioning of individual source-receiver pairs, the generated wavefields\nand FAS maps demonstrate spatial continuity, indicating our model’s capability to infer spatial\ncoherence from data.\nuniquelydefineasingleearthquakeevent, asthegeneratedwaveformsatstationscanexhibitvariability\nover the same source parameters that we used in the conditional variables. For example, this implies\nthat sources with identical locations and magnitudes may still differ in their source mechanisms. How-\never, our dataset mainly includes events from one fault system, where earthquake source mechanisms\nshow minimal variability. Thus, the variations of the wavefields among such different earthquakes are\nsmall. We demonstrate that median ground motions represent spatial continuity and have incorpo-\nratedfurtherclarificationsintheAbstract(Page1)andSectionIntroduction (Pages3-4)oftherevised\nmanuscript.\nComment 2: Since the model is not learning wave propagation, the title LEARNING PHYSICS\nis false. Its learning a mapping from one point to another, independent to each other, but no\ncoherency; and since the medium from the source to the surface is fixed, the high frequencies\n(here 10 Hz) are learning some shallow crustal amplification functions which one could similarly\ncapture by replacing the source-station coordinates with Vs30. Instead, the authors are trying to\nlearn a 3D problem with a 1D model which, not only is impossible, but also leads to erroneously\nspiked strong motion records like the ones they show in Figure 7.\nReply: As stated in our response to Comment 1, the model is not designed to reproduce full 3D\nwave propagation physics. Instead, it learns a stochastic mapping that captures the effective wave\npropagation characteristics between earthquake sources and stations. Importantly, the model demon-\nstrates its ability to capture spatial continuity, as evidenced by the generated FAS maps (see response\nto Comment 3 for further details). Regarding the use of the term “learning physics” in the title,\nwe agree with the reviewer that our model does not learn the 3D wave propagation equations. We\nclarify that this phrase is intended to convey the model’s ability to learn a stochastic mapping for wave\npropagation from data within a data-driven framework. Lest other readers have the same confusion,\nwe have revised the Abstract and Section INTRODUCTION to clarify this interpretation. Please refer\nto Pages 3-4 of the revised manuscript.\nSecond, as we clarify in our response to Comment 1, we train a single ML model that can generate\narbitrary pairs of source and sensors. To reproduce both dynamics (e.g., FAS) and kinematics (e.g., P\nand S wave arrival times and durations, Figure 4) for source-sensor pairs not included in the training\ndataset, learning site conditions is not sufficient. For example, site condition variations cannot explain\ndifferences in arrival times at the same source-site distances. Therefore, the model also needs to learn\nproperties of the subsurface velocity model and how seismic waves propagate through it.\nThird, not all generated waveforms exhibit spikes. We present selected examples in Section B.1\nWaveforms showing moderate performance to illustrate specific cases. VAE models are known to\ngenerate relatively smooth amplitude information, which might introduce numerical artifacts during\nphase retrieval. This may partially explain the occurrence of spikes. The quality of the generated\nwaveformscouldbeimprovedbyemployingalternativegenerativeframeworks, suchasdiffusionmodels\nin our follow-up work [1]. We have included more clarifications in Section DISCUSSION (Page 13)\nand B.1 Waveforms showing moderate performance (Page 32) of the revised manuscript.\nReferences:\n[1] Bi, Zhengfa, et al. “Advancing data-driven broadband seismic wavefield simulation with multi-\nconditional diffusion model.” arXiv preprint arXiv:2501.14348 (2025).\nComment 3: Another way to see what I am saying is to observe the incoherent pattern that\nemerges in the longer period range, at 2Hz for example, where ‘path’ and ‘site’ effects become\n3\nmuddled. There, Figure 14 shows that adjacent stations respond independently of each other for\nthe same ‘event’. In other words, there is no spatial continuity. The authors are just querying\nthe model at adjacent stations, whose waveforms, however, even for the same location of source\nand magnitude of earthquake, will not necessarily be spatially correlated like they are originating\nfrom the same event.\nReply: Thanks for your comment. We provide a detailed clarification on the spatial continuity in\nthis response. We address your concern from three perspectives: (i) analysis of median FAS maps\n2hz at different frequencies, (ii) evaluation of semi-variograms of the predicted ground motions, and (iii)\na sensitivity test using a modified earthquake scenario with a slightly perturbed source location. In\nour study, the spatial continuity is achieved by conditioning the model on continuous geographic\ncoordinates. \nFirst, we noticed that the FAS map at 2 Hz presented in the original Figure 14 exhibited some\nnoise and incoherent spatial patterns in certain regions. As the reviewer has suggested, we also\nupdated all relevant FAS maps in the manuscript (Figures 2, 13-16, 23) to display the median of 30\nrealizations, which now provides a more stable and representative summary of the model output.\nSpecifically, we show the revised FAS maps at 2 Hz and 5 Hz in FIG. 1 and 2, which have exhibited\nspatial coherence.\na CGM-GM b Ergodic GMM c Non-ergodic GMM\nFIG. 1. The comparison of FAS maps at 2 Hz from CGM-GM, ergodic GMM, and non-ergodic GMM. The\nCGM-GM result is based on the median FAS map from 30 realizations. The red star and blue triangles represent\nthe source and stations, respectively.\nSecond, in addition to the visual inspection, we compute the semi-variogram [1] of the prediction as\na function of separation distance (∆S) between the predicted stations, which is another direction for\nmeasuring spatial correlations between ground motions. As shown in FIG. 3, we observe that at very\nshort separation distances (e.g., [0,5] km), the semi-variogram value is very close to 0, meaning that\nthe correlation of the predicted values between nearby stations is very close to 1, which indicates that\nthereisspatialcontinuityinthepredictionsacrossstationlocations. Moreover,weprovideananalytical\nkernel function fitted to the empirical semivariogram, expressed as\n(cid:18) b·∆S2(cid:19)\nγ(∆S) = 1−exp − . (1)\nThe estimated values of b for the 2 Hz, 5 Hz, and 10 Hz cases are 4.858×10−3, 4.033×10−3, and\n3.86×10−3, respectively. These correspond to correlation lengths of ρ = 14.35, 15.75, and 16.10\n4\n5hz\na CGM-GM b Ergodic GMM c Non-ergodic GMM\nFIG. 2. The comparison of FAS maps at 5 Hz from CGM-GM, ergodic GMM, and non-ergodic GMM. The\nCGM-GM result is based on the median FAS map from 30 realizations.\nkm, respectively. The correlation length defines the spatial distance beyond which the correlation\nbetween data points significantly decreases. A correlation length of approximately 15 km is reasonable,\nrejecting the assumption of white noise. These findings imply the presence of spatial structure and\nprovidefurtherquantitativeevidencesupportingthespatialcontinuityofthegeneratedgroundmotions.\nFIG. 3. The empirical correlations at 2 Hz (left), 5 Hz (mid), and 10 Hz (right).\nTo further demonstrate spatial continuity between predictions across event locations, we slightly move\nthe coordinates of the event location and compute new FAS maps. Specifically, the original seismic\nevent is characterized by a magnitude of 3.84 and a depth of 7.94 km. The epicenter, denoted by a red\nstar, is located at a geographic position with a latitude of 37◦51.6′N and a longitude of 122◦15.6′W.\nWe modify the source location 1.4 km South while keeping all other physical parameters unchanged.\nThecomparisonbetweenthenewandoriginalFASmapsat10HzisshowninFIG.4. Theresultingmap\nexhibits a good similarity to the original, further demonstrating that our CGM-GM model maintains\nspatial continuity in its predictions across different event locations.\nWe have included an additional subsection in Supplementary Material to demonstrate the learned\nspatial continuity. Please refer to Section B.11 Spatial correlations on Pages 42-43.\nReference:\n[1] Jayaram, Nirmal, and Jack W. Baker. “Correlation model for spatially distributed ground-motion\nintensities.” Earthquake Engineering & Structural Dynamics 38.15 (2009): 1687-1708.\n5\na b\nNewsourcelocation Originalsourcelocation\nFIG. 4. The predicted FAS maps (10 Hz) with new and original source locations from our CGM-GM model. The\nresults are both based on the median FAS values from 30 realizations.\nComment 4: In that sense, the best comparison the authors can make is against the ergodic\nGMMs in Nor Cal. These treat each GM time series separately with some conditional variables.\nThenon-ergodicGMMsimplementwithineventvariabilitytermsbystudyingtheresidualsrelative\nto the ergodic terms of all the event terms at the same time across the strong motion network\nusing GPs. Why dont they show in Figure 2 comparison against the ergodic GMMs? This is\na much more fair comparison since the non-ergodic GMM model is a higher dimensional model\nthan the cGM-GM. And definitely statements like ‘outperforming a state-of the art non-ergodic\nGMM’ should be deleted since they are misinforming the audience that the model is capturing\neffects that it is not.\nReply: Thankyouforthesuggestion. WehaveaddedacomparisonbetweentheCGM-GMpredictions\nand the ASK14 model in terms of spectral acceleration (Sa) within the 2–15 Hz frequency range for\nevents with magnitudes greater than 3. To compute the ASK14 predictions, we extract the V\nS30\nand Z information from the USGS San Francisco Bay Region 3D Seismic Velocity Model (Version\n1.0\n21.1). It is worth noting that some deviations and potential adjustments to these parameters have\nbeen suggested in prior studies [1, 2]. As shown in FIG. 5, the results show that CGM-GM consistently\nreduces both the residuals and their standard deviation relative to the ergodic ASK14 predictions\nacross the entire frequency range. We have included an additional section to clarify this. Please refer\nto Section B.9 Residuals computed for the ASK14 model on Page 41 of the revised manuscript.\nMoreover,ourintentionwastohighlightthestrongpotentialofourapproachtocomplementingexisting\nmethods for generating ground motions. We deleted the sentence of “...outperforming a state-of-the-\nart non-ergodic GMM...” and clarified that we are aiming to complement physics-based simulations\nand non-ergodic empirical ground motion models. This is aspirational rather than claiming that it has\nalready done so. In line with the reviewer’s suggestion, we have revised the corresponding statement\nin the Abstract of the revised manuscript to clarify this point (Page 1).\nReferences:\n[1] Pinilla-Ramos, Camilo, et al. “Performance evaluation of the USGS velocity model for the San\nFrancisco Bay Area.” Earthquake Spectra 41.1 (2025): 457-494.\n[2] Lavrentiadis, Grigorios, et al. “Data-driven characterization of near-surface velocity in the San\nare defined as the logarithmic difference between the model predictions and the observed ground truth values.\nThe solid line and the shaded area denote the mean curves and the region of mean ± std.\nFrancisco bay area: A stationary and spatially varying approach.” Earthquake Spectra (2024):\n87552930251320666.\nComment 5: Despite the fact that cGM-GM uses event and station coordinates, because the\nstations ‘dont know’ that their motions come from the same event, it cannot be used to plot\nscenaria. In that sense, all the scenaria plots like Figure 2b and 2c should be replaced with some\nmedian amplitude of many realizations, to be compared with median amplitudes of ergodic GMM\nplots.\nReply: We agree with the reviewer that the generated ground motions from CGM-GM at a given\nstation cannot uniquely reflect a specific earthquake instance with the same magnitude and location,\ndue to variability in rupture processes, source mechanisms, and source-time functions. The detailed\nclarifications can be found in the response to Comment 1. Moreover, as suggested by the reviewer,\nwe have updated the relevant figure to show the median FAS computed from 30 realizations, which\nfurther improves spatial coherence while mitigating the effects of stochastic variability. Please refer to\nFigures (2, 13-16, 23) in the revised manuscript.\nComment 6: Page 7, line 3 of the text: [... implicit assumption of a homogeneous Earth...]\nshould be replaced with the [... implicit assumption of a radially homogeneous Earth...]. We\nshould also note here that this is a limitation of the baseline model that the authors selected;\na more fair comparison would be to select a baseline model conditional on magnitude, distance,\nand a site response proxy (the same they used in their empirical GMM maps) in which case they\nwould recover most of the capabilities that they are able to replicate with their cGM-GM model,\npossibly with a better time series synthesis capability.\nReply: Thanks for the suggestion. We have corrected the corresponding sentence to “...implicit\nassumption of a radially homogeneous Earth...” (Page 7). The non-ergodic and ergodic empirical\nmodels used in this study do not incorporate site response proxies, due to the unavailability of reliable\nsite condition information, as described in Section Embedding conditional variables (Page 17) and\nbaseline model.\n.\nComment 7: Topofpage8: againerroneousstatement: theagreementofthenon-ergodicmedian\nhere is purely phenomenological and if they plot the ergodic median the agreement will likely be\nbetter. Also, they should not be plotting scenaria, but median of scenaria, since their model is\nalso ergodic (in the sense that it is learning independent pairs of source-stations and not within\nevent variability like non-ergodic empirical models).\nReply: Weagreewiththereviewerthatpresentingmapsfromasinglerealizationmaynotbesufficient.\nWehaveupdatedallrelevantFASmapsintherevisedmanuscriptwiththemedianacross30realizations\ninstead of a single scenario. This strategy has effectively reduced pixel-level variability and provided a\nmore stable and representative illustration of the results. Additionally, in comparison with the ASK14\nmodel discussed in our response to Comment 4, we also compute the residuals for CGM-GM using\nthe median values from 30 realizations. The results show that CGM-GM consistently achieves lower\nresiduals and reduced standard deviation compared to the ergodic ASK14 predictions across the full\nfrequency range. Further details are provided in the response to Comment 4.\nComment 8: Page 10: bottom of first paragraph: Figure 4e shows bias in the arrival times, not\na good correlation as the authors claim. The heat map is clearly off center.\nReply: Thanks for pointing it out. First of all, we would like to emphasize that the duration plot\n(Figure 4e) is more scattered than arrival times (Figures 4c and 4d), which indicates the difficulties of\naccurately generating and measuring seismic durations. Second, the bias primarily occurs for durations\nwithin the 0–3 second range, corresponding to smaller-magnitude earthquakes. Due to lower signal-to-\nnoise(S/N)ratiosinsuchcases, itischallengingtopreciselymeasurethedurationsfromactualseismic\ndata since the relatively high noise level can artificially reduce the measured duration compared to the\nactual signal duration. Note that the generative model does not explicitly simulate noise, resulting in\ninherently higher S/N ratios. This potentially amplifies the discrepancy. Nonetheless, despite these\nlimitations, the overall duration trends between observed and generated data demonstrate a good\ngeneral correlation. We have included more clarifications regarding this aspect in Section P and S\narrival times and durations in the revised manuscript (Page 10).\nComment 9: Figure 5b and description of Figure on page 10: The SA trends at 0.1sec are\ncompletely off: the empirical model saturates at M8 due possibly to nonlinear site response and\nthe cGM-GM increases almost exponential (due to the fictitious spikes we described above). The\nfact that the lines cross at some point doesn’t mean there is agreement.\nReply: Thank you for your comment. We would like to clarify that our intention was not to suggest\nthat the CGM-GM performance at 0.1 sec is acceptable. Rather, we aim to illustrate the potential\ncapability of the CGM-GM model in capturing ground motion intensity. Specifically, by “potential”,\nwe refer to the model’s ability to prevent an excessive and rapid increase in predicted intensities at\nhigher magnitudes. To avoid confusion, we have clarified this point in the revised text and emphasized\nthe current challenge explicitly. Please see Section Generalization to larger-magnitude scenarios on\nPages 11-12 of the revised manuscript.\nplotted along with the 10Hz. Engineers don’t care about frequencies as high as 10Hz nearly as\nmuch as they care about low frequencies.\nReply: Thankyouforyourcomment. Weagreewiththereviewerandhaveaddedmapsofthemedian\npredictions from CGM-GM versus ergodic GMM at 2 Hz and 5 Hz. Please see Figures 13-14 on Page\n38. However, we would like to emphasize that ground-motion prediction at higher frequencies is of\nhigh importance to engineers for the design of critical infrastructures such as dams or nuclear power\nplants, for which the natural frequency of sensitive equipment lies in the higher frequency domain [1].\nReference:\n[1] Singh, Sugandha, and Abhinav Gupta. “Seismic response of electrical equipment subjected to\nhigh–frequency ground motions.” Nuclear Engineering and Design 374 (2021): 111046.\nIII. RESPONSE TO REVIEWER 3\nGeneral comment: The authors satisfactorily addressed almost all of the comments raised by\nthis reviewer in the first round, except Comment 2. The sense of the question was if CGM-GM\ncould be adopted zero-shot. As per the explanation given, the suspect is that the CGM-GM is\nnot fit for zero-shot. A comment on this would be highly appreciated.\nThe transparency of showing the goodness of fit in Fig.10 is highly appreciated.\nReply: Thank you very much for your positive feedback. Below, we provide clarification on your\nquestion regarding zero-shot generalization.\nCGM-GM is trained on a dataset of small events (M < 4). Within that magnitude range, the model\ncanbeusednaturallyina“zero-shot”manner. However,scalingthepretrainednetworkdirectlytolarge\nevents introduces challenges that exceed the information (e.g., wave physics) encoded in the original\ntraining set. Specifically, the earthquake source rupture processes change from point sources to finite\nrupture planes. Hence, we consider supervised fine-tuning for the generalization to larger events. This\nreason is presented in Section Generalization to larger-magnitude scenarios (Pages 10-12).\nAdditionally, we show both zero-shot generalization (pre-trained) and fine-tuning performance in Sec-\ntion B.10 Details of fine-tuning and GMM comparisons, especially Figure 21 (Pages 41-42). The\ncomparative results have demonstrated that the fine-tuning strategy can effectively improve the mean\npredictions and reduce the variations of Sa values under large-magnitude scenarios.\nWe hope this clarifies the intended operating regime of CGM-GM and the practical procedure for\nextending it to larger events.\n9"}, {"round": 3, "reviewer_comments": {"reviewer_2": "I would like to thank the authors for addressing my comments to some extent. There still are unresolved concerns in the revised version that I need to see addressed, however. I want to make sure that the authors properly highlight the contributions of the manuscript, which are important, without misleading the earthquake engineering audience who may be interested in using the results of their work.\nFor this purpose, I use sections of the text to point out inconsistencies between what the authors's model is and what they say it is in the text (line numbers would have been helpful in that sense but I'll use page numbers and approximate locations):\n1. p3 \"Our results demonstrate the excellent performance of our proposed model in learning the underlying wave propagation char- acteristics...\" There is no evidence, nor it is expected to be any evidence of the model to learn wave propagation: for that, the authors would have to show spatiotemporal evolution of the wavefield which they dont because this is not what the model is supposed to be doing. All references to how great a job the model does in learning wave propagation should be removed. They are false. The model is learning some form of frequency response median from irregularly spaced stations and that is great. Say that.\n2. p3. \"Our model learns a stochastic mapping that captures wave propagation characteristics between earthquake sources and stations...\" again, rephrase the misleading statement.\n3. p4. \"Nevertheless, we demonstrate that the model is capable of capturing spatial continuity in the median generated ground motions.\" No you dont. You demonstrate that the model is capable of capturing spatial continuity in the median generated FAS.\n4. p6. \"We utilize an existing earthquake event that occurred at 3:16 AM on October 21, 2011.\" This is misleading because, as mentioned earlier, the model is not designed to capture scenaria but median frequency responses. The authors should make it explicit that the model cannot capture event level spatiotemporal variability.\n5. p13. \"The primary objective of our proposed method is to generate spatially continuous and physically meaningful ground motions. Our model learns a stochastic mapping that captures wave propagation characteristics through a combined representation of source, site, and path effects.\" Again, rephrase the capturing of wave propagation. Call it frequency response or something that explicitly takes out the spatiotemporal continuity misinterpretation that runs across the manuscript.\n6. The abstract needs to be rephrased accordingly.\nOther comments:\n7. Why compare the scaling of large magnitude events to ASK14 and not to the non-ergodic GMM for California by Lavrentiadis et al (2023)? It seems to be a much more appropriate comparison given the flow of your work than switching metric and suddenly comparing to the Sa of ASK14.\n(Remarks on code availability)"}, "author_response": "Learning Physics for Unveiling Hidden Earthquake Ground\nMotions via Conditional Generative Modeling\nNCOMMS-24-45416B\nPu Ren, Rie Nakata, Maxime Lacour, Ilan Naiman, Nori Nakata, Jialin Song, Zhengfa Bi, Osman\nAsif Malik, Dmitriy Morozov, Omri Azencot, N. Benjamin Erichson, and Michael W. Mahoney\nDear Editor and Reviewer #2,\nWe sincerely appreciate your constructive comments and valuable suggestions. We are pleased that Reviewer #2 acknowledged the value of our research, and we have made the suggested changes\naccordingly. Doing interdisciplinary work (in this case, making contributions to machine learning,\nearthquakeengineering, andseismology)isdifficult, e.g., asterminologybetweendifferentfieldsdiffers.\nWe have modified our terminology, lest any readers be confused as to our claims and contributions.\nOur point-by-point responses are provided below marked in blue, and the corresponding manuscript\nrevisions are highlighted in red. We hope our manuscript will be found worthy of publication in this\nrevised form.\nI. RESPONSE TO REVIEWER 2\nGeneral comment: I would like to thank the authors for addressing my comments to some ex-\ntent. There still are unresolved concerns in the revised version that I need to see addressed, how-\never. I want to make sure that the authors properly highlight the contributions of the manuscript,\nwhich are important, without misleading the earthquake engineering audience who may be inter-\nested in using the results of their work.\nFor this purpose, I use sections of the text to point out inconsistencies between what the authors’s\nmodel is and what they say it is in the text (line numbers would have been helpful in that sense\nbut I’ll use page numbers and approximate locations):\nReply: We thank the reviewer for acknowledging the importance of our contributions. We agree\nthat the scope, capabilities, and limitations of our models should be precisely conveyed. Each of your\nconcerns is addressed in detail below.\nComment 1: p3 “Our results demonstrate the excellent performance of our proposed model\nin learning the underlying wave propagation char- acteristics...” There is no evidence, nor it\nis expected to be any evidence of the model to learn wave propagation: for that, the authors\nwould have to show spatiotemporal evolution of the wavefield which they dont because this is not\nwhat the model is supposed to be doing. All references to how great a job the model does in\nlearning wave propagation should be removed. They are false. The model is learning some form\nof frequency response median from irregularly spaced stations and that is great. Say that.\nReply: We have clarified our statement to avoid confusion about the model’s capability. The original\nsentenceclaimingthemodel“learnstheunderlyingwavepropagationcharacteristics”hasbeenreplaced\nwith:\n1\n“Our results demonstrate the performance of our proposed model in producing Fourier Amplitude\nSpectra (FAS) maps and their spatial heterogeneity as well as reproducing properties such as P and S\narrivals, and waveform durations.” (Page 3).\nWe demonstrate that the model is capable of capturing P and S arrivals and duration, as shown in\nthe cross-plots in Figure 4(c-e). In the previous manuscript, we referred to these features collectively\nas “wave propagation characteristics”. This wording might not be precise, and we have revised it for\naccuracy.\nAdditionally, we have revised other related statements regarding “learning wave propagation”. Specifi-\ncally, on Page 3, we have changed “This approach enables the representation of wave propagation and\nmeaningful spatial variations in ground motions” to:\n“This approach presents spatial continuity in median spectral components of the generated ground\nmotions.”.\nOn Page 13, we have revised “The most exciting part is that our CGM-GM can capture the underlying\nspatialheterogeneityandphysicalcharacteristics, asevidencedbythegenerationofrealisticFASmaps.\nWe conduct a comparative assessment of our CGM-GM against baseline models, including the CGM-\nbaseline and a state-of-the-art non-ergodic GMM. The results demonstrate that our method performs\ncomparably to the most advanced non-ergodic GMM. This validates that incorporating geospatial\ncoordinates as conditional variables effectively enables our model to learn spatial heterogeneities.” to:\n“The most exciting part is that our CGM-GM can capture the spatial heterogeneity of FAS without\nexplicit physics constraints. We conduct a comparative assessment of our CGM-GM against baseline\nmodels, including the CGM-baseline and state-of-the-art non-ergodic GMMs. The results demonstrate\nthat our method performs comparably to the most advanced non-ergodic GMMs.”.\nComment 2: p3. “Our model learns a stochastic mapping that captures wave propagation char-\nacteristics between earthquake sources and stations...” again, rephrase the misleading statement.\nReply: We have removed this sentence “Our model learns a stochastic mapping that captures wave\npropagation characteristics between earthquake sources and stations, governed by physical principles\nand Earth structures.”, to avoid any misleading implications. Please see Pages 3-4.\nComment 3: p4. “Nevertheless, we demonstrate that the model is capable of capturing spatial\ncontinuity in the median generated ground motions.” No you dont. You demonstrate that the\nmodel is capable of capturing spatial continuity in the median generated FAS.\nReply: The sentence has been revised from “Nevertheless, we demonstrate that the model is capable\nof capturing spatial continuity in the median generated ground motions.” to\n“Nevertheless, we demonstrate that the model is capable of capturing spatial continuity in the median\ngenerated FAS.” (Page 4).\nComment 4: p6. “We utilize an existing earthquake event that occurred at 3:16 AM on October\n21, 2011.” This is misleading because, as mentioned earlier, the model is not designed to capture\nscenaria but median frequency responses. The authors should make it explicit that the model\ncannot capture event level spatiotemporal variability.\nReply: Wehaverevisedthetexttoavoidmisleadinginterpretations. Theoriginalsentence“Weutilize\nan existing earthquake event that occurred at 3:16 AM on October 21, 2011.” has been changed to:\n2\n“We use the source parameters of an existing earthquake event that occurred at 3:16 AM on October\n21, 2011. Our model is not intended to reproduce ground motions of this specific event; rather, the\nsource parameter information is used as a representative case of an earthquake along the Hayward\nfault.” (Page 6).\nComment 5: p13. “The primary objective of our proposed method is to generate spatially\ncontinuous and physically meaningful ground motions. Our model learns a stochastic mapping\nthat captures wave propagation characteristics through a combined representation of source, site,\nand path effects.” Again, rephrase the capturing of wave propagation. Call it frequency response\nor something that explicitly takes out the spatiotemporal continuity misinterpretation that runs\nacross the manuscript.\nReply: Wehavemodifiedthetextfrom“Theprimaryobjectiveofourproposedmethodistogenerate\nspatiallycontinuousandphysicallymeaningfulgroundmotions. Ourmodellearnsastochasticmapping\nthat captures wave propagation characteristics through a combined representation of source, site, and\npath effects.” to:\n“Our model learns to capture the properties such as spatially continuous FAS, through a combined\nrepresentation of source, site, and path effects.” (Page 13).\nComment 6: The abstract needs to be rephrased accordingly.\nReply: We have revised the abstract to align with the clarified scope and contributions of our work.\nSpecifically, we have modified the text from “We propose a novel artificial intelligence (AI) simulator,\nConditional Generative Modeling for Ground Motion (CGM-GM), to synthesize high-frequency and\nspatiallycontinuousearthquakegroundmotionwaveforms. CGM-GMleveragesearthquakemagnitudes\nand geographic coordinates of earthquakes and sensors as inputs, learning a stochastic mapping that\ncaptures wave propagation characteristics governed by physical principles and Earth heterogeneities,\nwithout explicit physics constraints” to\n“We propose a novel artificial intelligence (AI) waveform generator, Conditional Generative Modeling\nfor Ground Motion (CGM-GM). CGM-GM leverages earthquake magnitudes and geographic coordi-\nnates of earthquakes and sensors as inputs, capturing spatially continuous Fourier amplitude spectra\n(FAS) as well as properties such as P and S arrivals, and waveform durations without explicit physics\nconstraints.”.\nComment 7: Why compare the scaling of large magnitude events to ASK14 and not to the non-\nergodicGMM for California by Lavrentiadis et al(2023)? It seems to be a much more appropriate\ncomparison given the flow of your work than switching metric and suddenly comparing to the Sa\nof ASK14.\nReply: We agree with the reviewer and have updated the manuscript to compare Fourier Amplitude\nSpectra (FAS) results with the non-ergodic GMM for California by Lavrentiadis et al (2021) [1], as\nLavrentiadis et al (2023) [2] was developed for spectral acceleration. Please see the updated Figures\n5 (Page 12) and 20 (Page 42) and related descriptions (Pages 10-12, 42-43) in the manuscript.\nReferences:\n[1] G.Lavrentiadis,N.A.Abrahamson,N.M.Kuehn,Anon-ergodiceffectiveamplitudeground-motion\nmodel for California, Bulletin of Earthquake Engineering (2021) 1–32doi:10.1007/ s10518-021-01206-\nw.\n[2] G. Lavrentiadis, N. Abrahamson, A non-ergodic spectral acceleration ground motion model for\nCalifornia developed with random vibration theory, Bulletin of Earthquake Engineering 21 (11) (2023)\n5265–5291.\n3"}, {"round": 4, "reviewer_comments": {"reviewer_2": "Again, I would like to thank the authors for taking the time to revise their manuscript and bring the manuscript content closer to what their work actually is about.\nIn their eagerness to achieve more, however, they --again-- have introduced incorrect representations of the material presented. I am baffled by this attitude, but I do believe that this group of young computer scientists with fresh ideas should gain visibility in the community which is why I continue to suggest revisions and do not recommend that the manuscript is rejected. Here are some more items for the authors to reconsider:\n1. The title and abstract need revision. The word 'physics' is misleading in the title and now the abstract is also misleading:\n[Predicting high-fidelity ground motions for future earthquakes is crucial for seismic hazard assessment and infrastructure\nresilience. Conventional empirical simulations suffer from sparse sensor distribution and geographically localized\nearthquake locations, while physics-based methods are computationally intensive and require accurate representations of\nEarth structures and earthquake sources. --> We propose a novel artificial intelligence (AI) waveform generator, Conditional\nGenerative Modeling for Ground Motion (CGM-GM). CGM-GM leverages earthquake magnitudes and geographic\ncoordinates of earthquakes and sensors as inputs, capturing spatially continuous Fourier amplitude spectra (FAS) as well as\nproperties such as P and S arrivals, and waveform durations without explicit physics constraints....]\nThe abstract has been edited by thereby lies the problem: the authors are not proposing a waveform generator. They are\nproposing an FAS median mapping tool that, if fed with a white noise phase, can produce random phase waveforms, but\nthey are not physically correct! This is not acceptable at least on my end.\n2. I appreciate that the authors went ahead and incorporated the LA21 NGMM for larger events in Figure 5. Where are the\ncontinuous changes in slope in the magnitude scaling coming from ? The implementation of the NGMM is clearly incorrect\nand I am not sure where the error is by looking at the figure.\n3. Furthermore, the trend of the cGM-GM and the NGMM is opposite. The proposed model doesn't saturate at large\nmagnitude events, on the contrary -- it explodes, yet another piece of evidence that it doesnt really learn physics. The\nauthors should think really carefully of what this figure represents.\n4. The same is evident in the high frequency regime of the small events. The model is learning a lot of noise, which is why it\nis predicting --at best-- median amplitudes up to 11Hz. Everything above that, including the trend, should be commented on\nin the text, since it shows that the model doesn't extrapolate well.\n(Remarks on code availability)\nN/A"}, "author_response": "Learning Physics for Unveiling Hidden Earthquake Ground\nMotions via Conditional Generative Modeling\nNCOMMS-24-45416C\nPu Ren, Rie Nakata, Maxime Lacour, Ilan Naiman, Nori Nakata, Jialin Song, Zhengfa Bi, Osman\nAsif Malik, Dmitriy Morozov, Omri Azencot, N. Benjamin Erichson, and Michael W. Mahoney\nDear Editor and Reviewer #2,\nThank you very much for your constructive comments and valuable suggestions. We particularly\nappreciate that Reviewer 2, who has concerns with our terminology, is helping to push forward inter-\ndisciplinary research by saying that “this group of young computer scientists with fresh ideas should\ngain visibility in the community.” Based on the most recent feedback, it is our understanding that the\nremaining concerns center on (i) the terminology involving “learning physics” and (ii) the implemen-\ntation of the non-ergodic GMM. To address them, we have made two main revisions as suggested:\n• Updated the title to remove “learning physics”;\n• Updated the abstract to clarify the model scope.\nWe have also made clarifications regarding NGMM implementation below. Please see our point-by-\npoint responses marked in blue, and the corresponding manuscript revisions are highlighted in red. We\nhope our manuscript will be found worthy of publication in this revised form.\nI. RESPONSE TO REVIEWER 2\nGeneral comment: Again, I would like to thank the authors for taking the time to revise their\nmanuscript and bring the manuscript content closer to what their work actually is about.\nIntheireagernesstoachievemore,however,they–again–haveintroducedincorrectrepresentations\nof the material presented. I am baffled by this attitude, but I do believe that this group of young\ncomputerscientistswithfreshideasshouldgainvisibilityinthecommunitywhichiswhyIcontinue\nto suggest revisions and do not recommend that the manuscript is rejected. Here are some more\nitems for the authors to reconsider:\nReply: We thank the reviewer for their continued effort in helping us make the manuscript more\nscientifically grounded. First of all, we would like to clarify that we are not a group of young computer\nscientists: the lead author has a background in civil and environmental engineering and now works in\nscientificmachinelearning; wehavetwoseismologists,eachwithwelloveradecadeofpost-dissertation\nresearchexperience,ontheteam;oneengineeringseismologistwhospecializesinstrongmotionworking\nwith leading experts in the field; and the last author is a physicist turned statistician and machine\nlearning expert, with decades of interdisciplinary experience.\nDoing interdisciplinary research is challenging since it could lead to differences in terminology across\nfields, whichmaycauseunintendedmisunderstandings. Wehavetakentheconcernsseriouslyandhave\nrevised the title and abstract, as well as clarified the NGMM implementations. We address each of\nyour concerns in detail below.\nComment 1: The title and abstract need revision. The word 'physics' is misleading in the title\nand now the abstract is also misleading:\n[Predictinghigh-fidelitygroundmotionsforfutureearthquakesiscrucialforseismichazardassess-\nment and infrastructure resilience. Conventional empirical simulations suffer from sparse sensor\ndistribution and geographically localized earthquake locations, while physics-based methods are\ncomputationallyintensiveandrequireaccuraterepresentationsofEarthstructuresandearthquake\nsources. → We propose a novel artificial intelligence (AI) waveform generator, Conditional Gen-\nerative Modeling for Ground Motion (CGM-GM). CGM-GM leverages earthquake magnitudes\nand geographic coordinates of earthquakes and sensors as inputs, capturing spatially continuous\nFourier amplitude spectra (FAS) as well as properties such as P and S arrivals, and waveform\ndurations without explicit physics constraints....]\nThe abstract has been edited by thereby lies the problem: the authors are not proposing a\nwaveform generator. They are proposing an FAS median mapping tool that, if fed with a white\nnoise phase, can produce random phase waveforms, but they are not physically correct! This is\nnot acceptable at least on my end.\nReply: Thank you for your comment. It appears that the current title suffers from the terminology\nchallenge, as it may give the impression that we claim to be learning complete earthquake waveform\npropagation physics, which is not our intent. Our method learns data-driven patterns that relate to\nphysical properties of ground motion, including FAS, P, S-wave arrivals, and durations (as we replied\nin the 2nd and 3rd revisions), and we are not replacing a physics-based simulation. However, we\nare willing to change the title to “Learning Earthquake Ground Motions via Conditional Generative\nModeling”. We are also open to any other suggestions that you may have.\nSecond, we are willing to further revise the abstract from\n“...Weproposeanovelartificialintelligence(AI)waveformgenerator, ConditionalGenerativeModeling\nfor Ground Motion (CGM-GM). CGM-GM leverages earthquake magnitudes and geographic coordi-\nnates of earthquakes and sensors as inputs, capturing spatially continuous Fourier amplitude spectra\n(FAS) as well as properties such as P and S arrivals, and waveform durations without explicit physics\nconstraints...”\nto\n“...We propose a novel artificial intelligence (AI) spectrogram generator, Conditional Generative Mod-\neling for Ground Motion (CGM-GM). CGM-GM leverages earthquake magnitudes and geographic co-\nordinates of earthquakes and sensors as inputs, when postprocessed with phase information, capturing\nspatially continuous Fourier amplitude spectra (FAS) as well as properties such as P and S arrivals,\nand waveform durations, without explicit physics constraints...”\nRespectfully, apartfromFASmedianinformation,ourmethodcanalsopredictpropertiessuchasPand\nS arrivals, and waveform durations, as shown in Figure 4 (Page 12 of our manuscript). The conditional\nVAE generates spectrograms, and the phase information is obtained from the corresponding predicted\namplitude data using a numerical method, the Griffin-Lim algorithm, as described in Section Network\nDesign (Pages 15-16) of our manuscript. This phase estimation method has also been used in the\nprior AI-enabled waveform generation workflow [1].\nReference:\n[1]Esfahani,RezaDD,etal. “TFCGAN:Nonstationaryground-motionsimulationinthetime–frequency\ndomain using conditional generative adversarial network (CGAN) and phase retrieval methods.” Bul-\nletin of the Seismological Society of America 113.1 (2023): 453-467.\nComment 2: I appreciate that the authors went ahead and incorporated the LA21 NGMM for\nlarger events in Figure 5. Where are the continuous changes in slope in the magnitude scaling\ncoming from ? The implementation of the NGMM is clearly incorrect and I am not sure where\nthe error is by looking at the figure.\nReply: The NGMM prediction was obtained by combining the ergodic BA19 model [2] and the non-\nergodic terms obtained by LA21 [3]. We use the values of the non-ergodic terms provided by Dr.\nLavrentiadis, the author and developer of the LA21, and a software written by Dr. Bayless, the author\nand developer of the BA19. These softwares were used and validated in multiple engineering projects\nled by Dr. Abrahamson, with whom Dr. Lacour works. We confirmed with Dr. Lavrentiadis that the\npredictions are valid and correct. We also confirmed that the continuous change does not exist for a\nmagnitude scaling plot of a single event-station pair, and appears when we plot all 73 stations. We\nconsider that these were caused by substantial variations in input parameters (e.g., Vs30) across the\nstations in the San Francisco Bay Area [4]. Therefore, this is not a problem with our implementation.\nReference:\n[2] Bayless, Jeff, and Norman A. Abrahamson. “Summary of the BA18 ground-motion model for\nFourier amplitude spectra for crustal earthquakes in California.” Bulletin of the Seismological Society\nof America 109.5 (2019): 2088-2105.\n[3] Lavrentiadis, Grigorios, Norman A. Abrahamson, and Nicolas M. Kuehn. “A non-ergodic effective\namplitude ground-motion model for California.” Bulletin of Earthquake Engineering 21.11 (2023):\n5233-5264.\n[4] Hirakawa, Evan, and Brad Aagaard. “Evaluation and updates for the USGS San Francisco bay\nregion 3d seismic velocity model in the east and North Bay portions.” Bulletin of the Seismological\nSociety of America 112.4 (2022): 2070-2096.\nComments 3 and 4: Furthermore, the trend of the cGM-GM and the NGMM is opposite. The\nproposed model doesn’t saturate at large magnitude events, on the contrary – it explodes, yet\nanother piece of evidence that it doesnt really learn physics. The authors should think really\ncarefully of what this figure represents.\nThe same is evident in the high frequency regime of the small events. The model is learning a lot\nof noise, which is why it is predicting –at best– median amplitudes up to 11Hz. Everything above\nthat, including the trend, should be commented on in the text, since it shows that the model\ndoesn’t extrapolate well.\nReply: Because the model is trained on small earthquake data within a limited frequency range\nwithout explicitly incorporating wave equations or subsurface information, it does not generalize to\nlarge-magnitude earthquakes (addressed in the 1st and 2nd revisions) nor perform perfectly on high-\nfrequency ranges (as described in Lines 17-18 of Page 8 in the manuscript)."}, {"round": 5, "reviewer_comments": {"reviewer_2": "I am satisfied with the changes and thank the authors for sticking with me through several revisions of their manuscript. I do not have any further comments on this manuscript.\n(Remarks on code availability)\nN/A"}, "author_response": ""}]

Decision letter
[Editorial Note: See attachment]