# Neural Networks as Universal Probes of Many-Body Localization in Quantum
Graphs

Cameron Beetar<EMAIL_ADDRESS>Jeff Murugan<EMAIL_ADDRESS>The
Laboratory for Quantum Gravity & Strings,
Department of Mathematics & Applied Mathematics,
University of Cape Town, Rondebosch, Cape Town, South Africa  Dario Rosa
<EMAIL_ADDRESS>Center for Theoretical Physics of Complex Systems,
Institute for Basic Science (IBS),
Daejeon - 34126, Korea

###### Abstract

We show that a neural network, trained on the entanglement spectra of a
nearest neighbor Heisenberg chain in a random transverse magnetic field, can
be used to efficiently study the ergodic/many-body localized properties of a
number of other quantum systems, without further re-training. We benchmark our
computational architecture against a $J_{1}\\!\\!-\\!\\!J_{2}$ \- model, which
extends the Heisenberg chain to include next-to-nearest neighbor interactions,
with excellent agreement with known results. When applied to Hamiltonians that
differ from the training model by the topology of the underlying graph, the
neural network is able to predict the critical disorders and, more generally,
shapes of the mobility edge that agree with our heuristic expectation. We take
this as proof of principle that machine learning algorithms furnish a powerful
new set of computational tools to explore the many-body physics frontier.

††preprint: APS/123-QED

_I. Introduction.—_ Science covets universal tools. Calculus, statistical data
analysis and perturbation theory all fall into this category - powerful
conceptional and computational frameworks applicable across a wide spectrum of
disciplines. The advent of computers, and scientific computing in particular,
was another paradigm shift which reshaped the research landscape since the
1980’s. Today it is evident that we stand on the cusp of another such
revolution, as machine learning algorithms begin to make an impact on nearly
every facet of modern scientific research, from medical diagnostics to
mathematical theorem proving. The seemingly universal applicability, and
remarkable adaptability, of such algorithms makes it an ideal computational
tool to address problems at the frontier of numerical (and analytic)
accessibility in physics.

Many-body quantum physics is replete with such problems. Among these, the
localization properties of disordered quantum matter stand out, as much for
their complexity as for their broad relevance. It is for these precise
reasons, armed with new tools from machine learning, that we return to this
problem in the current article. This is, of course, not the first time that
neural networks have been turned loose on the problem of many-body
localization (MBL). In fact, our strategy for attacking this problem is based
largely on [1] in which the thermalization to localization transition - the
so-called mobility edge \- was computed for a disordered, nearest-neighbour
Heisenberg spin chain, as a function of a random magnetic field $\bar{h}$. The
computation of the mobility edge even in a toy model like the Heisenberg chain
is computationally demanding, requiring sophisticated numerical tools, large
system sizes (up to 24 sites) and many disorder realizations [2, 3]. In
comparison, even a single hidden-layer neural network trained on the
entanglement spectrum of the Heisenberg chain is able to reproduce all the
qualitative features of the MBL phase diagram with significantly smaller
system sizes (a 16-site spin chain) and substantially fewer disorder
realizations [1].

In this article, we will push this analysis even further and ask whether a
neural network, trained on the entanglement spectra of the same (nearest-
neighbour) disordered Heisenberg spin chain as in [1], is able to predict the
mobility edge for a wide range of Hamiltonians defined on topologically
different networks, without additional re-training. The quantum networks that
we consider include (i) next-to-nearest neighbour interactions whose graph
exhibits clustering, (ii) star graphs with a very small interaction length and
(iii) a superposition of nearest neighbour and star interactions whose graph,
with both simple clustering as well as a small diameter, is an example of a
scale-free network which is a special case of a small-world network [4]. This
is an (albeit limited) example of what is called transfer-learning in the
machine-learning literature - think of a child, having learnt how to walk
adapting those mechanical movements to learn how to run, again, without having
to re-train on a different set of movements.

In what follows, we will set up the physics problem we would like to study -
the disordered spin chain defined on quantum networks - in sections II and
III. In section IV we describe briefly the neural network architecture that we
use to probe the quantum systems in question, with our results described in
section V. Further details of the neural network, as well as comparisons of
our results to previous work can be found in the supplementary material.

_II. The ETH/MBL transition.—_ Otherwise thermalizing systems can sometimes
exhibit properties of localization due to disorder [5, 6, 7]. The standard
example for this phenomenon, by now well-understood and accepted, is the
Anderson model — a non-interacting free fermion model hopping on a lattice and
subject to an on-site, random, potential [8]. Depending on the strength of the
disorder, the energy eigenstates can be extended (for small enough disorder)
or localized, i.e. having non-vanishing support only on few sites. The
disorder strength at which the transition from extended to localized behavior
happens is usually called critical disorder and it is in general highly
dependent on the dimensionality on which the model is defined [9].

Many-body localization addresses the fate of the localized phase in presence
of interactions [10, 11, 12]. In this case, the problem becomes immediately
much more involved. For example, from a computational point of view, one has
to deal with an exponentially large (in the number of particles, $N$) Hilbert
space; the latter being an intrinsic property of any quantum many-body system.
Consequently, any exact numerical study is constrained to deal with systems
having just few degrees of freedom [13, 14, 15, 2, 3], thus making the
extrapolation of the finite $N$ results to the thermodynamic limit a
formidable task, as well as a matter of ongoing debate [16, 17, 18, 19, 20,
21, 22, 23, 24, 25].

