Exploring percolation phase transition in the three-dimensional Ising model with machine learning

Exploring the phase diagram of Quantum Chromodynamics (QCD) and searching for the critical point (CP) of the phase transition from hadrons to Quark-Gluon Plasma (QGP) is a topic of intense research in relativistic heavy-ion collisions [1−4]. Past research has suggested that at vanishing or low baryon chemical potential $ \mu_B $ and high temperature T, the transition is a smooth crossover [5]. Predictions from effective field theory indicate that at low T and high $ \mu_B $, the phase transition may be first-order [6, 7], with the CP representing the end point of the first-order transition line. Physicists are currently dedicated to investigating the boundaries of the QCD phase diagram and determining the location of the CP, as these subjects hold considerable scientific interest [8, 9].

The QCD phase diagram is not yet well understood, either experimentally or theoretically. Lattice QCD calculations are currently limited to scenarios with zero or small baryon chemical potentials due to the sign problem [10, 11]. As a result, various phenomenological models, including spin models, have been proposed to study the universality aspect of critical phenomena [12−16]. Phase transitions arise from spontaneous symmetry breaking, and systems sharing the same symmetry fall into the same universality class, exhibiting identical critical exponents and potentially similar critical behaviors. It is argued that the QCD critical point belongs to the Z(2) universality class, which contains the 3D Ising model as well [1, 17−20]. By mapping the parametric equation-of-state of the Ising model, one can establish a connection between the phase diagram of the 3D Ising model in the (T, H) plane and that of QCD in the (T, $ \mu_{B} $) plane [21, 22].

Over the past few decades, there has been increasing interest in the exploration of the geometrical characteristics of phase transitions [23−25]. Geometry holds a unique significance in its ability to provide insights into the underlying mechanisms and behaviors of critical phenomena. By analyzing geometric properties, one can predict critical exponents, scaling laws, and universal properties of a critical system [25−29]. In particular, percolation theory serves as an optimal framework for this purpose. As a statistical mechanism, percolation theory is dedicated to investigate the behavior of interconnected clusters within various media and how this connectivity affects the macro-scale properties of the system. This theory has demonstrated its utility in detecting the distinctive features of the QCD phase transition. It has been used to investigate the initial phase of deconfinement transition in heavy-ion collisions [30−32]. The method has also been extended to explore the confinement-deconfinement transition of both $S U $(2) [33−35] and $S U $(3) [36, 37] lattice gauge theories. Through a suitable definition of clusters, the deconfinement transition of gauge theory can be characterized by the percolation of clusters within the dominant sector. Additionally, the QGP to hadron crossover transition has been characterized via the temperature dependence of the shear viscosity over entropy density ratio, utilizing the percolation framework within the Color String Percolation Model [38].

Current studies on the phase transitions of the Ising model have commonly employed magnetization to explore the thermal magnetization transition under alterations in temperature or external field. This thermal magnetization transition is characterized by a pivotal Curie temperature ($ T_c $). For temperatures above$ T_c $, the system assumes a disordered phase, and for temperatures below T_c, the system transitions into an ordered phase characterized by non-zero spontaneous magnetization. The geometrical percolation transition, on the other hand, is related to the connectivity of spin clusters [39]. Envision a lattice of spins where each spin interacts with its nearest neighbors. At a specific critical percolation threshold, denoted as $ T_p $, the emergence of a percolation cluster, in which the microscopic elements become connected and form a sample-spanning path across the system, is suggested to be an indicator of the occurrence of a continuous percolation phase transition [40].

