Site-Specific Ground Motion Generative Model for Crustal Earthquakes in Japan Based on Generative Adversarial Networks
Yuma Matsumoto Taro Yaoyama Sangwon Lee Takenori Hida Tatsuya Itoi

Affordance Blending Networks
    Hakan Aktas^1, Yukie Nagai^2, Minoru Asada^3,4, Erhan Oztop^3,5, Emre Ugur^1
^1Hakan Aktas and Emre Ugur are with Computer Engineering Department, Bogazici University
        hakanwhitestone@gmail.com
^2 Yukie Nagai is with IRCN, the University of Tokyo
^3Minoru Asada and Erhan Oztop are with SISREC, Osaka University
^4Minoru Asada is also affiliated with the International Professional University of Technology in Osaka ^5Erhan Oztop is also affiliated with Computer Engineering Department, Ozyegin University.


    
=======================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================




  We develop a site-specific ground-motion model (GMM) for crustal earthquakes in Japan that can directly model the probability distribution of ground motion acceleration time histories based on generative adversarial networks (GANs). The proposed model can generate ground motions conditioned on moment magnitude, rupture distance, and detailed site conditions defined by the average shear-wave velocity in the top 5 m, 10 m, and 20 m (V_S5, V_S10, V_S20) and the depth to shear-wave velocities of 1.0 km/s and 1.4 km/s (Z_1.0, Z_1.4). We construct the neural networks based on styleGAN2 and introduce a novel neural network architecture to generate ground motions considering the effect of source, path, and such detailed site conditions. 5% damped spectral acceleration of ground motions generated by the proposed GMM is consistent with empirical GMMs in terms of magnitude and distance scaling. The proposed GMM can also generate ground motions accounting for the shear-wave velocity profiles of surface soil with different magnitudes and distances, and represent characteristic that are not explained solely by V_S30.




§ INTRODUCTION

Ground-motion models (GMMs) are a critical component in probabilistic seismic hazard analysis (PSHA) <cit.>, <cit.>, serving as significant inputs for earthquake engineering. The GMM can evaluate the median as well as the variability of ground motions at a specific site, as functions of source characteristics, propagation path effects, and site conditions. In Japan, various empirical GMMs have been constructed for crustal earthquakes, subduction zone earthquakes, or both (e.g., <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>).
These GMMs were developed using the abundant strong motion data recorded by networks such as K-NET and KiK-net <cit.>.
Some GMMs have been constructed using a global database, while introducing a region-dependent term.
This approach allows for the development of GMMs with the vast amount of data obtained worldwide, while also taking regional characteristics into consideration.
For instance, in the Next Generation Attenuation-West2 (NGA-West2) project <cit.>, a dataset for crustal earthquakes was compiled <cit.>
and GMMs that can be used in Japan were developed <cit.>, <cit.>, <cit.>, <cit.>.
For subduction zone earthquakes, in NGA-Subduction (NGA-Sub) program <cit.>, GMMs were developed in the same manners <cit.>, <cit.>, <cit.>.

The aforementioned GMMs were constructed for peak ground acceleration (PGA), peak ground velocity (PGV), and spectral acceleration.
Consequently, hazard curves for these ground motion intensity measures are obtained by PSHA.
However, in recent years, it has become common to conduct non-linear dynamic response analysis using acceleration time histories of ground motion as input <cit.>
for detailed risk assessment of structures <cit.>.
One approach to synthesizing ground motion time histories using GMMs for dynamic response analysis involves the use of stochastic ground motion models (SGMMs) (e.g., <cit.>, <cit.>, <cit.>).
In the SGMMs, ground motions are described by a non-stationary stochastic model whose parameters are related to the source, path, and site conditions, and
uncertainties in the model parameters are accounted for to represent the variability of the ground motions under the given conditions.
Another approach involves the selection and scaling of ground motions (e.g., <cit.>, <cit.>, <cit.>).
Ground motions are selected from compiled observed record databases to match the target response spectrum under given conditions.
Although this method is widely used in earthquake engineering, it does not always ensure a sufficient number of records for a specific source, path, and site conditions.

In this study, we utilize a deep generative model to develop a GMM capable of directly modeling the ground motion time history data.
Deep generative models are probabilistic models that employ deep learning techniques.
A key concept of a deep generative model is its ability to capture the inherent distribution of the learned data, and to generate new data that follows this learned distribution.
In other words, instead of simply replicating the learned data, a generative model can generate a set of new data that is statistically similar to the original by capturing its underlying probability distribution <cit.>.
It is generally known that deep generative models, which consist of neural networks with many hidden layers, are capable of learning high-dimensional and complex probability distributions <cit.>.
By applying a deep generative model for ground motion time history data, it is expected to be possible to construct probabilistic models for such high-dimensional data,
which has been difficult with existing empirical GMMs and SGMMs.
In this study, we refer to such a deep generative model-based GMM as a ground motion generative model (GMGM).
The GMGM could become one option for the application of GMMs in earthquake engineering, such as in dynamic response analysis.

Several studies have examined the application of deep generative models for ground motions.
Esfahani et al. <cit.> utilized an autoencoder on the Fourier amplitude spectra (FAS) of ground motions to estimate
the minimum number of predictor variables that is required for a GMM. They also used the trained autoencoder to generate FAS for specific magnitudes and distances.
Among the studies that applied deep generative models for ground motion time histories, generative adversarial networks (GANs) <cit.> have been widely adopted.
Wang et al. <cit.> and Li et al. <cit.> examined data augmentation by applying GANs to the generation of ground motion time histories, targeting applications in earthquake detection problems.
Similarly, Li et al. <cit.> and Wang et al. <cit.> have conducted studies on data augmentation for ground motions using a technique
called conditional GANs (cGAN) <cit.> that can specify the generated data with some condition labels.
Gatti and Clouteau <cit.> proposed a method for generating ground motions up to high-frequency components by combining physics-based simulation methods and GANs. 
Grijalva et al. <cit.> applied GANs to the FAS of ground motions obtained in volcanic events.
Matinfar et al. <cit.> trained GANs on wavelet-transformed ground motions and developed a method to generate ground motions matching a target response spectrum.
Matsumoto et al. <cit.> trained GANs for ground motion time histories, and demonstrated that the trained GANs model
could adequately approximate the distribution of observed record database.
Addressing the characteristics of source, path, and site conditions, similar to GMMs,
Florez et al. <cit.> demonstrated the capability of cGAN in generating ground motion time histories conditioned on magnitude, distance, and
the average shear-wave velocity in the top 30 meters (V_S30). Following a similar approach, Esfahani et al. <cit.> developed a model named TFCGAN,
which learns the time-frequency domain amplitudes of ground motions conditioned on magnitude, distance, and V_S30. 
Additionally, they demonstrated how the time-frequency amplitudes produced by the trained cGAN could be used to retrieve the ground motion time histories using a phase retrieval technique.
Shi et al. <cit.> employed an extension of GANs, known as the generative adversarial neural operator <cit.>,
to construct a model capable of generating ground motion time histories conditioned on magnitude, distance, V_S30, and the style of faulting.

Although various approaches were proposed to apply GANs for ground motion data, most studies have typically described propagation path effects using the source distance and site conditions using V_S30.
However, such modeling is based on the ergodic assumption <cit.>, which may lead to an overestimation of variability in ground motions.
Recently, studies on non-ergodic GMMs have been actively conducted to eliminate this ergodic assumption <cit.>, and in Japan, site-specific GMMs (e.g., <cit.>) and non-ergodic GMMs (e.g., <cit.>) have also been developed.
It is equally important to address the ergodic assumption in GMGM.
While the fully non-ergodic GMMs take into account all of the source, path, and site effects,
we initially focus on the site effects only.
A GMGM that accounts for detailed site conditions could have a potential to eliminate parts of the ergodic assumption.

In this study, we propose a site-specific GMGM (SS-GMGM) for crustal earthquakes in Japan,
considering detailed site conditions and utilizing GANs on ground motion time history data.

The proposed SS-GMGM specifies site conditions using five condition labels: the average shear-wave velocities in the top 5 m, 10 m, 20 m (V_S5, V_S10, V_S20), and the depth to the layer with shear-wave velocities of 1.0 km/s and 1.4 km/s (Z_1.0, Z_1.4).
By combining these five site conditions with the moment magnitude (M_W)
and rupture distance (R_RUP), the SS-GMGM is trained on ground motion time history data conditioned on a seven-dimensional vector: [M_W, R_RUP, V_S5, V_S10, V_S20, Z_1.0, Z_1.4].
We also propose a novel neural network architecture that can generate ground motion time histories with this seven-dimensional vectors.
The quality and distribution of the generated ground motions from the trained SS-GMGM are evaluated, and the performance of the SS-GMGM is demonstrated by comparing it with existing empirical GMMs.
We also discuss how well the specified site conditions correlate with the generated ground motions.