Independently of the behavior in the thermodynamic limit, many finite $N$
studies of MBL systems have been already performed in the past. Among the
models which have been studied, the Heisenberg model in a random static
magnetic field has by now become the “gold standard” example of an MBL system
and, for this reason, this is the first model we discuss and it will
constitute the model that we will use to train our neural network. This step
will allow us to perform some sanity checks on our results against the amount
of results already at disposal in the literature, including the results of [1]
which uses the same neural network architecture.

The $N$-site Heisenberg model describes a set of $N$ spin-1/2 particles,
living on a chain with periodic boundary conditions, in presence of a random
static transverse magnetic field. The model is described by the Hamiltonian
[15]

$H=\frac{1}{4}\sum_{i,j}J_{ij}\bm{\sigma}_{i}\cdot\bm{\sigma}_{j}+\frac{1}{2}\sum_{i=1}^{N}h_{i}\sigma^{z}_{i},$
(1)

where the $\sigma^{z}_{i}$ corresponds to the $z$-component Pauli matrix at
the $i^{\text{th}}$ site, the $\bm{\sigma}_{i}\cdot\bm{\sigma}_{j}$ is the
scalar product between the Pauli matrices at the $i^{\text{th}}$ and
$j^{\text{th}}$ sites and the $h_{i}$ are the random static magnetic field
values, pseudo-randomly sampled from a uniform distribution
$h_{i}\in\left[-W,W\right]$, $W\in\mathbb{R}^{+}$ being the disorder strength.
In the Hamiltonian above, the interaction coupling matrix $J_{ij}$ can be
thought of as the adjacency matrix for a nearest-neighbour graph, on which the
model lives, and so it is defined to be $1$ when $i$ and $j$ are adjacent and
$0$ otherwise (i.e. $J_{ij}\to J\delta_{i,i+1}$). Since we have implicitly set
$\hbar=1$, we measure energy in units of $J\hbar$. As already mentioned, this
model has been extensively studied in recent years. In particular, while
exactly solvable for $W=0$, it turns out to satisfy ETH at small disorder,
while the eigenstates become localized at large values of the disorder. The
value of the disorder at which we measure the transition from extended to
localized behavior turns out to be energy-dependent, a feature that defines
the presence of a mobility edge [9, 2]. This means that, at intermediate
disorder strengths, regions of the spectrum in the extended regime coexist
with other regions (at different energies) in the localized regime. Given that
the model, at least at finite $N$, is by now well-understood, it constitutes a
perfect playground in which to train our neural network using a supervised
approach. Details of the training algorithm as well as a comparison to
existing predictions for the phase diagram can be found in the Appendix.

_III. Rewiring the Interactions.—_ Once the network is trained, and results
benchmarked against known results for the nearest-neighbour Heisenberg chain,
we turn our attention to the real problem of interest: systems with different
entanglement topologies as coded in the underlying interaction graph. In
particular, we will consider three classes of graphs, each encoding a specific
property that we would like to explore (see Figure 1):

The next-to-nearest neighbour Hamiltonian is specified by setting
$J_{ij}=J_{1}\,\delta_{i,i+1}+J_{2}\,\delta_{i,i+2}$ and introduces two
qualitatively new properties into the system111In the special case
$J_{1}=2J_{2}$ and $W=0$, this becomes the celebrated and exactly solvable
Majumdar-Ghosh model [32] for which the Hamiltonian is directly expressed in
terms of the quadratic Casimir of the 3-site spin algebra.; longer-range
interactions and clustering. This latter property is a feature of any graph
with a large fraction of triangles present and can be quantified by a
clustering coefficient which, at site $i$ is counted as
$C_{i}=2\Delta_{i}/(k_{i}(k_{1}-1))$, if the site had degree $k_{i}$ and has
$\Delta_{i}$ triangles attached to it. The corresponding graph-averaged
clustering coefficient is then $\overline{C}=\sum_{i}C_{i}/N$. In addition,
numerical studies for the ETH/MBL transition for this model have been already
obtained [16], thus making it as a control case to check the performance of
the NN on new models.

The $N$-site star graph in which one site, say $\alpha$, is singled out and
connected to each of the remaining $N-1$ sites with no other edges, so that
$J_{ij}=J_{\alpha}\delta_{\alpha i}$. The star graph is an example of a
complete bipartite graph with paths in the graph having either length 1 (for
sites connecting to $\alpha$) or 2 (for sites connected via $\alpha$). The
resulting average distance in the network
$\bar{l}=(1/N)\sum_{i}\bar{l}_{i}=2-4/N+3/N^{2}$ is exact, and
characteristically small.

The bicycle-wheel graph defined by the adjacency matrix
$J_{ij}=J_{1}\delta_{i\neq\alpha,i+1}+J_{3}\delta_{\alpha i}$ superposes the
nearest-neighbour Hamiltonian (1) with the star graph and is actually a simple
example of a broader class of small-world graphs [4] that exhibit both
clustering (due to the graph transitivity) as well as a very short average
path length (as a result of shortcuts through the hub). Here again, the graph,
as a representative of its class, is sufficiently simple that the network
characteristics can be computed exactly. Specifically, the average clustering
coefficient is $\overline{C}=2/3-2/3N+2/(N^{2}-N)$, while the average path
length reads $\bar{l}=2-6/N+5/N^{2}$.