The investigations of percolation transitions within the Ising model originally concentrate on the Geometrical Spin (GS) clusters, which comprise nearest neighbor spins with the same sign in the lattice. The results of the 2D Ising model show that percolation transition occurs exactly at the critical temperature $ T_c $ of the thermal magnetization transition [41]. However, a noticeable disparity between these two transitions is observed in three dimensions [42]. Consequently, the introduction of Fortuin–Kasteleyn (FK) clusters, wherein nearest neighboring spins of the same sign are considered to belong to the same cluster with a certain probability, is proposed [43, 44]. By introducing a parameter of the bond probability, the critical point and critical exponents of the FK clusters between the thermal and percolation transitions coincide. Nevertheless, precise solutions of the percolation temperature and associate theoretical values of critical exponents concerning the percolation phase transition of the GS clusters within the 3D Ising model remain as open questions.

In traditional statistical physics, the order parameter is commonly used to classify different states of matter and identify phase transitions within a system. However, the intricate nature of robust interactions introduces considerable challenges in ascertaining the order parameter associated with the QCD phase transition and subsequently measuring it through experiments. Recent advancements have proposed integrating machine learning (ML) [45, 46] as a promising approach to explore this complex physical issue [47, 48]. The advantage of ML methods lies in their strong adaptability and generalization capabilities. They can automatically learn features, handle large-scale data, recognize complex patterns and relationships, and support multi-task learning. This makes ML a powerful tool not only for enhancing efficiency and discovery but also for solving complex problems and achieving intelligent decision-making. Empirical evidence now substantiates that ML techniques can effectively discern the order parameter [49−54], learn percolation phase transitions in various models [55, 56] and identify thermal magnetization phase transitions [12, 57−61] in the Ising model through various ML methodologies. Furthermore, recent studies have demonstrated that the same neural network can be utilized to identify different phase transitions belonging to the same universality class [62]. This finding enhances the feasibility of using ML methods to investigate phase transitions across a wide range of scenarios.

In this study, we aim to investigate the percolation transition of GS clusters in the 3D Ising model using state-of-the-art ML techniques. Unlike conventional approaches, these methods do not require prior knowledge of order parameters or any additional information. We demonstrate the efficacy of both the unsupervised and semi-supervised ML methodologies in categorizing different phases of the percolation transition. Furthermore, we attempt to extract the temperature of the percolation threshold and a critical exponent associated with correlation length directly from the output layer of the neural networks. The subsequent sections of this paper are structured as follows: In Sec. II, we offer a brief introduction to the percolation transition in the 3D Ising model. Sec. III gives an overview of the network algorithms employed in this study. In Secs. IV and V, we present and discuss the results pertaining to the identification of percolation phase transitions within the Ising model, utilizing Principal Component Analysis (PCA) and the Domain Adversarial Neural Network (DANN) methods, respectively. Finally, in Sec. VI, we summarize our discoveries and provide an outlook for future research directions.

The Ising model [63] is a classic spin model utilized for investigating phase transitions. We consider a three-dimensional cubic lattice composed of $ N=L\times L\times L $ sites with periodic boundary conditions. Each site is occupied by a spin, $ s_i $. We assume that the spin of individual sites can have one of the two states, spin-up, $ s_i = +1 $, or spin-down, $ s_i = -1 $. In this analysis, we only study the 3D Ising model in zero external field, as described by the Hamiltonian

where J is the coupling constant between two spins $ s_i $ and $ s_j $. Here, we set $ J=1 $ as the energy unit. The sum is taken only over nearest-neighbor pairs in the grid, and it describes the interaction of the spins with each other.

The Wolff algorithm, a Monte-Carlo method, is frequently employed to generate equilibrium configurations of the Ising model under various conditions, such as changes in system size or temperature. In the present investigation, we explore a temperature range spanning from $ T=0.02 $ to $ T=7.92 $ with an external magnetic field $ H=0 $. We generate a set of 2000 independent spin configurations for each selected temperature with a given system size.