The structure of this paper is as follows.
Section <ref> describes the dataset used for training of the SS-GMGM. 
Section <ref> outlines the GANs method and the proposed neural network architectures.
The results of the proposed SS-GMGM are shown in section <ref>,
and section <ref> presents conclusions and future perspectives based on the detailed site conditions-specified training outcomes.
The program code we used for deep learning is available in a GitHub repository,
which can be accessed at <https://github.com/Mat-main-00/ss_gmgm>.



§ TRAINING DATASETS



 §.§ Data selection and correction

To compile the training dataset, we collected observed records from shallow crustal earthquakes in Japan.
The selection criteria for the earthquakes and records are outlined below.

  
  * Crustal earthquakes that occurred in the Eurasian plate between 1997 and 2016.
  
  * Moment magnitude M_W > 5.
  
  * Hypocentral depth less than 30 km.
  
  * Inclusion of both mainshocks and aftershocks that meet the above criteria.
  
  * Observed records at the K-NET stations.
  
  * The rupture distance R_RUP≤ 100 km.
  
  * Use the two horizontal components of ground motions assuming that they are independent, and ground motions for the rotation angle of 45 degrees are added to augment data.

Moment magnitude M_W was determined using the moment tensor solution from the F-net (Full Range Seismograph Network of Japan) database <cit.>,
and the lower threshold of M_W was set to 5.0 referring to the current practices in PSHA in Japan <cit.>.
Rupture distance R_RUP was calculated as the shortest distance from the rupture area to the station. When M_W is enough small,
and the earthquake can be considered as a point source,
R_RUP was calculated as the hypocenter distance.
The final training dataset consists of 21,696 records from 62 earthquakes.
Each record is a horizontal one-component acceleration time history with 100 Hz sampling.
Details of the selected earthquakes are further described in Matsumoto et al. <cit.>.
The locations of the earthquake epicenters and stations are shown in Figure <ref>,
and the magnitude-distance distribution is shown in Figure <ref>.









 §.§ Parameters for site conditions

The K-NET database provides P-S logging results at one-meter intervals down to a depth of 20 meters from the surface at each station.
Although the SS-GMGM could directly utilize this data as site conditions, we chose to represent surface soil conditions with V_S5, V_S10, and V_S20 considering practical applications.
Therefore, the surface soil is modeled as a three-layer structure in the SS-GMGM.
In many existing GMMs, the surface soil conditions are expressed only by V_S30.
For comparative purposes, we also calculated V_S30 using the following empirical formula proposed by Kanno et al. <cit.>:

    V_S30 = 1.13V_S20 + 19.5

It is important to note that V_S30 is not used as a model parameter in the SS-GMGM but is only used for comparison.
To account for amplification by deep sedimentary layers, Morikawa and Fujiwara <cit.> used Z_1.4 in their GMM.
Abrahamson et al. <cit.>, Boore et al. <cit.>, and Chiou and Youngs <cit.> used Z_1.0 as the parameter to represent that amplification.
Referring to these studies, we use both Z_1.0 and Z_1.4 as the site condition parameters of the SS-GMGM.
We obtained the values of Z_1.0 and Z_1.4 at the each station included in the complied dataset from the deep subsurface structural model V3.2 provided in Japan Seismic Hazard Information Station (J-SHIS) <cit.>.
It should be noted that the values of Z_1.0 used in the above three GMMs refer to the depth from the ground surface,
whereas the values in J-SHIS database represent the depth from the engineering bedrock, thus the condition settings are not strictly same.
Furthermore, Z_1.0 and Z_1.4 are selected primarily for the purpose of comparison with existing GMMs, and may not necessarily be optimal as the parameters representing the characteristics of ground motions in our dataset.
A future task will be to figure out how to appropriately incorporate the characteristics of deep sedimentary layers into the GMGMs.
Through these procedures, we derived a seven-dimensional vector [M_W, R_RUP, V_S5, V_S10, V_S20, Z_1.0, Z_1.4] as a condition label for each observed record.



 §.§ Data pre-processing

In the compiled dataset, each ground motion was aligned with respect to the P-wave arrival time,
which wes manually determined by visually inspecting the waveform.
Oscillations caused by different events, such as aftershocks occurring shortly after the mainshock, were identified and removed based on visual inspection of the waveform.
The data length was then set to 7,992 (79.92 seconds) by truncating the end of each record to ensure that regions of large amplitude were retained for almost all observed records.
When the length of the original record was less than 7,992, zeros were appended to the end to equalize the data length across the dataset.
A cosine taper was applied to the final 100 steps, and 100 zeros were appended to the start and finish of the records in order to lessen edge effects and stabilize the learning process.
The sampling frequency was not changed from 100 Hz, and no band-pass filter was applied.
Consequently, ground motion acceleration time histories with a duration of 81.92 seconds (8,192 samples), starting 1.0 second before the onset of the P-wave,
were obtained.

Training of GANs models can be unstable <cit.>.
To improve the stability of the learning process, normalization techniques are commonly employed within the neural network architecture <cit.>. 
Therefore, in this study, the waveforms and their amplitudes are learned separately to improve model performance.
Each observed record amplitude was normalized by its PGA, and the PGA value was appended to the corresponding condition label to form an eight-dimensional vector.
Then, normalization was also performed on each element of the condition label vectors, ensuring a mean of zero and a standard deviation of 0.1.




§ METHODS



 §.§ GANs

GANs consist of two deep neural networks: a generator G and a discriminator D.
The generator G receives a noise vector z∼𝒩(0, I) as input to generate new data G(z) (referred to as generated data).
The discriminator D takes both observed data x and generated data G(z) as input and estimate the probability that the input is observed data.
GANs are trained through an iterative process where the generator and the discriminator are alternately updated. 
The discriminator is trained to correctly distinguish between the observed and generated data, whereas the generator attempts to produce data that the discriminator will mistakenly identify as observed data.
When appropriately trained through such procedures, the discriminator is known to accurately capture the distribution of the learned dataset, enabling the generator to generate realistic new data that follows this learned distribution <cit.>.



In our previous study <cit.>, we used a GANs model known as Wasserstein GAN with gradient penalty (WGAN-GP) <cit.>. However, in this study,
we propose a model based on styleGAN2 <cit.>, which has achieved higher quality data generation, to construct the SS-GMGM.
A key feature of styleGAN2 is the architecture of its generator.
The generator G of styleGAN2 is composed of two neural networks: the mapping network f and the synthesis network g.
The mapping network takes the noise vector z as input and output a vector w, referred to as the intermediate latent variable.
The synthesis network then takes this intermediate latent variable as input and output generated data g(w) (equivalent to G(z)).
In conventional GANs, the generator learns the correspondence between the noise vector z that follows a normal distribution and the observed data x.
However, observed data x are typically not distributed according to a normal distribution.
By transforming the noise vector z into an intermediate latent variable w, the input of the synthesis network  is not sampled according to any fixed distribution, but its sampling density is induced by a learned mapping f(z). This would make it easier for the synthesis network to learn the relationship between the intermediate latent variables w and observed data x <cit.>.

Following the method of Karras et al. <cit.>, we set the objective function for the generator training as follows:

    min_θ𝔼_z[log(1 - D(G(z; θ)))]
      + 𝔼_w, u∼𝒩(0, I)(J^T_wu_2 - a)^2

where θ is the parameters of the generator, J_w
is the Jacobian matrix J_w = ∂ g(w)/∂w, and a is the constant.
The first term is the logistic loss of the normal GANs <cit.>, and the second term is the regularization term.
Similarly, we set the objective function for the discriminator training as follows:

    min_ψ -𝔼_x[log D(x; ψ)]
      - 𝔼_z[log(1 - D(G(z); ψ))] + γ/2𝔼_x[∇ D(x;ψ)
      ^2]

where ψ is the parameters of the discriminator, and γ is the constant.
The first and second terms are the logistic losses of the normal GANs <cit.>, and the third term is the regularization term called R_1 regularization <cit.>.





 §.§ Proposed neural network architecture of the SS-GMGM