In each case, we will retain the random static transverse magnetic field term
$\sim\sum h_{i}\sigma_{i}^{z}$ that disorders the system. In this way, we are
able to investigate the competing effects of disorder which tends to promote
many-body localization and interactions which favour thermalization while
keeping the system size $N$ within a computationally managable range.

Figure 1: Interaction graphs for an 8-site system. Top: Nearest-neighbour
interactions on which the neural network is trained. Bottom, from left to
right: superposition of star and nearest-neighbour interactions; a star graph
and next to nearest-neighbour circulate graph. The central nodes in the lower-
left images corresponds to the $\alpha$-site in the defining adjacency
matrices.

_IV. Binary Classification Neural Network.—_ Distinguishing between ETH and
MBL regimes in these many-body systems is essentially a classification
problem. As such, the mobility edge that forms the boundary between these two
phases corresponds to a band of classification uncertainty in the algorithm.
To see how this works, we will adopt a binary classification neural network
scheme similar to that in [1]. In short, the neural network maps an
$n$-dimensional vector, consisting of the entanglement spectrum at a given
eigenenergy and a specified disorder, to a 2-dimensional vector that
quantifies how confident the network is that the entanglement spectrum
corresponds to an ETH or localizing phase. We outline here some of the key
features of our neural network, and refer the reader to the Supplementary
Material for further elaboration.

A defining feature of artificial neural networks (ANN) are their layers.
Mathematically speaking, a layer consists of a linear map and a nonlinear
activation function. The layer is called hidden if it is not the final layer
of the network. For the purposes of this work, we use a single hidden layer
ANN (see figure 2 for a visual representation).

Figure 2: A visualisation of the neural network, and its layers, associated
with $N=16$ systems. Numbers above the directed edges correspond to the number
of neurons between each transformation. Icons labelled 1 and 3 correspond to
linear transformations while icons labelled 2 and 4 correspond to nonlinear
activation functions.

In this context, the dimension of the output of a given layer corresponds to
the number of neurons in that layer. Specifically, if there are $n^{(i)}$
neurons at layer $i$ in the network, the linear map from the $i^{\text{th}}$
to the $(i+1)^{\text{th}}$ layer is specified by the transformation
$\bm{v}^{(i+1)}=M^{(i+1,i)}\bm{v}^{(i)}+\bm{b}^{(i+1)},$ with weight matrix
$M\in\mathbb{R}^{n^{(i+1)}\times n^{(i)}}$ and bias vector
$\bm{b}\in\mathbb{R}^{n^{(i+1)}}$.

The two activation functions used in our network are an element-wise SoftSign
function in the hidden layer $\text{SS}_{i}(\bm{v})\equiv v_{i}/(1+|v_{i}|)$,
where $v_{i}$ is the $i^{\text{th}}$ component of the Real vector $\bm{v}$,
and the Softmax function $\text{SM}_{i}(\bm{v})\equiv
e^{-v_{i}}/\sum_{j}e^{-v_{j}}$, in the output layer. In this way, the neural
network can be summarised as the map

$\begin{split}&f_{\text{net}}\text{:
}\mathbb{R}^{n^{(1)}}\rightarrow\mathbb{R}^{n^{(3)}};\quad\bm{v}^{(1)}\mapsto\bm{v}^{(3)},\\\
&f_{\text{net}}\left(\bm{v}^{(1)}\right)=\text{SM}\left[M^{(3,2)}\text{SS}\left(M^{(2,1)}\bm{v}^{(1)}+\bm{b}^{(2)}\right)+\bm{b}^{(3)}\right]\,.\end{split}$

The choice $M^{(3,2)}\in\mathbb{R}^{2\times n^{(2)}}$ ensures that the outputs
are 2-dimensional vectors. The Softmax function, in turn, converts this into a
vector $(v^{(3)}_{1},v^{(3)}_{2})^{\mathrm{T}}$ with
$v^{(3)}_{1}+v^{(3)}_{2}=1$. Here, $v^{(3)}_{1}$ corresponds to the network
confidence that the input is in an ETH phase, while $v^{(3)}_{2}$ gives the
corresponding confidence level for the MBL phase.

In training the network on the nearest-neighbour model, we assume that systems
with disorder sampled from $[-0.25,0.25]$ obey the ETH and map to $(1,0)$,
while those with disorder sampled from $[-12,12]$ are in the MBL phase and map
to $(0,1)$. These values are known to be well-inside the ETH (MBL) regions and
we train our neural network at these disorder samplings, only. This scheme,
arising from a training method that includes cross-validation, dropout
regularisation, weight-decay, and lack-of-confidence penalties (see
Supplementary Material), produces a neural network that can effectively
classify entanglement spectra as being in the ETH or MBL phase, while
systematically avoiding both over-fitting and classification agnosticism.

_V. Results.—_ To produce phase diagrams for a given system, we consider 48
unique disorder realisations at $W=1.0$ and scale the disorder at each site,
for each realisation, by multiplying their values by the disorder value of
interest. The entanglement spectra of each realisation and each eigenenergy is
classified by the neural network at interval disorder values and these results
are averaged either by grouped eigenenergy index, or by a shifted-normalised
eigenenergy binning. However, as we focussed on Hamiltonians with interaction
graphs distinct from nearest-neighbour topology, we adopt the former since it
is robust against degenerate eigenenergies. The resulting phase diagrams are
presented in figure 3(a). We note that the shape of the mobility edge in the
nearest-neighbour case is in excellent agreement with existing literature for
this case. This is particularly striking when considering the shifted-
normalised eigenenergy binning [1].