In the Ising model, the GS clusters are identified as groups of nearest-neighboring sites with the same spin direction. As the temperature of the system changes, the spins tend to align with their neighbors due to interactions between them. The percolation phase transition infers a significant change in the way clusters of aligned spins form and extend within the whole lattice at a certain critical threshold $ T_p $. Beyond the percolation temperature, small clusters of aligned spins are isolated and do not span the entire lattice. However, as the temperature crosses $ T_p $, these smaller clusters start to coalesce and connect, leading to the emergence of a percolating cluster that wraps around the whole lattice. This percolating cluster signifies a sudden change in the system behavior, as the alignment of spins becomes correlated over long distances.

To facilitate an intuitive comprehension, Fig. 1 (a) shows the spin configurations within the 3D Ising model for a system size of $ 32\times 32\times 32 $ at a temperature of $ T=7.52 $. In this depiction, red lattices represent spin-up orientations, while blue ones denote spin-down. At high temperature, owing to the stochastic distribution of the spin states, no percolation clusters are discernible. Moving to Fig. 1 (b), which illustrates the configuration at a low temperature of $ T=4.42 $, we observe a predominant cluster colored in red. This cluster spans the whole lattice, signifying the presence of percolation clusters. The percolation transition is of great interest because it often leads to emergent behaviors and critical phenomena, where small changes in a model parameter can lead to drastic alterations in the overall behavior of the system.

Figure 1.  (color online) The spin configurations of the 3D Ising model with a lattice size of $ L = 32 $ are depicted (a) in the absence of percolation clusters at $ T = 7.52 $ and (b) in the presence of a percolation cluster at $ T = 4.42 $.

Supervised learning, unsupervised learning, and semi-supervised learning are fundamental paradigms in machine learning [64], each with different approaches and applications. Supervised learning involves training a model on labeled data, where input-output pairs are provided, enabling the model to learn patterns and make accurate predictions on new, unseen data. Unsupervised learning, on the other hand, deals with unlabeled data, aiming to uncover inherent structures, clusters, or relationships within the data. Semi-supervised learning merges elements of both, incorporating labeled and unlabeled data to enhance model performance. By leveraging the small amount of labeled data alongside the larger pool of unlabeled data, semi-supervised learning strikes a balance between efficiency and accuracy, making it valuable when acquiring fully labeled datasets is expensive or time-consuming.

To enhance the versatility and future applicability of our approach, we used both the unsupervised PCA method and the semi-supervised DANN methodology. It allows the extensino of the possible utility of our methods to the analysis of experimental data, ensuring a broader scope and improved adaptability for future applications in high energy physics.

PCA is among the most widely utilized multivariate techniques [65], with its origins tracing back to Pearson's pioneering work in statistics [66]. Pearson's formulation involves identifying optimal lines and planes to closely align with point distributions in space. PCA aims to simplify complex datasets by identifying a new set of orthogonal axes, called principal components, that capture the most significant variations in the data. These components are ordered in terms of the amount of variance they explain, allowing for the reduction of high-dimensional data into a lower-dimensional space while retaining as much relevant information as possible. Functioning as a classic example of unsupervised learning, PCAs are used in tasks such as data clustering and dimensionality reduction [49, 50, 67].

The sketch in Fig. 2 illustrates the network architecture employed by the PCA algorithm in this investigation. It can be segmented into two main components: the fit and transform stages. To begin with, the network takes the spin configurations of the data denoted as $ X_{(m,n)} $and the predetermined number of principal component features denoted as K as inputs. In the fit stage, PCA computes the mean and covariance matrix for the dataset $ X_{(m,n)} $. Subsequently, the Singular Value Decomposition (SVD) decomposes the covariance matrix into eigenvalues and their corresponding eigenvectors, which are denoted as V. Figure 2 presents the mathematical expression of SVD as $ U \ast \Sigma \ast V{}^{T}_{\left( n,n\right)} $. In this expression, Σ contains the singular values of the original matrix, while U and $ V_{\left( n,n\right)} $ represent the left and right singular vectors of the original matrix, respectively. $ V{}^{T}_{\left( n,n\right)} $ denotes the transpose of $ V_{\left( n,n\right)} $. Subsequently, based on the predefined number of features K, we retain the first K rows of the feature matrix V, represented as $ V_{(K,n)} = V[0:K] $. Shifting to the transform stage, the network executes a projection of the original data $ X_{(m,n)} $ onto the selected principal components $ V{}^{T}_{\left( K,n\right)} $. As a result, the data is transformed into a reduced-dimensional representation, denoted as $ X^{\prime}_{m,K} $. This accomplishment effectively realizes the objectives of both dimensionality reduction and feature extraction.

