Deep learning of topological phase transitions from entanglement aspects: An unsupervised way

An autoencoder (AE) is a type of neural network that aims at compressing input data into more efficient representation in an unsupervised manner. It mainly consists of two parts, an encoder and a decoder. The encoder $f$ takes $d$-dimensional input data $x$ and outputs a $l$-dimensional latent variable $z=f(x)$; the decoder $g$ then maps $z$ back to $x^{\prime}=g(z)$ in $d$-dimension. The learnable parameters of the model are trained by performing gradient descent updates in order to minimize the reconstruction loss, $L(x,x^{\prime})$, usually chosen to be the mean squared error. Since typically $l<d$ after encoding the data, using AE is often viewed as a nonlinear dimension reduction method.

A typical model architecture we used in this study is schematically shown in Fig. 3(b). In the encoder part, it is made of a convolutional module (ResNet-like structure He16 ) followed by a linear module consisting of a few fully connected hidden layers. In the decoder part, it is basically made of similar layer structures but arranged in a reversed order with respect to the encoder. Note that, however, the convolutional layers here are replaced by transposed convolutional ones. The activation functions of the intermediate layers are always rectified linear units (ReLUs) Jarrett09 , except for the final layer where the sigmoid function is used.

Although the traditional AE can learn a function to encode two input data points into distinct latent variables $z_{1}$ and $z_{2}$, one may have no idea what would be the decoded result when giving the input $(z_{1}+z_{2})/2$. To overcome this arbitrariness in AEs, we also consider the variational autoencoders (VAE) which can learn a latent variable model $g(x,z)$ with a joint distribution of a latent variable $z$ and input $x$ Kingma13 . In sharp contrast with traditional AEs, $z$ here is basically drawn from some prior probability distribution $p(z)$, which is almost always chosen to be a multivariate Gaussian, and thus leads to certain controllability. In addition, the weights of the VAE are now trained by simultaneously optimizing two loss functions, a reconstruction loss and the Kullback-Leibler (KL) divergence between the learned latent distribution and a prior unit Gaussian. Such an additional KL divergence loss can be viewed as a regularization term in a traditional AE. The VAE has a similar model architecture compared to that in Fig. 3(b) by simply getting rid of convolutional modules at the head and tail. Furthermore, as shown in Fig. 3(c) the middle most layer would output the multiple means and the variances, depending on the number for encoded features we need, in the latent space.

Principal component analysis (PCA) is a standard yet simple method for dimensional reduction and data visualization Pearson01 ; Jolliffe02 . It is an orthogonal, linear transformation of the input features to a sorted set of new variables by their variance. Such method is motivated by the experience that in many cases, the most relevant information can be revealed in the direction with largest variance for a given signal, while directions with small variance usually indicate noises and may be neglected.

Concretely, let us consider $N$ $p$-dimensional feature vectors, ${\boldsymbol{X}}=\{{\boldsymbol{x}}_{1},{\boldsymbol{x}}_{2},\cdots,{\boldsymbol{x}}_{N}\}$. One can assume that the mean of all vectors, $\sum_{i}{\boldsymbol{x}}_{i}=0$, without loss of generality, and then ${\boldsymbol{X}}$ is called a (zero-mean) centered matrix. By definition, the transformation weight vector for producing the first principle component ${\boldsymbol{w}}_{1}$ can be found by

The next ordered weight vectors are then obtained by repeating Eq. (8) after subtracting out the calculated principle components from ${\boldsymbol{X}}$. However, in practice, one can prove that this procedure for getting ${\boldsymbol{w}}_{i}$ is equivalent to find out the eigenvectors of the $N\times N$ symmetric matrix ${\boldsymbol{X}}^{T}{\boldsymbol{X}}$, Mehta19 i.e.,

where we have assumed that the obtained eigenvalues are sorted such that $\lambda_{1}\geq\lambda_{2}\geq\cdots\geq\lambda_{N}\geq 0$, representing variances for the input feature vectors. It is also useful to define the relative variance, $\tilde{\lambda}_{i}=\lambda_{i}/\sum_{i}{\lambda_{i}}$ in order to count accumulated variance percentage.