We find that the ANN is also able to probe the phase boundary. To do so, we
make the assumption that all data having MBL-classification confidence $p\geq
0.9$ ($p\leq 0.1$) are MBL (ETH) spectra, and are assigned a value of 0. On
the other hand, spectra with $0.1<p<0.9$ are transition spectra, and are
assigned a value of 1 [1]. With these criteria and the same averaging over 48
realisations as was used in producing the phase diagram, we can effectively
isolate the mobility edge between the ETH and MBL phases (see Supplementary
Material).

(a)

(b)

(c)

(d)

Figure 3: Phase diagrams for 48 disorder realisations of $N=12$ site systems
with different graph topologies and disorder steps $\Delta W=0.01$. All data
is generated from the same 48 disorder realisations. In all diagrams, the data
is grouped into bins of 8 consecutive energy eigenvalues and averaged inside
each bin. We choose this binning over the $\epsilon$-binning as we do not
claim to know anything about the distribution of the energy eigenvalues of all
the topologies. The phase diagrams display the ETH/MBL transition for (a) the
nearest-neighbour system (b) the next-to-nearest-neighbour system, (c) the
star graph, and (d) a superposition of the nearest-neighbour and star graphs.
Dashed lines in (a) and (b) approximately correspond to regions exhibiting a
sharp drop in the ratio of the system size, $N$, to the BKT correlation length
$\xi_{\text{BKT}}$ at $N/\xi_{\text{BKT}}\approx 0$ [16].

We are now in a position to return to the question posed in the Introduction
and turn our network, trained on the disordered nearest-neighbour chain, loose
on more complex interaction graphs to see if is able to classify entanglement
spectra in systems for which the neural network has had no training. As
motivated in section III, we consider three prototypical models, specified by
their adjacency matrices, as representatives of broader classes of quantum
networks. Applying precisely the same procedure as was done to produce figure
3(a), the network produces phase diagrams 3(b)-3(d). At this point, there are
several key observations to make:

—As the connectivity between sites is increased, we observe that the ETH
region is ’pushed out’ to larger disorder values. Heuristically, this is to be
expected, as the coupling to the part of the system acting as the thermal bath
is increased with increased connectivity. Along these lines, there are studies
questioning and debating the existence of the MBL phase for long-range
interactions and so our results sound very plausible [27, 28, 29]. There are
also direct comparisons, [16], of the Heisenberg and $J_{1}\\!\\!-\\!\\!J_{2}$
models whose breakdowns of the ETH appear to agree with the results found here
(see dashed regions in figures 3(a) and 3(b)), and these provide justification
of the scaling of the ETH region found by the neural network in [1].

—Finally, in the case of the star graph (figure 3(c)) and the superposition of
nearest-neighbour and star networks (figure 3(d)) we see further expansion of
the ETH region, as expected. Interestingly, the transition region is
essentially the same in the two graphs, thus suggesting that connectivity is
much more prominent than clustering in determining the MBL transition.
Additionally, in the system with star graph connectivity, eigenenergy
degeneracy that we had not anticipated was found at very-low disorder values
and, as there is no ‘canonical’ ordering of these degenerate eigenvalues, the
entanglement spectrum is ill-defined. This manifests in the appearance of
small ‘uncertainty bubbles’ in the very low disorder region of figure 3(c).
This degeneracy is removed by nearest-neighbour couplings as evidenced by the
lack of these uncertainty bubbles in figure 3(d).

_VI. Conclusions and Discussion —_ We have shown that a neural network,
trained to recognize the ETH/MBL transition by means of the entanglement
spectra in a prototype example, can efficiently detect the same kind of
transition on new models without further retraining.

This result is particularly relevant in promoting neural networks to the
status of effective tools in investigations of MBL/ETH transitions. Indeed,
confronted by a new model, one often does not know a disorder strength for
which the model can be safely considered to be in the ETH (MBL) phase,
hindering both training as well as subsequent results obtained by the neural
network.

Together with the results of [1], showing that a simple neural network detects
the qualitative features of an MBL phase diagram of a model used to train the
network; we have shown here that ANNs should be considered a powerful tool to
perform preliminary investigations of the phase diagram of a given quantum
many-body system. These investigations can (and should [30]) be confirmed, and
made more precise, at a later stage by more conventional (and computationally
demanding) numerical tools [2, 3].

More speculatively, one can ask what is the ultimate reason for the excellent
performances of neural networks on models different from the training model.
We believe that this is yet another manifestation of random matrix theory
universality: in the ETH phase, typical eigenstates are supposed to be almost
random while in the MBL phase the eigenstates should exhibit highly model-
dependent features. Hence, we believe that the neural network learns how to
distinguish universal from non-universal features and this should be the key
to its performance while changing models.

Our work can be extended in several ways. Allowing that it is an excellent
tool for preliminary investigations, a natural question is whether it is
possible to feed the neural network with simpler data than the entanglement
spectrum, while retaining its efficiency in new models. For example, how does
the network perform working with the energy spectrum? [31]. At the physics
level, our results show that the connectivity of the graph strongly affects
the transition region. It would be important to make more robust and precise
these findings, through more conventional techniques. We plan to address all
these points in the near future.