In this section, we introduce the overview of our newly proposed neural network architecture of the SS-GMGM.
For more detailed information on the neural networks architectures and their parameter settings, please refer to our GitHub repository <https://github.com/Mat-main-00/ss_gmgm>.

The overall architecture of the SS-GMGM is shown in Figure <ref>.
The generator takes the noise vector as input and generates a ground motion along with a corresponding condition label.
The discriminator, receiving pairs of ground motion and condition label as input, outputs the probability that the inputs are observed records.
Since the conventional styleGAN2 model cannot handle condition labels,
we made several modifications to the network configurations of Karras et al. <cit.>.
An overview of the proposed neural network architecture of the generator is shown in Figure <ref> (a).
A new label generation network h was added to generate the condition labels using the feature maps in the final layer of the synthesis network as input.
The architectures of the mapping network and synthesis network are almost same as the model configuration by Karras et al. <cit.>,
however, the parameters of the neural networks were modified to fit the data shape of the ground motions.
The mapping network consists of eight fully connected layers, using leaky ReLU <cit.> as the activation function.
Skip connection <cit.> is used to construct the synthesis network, and both ELU <cit.> and leaky ReLU are used as the activation function.
The label generation network consists of eight fully connected layers, with leaky ReLU employed as the activation function.

An overview of neural network architecture of the discriminator is shown in Figure <ref> (b).
The discriminator receives the ground motion data and label data in different neural networks.
The neural network for the ground motion data was constructed in the same way as the generator, modifying the neural network configurations and parameters of Karras et al. <cit.>. The residual network <cit.> is utilized, and leaky ReLU is applied as the activation function.
We introduced a three fully connected layers for label data referring to a configuration called projection discriminator <cit.>.
The outputs of this introduced network are combined with the feature maps of the network for ground motion data by taking inner product.
Through this process, the discriminator is able to comprehend the information from the condition labels in an appropriate manner <cit.>.




The hyperparameters were determined according to Karras et al. <cit.>.
The learning rate was set to 0.002, the batch size was 64, and the dimensions of z and w were both set to 512.
Adam <cit.> was used as the optimization method.
The deep learning and the construction of neural networks were carried out using the Python library PyTorch <cit.>.
For other hyperparameter settings, please refer to our GitHub repository.






§ ANALYSIS RESULTS AND MODEL EVALUATION

The SS-GMGM was trained using the compiled dataset, and the trained model generated 100,000 different ground motions and corresponding condition labels.
Several post-processing steps were performed on the each generated ground motion.
Initially, the offset was removed so that mean acceleration becomes zero. Then a fourth-order Butterworth filter was applied
with 0.1 Hz low-frequency cutoff and 20 Hz high-frequency cutoff.
The accelerations in the first and last two seconds were set to zero, similar to the conditions of the training data.
Finally, the amplitudes of each generated ground motion were recovered by multiplying the values of PGA in the corresponding generated condition labels.
In this section, we evaluate the performance of the SS-GMGM using these generated ground motions and condition labels.



 §.§ Generated ground motion waveforms

We first evaluate the performance of the SS-GMGM by visually checking the generated data.
Figures <ref> and <ref> show the ground motion waveforms of observed records and generated ground motions for different 
magnitude, distance, and V_S30 scenarios, respectively.
V_S30 for each generated ground motion was calculated from the generated V_S20 using equation <ref>.
The SS-GMGM appropriately captures waveform characteristics such as the onset of P-waves and S-waves, as well as the envelope shapes.
The amplitude scaling with magnitude and distance is generally well represented, and the relationships between distance and duration, as well as between distance and P-S time, are also appropriate.
Moreover, focusing on the values of V_S30, lower frequency ground motions are generated in conditions of soft soil where V_S30 < 300 m/s, compared to those of stiff soil where V_S30 > 300 m/s.
In the data of row 3 and column 2 in Figures <ref> and <ref>,
a tendency for subsequent vibrations to have lower frequency components in soft soil condition is also captured.
The more detailed examinations of the correspondence between ground motions and site conditions are conducted in subsection <ref>.



 §.§ Overall characteristics of generated ground motions


A key characteristic of the SS-GMGM is its ability to directly learn the probability distribution of ground motion time histories.
Therefore, we examine the distribution of the generated ground motions by comparing it with that of training dataset.
The evaluation is conducted by examining the following five
indices used in Rezaeian and Der Kiureghian <cit.>, <cit.> to represent the characteristic of ground motions:

  
  * Arias intensity, I_A.
  
  * Significant/Relative duration of ground motion, D_5-95.
  
  * Significant/Relative duration of ground motion, D_5-45.
  
  * Zero-level crossing rate, ν.
  
  * Mean of the rates of negative maxima and positive minima, η.

Arias intensity <cit.> is a measure of total energy contained in the ground motion and is defined as:

    I_A = π/2g∫_0^t_da^2(t)dt

where g is the gravitational acceleration, t_d is the total duration of ground motion, and a(t) is the acceleration at time t.
We set g = 980.665 cm/s^2.
Significant duration <cit.> is widely used as an index for assessing the ground motion duration <cit.>.
The value of D_5-95 is defined as the time interval required for the cumulative power of the ground motion to reach from 5% to 95% of the total cumulative power,
and generally corresponds to the duration of strong motion <cit.>.
D_5-45 is also defined as the time interval from 5% to 45%, and corresponds to the time at the middle of the strong-shaking phase <cit.>.
The Zero-level crossing rate ν is adopted to characterize the dominant frequency in ground motion <cit.>,
and we calculated the value of ν as the mean of zero-level up-crossing rate and zero-level down-crossing rate. 
Negative maxima and positive minima are used to characterize the bandwidth of ground motions <cit.>.
It is known that a ground motion with wider bandwidth tends to have larger η. We defined η as the mean of the rates of negative maxima and positive minima.





We calculated the values of these five indices for both observed records and generated ground motions.
Note that the values of ν and η were calculated for the vibrations in the time interval corresponding to D_5-95,
and observed records are filtered with the same Butterworth filter used for post-processing of generated ground motions.
Figure <ref> shows the comparison results of cumulative distribution function (CDF) of each index.
Although the generated ground motions tend to have a little smaller D_5-45 values,
the distributions of I_A and D_5-95 are consistent with that of observed records.
This means that the SS-GMGM captures the temporal characteristic of observed records with reasonable accuracy.
The distribution of ν also shows a good agreement.
The distribution of η indicate that ground motions in the tails of the distribution of observed records are not generated extensively,
however, the overall distribution of generated ground motions is generally matched with that of observed records.
The SS-GMGM captures the overall characteristics of the observed records in terms of both temporal and frequency characteristics,
and is an appropriate probabilistic model of ground motion time histories that
approximates the distribution of learned dataset. 






 §.§ Characteristics of near-field ground motions

Near-field ground motions containing strong velocity pulse have caused destructive damage to structures (e.g., <cit.>, <cit.>),
and are one of the critical ground motions that should be considered in earthquake engineering.
Following the numerical integration of the generated ground motions, the method of Baker <cit.> was applied to classify and extract the velocity pulses.
In the numerical integration of the generated data, a fourth-order Butterworth filter (cutoff frequency 0.2Hz) was used to remove the low frequency components to prevent the noise from them.
It should be noted that while the method of Baker <cit.> was designed to be applied to the fault normal components of ground motions,
the generated ground motions by the SS-GMGM do not contain information on the positional relationship with the fault plane. Therefore, it was applied directly to the generated ground motions.


First, following the evaluation in our previous study <cit.>,
we compare the distribution of following pulse indicator values proposed by Baker <cit.>:

    Pulse indicator = 1/1 + exp(-23.3 + 14.6v_r + 20.5E_r)

where v_r is the PGV of the residual data, which is obtained by subtracting the extracted pulse from original velocity data,
and E_r is the total accumulate power of residual data divided by that of original data.
Pulse-like ground motion tends to have high pulse indicator value.
Figure <ref> shows the CDF of the pulse indicator of observed records and generated ground motions.
Although the proportion of pulse-like ground motion in the generated data is slightly smaller than in the observed records dataset,
the overall trend of the distribution is generally consistent.
Figure <ref> shows the examples of the original velocity waveforms which were classified as pulse-like ground motion and corresponding extracted pulse waveforms.
A clear velocity pulse is found in generated ground motion with the near-filed setting.
The proposed SS-GMGM is capable of generating ground motions with engineering-significant characteristic.






 §.§ Statistical evaluation of the FAS