K-means clustering is a simple and easily understandable clustering algorithm without any supervision MacQueen67 ; Lloyd82 . Given a prior knowledge about the number of clusters $K$, the basic idea is to find the best cluster means such that the variance within each cluster is minimized. To put it more precisely, consider a set of $N$ $p$-dimensional data points without labels, ${\boldsymbol{X}}=\{{\boldsymbol{x}}_{1},{\boldsymbol{x}}_{2},\cdots,{\boldsymbol{x}}_{N}\}$ and call the set, ${\boldsymbol{C}}=\{{\boldsymbol{\mu}}_{1},{\boldsymbol{\mu}}_{2},\cdots,{\boldsymbol{\mu}}_{K}\}$ (${\boldsymbol{\mu}}$ is also $p$-dimensional), as the $K$ centers for the whole data. The objective of K-means method is then to assign each ${\boldsymbol{x}}_{i}$ to an appropriate cluster such that the loss function,

is minimized. Note that the assignment $a_{ik}$ is 1 if ${\boldsymbol{x}}_{i}$ is assigned to cluster $k$ while 0 otherwise and $\sum_{k}a_{ik}=1$ for every $i$. The implementation is usually done by iteration until certain convergence with chosen tolerance level has been achieved.

Although K-means clustering is quite intuitive, one still needs to provide the number of clusters $n$ as a priori knowledge. To obtain a more reasonable estimation of this number and somehow eliminate the effect of distance function chosen in K-means method, we employ Silhouette analysis (SA) to justify it Rousseeuw87 ; deAmorim15 . For a given set of clusters $\{\mathcal{C}_{i}\}$, SA assigns a value to each data point $x\in\mathcal{C}_{i}$ by

where $b(x)=\min_{j\neq i}b_{j}(x)$ with $b_{j}(x)$ being the mean distance of $x$ to all points in the cluster $j$ and $a(x)$ is the mean distance between $x$ and all other data points in the same cluster $i$. In other words, the Silhouette value $s(x)$, bounded between $\pm 1$, is a measure of how similar $x$ is to its own cluster (cohesion) compared to the other clusters (separation). Considering the mean of Silhouette values, $s$-score, as a function of $n$ (after K-means clustering), the best estimation of $n$ is simply the one gives maximum of $s$-score.

In addition, once the best $n$ is determined, the data point with the highest Silhouette value within each cluster may then be taken as a confident seed to build a “labeled” training set to train a neural network via supervised learning approach. It in turn could be used to make sharper phase boundaries in a phase diagram, originally obtained in an unsupervised way over interested parameter space.

Following our ML pipeline mentioned in Section III, the first step is to train a neural network to encode our type-I input data to effective representations in the latent space. However, in order to determine the minimal dimension $d_{z}$ of the latent space, we train a series of AEs with same model architecture appendix except for the number of hidden neurons in the middlemost layer ($n_{mid}$, from 2 to 10). For each $n_{mid}$, we record the necessary dimension of the (converged) latent representations to keep at least $99\%$ variance of them by PCA. Due to the unconstrained nature of the latent representations in AEs, we repeat 100 times of the same training procedure with the same initial weight distribution. As shown in Fig. 4(a), we observe that the minimal dimension $d_{z}$ would be 4, because such number becomes dominant in the discrete distribution of $d_{z}$ when $n_{mid}$ increases.

Next, SA is utilized to estimate the best number of clusters $n$ for the latent representations of all input data via K-means method. Note that these representations are provided from previous trained AEs with $n_{mid}=4$ (as suggested above) hidden neurons in the middlemost layer. The box plot, Fig. 4(b), clearly shows that the mean of Silhouette value achieves the highest one when $n=3$. By projecting the 4D latent representations into a 2D space spanned by the first two principle components (features), we obtain Fig. 4(c). This feature plot gives us an insight about how the system could be reasonably divided into three clusters (phases).

Again, due to the unconstrained nature of the latent representations in AEs, the phase transition points found are statistically at mean values -1.012 and 0.980 with standard deviations 0.093 and 0.081, respectively, after collecting clustering results from 100 sets of latent representations via different trained AEs. This is quite close to the theoretically calculated values -1 and 1, but having relatively large deviation. To make phase boundaries sharper, we further train a CNN classifier appendix in supervised learning manner to predict the whole phase diagram. To prepare a “labeled” dataset, we first pick up three seeds in the input “images”, each of which gets the highest Silhouette value in each cluster and is thus believed to be inside each phase with strong confidence. In our example, the three seeds are located at $\mu/2t=-2.988,0.182,3.019$ (from latent representations traced back to original input data points), and for each cluster (phase) we expand symmetrically around the seed by a window width $0.1t$ to obtain 2000 points with equal spacing. These 6000 data points are finally formed our training dataset, while the original 20,001 ones become our test set without labels. As shown in Fig. 4(d), the phase boundaries obtained by trained CNN classifiers are clearly sharper at mean values -1.015 and 1.038 with smaller deviations 0.018 and 0.023 (after repeating same training procedure 100 times).

