Goodfellow et al. proposed a novel architecture based on convolutional neural networks in 2014, which introduced a fundamental departure from traditional neural networks by comprising two components: a generator and a discriminator. This framework was named Generative Adversarial Networks (GAN) (Wu et al., 2021).

As illustrated in Figure 2 the input data includes “ingredients” required to synthesize realistic data: original feature values (labels), random noise, and raw data. Through the generator, new labeled data is produced and passes the discriminator’s evaluation. After iterative training, the network ultimately generates data indistinguishable from real samples (Mao et al., 2017).

Traditional GANs employ KL divergence or JS divergence for training, but these methods suffer from issues like gradient vanishing and mode collapse. To address these problems, the Wasserstein Generative Adversarial Network (WGAN) was proposed. The core innovation of WGAN lies in replacing traditional divergence metrics with the Wasserstein distance to measure the distance between probability distributions (Equation 1). W P r , P g = inf γ ∼ ∏ P r , P g ⁡ E x , y ∼ γ ∥ x − y ∥ (1)

In the above formula, ∏ P r , P g represents the set of all joint distributions of P r and P g . For each possible joint distribution γ ,the distance ∥ x − y ∥ between the real sample x and the fake sample y can be calculated. Then, the expected value E x , y ∼ γ ∥ x − y ∥ of the samples and the distance for the joint distribution γ can be computed. The lower bound inf γ ∼ ∏ P r , P g ⁡ E x , y ∼ γ ∥ x − y ∥ of the entire probability is the Wasserstein distance. It can be understood that E x , y ∼ γ ∥ x − y ∥ represents the energy consumption required to push P r to P g under the overall situation of γ and W P r , P g is the minimum consumption for the shortest distance.

Since the solution of inf γ ∼ ∏ P r , P g cannot be obtained directly, it is changed to a Lipschitz - continuous form. That is, a restriction condition is added to the function to obtain an optimal solution that can be used as a substitute. As shown in Formula 2. W P r , P g = 1 K sup ∥ f ∥ L ≤ K ⁡ E x ∼ P r f x − E x ∼ P g f x (2)

In the above - mentioned (Formula 3), K ≥ 0 and any two elements must satisfy. f x 1 − f x 2 ≤ K x 1 − x 2 (3) which is indicated by the arrow (3). K is called the Lipschitz constant. Then, Formula 5 can be understood as taking the upper bound of E x ∼ P r f x − E x ∼ P g f x for all eligible under the condition of the constant K ,and then shrinking it by K times. In generative adversarial networks, the parameter ω is needed to complete the update of the function f ω , as shown in the following Formula 4. K · W P r , P g ≈ max ω : f ω L ≤ K ⁡ E x ∼ P r f ω x − E x ∼ P g f ω x (4) L = E x ∼ P r f ω x − E x ∼ P g f ω x (5)

The results obtained in the above formula are approximate to the Wasserstein distance between the real samples and the generated samples. By minimizing the Wasserstein distance, that is, minimizing L , the problem of gradient vanishing can be solved. Since the activation function in the last layer of the discriminator network in the original GAN is designed for a binary classification task, but now it is improved using the Wasserstein distance, there is no need to continue using the Sigmoid function as the activation function in the last layer. From Equations 6, 7, we can obtain the latest loss functions of the two networks. L G = − E x ∼ P g f ω x (6) L D = E x ∼ P r f ω x − E x ∼ P g f ω x (7)

In summary, the key improvements of WGAN can be categorized into four principal aspects: (1) Removing the Sigmoid function in the last layer of the discriminator network; (2) Adopting the improved forms of the loss functions for both the generator and discriminator networks as mentioned above, instead of using the log - form; (3) Ensuring that the absolute value of the updated parameters does not exceed a fixed constant; (4) Selecting appropriate optimization algorithms. Sometimes, using momentum gradient descent and adaptive moment estimation algorithms may cause instability in the loss function, leading to substantial variations in the generated samples.