In this subsection, we conduct statistical performance evaluation of the SS-GMGM.
Following the evaluation method of Florez et al. <cit.> and Esfahani et al. <cit.>,
we compare the distributions of observed and generated FAS with different M_W, R_RUP, and V_S30 scenarios.
Figure <ref> compares the logarithmic means and logarithmic standard deviations (std) of the FAS of both observed records and generated ground motions.
The bins of the condition labels are determined so that enough number of observed records are contained.
The generated FAS are generally consistent with observed FAS in frequency range [1, 20] Hz.
In comparisons under relatively small M_W conditions, the generated ground motions contain many low frequency components regardless of R_RUP or V_S30,
resulting in an overestimation particularly in the region below 0.5 Hz.
On the other hand, in the region with larger magnitudes range 6.7 ≤ M_W ≤ 6.9 (Figure <ref> (a)),
the generated data agree with the observed records in wider frequency range [0.1, 20] Hz.
This trend may be attributed to the influence of the generative processes of the SS-GMGM.
The constraint that the acceleration returns to zero under normal conditions is not included in the SS-GMGM. As a result, regions with small amplitudes periodically experience slight deviations from the zero line, which could result in low frequency noise.
As small magnitude earthquakes do not cause oscillations with low frequency components, the effect of low frequency noise is notably evident in the range of 5.5 ≤ M_W ≤ 5.7 (Figure <ref> (c)).
The frequency bands for which the distribution of the generated FAS agrees with the FAS of observed records are generally comparable to those examined in Florez et al. <cit.> and Esfahani et al.<cit.>.





 §.§ Comparison with empirical GMMs

The generated results of the SS-GMGM are compared with the prediction of the following two empirical GMMs:

  
  * Morikawa and Fujiwara <cit.> (MF13) GMM 
  
  * Abrahamson et al. <cit.> (ASK14) GMM 


Morikawa and Fujiwara <cit.> proposed two models which differ in amplitude scaling with reference to magnitudes.
We utilize a model with a quadratic magnitude term and perform predictions using the following formula:

    log_10S_a    = a(M_W^' - M_W_1)^2 + bR_RUP + c - log_10(R_RUP + d· 10^eM_W^') + G_d + G_s 
    
      M_W^'   = min(M_W, M_W_0) 
    
      G_z    = p_zlog_10(max(Z_min, Z_1.4)/Z_0) 
    
      G_s    = p_slog_10(min(V_Smax, V_S30)/V_0)

where S_a is the spectral acceleration value for specific period, a, b, c, d, e, M_W_0, M_W_1, p_z, p_s, Z_min, and V_Smax are coefficients.
G_z is the correction term for amplification by deep sedimentary layers, and G_s is the correction term for amplification by shallow soft soils.

The ASK14 GMM contains some explanatory variables to specify the source conditions.
Since the SS-GMGM cannot take into account such detailed source conditions,
the following equation is used for the prediction, considering only terms related to site conditions and the region-specific term:

    log S_a    = f_1(M_W, R_RUP) + f_5(Sa_1180, V_S30) + f_10(Z_1.0, V_S30) + Regional(V_S30, R_RUP)

where Sa_1180 is the median spectral acceleration on hard rock.
For detailed formulation of f_1, f_5, f_10, and Regional(·), please refer to Abrahamson et al. <cit.>.






  §.§.§ Distribution

Figure <ref> compares the spectral acceleration values (5% damping) at periods T = 0.2 and 1.0 seconds between the generated ground motions and prediction results from the MF13 GMM and ASK14 GMM. The distance scaling of the SS-GMGM is generally consistent with the MF13 GMM,
and the variability of generated data is also agree with the prediction of the MF13 GMM.
Compared to the median spectral acceleration values of the ASK14 GMM at T = 0.2 seconds, the values from the SS-GMGM tend to be slightly higher, especially under soft soil conditions. 
However, the distributions of the spectral acceleration values at T = 1.0 second are generally in agreement.






  §.§.§ Residual analysis

The median values of the MF13 GMM and ASK14 GMM corresponding to the each generated condition label
are predicted, and the residual values (log_10[gen/pre], gen is the generated data and pre is the median of the GMMs) are calculated.
To align the condition of the dataset, only generate data with M_W ≥ 5.5 is used for the MF13 GMM predictions.
Figure <ref> shows the residuals with reference to M_W, R_RUP, V_S30, and Z_1.4 compared with the MF13 GMM.
The residuals are generally centered around zero regardless of M_W, V_S30, and Z_1.4 values,
indicating that the SS-GMGM has appropriately learned the magnitude scaling as well as the amplification characteristics due to shallow soil and deep sedimentary layers.
In terms of distance scaling, there is a tendency to slightly overestimate in regions where R_RUP is large,
and to slightly underestimate in regions with short distances at T = 1.0 seconds.
However, across a wide range, the residuals are generally centered around zero, showing that the distance scaling of the SS-GMGM is consistent with reasonable accuracy.
Figure <ref> shows the residual plots of the SS-GMGM against the ASK14 GMM prediction for each value of M_W, R_RUP, V_S30, and Z_1.0.
The generated spectral accelerations tend to be slightly larger at T = 0.2 seconds, as also shown in Figure <ref>.
However, the mean values of the residuals remain relatively constant regardless of the values of M_W, R_RUP, V_S30, and Z_1.0.
This confirms that the SS-GMGM has adequately learned the scaling of magnitude and distance, as well as of the amplifications by the shallow soil and deep sedimentary layers when compared with the ASK14 GMM.



 §.§ Evaluation as a site-specific model

The performance of the SS-GMGM is assessed in this subsection when site conditions are given in detail by looking at the relationship between generated ground motions and condition labels.
First, the shear-wave velocity profile modeled in three layers for the top 20 meters is back-calculated from the generated values of V_S5, V_S10, and V_S20
as follows:

    v_0-5   = V_S5
    
      v_5-10   = V_S5V_S10/2V_S5 - V_S10
    
      v_10-20   = 2v_5-10V_S5V_S20/4v_5-10V_S5 - v_5-10V_S20 - V_S5V_S20

where v_0-5, v_5-10, and v_10-20 are the average shear-wave velocity of the first-layer (0-5 m), second-layer (5-10 m), and third-layer (10-20 m) of the surface soil, respectively.
Subsequently, generated data with nearly identical values of M_W, R_RUP, and V_S30, which are used as condition labels of some GMGM (<cit.>, <cit.>), are sampled.
Figure <ref> presents the time-history waveforms, FAS, and corresponding velocity profiles for three generated samples,
and Table <ref> represents the corresponding condition label values.
Note that the FAS was smoothed using a Parzen window with a bandwidth of 0.2 Hz.
Even when the values of M_W, R_RUP, and V_S30 are nearly the same, it is evident that the characteristics of the ground motions differ significantly.
In Figure <ref>, the generated ground motion in (i) has large amplitude, and has a peak near 8 Hz in FAS,
whereas the generated ground motion in (ii) has relatively smaller amplitude with a peak near 5 Hz in FAS.
Examining the corresponding shear-wave velocity profiles, the data for case (i) reveals a velocity profile where the first-layer has a shear-wave velocity of v_0-5 = 134 m/s,
and the second-layer has a shear-wave velocity of v_0-5 = 512 m/s.
This indicates a soft upper layer of 5 meters, underlain by a harder layer beyond 5 meters. Similarly, in case (ii),
the top 5 meters also consists of a soft layer, but v_0-5 is slightly larger than in (i), and the difference in shear-wave velocities between the first and second layer is smaller.
The high frequency vibrations are amplified due to multiple reflections at the soft surface soil,
suggesting that the generated data in (i) shows larger amplitude and higher predominant frequency compared to (ii) due to this amplification.





Furthermore, the data in (iii), while having similar amplitude and shear-wave velocity profile compared to (ii), contain many low frequency components.
The Z_1.0 and Z_1.4 values for (iii) are relatively high,
indicating that the generated ground motions are also consistent with the effect of deep sedimentary layers.
The SS-GMGM is capable of representing characteristic of ground motions that cannot be explained solely by V_S30 value.











 §.§ Distribution of generated ground motions with specified magnitude, distance, and detailed site conditions

We examine the distribution of generated ground motions when magnitude, distance, and detailed site conditions are specified.
Initially, generated data fitting the following criteria were sampled, and 131 data are obtained:

  
  * 5.9 ≤ M_W ≤ 6.1
  
  * 25 km≤ R_RUP≤ 35 km
  
  * 404 m/s≤ V_S30≤ 426 m/s