_Acknowledgments.—_ DR acknowledges the support by the Institute for Basic
Science in Korea (IBS-R024-Y2 and IBS-R024-D1). JM is supported in part by a
Simons Associateship at the ICTP, Trieste. CB acknowledges support from the
South African Research Chairs Initiative of the NRF. Computations were
performed using facilities provided by the University of Cape Town’s ICTS High
Performance Computing team. All graphics were produced by Wolfram Mathematica
version 12.0.0.

## References

  * Schindler  _et al._ [2017] F. Schindler, N. Regnault, and T. Neupert, Probing many-body localization with neural networks, Phys. Rev. B 95, 245134 (2017).
  * Luitz  _et al._ [2015] D. J. Luitz, N. Laflorencie, and F. Alet, Many-body localization edge in the random-field heisenberg chain, Phys. Rev. B 91, 081103 (2015).
  * Sierant  _et al._ [2020] P. Sierant, M. Lewenstein, and J. Zakrzewski, Polynomially filtered exact diagonalization approach to many-body localization, Phys. Rev. Lett. 125, 156601 (2020).
  * Watts and Strogatz [1998] D. J. Watts and S. H. Strogatz, Collective dynamics of ‘small-world’ networks, Nature 393, 440 (1998).
  * Nandkishore and Huse [2015] R. Nandkishore and D. A. Huse, Many-body localization and thermalization in quantum statistical mechanics, Annual Review of Condensed Matter Physics 6, 15 (2015).
  * Alet and Laflorencie [2018] F. Alet and N. Laflorencie, Many-body localization: An introduction and selected topics, Comptes Rendus Physique 19, 498 (2018).
  * Abanin  _et al._ [2019] D. A. Abanin, E. Altman, I. Bloch, and M. Serbyn, Colloquium: Many-body localization, thermalization, and entanglement, Rev. Mod. Phys. 91, 021001 (2019).
  * Anderson [1958] P. W. Anderson, Absence of diffusion in certain random lattices, Phys. Rev. 109, 1492 (1958).
  * Abrahams [2010] E. Abrahams,  _50 Years of Anderson Localization_ (WORLD SCIENTIFIC, 2010).
  * Gornyi  _et al._ [2005] I. V. Gornyi, A. D. Mirlin, and D. G. Polyakov, Interacting electrons in disordered wires: Anderson localization and low-$t$ transport, Phys. Rev. Lett. 95, 206603 (2005).
  * Basko  _et al._ [2006] D. Basko, I. Aleiner, and B. Altshuler, Metal–insulator transition in a weakly interacting many-electron system with localized single-particle states, Annals of Physics 321, 1126 (2006).
  * Basko  _et al._ [2007] D. M. Basko, I. L. Aleiner, and B. L. Altshuler, Possible experimental manifestations of the many-body localization, Phys. Rev. B 76, 052203 (2007).
  * Oganesyan and Huse [2007] V. Oganesyan and D. A. Huse, Localization of interacting fermions at high temperature, Phys. Rev. B 75, 155111 (2007).
  * Žnidarič  _et al._ [2008] M. Žnidarič, T. c. v. Prosen, and P. Prelovšek, Many-body localization in the heisenberg $xxz$ magnet in a random field, Phys. Rev. B 77, 064426 (2008).
  * Pal and Huse [2010] A. Pal and D. A. Huse, Many-body localization phase transition, Phys. Rev. B 82, 174411 (2010).
  * Šuntajs  _et al._ [2020a] J. Šuntajs, J. Bonča, T. c. v. Prosen, and L. Vidmar, Quantum chaos challenges many-body localization, Phys. Rev. E 102, 062144 (2020a).
  * Abanin  _et al._ [2021] D. Abanin, J. Bardarson, G. De Tomasi, S. Gopalakrishnan, V. Khemani, S. Parameswaran, F. Pollmann, A. Potter, M. Serbyn, and R. Vasseur, Distinguishing localization from chaos: Challenges in finite-size systems, Annals of Physics 427, 168415 (2021).
  * Šuntajs  _et al._ [2020b] J. Šuntajs, J. Bonča, T. c. v. Prosen, and L. Vidmar, Ergodicity breaking transition in finite disordered spin chains, Phys. Rev. B 102, 064207 (2020b).
  * Corps  _et al._ [2020] A. L. Corps, R. A. Molina, and A. Relaño, Thouless energy challenges thermalization on the ergodic side of the many-body localization transition, Phys. Rev. B 102, 014201 (2020).
  * Panda  _et al._ [2020] R. K. Panda, A. Scardicchio, M. Schulz, S. R. Taylor, and M. Žnidarič, Can we study the many-body localisation transition?, EPL (Europhysics Letters) 128, 67003 (2020).
  * Luitz and Lev [2020] D. J. Luitz and Y. B. Lev, Absence of slow particle transport in the many-body localized phase, Phys. Rev. B 102, 100202 (2020).
  * Morningstar  _et al._ [2020] A. Morningstar, D. A. Huse, and J. Z. Imbrie, Many-body localization near the critical point, Phys. Rev. B 102, 125134 (2020).
  * Morningstar  _et al._ [2021] A. Morningstar, L. Colmenarez, V. Khemani, D. J. Luitz, and D. A. Huse, Avalanches and many-body resonances in many-body localized systems (2021), arXiv:2107.05642 [cond-mat.dis-nn] .
  * De Roeck  _et al._ [2016] W. De Roeck, F. Huveneers, M. Müller, and M. Schiulaz, Absence of many-body mobility edges, Phys. Rev. B 93, 014203 (2016).
  * Brighi  _et al._ [2020] P. Brighi, D. A. Abanin, and M. Serbyn, Stability of mobility edges in disordered interacting systems, Phys. Rev. B 102, 060202 (2020).
  * Note [1] In the special case $J_{1}=2J_{2}$ and $W=0$, this becomes the celebrated and exactly solvable Majumdar-Ghosh model [32] for which the Hamiltonian is directly expressed in terms of the quadratic Casimir of the 3-site spin algebra.
  * Nag and Garg [2019] S. Nag and A. Garg, Many-body localization in the presence of long-range interactions and long-range hopping, Phys. Rev. B 99, 224203 (2019).
  * Nandkishore and Sondhi [2017] R. M. Nandkishore and S. L. Sondhi, Many-body localization with long-range interactions, Phys. Rev. X 7, 041021 (2017).
  * Yao  _et al._ [2014] N. Y. Yao, C. R. Laumann, S. Gopalakrishnan, M. Knap, M. Müller, E. A. Demler, and M. D. Lukin, Many-body localization in dipolar systems, Phys. Rev. Lett. 113, 243002 (2014).
  * Théveniaut and Alet [2019] H. Théveniaut and F. Alet, Neural network setups for a precise detection of the many-body localization transition: Finite-size scaling and limitations, Phys. Rev. B 100, 224202 (2019).
  * Kausar  _et al._ [2020] R. Kausar, W.-J. Rao, and X. Wan, Learning what a machine learns in a many-body localization transition, Journal of Physics: Condensed Matter 32, 415605 (2020).
  * Majumdar and Ghosh [1969] C. K. Majumdar and D. K. Ghosh, On Next-Nearest-Neighbor Interaction in Linear Chain. I, Journal of Mathematical Physics 10, 1388 (1969).
  * Li and Haldane [2008] H. Li and F. D. M. Haldane, Entanglement spectrum as a generalization of entanglement entropy: Identification of topological order in non-abelian fractional quantum hall effect states, Phys. Rev. Lett. 101, 010504 (2008).
  * Mehta  _et al._ [2019] P. Mehta, M. Bukov, C.-H. Wang, A. G. Day, C. Richardson, C. K. Fisher, and D. J. Schwab, A high-bias, low-variance introduction to machine learning for physicists, Physics Reports 810, 1 (2019).
  * Shannon [1948a] C. E. Shannon, A mathematical theory of communication, The Bell System Technical Journal 27, 379 (1948a).
  * Shannon [1948b] C. E. Shannon, A mathematical theory of communication, The Bell System Technical Journal 27, 623 (1948b).