The utilization of highly geologically plausible and high-fidelity labeled data can significantly enhance the predictive performance of neural networks. This study employs a WGAN (Wasserstein Generative Adversarial Network) to generate high-quality labeled data. The process first requires obtaining prediction results from actual seismic data processed through single-attribute analysis within the target area, which serves as the “raw material” for training the WGAN. For validation purposes the North Sea F3 Block in the Netherlands, a publicly available seismic dataset collected in 1987, is widely referenced in studies on AI-based seismic fault detection. This dataset contains numerous faults with distinct characteristics, making it an excellent benchmark for validating the model’s training effectiveness.

First, the seismic data is processed using the single-attribute method (Figures 3–6), and its prediction results are input into the WGAN network. Through multiple iterations, this generates labeled data that more closely resembles real geological samples.

After processing the single-attribute-generated results through slicing and grayscale/binary image processing, geologically realistic labeled samples are obtained. These pseudo-labeled samples are then input into the generator network. The generated pseudo-labeled data and real data are subsequently assessed for authenticity by the discriminator network. Through iterative refinement, the gap between the generator’s output and the real data gradually diminishes until high-quality labeled data meeting research requirements is achieved, thereby completing the expansion of the labeled dataset.

The WGAN-generated label samples successfully mitigate the issue of sample scarcity in previous methods. As demonstrated in the figure below, which depicts label samples obtained from the trained network, the generated labels exhibit enhanced realism and diversity. Unlike prior approaches that were constrained to detecting only prominent fractures and simplistic/linear labels, these improved results demonstrate finer structural details and more natural variations, thereby achieving closely alignmen with real-world geological features (Figure 7).

The study area is located in the middle-northern segment of the Western Sichuan Depression, within the Xujiahe Formation, extending in a near-east-west direction. Its boundaries are defined by the Xiaquan Structural Belt, Fenggu Structural Belt, Zitong Sag, Chengdu Sag, and its eastern margin (Wang et al., 2024). This region has undergone multiple tectonic events since the Late Triassic, forming a complex ancient large uplifted region (Deng et al., 2022; Zhao et al., 2024). The Xujiahe Formation in the Western Sichuan Depression has experienced successive tectonic movements during the Indosinian, Yanshanian, and Himalayan orogenic periods, ultimately shaping a regional structural framework characterized by high relief in the northwest and low-lying topography in the southeast (Zhao et al., 2017).

Well logging interpretation is essential for fracture prediction, as it provides insights into fracture development. This study primarily utilizes well logging data to conduct comparative interpretations with seismic attributes and deep learning methods, thereby verifying the accuracy of our fracture prediction results. The presence of fractures induces anomalous responses in conventional well logging curves (Figure 12), such as distinctive patterns in acoustic curves at fracture intervals. While imaging well logging offers superior fracture identification capabilities, its high cost limits data availability in the study area. Consequently, this research combines conventional well logging data with supplementary imaging well logging data to enhance fracture identification.

Through conventional well logging and imaging well logging methods, comprehensive interpretations have been conducted to describe fracture feature variations across different regions. Below are detailed analyses of two representative wells. Figure 13 illustrates the rose diagram of fracture orientations and inclination angle distribution histogram for Well 1. The results indicate that fractures in Well one and Well three predominantly exhibit SE-NW (120°–130°) and SWW-NEE (70°–90°) orientations. Fracture inclinations are primarily 0°–10°, with moderate-to-high-angle fractures (10°–30°) also well-developed. Additionally, fracture aperture increases progressively with depth. Figure 13 displays the rose diagram and inclination angle distribution histogram for Well 2. Fractures in this well are dominated by low-angle orientations (10°–20°), while moderate-angle fractures are underdeveloped, and high-angle fractures are absent. In the Xu2 Member, fracture linear density initially increases with depth and then fluctuates downward.