The black dashed line with circles in Figure <ref> (a) shows the mean and standard deviation of the acceleration spectra for all 131 sampled generated ground motions.
Among the 131 samples, we further extracted two clusters that share similar shear-wave velocity profiles (Cases 1 and 2).
The red dashed line with squares (Case 1) and blue dotted line with triangles (Case 2) in Figure <ref> (a)
shows the mean acceleration spectra of the extracted clusters, and their corresponding shear-wave velocity profiles of each sampled data
are shown in Figure <ref> (b).
Note that the standard deviation of Case 1 and Case 2 are calculated together to ensure a sufficient sample size as follows:

    σ = √((n_1 - 1)σ̂_1^2 + (n_2 - 1)σ̂_2^2/n_1 + n_2 - 2)

where n_1 and n_2 are the number of data in Case 1 and 2, respectively, and σ̂_1^2 and σ̂_2^2 are the unbiased variance of data in Case 1 and 2, respectively.
The mean acceleration spectra varies depending on the soil profile.
For example, Case 1 has relatively soft soil condition of top 5 meters, and the spectral acceleration of Case 1 at short period range had larger values.
The standard deviations at each period have decreased by specifying the shear-wave velocity profile,
and the averaged standard deviation for each period of Case 1 and Case 2 results in about 0.35.
According to the study by Morikawa et al. <cit.>, the standard deviation (total of the within-event and between-events) of the 5% damped spectral acceleration at specific station
generally ranges between 0.35 and 0.45.
Hikita and Tomozawa <cit.> studied the variability of spectral acceleration of single-path ground motions,
and the standard deviations approximately ranged between 0.3 and 0.5 in their examination.
Although a strict comparison is challenging due to different settings, the standard deviation obtained in this study is considered to be generally reasonable.



 §.§ Data interpolation by SS-GMGM

GANs are known to generate new data samples by interpolating the distribution of a learned dataset <cit.>.
In this subsection, we evaluate the data interpolation of the SS-GMGM for generated ground motions and condition labels.

First, Figure <ref> compares the distribution of M_W and V_S30 between the observed records dataset and generated data.
While no generated data correspond to the slightly included area of V_S30 > 1000 m/s in the training dataset,
new data have been generated within the region where observed records are distributed.
For instance, there are no observed records corresponding to M_W = 6.5 in the dataset, 
yet the SS-GMGM also generates data for such regions.
Looking at the residual plots in Figure <ref>, the generated results for M_W = 6.5 are consistent with the predictions of the MF13 GMM, indicating that the generated ground motions are reasonable.
The SS-GMGM can generate data for combinations of condition labels that are not included in the observed record database.

Next, we examine the case specifying M_W, R_RUP, and the shear-wave velocity profile.
Here, we use the generated results of Case 2 shown in Figure <ref> as an example.
Table <ref> lists the values of the condition labels of six generated ground motions in Case 2.
For these generated data, the values for Z_1.0 and Z_1.4, which were not specified during data sampling,
are almost identical.
Similar trends were occasionally observed when different combinations of the M_W, R_RUP, and shear-wave velocity profile were specified.
Table <ref> also shows the site conditions of YMG007 station, a K-NET station included in our training dataset.
It is important to note that since only shear-wave velocity profile up to 10 meters has been obtained at YMG007 station in K-NET database,
it is assumed that the shear-wave velocity from 10 meters to 20 meters is identical to the upper layer.
The shear-wave velocity profile and values of Z_1.0 and Z_1.4 at YMG007 station are very similar to those of the generated data in Case 2.
This suggests that while the SS-GMGM can interpolate data when only a few condition labels are considered, as shown in Figure <ref>,
such interpolation becomes more challenging as the number of considered condition labels increases.
This will be discussed further in section <ref>.







Finally, we examine the case where only site conditions were specified.
At YMG007 station, two observed records from different earthquakes are obtained in our training dataset.
The relationships between these earthquakes and the YMG007 station are shown in
Figure <ref>.
Then, generated data having the same site conditions as YMG007 station are sampled, corresponding to the following criteria:


  * 225 m/s≤ v_0-5≤ 275 m/s

  * 355 m/s≤ v_5-10≤ 405 m/s

  * 395 m/s≤ v_10-20≤ 445 m/s

  * 3 m≤ Z_1.0≤ 5 m

  * 9 m≤ Z_1.4≤ 11 m

Figure <ref> shows the magnitude-distance distribution of the two observed records at YMG007 station and the sampled generated data.
Although no data have been generated for region where R_RUP is approximately 100 km,
it is evident that data corresponding to combinations of magnitude and distance not included in the observed records have been generated even when specifying detailed site conditions.
Therefore, even when the number of considered condition labels is increased, the SS-GMGM is capable of generating data not included in the observed records dataset, provided that the number of condition labels used for interpolation remains small.
For the data for 42 km≤ R_RUP≤ 44 km in Figure <ref>,
time history waveforms and their corresponding acceleration spectra and shear-wave velocity profiles are shown in Figure <ref>.
The frequency characteristics of generated ground motions are generally similar to the observed record,
and the magnitude scaling of the amplitude is generally appropriate.







§ CONCLUSION AND DISCUSSION

In this study, we developed a GMM (SS-GMGM) for crustal earthquakes in Japan that can directly generate ground motion time histories with specifying detailed site conditions.
The proposed SS-GMGM were developed based on a type of deep generative model, called styleGAN2, and a novel neural network architecture that could serve as
a generative model with detailed condition labels were also proposed.
The neural network architecture of the SS-GMGM is capable of account for site conditions with five values, [V_S5, V_S10, V_S20, Z_1.0, Z_1.4], in addition to magnitude and distance information.
We demonstrated that the characteristics of the generated ground motions are consistent with these condition labels,
and ground motions with engineering-significant characteristic can be generated.
Furthermore, the amplifications by shallow soil and deep sedimentary layer were shown to be accurately represented, and the SSGMGM's magnitude and distance scaling were shown to match those of empirical GMMs.
Additionally, by modeling the surface soil in three layers, the SS-GMGM could express differences in ground motion characteristics that could not be explained solely by V_S30,
and generated new data samples with different magnitude-distance settings even when detailed site conditions were specified.




The SS-GMGM constructed in this study does not always produce high-quality ground motions, in some cases, significant noises occur at the end of the data.
Figure <ref> shows an example of such generated result.
This type of data tends to be generated, although infrequently, under combinations of condition labels where there are fewer data points in training dataset.
The cause of this tendency might be attributed to the influence of the vibrations included in the observed records.
Subsequent vibrations believed to be caused by aftershocks have been eliminated, as stated in the data preprocessing subsection. On the other hand,
aftershocks left some vibrations that were difficult to pinpoint.


Therefore, the observed records in our training dataset might include vibrations like those in the time interval from 40 to 60 seconds in Figure <ref>,
and it is conceivable that the SS-GMGM has learned such tendencies.
The vibrations seen after about 75 seconds in Figure <ref> are not mere noise but resemble a ground motion.
This could be due to the characteristics of convolutional neural networks (CNNs).
The parts of neural networks of the SS-GMGM have constructed by CNNs,
which are known to have a property called equivariance to translation due to parameter sharing <cit.>.
It is important to improve the performance of the SS-GMGM by carefully examining the data preprocessing procedures and the
configuration of the neural networks suitable for generating ground motion data.

GANs generate new data samples by interpolating the distribution of learned dataset.
Although interpolation could be performed when a few condition labels are specified,
the site conditions specified by the five-dimensional vector [V_S5, V_S10, V_S20, Z_1.0, Z_1.4]
were not well-interpolated.
This issue might stem from the curse of dimensionality.
In the compiled dataset, there are 728 observation stations.
The distribution of soil conditions defined by five-dimensional vectors is relatively high-dimensional considering this number of stations.
Data points may not be sufficiently dense for effective interpolation when many conditions are specified.
This challenge is one of the difficulties in data-driven approaches.
As more condition labels are considered in detail for the source, path, and site conditions,
the lack of data becomes a more significant issue.
Potential solutions to this problem about site conditions include increasing the number of observation stations used in compiled dataset.
The performance of the SS-GMGM may be improved by including observed records at KiK-net stations, as they were not used in this investigation. Using simulation-based techniques to estimate ground motion amplifications and predict ground motions at the bedrock using GMGMs could be another way to solve this problem.
This approach would require detailed site condition information, but would allow the avoidance of issues mentioned above in constructing the SS-GMGM.