##  Appendix A Supplementary Material for “Neural Networks as Universal Probes
of Many-Body Localization in Quantum Graphs”

Appreciating that the machine learning tools that are at the core of this
letter may not be entirely familiar to the bulk of its readership, and in the
interest of self-containment, here we will provide an extended discussion of
our neural network - its construction, architecture and training, as well as
some of the terminology used - together with some of the preliminary results
and checks that we undertook in order to increase our confidence in the new
results reported in the letter.

###  A.1 I. Extended Neural Network Discussion

####  A.1.1 Input Data

An $N$-site system of spin-1/2 particles has $2^{N}$ eigenstates. We will
consider only those eigenstates with no net spin; i.e. $S_{\text{tot}}^{z}=0$.
Consequently, we will enforce only even values for $N$ so that the number of
energy eigenstates under consideration is $N\choose N/2$. When partitioning
the system, we only ever split the system into equal halves
($N_{A}=N/2=N_{B}$), and trace out the $B$-subsystem. The constraint of zero
net spin then enforces that $N$ be divisible by 4. Following [1], we also
train the network on entanglement data. The entanglement (or modular)
Hamiltonian [33], $H_{e}$, is defined through the reduced density matrix,
$\rho_{A}$, obtained by tracing out the $B$-subsystem, as

$e^{-H_{e}}\equiv\rho_{A}=\text{Tr}_{B}\ket{\Psi}\bra{\Psi}.$ (1)

In this context, each energy eigenvalue is associated to an entanglement
spectrum of size $N/2\choose N/4$. Additionally, we will use only the smallest
half of these values, since larger values give vanishingly small contributions
to the negative-exponential. Thus, as input for our network we always take an
${N/2\choose N/4}/2$ length vector that corresponds to the smallest half of
the entanglement spectrum associated with a given eigenenergy (which, itself,
will vary with the strength of the disorder). Hereafter we use the terms
“input vector” and “entanglement spectrum” interchangeably. The final
limitation on the data used in training is that we consider only the central
80% of eigenergies for training; we remove the lower and upper 10% of
eigenergies as they can have behaviour that strongly deviates from that of the
central region, as is evidenced by figure 3(a) in the main text, where we see
that, for low disorder, the low energy states can exhibit MBL-like behavior,
[1].

####  A.1.2 Training Network Architecture

The training network architecture is a bit more complicated than that of the
extracted network map described in the letter. To elaborate on this, we will
make use of figure 1 which describes a training methodology built on four key
attributes: a confidence penalty, dropout regularization, weight decay, and
cross validation.

Figure 1: The full network structure used for training for $N=16$ input sizes.
Notice that there are two inputs (far left), two outputs (far right, middle),
and two cost functions (far right, top and bottom). Of the larger nodes, nodes
4-7 constitute the part of the network we extract once training is complete.