Figure 2.  (color online) The PCA architecture used in our analysis.

DANN [68] is a specific type of transfer learning [69, 70] that emphasizes the mapping relationship between labeled source domain data and unlabeled target domain data, rather than the clustering and dimensionality reduction features provided by PCA. DANN uniquely integrates deep feature learning and domain adaptation, enabling the classification decision to be based on both discriminative and domain-invariant features for accurate data classification. By leveraging DANN, the classification process can effectively utilize features that are both informative and immune to domain variations.

The overall structure of DANN is shown in Fig. 3, comprising three main components: a feature extractor (green), a label predictor (blue), and a domain classifier (red). The feature extractor captures informative features from the input data, transforming them into the feature vector f. It is composed of a convolutional layer, a max-pooling layer, and a fully connected hidden layer. In the convolutional layer, there are 16 filters, each measuring $ 4 \times 4 \times 4 $ with a stride of 1, applied to the input data to generate feature maps. The max-pooling layer uses a filter of size $ 2 \times 2 \times 2 $ with strides of 2 to reduce the dimensionality of the feature maps. These feature maps are then flattened and passed to a fully connected (FC) layer containing 50 neurons. Additionally, Batch Normalization (BN) and ReLU activation are applied to prevent overfitting and speed up the training process. Once the feature extractor has extracted feature vectors f from both the source and target domain data, the network forwards these vectors f to the label predictor and domain classifier. The label predictor consists of a fully connected layer with 2 neurons, applying BN and softmax activation. Its output is expressed as a vector ($ P_0 $, $ P_1 $), indicating the probabilities of the input data belonging to different event types. The domain classifier comprises a fully connected layer with 1 neuron, utilizing BN and sigmoid activation. It outputs the vector $ P_d $, determining whether the feature vector f originates from the source domain or the target domain; if it is from the source domain, $ P_d $ = 0; if from the target domain, $ P_d $ = 1. In our neural network architecture, the feature extractor and label predictor form a conventional feedforward neural network. Additionally, the feature extractor constitutes an adversarial network by connecting a gradient reversal layer to the domain classifier. The loss function for the domain adversarial network comprises two components: one for the label predictor and the other for the domain classifier. The loss function of the network is defined as: $ L(\theta_f, \theta_{y}, \theta_{d})= L_{y}( \theta_{f} , \theta_{y} )-L_{d}( \theta_{f} , \theta_{d} ), $ where y, f, and d refer to the labels of the three components of DANN, and θ represents the internal parameters of the network. The learning rate λ is set to 0.0005. During the iterative training process, the $ Adam $ optimizer in TensorFlow 2.4.1 is utilized to minimize the loss functions for each component of the network. Once training is complete, the trained model can be employed to predict unlabeled samples and generate prediction results. This approach utilizes both labeled and unlabeled samples to train different parts of the model, facilitating comprehensive training and optimization of the model.

Figure 3.  (color online) The DANN architecture used in the analysis.

The geometrical percolation transition typically involves several key quantities, such as percolation strength, the largest cluster size, average cluster size, as discussed in reference [71]. The cluster size, denoted as S, is determined by the count of spin sites it encompasses. To ensure sufficient data for subsequent calculations, our primary focus was on the largest cluster size, referred to as $S_{\max}$, which is influenced by both temperature and volume of the system.