Based on calculating the linear density of fractures in the target stratigraphic interval (Xu2 Member of the Xujiahe Formation), high-value zones of fracture density are identified within the three-stage tectonic stress integration region. The well logging results will serve as the primary validation basis for subsequent fracture prediction outcomes derived from multi-attribute prediction and the multi-attribute fusion generative adversarial network (WGAN) method. Through integration with well logging interpretations and fracture development characteristics, the final fracture characterization results will achieve greater credibility.

Seismic data with distinct attribute characteristics play a crucial role in fracture identification. The Xujiahe Formation primarily consists of tight sandstone deposits, and its seismic data exhibit high quality, with notable resolution and frequency attributes. As shown in Figure 14, the seismic bandwidth is 63 Hz, and spectral analysis reveals that the dominant frequencies are primarily distributed between 20 Hz and 80 Hz.

The identification results of fractures vary under different frequency characteristics. In this study, we employ Matching Pursuit Spectral Decomposition (MPSD) processing within the effective frequency band to obtain distinct frequency-component volumes. Subsequently, structure-oriented filtering is applied to these frequency volumes to enhance the detectability of fault dip strength and azimuth attributes (Figure 15). This integrated approach yields results that more accurately characterize fracture features, as the frequency-dependent processing better captures fracture-related seismic responses while structural filtering improves the continuity and interpretability of fracture systems. The methodology demonstrates superior performance in identifying multi-scale fractures compared to conventional full-spectrum approaches, particularly in complex tight sandstone reservoirs where fracture manifestations vary significantly across frequency bands.

As illustrated in Figure 16, this study selects frequency-component volumes at 30 Hz, 40 Hz, 50 Hz, and 60 Hz for coherence attribute processing, yielding distinct fracture identification outcomes.

As shown in the figure above, different frequency bands exhibit distinct responses to fractures of varying scales. Coherence processing of low-frequency components can effectively characterize medium-to-large scale fractures, while coherence processing of high-frequency components can delineate small-scale fractures and fractures zones (Bahorich and Farmer, 1995). Through coherence slices at different scales, it can be observed that low-frequency components mainly reflect the macro-scale distribution of fracture development, whereas high-frequency components can provide more detailed characterization of small-scale fractures (Gersztenkorn and Marfurt, 2002).

To evaluate the method proposed in this study, we performed additional fracture identification using four attributes (Huang et al., 2025). The specific results are shown below (All the result values in the figure are obtained after standardized processing using this method (Niu et al., 2025; Ren et al., 2024). The larger the value, the greater the probability of fracture distribution (Li et al., 2024; Zhang et al., 2025).

Fracture prediction based on maximum likelihood attributes demonstrates that both standard and refined maximum likelihood approaches achieve comprehensive delineation of large fractures. Particularly noteworthy is that the refined maximum likelihood method enhances the identification of small fractures (Figure 17).

The fundamental principle of fracture prediction based on coherence attributes lies in delineating fractures by analyzing lateral variations in seismic waveform continuity. This method enables clear identification of large-scale fractures within the seismic dataset (Figure 18; Kurt et al., 2012).

Compared with coherence attributes and maximum likelihood attributes, the fracture layering description derived from gradient structure tensor attributes provides clearer delineation of fracture angles and large-scale structural features (Chen et al., 2012). However, it exhibits limitations in characterizing small fractures (Wang et al., 2018). In essence, this approach yields more accurate information on the overall fracture trends within the Xujiahe Formation (Figure 19; Chopra and Marfurt, 2007; Lou et al., 2022).

Fracture prediction based on curvature attributes exhibits a strong response to linear features, enabling clear characterization of fault structures and fracture delineation. Among the commonly used methods, maximum positive curvature and minimum negative curvature attributes yield the most effective results (Roberts, 2001). The curvature calculations of both approaches enhance the description of stratal bending, produce stronger responses at fracture locations, and provide richer geometric information about the formations (Figure 20; Suo et al., 2012).