We will explicitly consider the neural network used to classify an $N=16$ site
system. The only difference between this network and the $N=12$ site network
is in the input sizes - the same number of neurons are used at the hidden
layer and thereafter.

Confidence Optimization: For a $N=16$ system, the lower half of the
entanglement spectrum (the network input) has size ${8\choose 4}/2=35$,
whereas the training network in 1 has two of these as inputs. This is because,
at every training step, two input vectors are given - the first being taken
from the disorder region of assumed phase (either $0.25\rightarrow\text{ ETH}$
or $12.0\rightarrow\text{ MBL}$) which is accompanied by its target vector
(being $(1,0)$ and $(0,1)$, respectively) which is called supervised input,
and the second being a vector from a region of intermediate disorder
$W\in(0.25,12.0)$ where no target vector is specified which is called
unsupervised input. The first two larger nodes in 1 simply reshape this data
into a form the network can process. At the remaining large nodes are
identical copies of the network through which both inputs are modified. For
the supervised input we use the conventional cross-entropy cost function [34],

$\text{Sup.
Cost}(\hat{f},f_{\text{net}})=-\sum_{\bm{v}\in\text{TD}}\sum_{i=1}^{2}f_{\text{net},i}(\bm{v}_{i})\ln\hat{f}(\bm{v}_{i}),$
(2)

where $f_{\text{net}}$ is the network map, $\hat{f}$ is the target map and the
first sum is over the input vectors in the $W=0.25$ or $W=12.0$ regions
(‘Targeted Data’). The unsupervised input, not having a target vector to take
the cross-entropy with is used to enforce the confidence penalty by noticing
that the function $C(x)=x\log x$, with domain $x\in[0,1]$ has global minima at
values of 0 and 1. Thus we use the Shannon entropy as a cost function for the
unsupervised data [35, 36]

$\text{Unsup.
Cost}(f_{\text{net}})=-\sum_{\bm{v}\in\text{ID}}\sum_{i=1}^{2}f_{\text{net},i}(\bm{v}_{i})\ln
f_{\text{net},i}(\bm{v}_{i}),$ (3)

where the first sum is over the ‘Intermediate disorder Data’. The inclusion of
this cost function penalizes the neural network for finding values farther
from 0 and 1, pushing the network towards being as confident as possible in
the intermediate disorder region. One could choose to moderate the influence
of this penalty by multiplying through (3) by a hyperparameter, but
investigation finds that this value is optimal for values near 1, and so we do
not explicityly investigate the effects of this hyperparameter here [1].
Considering a cost function that combines (2) and (3), such a cost function
would push the netword towards classification of (assumed) definite-phase
states and unknown-phase states with minimal agnosticism.

Dropout Regularization: To avoid overweighting particular neurons, we make use
of dropout regularization; at every training step some of the neurons are not
trained at all. We do not enforce that a specific number of neurons not be
evaluated at each step, rather we enforce that every neuron has a probability
of 0.5 of being dropped at every training step.

Weight Decay: At any point during training the weights within a linear map may
be sporadically transformed into nonzero values to produce anomalous
overfittings. To avoid this, we implement weight decay of the transformation
matrix in the hidden layer. To be precise, we make use of a stochastic
gradient descent training algorithm that updates the weights and biases of the
linear layers at each training step by randomly selecting a subset of the
training data and minimizing the error on that subset of the training data,
that is we minimize the error (cost) on

$A_{i+1}=A_{i}-\lambda\frac{\partial f_{\text{cost}}}{\partial A_{i}},$ (4)

where $A_{i}$ are the weights and biases at training step $i$, $\lambda$ is
the learning rate hyperparameter, and $f_{\text{cost}}$ is the cost function,
for a subset of the total training data. The effect of including weight decay
is to add a decay term $-\mu A_{n}$ to (4), where $\mu$ is the weight decay
hyperparameter. The inclusion of this term in (4) sees a corresponding change
to the cost function, where we must also now include a $+\mu|A|^{2}$ (for the
$l_{2}$-norm $|\cdot|$). We choose to apply this weight-decay procedure only
to the hidden layer weight matrix. We used values of $\lambda=10^{-4}$ and
$\mu=1$ for training, as is done in [1].

Cross-Validation: Finally, before training starts all the training data is
randomly partitioned into two sets if equal cardinality. One set is used for
training, while the other (validation) set is used to ensure overfitting is
minimised by comparing how the updated network affects the error measured on a
validation set - a general neural network should be attempt to minimize the
error on both the training data and the validation data.

Including all considerations thus far, the final cost function is

$\begin{split}\text{Cost}&(\hat{f},f_{\text{net}})=-\sum_{\bm{v}\in\text{TD}}\sum_{i=1}^{2}f_{\text{net},i}(\bm{v}_{i})\ln\hat{f}(\bm{v}_{i})\\\
&-\sum_{\bm{v}\in\text{ID}}\sum_{i=1}^{2}f_{\text{net},i}(\bm{v}_{i})\ln
f_{\text{net},i}(\bm{v}_{i})+\mu|M|^{2},\end{split}$ (5)

where $M$ is the weight matrix of the hidden layer, and this is the cost
function used while training.