It has been verified that $ S_{\rm max} $ plays a crucial role in understanding the behavior of systems undergoing a phase transition from a disconnected state to a connected state [14]. The appearance of a percolating cluster serves as a critical indicator of a phase transition within the system. Along with other extensive variables such as susceptibility and correlation length, the size of the largest cluster obeys scaling laws near the critical threshold. Understanding the size and structure of the largest cluster is essential for gaining insights into the critical behavior of percolating systems and characterizing overall connectivity properties.

In the thermal magnetization phase transition of the Ising model, it has been found that the first principal component ($ K_1 $), defined as the direction that maximizes variance in the data, captures the highest level of variability in the dataset [52]. The results derived from $ K_1 $ confirm its linear relationship with magnetization within a finite system size. This implies that the first principal component of the input data can effectively capture and learn about magnetization, which is the characteristic order parameter for the thermal magnetization transition in the 3D Ising model. To assess the capacity of the unsupervised ML method for learning and identifying important features associated with the largest cluster size from the input data in the percolation transition, we used PCA with one principal component.

To implement the PCA network, we employed the PCA class methods from the scikit-learn library 1.0.1 in Python 3.7.11, focusing on data processing with the relevant functions. The key hyperparameters of our PCA network include the dimensionality reduction parameter $ K_1 $ and the singular value decomposition method set as $ svd\_solver = full $. Upon configuring these hyperparameters, we employed the $ pca.fit() $ function to conduct principal component analysis on the input data. This function fits the PCA model to the data, enabling the computation of principal components and explained variance ratios. It involves learning the transformation parameters from the data necessary for dimensionality reduction. Subsequently, the we employed the $ pca.transform() $ function to reduce the dimensionality of data by projecting it onto the new feature space defined by the principal components. This process effectively reduces the dimensionality of data while preserving the essential information captured by the principal components.

We generated event samples for a given system size of $ L = 28 $ and cover a temperature range of $ T=3.52 $ to $ T=5.52 $. The spin configurations of the largest cluster, where the cluster spins only are retained and all other spins are set to 0, were utilized as input data for each sample of the 3D Ising model for the PCA network architecture as depicted in Fig. 2. Subsequently, we computed the average sizes of the 2000 largest clusters at each temperature point and conduct a Pearson correlation analysis with respect to the first principal component obtained from PCA. The correlation between the average size of the largest clusters $ \langle S_{\rm max}\rangle $ and the first principal component $ K_1 $ of the PCA is graphically illustrated in Fig. 4. Our findings reveal a linear relationship between $ \langle S_{\rm max}\rangle $ and $ K_1 $. The computed Pearson correlation coefficient of 1.0 suggests a strong positive correlation between $ \langle S_{\rm max}\rangle $ and $ K_1 $.

Figure 4.  (color online) The average size of the largest clusters as a function of the output of the first principal component in PCA with $ L = 28. $

The first principal component represents the newly defined coordinate axis within PCA. This axis encapsulates the highest amount of information and exhibits the greatest power in distinguishing patterns in the data. It plays the most significant role in accounting for variations in the data and excels at elucidating changes in data. The observed strong positive correlation between the average size of the largest clusters $S_{\rm max}$ and $ K_1 $ provides further confirmation that the largest cluster carries a substantial amount of critical system information. In addition, this information can be effectively acquired from the first principal component of PCA in the geometric percolation phase transition.

Next, we explored the capability of PCA in identifying different phases and investigating the critical threshold associated with the percolation transition within the 3D Ising model. As detailed in Sec. II, the Ising model undergoes a geometric percolation phase transition at a critical threshold $ T_p $. This transition separates a disconnected state at high temperatures from a fully connected state at low temperatures. To achieve this, we establish a PCA network to conduct unsupervised learning directly on samples of spin configurations of the largest clusters.