Previous studies applying GANs to GMM (<cit.>, <cit.>, <cit.>) have constructed models based on cGAN,
specifying condition labels for generating ground motions.
In contrast, this study generates ground motions and condition labels from normal random inputs without specifying such labels.
The advantages and disadvantages of our method compared to cGAN based method are as follows. 
An advantage is its suitability for generating ground motions when detailed condition labels are specified.
As mentioned earlier, GANs approximate the distribution of the learned dataset, making them unsuitable for generating data that would be a complete extrapolation from the dataset range.
For example, in the dataset used in this study, the moment magnitudes of the records range from 5.0 to 7.1, making it difficult for the SS-GMGM to generate ground motions for magnitudes like 4 or 8.
The distribution of condition labels within the seven-dimensional space is sparse, considering the number of available observed records.
Since cGAN requires the specification of condition labels to generate data, such sparsity could make it challenging to appropriately set these labels without exceeding the applicable range.
Consequently, models that specify more detailed condition labels may generate ground motions that ideally should not be generated.
Such problems can be avoided by generating condition labels with ground motions at the same time.
On the other hand, a disadvantage of our method is the time required to generate ground motions that matches specific condition labels.
The frequency of generated data depends on its occurrence in the training dataset,
which means data for less common conditions is less likely to be generated.
To utilize the SS-GMGM in the filed of earthquake engineering, it is necessary to construct a framework that incorporates rare-event simulation techniques <cit.>, which are studied in fields such as structural reliability <cit.>.

Additionally, while this study constructed GMGMs for the horizontal one component of ground motions occurred in crustal earthquakes,
it is also important from the viewpoint of earthquake engineering to predict ground motions for subduction zone earthquakes and
to simultaneously predict three components.
Expanding the application range of the SS-GMGM is essential and is a future task.