The training method (stochastic gradient descent) and associated additional
techniques (cross-validation etc,) described above are in no way unique to
this work and there exist numerous resources. We point the reader to [34] and
references therein. Furthermore, the reader may note a strong similarity with
this work and the work done in [1]; this is intentional. We, however, make the
use of the SoftSign activation function in the hidden layer instead of the
Rectified Linear Unit activation function, as empirical investigation showed
that the SoftSign function produced the best classification confidence and a
minimal transition region (thinnest red region in, for instance, figure 2(b)).

####  A.1.3 Training Results

The results of the trained network are presented in table 1. These results
characterize the ability of the networks to recognize ETH/MBL type systems,
with greater accuracy achieved when greater system size is used.

Table 1: The trained network’s classification results. System size refers to the number of sites the system has (i.e.the value of $N$). The data type specifies for which sort of training data we are calculating statistics; The mean is the average confidence of the network in the data being of the type specified by data type, and the standard deviation is the standard deviation of the mean. The data set sizes before partitioning for cross-validation were $37050$ for $N=12$ and $514850$ for $N=16$ (excluding intermediate disorder data).  System Size | Data Type | Mean (%) | Std Dev. (%)
---|---|---|---
12 | ETH | $99.35$ | $1.3\times 10^{-1}$
12 | MBL | $99.99999$ | $3.0\times 10^{-5}$
16 | ETH | $99.933$ | $1.4\times 10^{-2}$
16 | MBL | $100.00$ | $<10^{-10}$

(a)

(b)

Figure 2: Phase diagrams and transition region probing (using the method
described in the supplementary text) for different system sizes averaged over
48 disorder realisations for nearest-neighbour systems. The $N=12$ site system
has values computed in disorder steps of $\Delta W=0.01$; the $N=16$ site
system has $\Delta W=0.125$. All figures run from $W=0.25$ to $W=4.0$. (a) The
phase diagrams for an $N=12$ site system (left) and an $N=16$ site system
(right). (b) The associated transition region probes associated with the phase
diagrams in (a).

In addition to the results specified in table 1, we also mention the least
confident classification value for $N=12$: 96.9428% in the ETH regime and
99.9996% in the MBL regime. Similarly, these worst confidence values for
$N=16$ were 99.78% and 100.00%, respectively. Assuming all states with
classification confidence $90\%<$ are in that state, the network classifies
100% of the training and validation data successfully post-training. These
results also show that the networks are particularly good at classifying MBL
phases.

###  A.2 II. Extended Results

In the paper, we make use of the phase diagrams found by averaging over 48
unique realisations of an $N=12$ site system, and further grouping and
averaging over 8 consecutive eigenvalues (see figure 3). While this
methodology is essential for the probing of new topologies whose eigenenergy
structure and spacing is not assumed to be known, an alternative methodology
that may be used is grouping by binned and shifted-normalised eigenergies,

$\epsilon=(E-E_{\text{min}})/(E_{\text{max}}-E_{\text{min}}),$ (6)

where $E$ is a given eigenenergy, and $E_{\text{min(max)}}$ is the
minimum(maximum) eigenenergy for a given disorder realisation at a given
disorder value. The result of using this grouping methodology is presented in
figure 2(a).

As we are interested in the phase transition from ETH to MBL, it is
instructive to highlight the region in which the network struggles to
confidently classify entanglement spectra. To do this, we produce and classify
48 unique disorder realisations over the range $W\in[0.25,12]$, as before.
However, once the network has classified the data, we assume all output
2-vectors having first entry $v_{1}<0.1$ are in the MBL phase and are set to
0, all $v_{1}>0.9$ are in the ETH phase and are also set to 0, and the
remaining data ($0.1\leq v_{1}\leq 0.9$) are in the transition region and are
set to 1. In this way we produce a heat map of the transition region. The
process of averaging over realisations eigenergies (either by eigenenergy-
grouping or $\epsilon$-binning) gives rise to yellow-tinted regions, as shown
in figure 2(b).

In particular, the shape of the $N=16$ phase diagram should be compared to [1]
and [6] after a purely superficial change in the aspect ratio of the image;
where there is notable similarity in the both cases up to a re-scaling of the
critical disorder values.

We can probe the transition region in all topologies using either averaging
methodology. The transition region probe for the phase diagrams in the paper
(see figure 3 in the paper) is presented in figure 3. The small ‘bubbles’ of
high uncertainty at low disorder in figure 3(c) may be of interest to the
reader. These bubbles of uncertainty occur as a result of an energy degeneracy
in the star system. As there is not a canonical way of ordering these
degenerate eigenvalues, the notion of an entanglement spectrum becomes ill-
defined and hence it does not make sense to try classify this region.
Nonetheless, we include it so that the various topologies may be compared
directly.

(a)

(b)

(c)

(d)

Figure 3: The associated transition regions to the phase diagrams in figure 3
in the paper for 48 realisations of $N=12$ site systems having different
topologies. In all diagrams, we group data into bins of 8 consecutive energy
eigenvalues and average the data inside each bin.(a) Results for the nearest-
neighbour system. (b) Results for the next-to-nearest-neighbour system. (c)
Results for the star system. (d) Results for the bicycle-wheel system
(superposition of nearest-neighbour and star systems). The dashed lines in (a)
and (b) approximately correspond to regions in which the there is a sharp drop
in the ratio of the system size, $N$, to the Berezinskii-Kosterlitz-Thouless
correlation length $\xi_{\text{BKT}}$ at $N/\xi_{\text{BKT}}\approx 0$ [16].