Alternatively, we next consider taking BCM generated quantities as our training inputs. There are two potential formats. We first train an AE to encode the type-II inputs, whose format looks simpler. Note that here we choose the middlemost layer to have 2 hidden neurons according to similar experiments done in type-I case [(see Fig. 4(a)]; moreover, all convolution-related modules are taken away in the model architecture appendix . As shown in Fig. 5(a), SA shows that 2 clusters, among other chosen number of distinct clusters, are the best result obtained via K-means method on the latent representations. It indicates that, when combining with the projected plot of the latent features, as seen in Fig. 5(b), this approach can only distinguish a topological phase (phase I) from non-topological ones (phases III and IV). The reason that phases III and IV can not be further distinguished can be attributed to rather limited information provided by type-II inputs, as indicated in our previous work tsai20 .

To gain more information, we use type-III inputs to train autoencoders with $n_{mid}$ hidden neurons in the middlemost layer to get the corresponding latent representations. As shown in Fig. 6(a), the necessary dimension of the latent space is proportional to $n_{mid}$ and thus we choose the dominant one at large $n_{mid}$, namely, $d_{z}=8$. Note that to avoid the arbitrariness of phase when computing eigenvectors, we have preprocessed the inputs by squaring each entry of the input matrix (“image”). After K-means clustering the latent representations obtained from AE, Fig. 6(b) depicts the results of SA and the input data are suggested to be separated into 3 clusters. Fig. 6(c) shows how the encoded representations in the 3D latent space after projecting to 2D may be divided into 3 clusters. In fact, by plotting the phase diagram as a function of $\mu$, the phase transition boundaries are somewhat shifted from theoretical values with relatively large deviation. This, however, can be improved if following the same supervised learning strategy as mentioned in the type-I input case typeiii . Note that the labeled training dataset used here is based on 3 seeds at $\mu/2t=-3.362,-0.026,3.359$, each of which gets the highest Silhouette value in the corresponding cluster. Even though these vales are taken in the subsystem A with $L=10$, Fig. 7 shows that the locations of the most confident Silhouette values do not change much as $L$ increases. Therefore, we finally take the subsystem of size $L=100$ to refine the transition boundaries and to reduce the possible finite-size effect all together. As shown in Fig. 6(d), the phase transition points found then are statistically at mean values -0.917 and 0.978 with standard deviations 0.149 and 0.0719, respectively, using well-trained CNN classifiers appendix .

In VAE approach, we try to impose certain constraints (multivariate Gaussian here) on the encoded representations and let a VAE learn the parameters of a probability distribution modeling the input data. Following the same ML pipeline, the first step is to train a series of VAEs with $n_{mid}=2$ to 10 pairs of hidden neurons in the middlemost layer appendix ; in each pair one neuron outputs the mean of encoded representation distribution and the other produces its corresponding variance. After repeating the same training procedure 100 times for a given $n_{mid}$, the number of necessary dimension for the latent representations obtained after PCA to keep at least $99\%$ variance of the encoded data is simply one (pair). In sharp contrast with AE, VAE provides a more stable result, as clearly shown in Fig. 8(a).

Once the minimal dimension of the encoded representations is set, we employs SA to estimate the best number of clusters $n$ for them via K-means method. From Fig. 8(b), it shows that the Silhouette value achieves the highest one when $n=3$. Fig. 8(c) depicts the encoded representations as a function of $\mu/2t$ and different colors indicate distinct clusters (phases). By plotting the phase diagram, we find that transition points are at mean values -0.992 and 0.9597 with standard deviations 0.047 and 0.026, respectively, after collecting results from 100 trained models with the same initial weight distribution. This small deviation is significantly different from AE approach, where the position of transition points found is more unstable. The supervised learning in the last step is not necessary in this case.

On the other hand, we take BCMs as our training input. Firstly, we train a VAE to encode the type-II inputs. Note that here we enforce the middlemost layer to output one set of mean and variance because similar PCA tests have been done to determine the minimal dimension of the encoding representations. Next, as shown in Fig. 9(a), SA shows that 2 clusters are the best result obtained through K-means method on the encoded representations [see Fig. 9(b)]. It turns out that this approach can distinguish topological phase (phase I) from non-topological ones (phases III and IV), but phases III and IV can not be further distinguished. This result is again consistent with previous work tsai20 , indicating too compressed information provided by type-II inputs.

Thus, we next examine to train another VAEs by feeding in type-III inputs, which have been preprocessed by squaring each entry in matrices. We find the minimal dimension of the latent dimension to be 1 based on PCA test, as shown in Fig. 10(a). After K-means clustering the latent representations, Fig. 10(b) depicts the results of SA and the input data are suggested to be separated into 3 clusters. Fig. 10(c) shows how the encoded representations in the latent space can be divided into 3 clusters. In fact, by drawing the output probability for each phase as a function of $\mu$, the phase transition points are somewhat deviated from theoretical values. This, however, can be largely improved if following the same supervised learning strategy as in the MCM case typeiii . This way not only shifts transition points back to mean values -0.995 and 1.003, but also reduces standard deviation from 0.335, 0.244 to 0.107, 0.067, respectively (statistically over 100 same-architecture CNN classifiers). The “labeled” training dataset used here is based on 3 seeds at $\mu/2t=-3.375,-0.063,2.892$, each of which gets the highest Silhouette value in the corresponding cluster. Note that, similar to the AE case, we again take the subsystem of size $L=100$ to refine the transition boundaries and to reduce the possible finite-size effect as well.

Similar to the 1D $p$-wave superconductor case, we first train a neural network to encode our type-I input data to latent representations. We determine the minimal dimension $d_{z}$ of the latent space by training a series of AEs with same model architecture appendix except for the number of hidden neurons in the middlemost layer ($n_{mid}$, from 2 to 10). In Fig. 11(a), we show the discrete distribution of $d_{z}$, which indicates the necessary dimension of the latent representations to keep at least $99\%$ variance of them for each $n_{mid}$ by PCA, after repeating 100 times of the same training procedure. Clearly, $d_{z}$ is suggested to be 4.

Once $d_{z}$ is known, we take the latent representations of all input data from previously trained AEs with $n_{mid}=d_{z}$ and do SA to estimate the optimal number of clusters $n$ via K-means method. As shown in Fig. 11(b), it points out that the mean of Silhouette value reaches the highest one when $n=2$. In addition, by projecting the 4D latent representations into a 2D space spanned by the first two principle components (features), we have Fig. 11(c). By noticing the density change of the feature points, the plot suggests how the system could be consistently divided into two clusters (phases).

Due to the unconstrained nature of the latent representations in AEs, the phase transition point found is statistically at mean value $v=0.987$ with standard deviations 0.207, after collecting clustering results from 100 sets of latent representations via different trained AEs ($n_{mid}=4$). Although such critical value is very close to the theoretical value 1, it still gets non-ignorable deviation. To reduce this variance, we again train a CNN classifier appendix in supervised learning manner to predict the whole phase diagram. We prepare a “labeled” dataset by picking up two seeds having highest Silhouette value in each cluster. They correspond to $v/w=0.332,5.799$ of the original input “images”, and for each point we expand symmetrically around its location by a window width $0.1w$ to obtain 2000 points with equal spacing. These 4000 data points are collected to be our training dataset, whereas the original 10,001 ones are formed our test set without labels. It turns out that the phase boundary predicted by trained CNN classifiers for the test set is clearly sharper at mean value 1.024 with much smaller deviation 0.041 (after repeating same training procedure 100 times), as shown in Fig. 11(d).

Following the insights obtained from the 1D $p$-wave SC case, we employ VAE approach to impose some constraints on the encoded representations for a more stable solution. As a first step, we train a series of VAEs with $n_{mid}=2$ to 10 pairs of hidden neurons in the middlemost layer appendix . Note that, similarly, one of the neurons in each pair outputs the mean of encoded representation distribution and the other produces its corresponding variance. After repeating the same training procedure 100 times for a given $n_{mid}$, the number of necessary dimension for the latent representations obtained after PCA to keep at least $99\%$ variance of them is simply one (pair), as clearly shown in Fig. 12(a). This again proves the stability of obtaining robust latent representations via VAE method.

Since the minimal dimension of the encoded representations is one (pair), we then employs SA to estimate the best number of clusters $n$ for them via K-means method. Fig. 12(b) clearly shows that the Silhouette value achieves the highest one when $n=2$. Moreover, Fig. 12(c) depicts the encoded representations as a function of $v/w$ and different colors indicate distinct clusters (phases). Finally, from our repeated clustering results, the phase transition point is statistically found to be at mean value 1.051 with standard deviations 0.069. This deviation is relatively smaller than the one obtained from the AE approach. Therefore, the supervised learning in the last step in Fig. 3(a) is neglected here.