§ DATA AND RESOURCES

  The strong motion records and the shear-wave velocity values of the surface soil
  used in this study can be downloaded through the website of the National Research Institute for Earth Science and Disaster Resilience
  (NIED; <https://www.kyoshin.bosai.go.jp/kyoshin/>). We only used the records of K-NET.
  The moment magnitude was obtained from the NIED F-net database <https://www.fnet.bosai.go.jp/>,
  and the values of Z_1.0 and Z_1.4 were obtained from the NIED J-SHIS website <https://www.j-shis.bosai.go.jp/>.
  The program code used in deep learning is available in the GitHub repository <https://github.com/Mat-main-00/ss_gmgm>.



§ DECLARATION OF COMPETING INTERESTS

The authors acknowledge that there are no conflicts of interest recorded.



§ ACKNOWLEDGMENTS

This study was supported by JSPS KAKENHI Grant Numbers JP23H00219 and JP22J23006.

unsrt


10

Cornel1971
Carl Allin Cornell.
Probabilistic analysis of damage to structures under seismic loads.
In D. A. Howells, P. Haigh, and C. Taylor, editors, Dynamic Waves in Civil Engineering, pages 473–488, 1971. Wiley, New York.

Baker2021seismic
Jack Baker, Brendon Bradley, and Peter Stafford.
Seismic Hazard and Risk Analysis.
Cambridge University Press, 2021.

FukushimaNew1990
Yoshimitsu Fukushima and Teiji Tanaka.
A new attenuation relation for peak horizontal acceleration of strong earthquake ground motion in Japan.
Bull. Seismol. Soc. Am., 80(4):757–783, 1990.

Molas1995
Gilbert L. Molas and Fumio Yamazaki.
Attenuation of earthquake ground motion in Japan including deep focus events.
Bull. Seismol. Soc. Am., 85(5):1343–1358, 1995.

Si1999
Hongjun Si and Saburoh Midorikawa.
New attenuation relationships for peak ground acceleration and velocity considering effects of fault type and site condition.
J. Struct. Constr. Eng., 523:63–70, 1999.
(in Japanese with English abstract).

Nishimura2003
Toshimitsu Nishimura and Masanori Horike.
The attenuation relationships of peak ground acceleration for the horizontal and the vertical components inferred from the Kyoshin network data.
J. Struct. Constr. Eng., 571:63–70, 2003.
(in Japanese with English abstract).

Kanno2006
Tatsuo Kanno, Akira Narita, Nobuyuki Morikawa, Hiroyuki Fujiwara, and Yoshimitsu Fukushima.
A new attenuation relation for strong ground motion in Japan based on recorded data.
Bull. Seismol. Soc. Am., 96(3):879–897, 2006.

Zhao2006Attenuation
John X. Zhao, Jian Zhang, Akihiro Asano, Yuki Ohno, Taishi Oouchi, Toshimasa Takahashi, Hiroshi Ogawa, Kojiro Irikura, Hong K. Thio, Paul G. Somerville, Yasuhiro Fukushima, and Yoshimitsu Fukushima.
Attenuation relations of strong ground motion in Japan using site classification based on predominant period.
Bull. Seismol. Soc. Am., 96(3):898–913, 2006.

Morikawa2013
Nobuyuki Morikawa and Hiroyuki Fujiwara.
A new ground motion prediction equation for Japan applicable up to M9 mega-earthquake.
J. Disaster Res., 8(5):878–888, 2013.

Ghofrani2014
Hadi Ghofrani and Gail M. Atkinson.
Ground-motion prediction equations for interface earthquakes of M7 to M9 based on empirical data from Japan.
Bull. Earthq. Eng., 12(2):549–571, 2014.

Zhao2016cr
John X. Zhao, Shuanglin Zhou, Jun Zhou, Chen Zhao, Heng Zhang, Yingbin Zhang, Pingjun Gao, Xiaowen Lan, David Rhoades, Yoshimitsu Fukushima, Paul G. Somerville, and Kojiro Irikura.
Ground‐motion prediction equations for shallow crustal and upper‐mantle earthquakes in Japan using site class and simple geometric attenuation functions.
Bull. Seismol. Soc. Am., 106(4):1552–1569, 2016.

Si2022Development
Hongjun Si, Saburoh Midorikawa, and Tadahiro Kishida.
Development of nga-sub ground-motion prediction equation of 5%-damped pseudo-spectral acceleration based on database of subduction earthquakes in Japan.
Earthq. Spectra, 38(4):2682–2706, 2022.

KNET
National Research Institute for Earth Science and Disaster Resilience.
NIED K-NET, KiK-net.
National Research Institute for Earth Science and Disaster Resilience, doi:10.17598/NIED.0004, 2019.

Bozorgnia2014NGA
Yousef Bozorgnia, Norman A. Abrahamson, Linda Al Atik, Timothy D. Ancheta, Gail M. Atkinson, Jack W. Baker, Annemarie Baltay, David M. Boore, Kenneth W. Campbell, Brian S.-J. Chiou, Robert Darragh, Steve Day, Jennifer Donahue, Robert W. Graves, Nick Gregor, Thomas Hanks, I. M. Idriss, Ronnie Kamai, Tadahiro Kishida, Albert Kottke, Stephen A. Mahin, Sanaz Rezaeian, Badie Rowshandel, Emel Seyhan, Shrey Shahi, Tom Shantz, Walter Silva, Paul Spudich, Jonathan P. Stewart, Jennie Watson-Lamprey, Kathryn Wooddell, and Robert Youngs.
NGA-West2 research project.
Earthq. Spectra, 30(3):973–987, 2014.

Ancheta2014nga
Timothy D. Ancheta, Robert B. Darragh, Jonathan P. Stewart, Emel Seyhan, Walter J. Silva, Brian S.-J. Chiou, Katie E. Wooddell, Robert W. Graves, Albert R. Kottke, David M. Boore, Tadahiro Kishida, and Jennifer L. Donahue.
NGA-West2 database.
Earthq. Spectra, 30(3):989–1005, 2014.

abrahamson2014summary
Norman A. Abrahamson, Walter J. Silva, and Ronnie Kamai.
Summary of the ASK14 ground motion relation for active crustal regions.
Earthq. Spectra, 30(3):1025–1055, 2014.

boore2014nga
David M. Boore, Jonathan P. Stewart, Emel Seyhan, and Gail M. Atkinson.
NGA-West2 equations for predicting PGA, PGV, and 5% damped PSA for shallow crustal earthquakes.
Earthq. Spectra, 30(3):1057–1085, 2014.

campbell2014nga
Kenneth W. Campbell and Yousef Bozorgnia.
NGA-West2 ground motion model for the average horizontal components of PGA, PGV, and 5% damped linear acceleration response spectra.
Earthq. Spectra, 30(3):1087–1115, 2014.

chiou2014update
Brian S.-J. Chiou and Robert R. Youngs.
Update of the Chiou and Youngs NGA model for the average horizontal component of peak ground motion and response spectra.
Earthq. Spectra, 30(3):1117–1153, 2014.

Bozorgnia2022NGA
Yousef Bozorgnia, Norman A. Abrahamson, Sean K Ahdi, Timothy D Ancheta, Linda Al Atik, Ralph J. Archuleta, Gail M. Atkinson, David M. Boore, Kenneth W. Campbell, Brian S-J Chiou, Victor Contreras, Robert B Darragh, Sahar Derakhshan, Jennifer L Donahue, Nick Gregor, Zeynep Gulerce, IM Idriss, Chen Ji, Tadahiro Kishida, Albert R. Kottke, Nicolas Kuehn, DongYoup Kwak, Annie O-L Kwok, P Lin, Jorge Macedo, Silvia Mazzoni, Saburoh Midorikawa, Sifat Muin, Grace A. Parker, Sanaz Rezaeian, Hongjun Si, Walter J. Silva, Jonathan P. Stewart, Melanie Walling, Katie Wooddell, and Robert R. Youngs.
NGA-Subduction research program.
Earthq. Spectra, 38(2):783–798, 2022.

abrahamson2022summary
Norman A. Abrahamson and Zeynep Gulerce.
Summary of the Abrahamson and Gulerce NGA-Sub ground-motion model for subduction earthquakes.
Earthq. Spectra, 38(4):2638–2681, 2022.

Parker2022NGA
Grace A. Parker, Jonathan P. Stewart, David M. Boore, Gail M. Atkinson, and Behzad Hassani.
NGA-subduction global ground motion models with regional adjustment factors.
Earthq. Spectra, 38(1):456–493, 2022.

Kuehn2023regionalized
Nicolas M. Kuehn, Yousef Bozorgnia, Kenneth W. Campbell, and Nicholas Gregor.
A regionalized partially nonergodic ground-motion model for subduction earthquakes using the NGA-Sub database.
Earthq. Spectra, 39(3):1625–1657, 2023.

Vamvatsikos2002
Dimitrios Vamvatsikos and Carl Allin Cornell.
Incremental dynamic analysis.
Earthq. Eng. Struct. Dynam., 31(3):491–514, 2002.

FEMAp58
Federal Emergency Management Agency (FEMA).
FEMA P-58-1: Seismic Performance Assessment of Buildings, Volume 1 - Methodology, Second Edition.
Federal Emergency Management Agency, Washington, DC.,
2018.

Kameda1994probabilistic
Hiroyuki Kameda.
Probabilistic seismic hazard and stochastic ground motions.
Eng. Struct., 16(7):547–557, 1994.

Rezaeian2008
Sanaz Rezaeian and Armen Der Kiureghian.
A stochastic ground motion model with separable temporal and spectral nonstationarities.
Earthq. Eng. Struct. Dynam., 37(13):1565–1584, 2008.

Rezaeian2010
Sanaz Rezaeian and Armen Der Kiureghian.
Simulation of synthetic ground motions for specified earthquake and site characteristics.
Earthq. Eng. Struct. Dynam., 39(10):1155–1180, 2010.

Naeim2004
Farzad Naeim, Arzhang Alimoradi, and Shahram Pezeshk.
Selection and scaling of ground motion time histories for structural design using genetic algorithms.
Earthq. Spectra, 20(2):413–426, 2004.

Bommer2004Use
Julian J. Bommer and Ana Beatriz Acevedo.
The use of real earthquake accelerograms as input to dynamic analysis.
J. Earthq. Eng., 8:43–91, 2004.

Baker2006
Jack W. Baker and Carl Allin Cornell.
Spectral shape, epsilon and record selection.
Earthq. Eng. Struct. Dynam., 35(9):1077–1095, 2006.

wang2024pivoting
Yuan Wang, Lokesh Kumar Sambasivan, Mingang Fu, and Prakhar Mehrotra.
Pivoting retail supply chain with deep generative techniques: Taxonomy, survey and insights, 2024.
available at https://arxiv.org/abs/2403.00861 (last accessed April 2024).

Ruthotto2021
Lars Ruthotto and Eldad Haber.
An introduction to deep generative modeling.
GAMM-Mitteilungen, 44(2):e202100008, 2021.

Esfahani2021Exploring
Reza Dokht Dolatabadi Esfahani, Kristin Vogel, Fabrice Cotton, Matthias Ohrnberger, Frank Scherbaum, and Marius Kriegerowski.
Exploring the dimensionality of ground‐motion data by applying autoencoder techniques.
Bull. Seismol. Soc. Am., 111(3):1563–1576, 2021.

Goodfellow2014
Ian J. Goodfellow, Jean Pouget-Abadie, Mehdi Mirza, Bing Xu, David Warde-Farley, Sherjil Ozair, Aaron Courville, and Yoshua Bengio.
Generative adversarial networks, 2014.
available at https://arxiv.org/abs/1406.2661 (last accessed April 2024).

Wang2019segam
Tiantong Wang, Zhongping Zhang, and Youzuo Li.
EarthquakeGen: Earthquake generator using generative adversarial networks.
SEG Int. Expo. Annu. Meet., page D043S119R004, 2019.

li2020seismic2
Yuanming Li, Bonhwa Ku, Gwantae Kim, Jae-Kwang Ahn, and Hanseok Ko.
Seismic signal synthesis by generative adversarial network with gated convolutional neural network structure.
In IGARSS 2020, pages 3857–3860, 2020.

li2020seismic1
Yuanming Li, Bonhwa Ku, Shou Zhang, Jae-Kwang Ahn, and Hanseok Ko.
Seismic data augmentation based on conditional generative adversarial networks.
Sensors, 20(23), 2020.

wang2021seismogen
Tiantong Wang, Daniel Trugman, and Youzuo Lin.
Seismogen: Seismic waveform synthesis using gan with application to seismic data augmentation.
J. Geophys. Res. Solid Earth, 126(4):e2020JB020077, 2021.

mirza2014conditional
Mehdi Mirza and Simon Osindero.
Conditional generative adversarial nets, 2014.
available at https://arxiv.org/abs/1411.1784 (last accessed April 2024).

gatti2020towards
Filippo Gatti and Didier Clouteau.
Towards blending physics-based numerical simulations and seismic databases using generative adversarial network.
Comput. Meth. Appl. Mech. Eng., 372:113421, 2020.

grijalva2021eseismic
Felipe Grijalva, Washington Ramos, Noel Pérez, Diego Benítez, Román Lara-Cueva, and Mario Ruiz.
ESeismic-GAN: A generative model for seismic events from cotopaxi volcano.
IEEE J. Sel. Top. Appl. Earth Obs. Remote Sens., 14:7111–7120, 2021.

matinfar2023deep
Mehrshad Matinfar, Naser Khaji, and Goodarz Ahmadi.
Deep convolutional generative adversarial networks for the generation of numerous artificial spectrum-compatible earthquake accelerograms using a limited number of ground motion records.
Comput.-Aided Civ. Infrastruct. Eng., 38(2):225–240, 2023.

matsumoto2023fundamental
Yuma Matsumoto, Taro Yaoyama, Sangwon Lee, Takenori Hida, and Tatsuya Itoi.
Fundamental study on probabilistic generative modeling of earthquake ground motion time histories using generative adversarial networks.
Jpn. Archit. Rev., 6(1):e12392, 2023.

Florez2022
Manuel A. Florez, Michaelangelo Caporale, Pakpoom Buabthong, Zachary E. Ross, Domniki Asimaki, and Men‐Andrin Meier.
Data‐driven synthesis of broadband earthquake ground motions using artificial intelligence.
Bull. Seismol. Soc. Am., 112(4):1979–1996, 2022.

esfahani2023tfcgan
Reza D. D. Esfahani, Fabrice Cotton, Matthias Ohrnberger, and Frank Scherbaum.
TFCGAN: Nonstationary ground‐motion simulation in the time–frequency domain using conditional generative adversarial network (cGAN) and phase retrieval methods.
Bull. Seismol. Soc. Am., 113(1):453–467, 2022.

shi2024broadband
Yaozhong Shi, Grigorios Lavrentiadis, Domniki Asimaki, Zachary E. Ross, and Kamyar Azizzadenesheli.
Broadband ground motion synthesis via generative adversarial neural operators: Development and validation, 2024.
available at https://arxiv.org/abs/2309.03447 (last accessed April 2024).

rahman2022generative
Md Ashiqur Rahman, Manuel A. Florez, Anima Anandkumar, Zachary E Ross, and Kamyar Azizzadenesheli.
Generative adversarial neural operators.
TMLR, 2022.

Anderson1999
John G. Anderson and James N. Brune.
Probabilistic seismic hazard analysis without the ergodic assumption.
Seismol. Res. Lett., 70(1):19–28, 1999.

lavrentiadis2023overview
Grigorios Lavrentiadis, Norman A. Abrahamson, Kuehn M. Nicolas, Yousef Bozorgnia, Christine A. Goulet, Anže Babič, Jorge Macedo, Matjaž Dolšek, Nicholas Gregor, Albert R. Kottke, et al.
Overview and introduction to development of non-ergodic earthquake ground-motion models.
Bull. Earthq. Eng., 21(11):5121–5150, 2023.

akaba2024
Hinata Akaba, Taro Yaoyama, and Tatsuya Itoi.
Construction of site-specific ground motion model based on multi-output gaussian process.
J. Jpn. Assoc. Earthq. Eng., 24(1):1_73–1_91, 2024.
(in Japanese with English abstract).

Sung2024
Chih‐Hsuan Sung, Hiroe Miyake, Norman Abrahamson, and Nobuyuki Morikawa.
Nonergodic ground‐motion models for subduction zone and crustal earthquakes in Japan.
Bull. Seismol. Soc. Am., 2024.

Fnet
National Research Institute for Earth Science and Disaster Resilience.
NIED F-net, 2024.
available at https://www.fnet.bosai.go.jp/ (last accessed April 2024).

FUJIWARA2023
Hiroyuki Fujiwara, Morikawa Nobuyuki, Takahiro Maeda, Asako Iwaki, Shigeki Senna, Shinichi Kawai, Hiroki Azuma, Ken Xiansheng Hao, Masajiro Imoto, Kazue Wakamatsu, Junichi Miyakoshi, Takeshi Morii, Naomi Shimazu, Mari Takahashi, and Masaki Akatsuka.
Improved seismic hazard assessment after the 2011 great east Japan earthquake (part 2).
Technical Note of the National Research Institute for Earth Science and Disaster Resilience, pages 1–388, 2023.

jshis
National Research Institute for Earth Science and Disaster Resilience.
NIED J-SHIS.
National Research Institute for Earth Science and Disaster Resilience, https://doi.org/10.17598/nied.0010, 2019.

Arjovsky2017
Martin Arjovsky, Soumith Chintala, and Léon Bottou.
Wasserstein generative adversarial networks.
In Proc. of the 34th Int. Conf. on Machine Learning, pages 214–223, 2017.

radford2016unsupervised
Alec Radford, Luke Metz, and Soumith Chintala.
Unsupervised representation learning with deep convolutional generative adversarial networks, 2016.
available at https://arxiv.org/abs/1511.06434 (last accessed April 2024).

Gulrajani2017
Ishaan Gulrajani, Faruk Ahmed, Martin Arjovsky, Vincent Dumoulin, and Aaron C Courville.
Improved training of wasserstein gans.
In Adv. Neural Inf. Process. Syst., volume 30, 2017.

karras2020analyzing
Tero Karras, Samuli Laine, Miika Aittala, Janne Hellsten, Jaakko Lehtinen, and Timo Aila.
Analyzing and improving the image quality of stylegan.
In Proc. IEEE/CVF Comput. Soc. Conf. Comput. Vis. Pattern Recognit., pages 8110–8119, 2020.

Karras2021
Tero Karras, Samuli Laine, and Timo Aila.
A style-based generator architecture for generative adversarial networks.
IEEE Trans. Pattern Anal. Mach. Intell., 43(12):4217–4228, 2021.

Mescheder2018
Lars Mescheder, Andreas Geiger, and Sebastian Nowozin.
Which training methods for GANs do actually converge?
In ICML, volume 80, pages 3481–3490, 2018.

maas2013rectifier
Andrew L Maas, Awni Y Hannun, Andrew Y Ng, et al.
Rectifier nonlinearities improve neural network acoustic models.
In ICML, volume 30, page 3, 2013.

karras2018progressive
Tero Karras, Timo Aila, Samuli Laine, and Jaakko Lehtinen.
Progressive growing of gans for improved quality, stability, and variation, 2018.
available at https://arxiv.org/abs/1710.10196 (last accessed April 2024).

clevert2016fast
Djork-Arné Clevert, Thomas Unterthiner, and Sepp Hochreiter.
Fast and accurate deep network learning by exponential linear units (ELUs), 2016.
available at https://arxiv.org/abs/1511.07289 (last accessed April 2024).

He2016CVPR
Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun.
Deep residual learning for image recognition.
In Proc. IEEE Comput. Soc. Conf. Comput. Vis. Pattern Recognit., pages 770–778, 2016.

miyato2018cgans
Takeru Miyato and Masanori Koyama.
cGANs with projection discriminator, 2018.
available at https://arxiv.org/abs/1802.05637 (last accessed April 2024).

kingma2017adam
Diederik P. Kingma and Jimmy Ba.
Adam: A method for stochastic optimization, 2017.
available at https://arxiv.org/abs/1412.6980 (last accessed April 2024).

Paszke2019
Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, Alban Desmaison, Andreas Kopf, Edward Yang, Zachary DeVito, Martin Raison, Alykhan Tejani, Sasank Chilamkurthy, Benoit Steiner, Lu Fang, Junjie Bai, and Soumith Chintala.
Pytorch: An imperative style, high-performance deep learning library.
In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, and R. Garnett, editors, Adv. Neural Inf. Process. Syst., volume 32, 2019.

arias1970measure
Arturo Arias.
A measure of earthquake intensity.
In Robert J. Hansen, editor, Seismic Design for Nuclear Power Plants, pages 438–483. MIT Press, 1970.

Bommer1999effective
Julian J. Bommer and Alejandro Martínez-Pereira.
The effective duration of earthquake strong motion.
J. Earthq. Eng., 3(2):127–172, 1999.

Hancock2006
Jonathan Hancock and Julian J. Bommer.
A state-of-knowledge review of the influence of strong-motion duration on structural damage.
Earthq. Spectra, 22(3):827–845, 2006.

Trifunac1975
M. D. Trifunac and A. G. Brady.
A study on the duration of strong earthquake ground motion.
Bull. Seismol. Soc. Am., 65(3):581–626, 1975.

Pitarka1998
Arben Pitarka, Kojiro Irikura, Tomotaka Iwata, and Haruko Sekiguchi.
Three-dimensional simulation of the near-fault ground motion for the 1995 Hyogo-Ken Nanbu (Kobe), Japan, earthquake.
Bull. Seismol. Soc. Am., 88(2):428–440, 1998.

Lin2018
Yen‐Yu Lin, Te‐Yang Yeh, Kuo‐Fong Ma, Teh‐Ru Alex Song, Shiann‐Jong Lee, Bor‐Shouh Huang, and Yih‐Min Wu.
Source characteristics of the 2016 Meinong (ML 6.6), Taiwan, earthquake, revealed from dense seismic arrays: Double sources and pulse‐like velocity ground motion.
Bull. Seismol. Soc. Am., 108(1):188–199, 2018.

Baker2007
Jack W. Baker.
Quantitative classification of near-fault ground motions using wavelet analysis.
Bull. Seismol. Soc. Am., 97(5):1486–1501, 2007.

morikawa2008strong
Nobuyuki Morikawa, Tatsuo Kanno, Akira Narita, Hiroyuki Fujiwara, Toshihiko Okumura, Yoshimitsu Fukushima, and Aybars Guerpinar.
Strong motion uncertainty determined from observed records by dense network in Japan.
J. Seismol., 12:529–546, 2008.

Hikita2023
Tomoki Hikita and Yusuke Tomozawa.
Variability of response spectra of single-path ground-motion record pairs.
Earthq. Spectra, 39(4):2470–2491, 2023.

luzi2021double
Lorenzo Luzi, Yehuda Dar, and Richard Baraniuk.
Double descent and other interpolation phenomena in gans, 2021.
available at https://arxiv.org/abs/2106.04003 (last accessed April 2024).

Goodfellow-et-al-2016
Ian Goodfellow, Yoshua Bengio, and Aaron Courville.
Deep Learning.
MIT Press, 2016.

bucklew2004introduction
James Antonio Bucklew and J Bucklew.
Introduction to rare event simulation.
Springer New York, 2004.

AU2001263
Siu-Kui Au and James L. Beck.
Estimation of small failure probabilities in high dimensions by subset simulation.
Probabilistic Eng. Mech., 16(4):263–277, 2001.