The numerical results obtained from five different system sizes and spanning the temperature range $ T\in \left[ 3.52,5.52\right] $ are illustrated in Fig. 5 (a). The first principal component $ K_1 $ initially exhibits a gradual decrease as the temperature increases for all sizes considered. This is followed by a sharp drop as the temperature approaches a specific point, ultimately reaching equilibrium at higher temperatures. Notably, the $ K_1 $ values for different sizes intersect at a particular temperature, corresponding to the critical threshold of the percolation phase transition. This demonstrates the effectiveness of PCA as a methodology for classifying the two different phases of the percolation transition in the 3D Ising model.

Figure 5.  (color online)(a) The first principal component of PCA as a function of temperature for five different system sizes. (b) Finite-size scaling analysis to determine percolation critical threshold based on the first principal component.

We perfomed five independent iterations of the PCA network, each with a different size of $ L = $ 12, 16, 20, 24, 28. The pseudo-critical threshold T for each size corresponds to the point where $ K_{1} =0 $. Subsequently, we employed finite-size scaling analysis to estimate the percolation transition temperature in the infinite L limit, using information derived from the first principal component [52, 72]. Sebsequently, using the methods employed in Refs. [73, 74], we obatined the percolation temperature $ T_p $ by extrapolating through the fit of $ |T - T_p| \sim L^{-1/\nu} $ as the limit $ 1/L $ approached zero, as illustrated in Fig. 5 (b). Statistical errors were estimated using the standard deviation and were smaller than the size of the data points. The percolation temperature determined from the fitting process was $ T_p = 4.475\pm 0.004 $. This critical threshold for the percolation transition in the 3D Ising model is slightly lower than the critical temperature of the thermal magnetization phase transition, which is $ T_c = 4.5115\pm 0.0001 $[75, 76]. This finding is qualitatively consistent with the results in Refs. [42, 77].

DANN is a deep learning technique mainly applied in the domain adaptation field. The primary advantage of DANN is its ability to adapt a machine learning model from one domain to another, particularly when the source and target domains have different data distributions. This neural network is trained in such a way that the feature representations of the two domains become indistinguishable to the domain classifier.

In our analysis, we labeled percolation phenomenon at extremely low temperatures as phase $'1'$ and those with no percolation phenomenon at exceedingly high temperatures as phase $ '0' $ in the Ising model. Consequently, we designated the spin data of the largest clusters at low and high temperatures as the source domain data during DANN network training, with the unlabeled one in the intermediate temperature range regarded as the target domain data. Detailed network architectures and the training process of DANN are illustrated in Fig. 3. To establish the optimal temperature range for each scale, we employed a technique involving the fixation of either low-temperature or high-temperature labels while decreasing the high-temperature labels or increasing the low-temperature labels, respectively. The specific selection of the temperature range for target domain data at various system sizes is outlined in Table 1. The source domain data can then be annotated with suitable labels according to the target domain range. In this process, DANN learns the mapping relationship between the source domain data and the target domain data. The ultimate label predictor proves effective at accurately forecasting the percolation classification of unlabeled data in the target domain through the utilization of domain adaptation and back-propagation techniques.

After training the DANN on the optimal domain, we assessed samples at various temperatures, and then used DANN to predict the probability of each configuration belonging to phase '1' or phase '0'. The numerical results obtained at various system sizes at a vanishing magnetic field are illustrated in Fig. 6 (a). The average outputs of different sizes cross at a specific temperature corresponding to the percolation temperature. It infers that the DANN can successfully classify the two different phases in the percolation transition in the Ising model.

Figure 6.  (color online) (a) The output layer averaged over test sets of DANN as a function of temperature for five system sizes. (b) Finite-size scaling analysis to determine the percolation critical threshold. (c) Data collapse of the average output layer as a function of $ tL^{1/v} $, where $ t = (T-T_{p})/T_{p} $ is the reduced temperature.

To determine the percolation critical threshold, we perfomed DANN training at five sizes and identified the pseudo-critical threshold T as the intersection point of the two curves corresponding to phase $ '0' $ and phase $ '1' $ for each size. The percolation threshold $ T_p $ was obtained by extrapolating the fit of $ |T - T_p| \sim L^{-1/\nu} $ as the limit $ 1/L $ approached zero, as illustrated in Fig. 6 (b). The obtained critical threshold for the percolation transition, $T_p = 4.466\pm 0.005$, is in agreement with the result from the PCA network within errors. This temperature is slightly lower than the critical temperature associated with the thermal magnetization phase transition in the 3D Ising model.

The critical exponents associated with the 3D Ising model play a crucial role in characterizing the behavior of various thermodynamic and correlation functions as the system approaches the critical point. They are considered universal in nature, $ i. e. $, they do not change with the specific characteristics of the physical system. Among these exponents, the correlation length exponent ν holds significance as it reveals how the correlation length undergoes a pronounced divergence near the critical point. This divergence is described by the relationship between the correlation length ξ and temperature, expressed as $ \xi\sim|T - T_p|^{-\nu} $. For a finite system, the correlation length is expected to scale proportionally to the size of the system. Consequently, one can establish a connection between the temperature and system size characterized by the expression $ |T - T_p| \sim L^{-1/\nu} $.

To extract the critical exponent ν, we employed the data collapse technique [78, 79] and the results are shown in Fig. 6 (c). The figure demonstrates that with different system sizes, the system has a compelling convergence with $ \nu = 0.72 \pm 0.03 $. This value quantitatively agrees with the previously measured cluster size heterogeneity result of $ \nu \approx 0.75 $ for geometric clusters in the 3D Ising model [80]. Earlier studies on FK clusters [34, 81] demonstrated that the calculated correlation length exponent ν is consistent with that of the thermal magnetization phase transition, reported as $ \nu = 0.629 912(86) $ using the MC method [75]. Our findings indicate that the extracted correlation length exponent for GS clusters surpasses that of the thermal magnetization phase transition in the 3D Ising model.

In summary, we employed both unsupervised and semi-supervised machine learning techniques to investigate the percolation phase transitions of the geometrical spin clusters within the 3D Ising model. We obtained a linear correlation between the average size of the largest clusters and primary component of PCA, suggesting that the largest clusters contain sufficient information pertaining to the percolation transition. We used two distinct machine learning approaches by utilizing the spin configurations from the largest clusters as input data, and applied the finite-size scaling method to estimate the critical threshold of the percolation phase transition. Both of these machine learning methodologies effectively classify disconnected states at higher temperatures and fully connected states at lower temperatures. As a result, we obtained the critical threshold for the percolation transition as $ T_p = 4.475\pm 0.004 $ and $ T_p = 4.466\pm 0.005 $ using PCA and DANN, respectively. These temperatures are slightly lower than the critical temperature of the thermal magnetization phase transition. The extracted correlation length exponent was $ \nu = 0.72 \pm 0.03 $, which is consistent with the cluster size heterogeneity of geometric clusters. This value is greater than the critical exponent associated with the thermal magnetization phase transition in the 3D Ising model.

Percolation theory has become a powerful tool for investigating phase transitions in various physical systems. The machine learning methods developed in this study have the remarkable ability to identify phase transitions and extract critical temperature and critical exponents with minimal or even no training data. This capability allows for convenient application in the investigation of criticality in high energy experimental data or other Monte Carlo models, which could aid the exploration of the underlying physical mechanisms governing the QCD phase transition in heavy-ion collisions.

We are grateful to Prof. Mingmei Xu, Lizhu Chen and Dr. Feiyi Liu for fruitful discussions and comments. We further thank Prof. Hengtong Ding for providing us with computing resources. The numerical calculation have been performed on the GPU cluster in the Nuclear Science Computing Center at Central China Normal University (NSC3).