Monte Carlo Methods II: Ising Model and Hopfield Network#
Introduction#
In 2024, the Nobel Prize in Physics was awarded to John J. Hopfield and Geoffrey E. Hinton for their pioneering work on neural networks. Their contributions have played a significant role in the development of modern artificial intelligence and computational models inspired by biological learning mechanisms.
John J. Hopfield introduced the Hopfield network, a type of recurrent neural network that models associative memory and collective computation in neural systems.
Geoffrey E. Hinton made significant advancements in deep learning, particularly in training deep neural networks through the backpropagation algorithm.
Their research has had profound implications in artificial intelligence, neuroscience, and computational physics.
Why Do Neural Networks Matter?#
Neural networks are now an essential part of many everyday technologies:
Facial Recognition: Used to unlock smartphones and enhance security.
Chatbots: Power customer service agents that interact with users.
Voice Assistants: Siri, Alexa, and Google Assistant rely on neural networks for speech recognition.
Self-Driving Cars: Autonomous vehicles process sensor data using deep learning models.
Importance of This Discovery#
Neural networks have revolutionized multiple domains, including:
Language Translation: Real-time translation services are powered by deep learning models.
Medical Diagnosis: AI-assisted analysis of medical images improves diagnostic accuracy.
Scientific Research: AI-driven simulations enhance modeling in physics, chemistry, and biology.
The Ising Model#
To understand how Hopfield networks operate, we first introduce the Ising model.
The Ising model is a mathematical model in statistical mechanics, introduced by Wilhelm Lenz in 1920 and solved for the one-dimensional case by his student Ernst Ising in 1925. The model was originally used to explain ferromagnetism, where magnetic materials exhibit spontaneous magnetization due to interactions between neighboring atomic spins.
Wilhelm Lenz conceived the model as a simplified representation of magnetic interactions in a lattice, where spins can either point “up” (+1) or “down” (-1).
Ernst Ising solved the one-dimensional version of the model in his doctoral thesis, showing that it did not exhibit phase transitions—a result that was surprising at the time.
Lars Onsager solved the two-dimensional version of the model in 1944, demonstrating that it undergoes a phase transition at a critical temperature, where spontaneous magnetization occurs.
Since then, the Ising model has become one of the most widely studied models in statistical physics and beyond. Its applications extend not only to physics (such as in magnetism and lattice gases) but also to fields like biology (neural networks), computer science (optimization problems), and even sociology (modeling opinion dynamics).
The Metropolis algorithm (1953) was developed for Monte Carlo simulations of systems like the Ising model, enabling the study of large, complex systems by simulating their thermal fluctuations and statistical properties. This method revolutionized computational physics and remains a powerful tool in many areas of research today.
Basic Idea#
Consider a grid (lattice) where each point represents a spin.
Each spin can take one of two states: up (+1) or down (-1).
Spins interact with their neighbors, aiming to align in a way that minimizes system energy.
The system’s energy parameter determines the likelihood of spins flipping, making temperature a crucial factor in system behavior.
Let’s start programming as we introduce the different concepts. Here, we create a python class to keep track of the state and temperature of an Ising Model.
# Subclass a numpy array to store the state of an Ising Model
import numpy as np
from numpy.random import randint
class Ising(np.ndarray):
def __new__(cls, shape=(64,64), kT=1, seed=None):
if seed is not None:
np.random.seed(seed)
obj = np.random.choice([-1,1], size=shape).view(cls)
obj.kT = kT
return obj
def __array_finalize__(self, obj):
if obj is None: return
self.kT = getattr(obj, 'kT', None)
With the class Ising
, we can now instanize a Ising model by creating a 2D grid of spin state:
# HANDSON: instanize a Ising model by creating a 2D grid of spin state.
We can visualize the spin state by plotting ising
as an image with plt.imshow()
:
# HANDSON: visualize the spin state by plotting the state as an image
Interaction and Energy Minimization#
The system’s total energy is determined by how well-aligned the spins are with their neighbors.
At high temperatures, spins flip frequently due to thermal fluctuations, leading to disorder.
At low temperatures, spins align into ordered configurations, exhibiting phase transitions similar to those in magnetic materials.
The Hamiltonian (energy function) for the Ising model is:
(458)#\[\begin{align} H = - \sum_{i,j} J_{ij} s_i s_j \end{align}\]where:
\(s_i\) represents the spin at site \(i\).
\(J_{ij}\) is the coupling strength between neighboring spins.
The system evolves dynamically to minimize \(H\), leading to patterns of stable configurations.
Energy Change Due to Spin Flip#
When a single spin \(s_i\) flips, the change in energy is calculated by comparing the energy before and after the flip. The energy difference \(\Delta E\) due to flipping the spin at site \(i\) is:
where the sum is over the nearest neighbors \(j\) of site \(i\). The factor of 2 arises because flipping the spin at site \(i\) changes its contribution to the energy from \(-s_i s_j\) to \(+s_i s_j\).
Choose \(J_{ij}\) so that only neighborhood cells have interaction 1, we can implement the following function:
# HANDSON: please implement the following function according to the above equation
def energy_change(ising, i, j):
I, J = ising.shape # Get lattice dimensions
spin = ising[i, j] # Current spin value
# Neighboring spins (periodic boundary conditions)
neighbors = ( ising[(i+1)%I, j]
+ ising[(i-1)%I, j]
+ ising[i, (j+1)%J]
+ ising[i, (j-1)%J])
# Energy difference due to spin flip
return 2 * spin * neighbors
We may now obtain the energy change for flipping a spin at different cells:
# HANDSON: choose different i and j and compute the energy change by calling `energy_change(state, i, j)`
Metropolis Algorithm#
To simulate the evolution of the system at a given temperature, we use the Metropolis algorithm. This algorithm probabilistically accepts or rejects a spin flip based on the energy change \(\Delta E\) and the temperature \(T\). The probability of accepting a spin flip is given by the Boltzmann factor:
where:
\(\Delta E\) is the energy change caused by the flip.
\(k\) is the Boltzmann constant. \(T\) is the temperature.
This allows the system to “explore” higher energy states at higher temperatures (thermal fluctuations), while favoring low-energy configurations as the system cools down.
We can implement a function flip(ising, i, j)
to determent if we need to flip the spin ising
at grid i
, j
:
# HANDSON: please implement the following function according to the above equation
def flip(ising, i, j):
...
Note the np.random.rand()
function in flip()
.
Choose a cell with positive energy change and run flip multiple times.
What do you see?
# HANDSON: try repeat flip() many times
What happen if you change the temperature?
# HANDSON: try change `kT` and then repeat flip()
With the above helper functions, we are ready to implement the Ising Model using the metropolis algorithm:
# HANDSON: please implement the following function according to the metropolis algorithm
def run(ising, N):
...
That’s it! We can how compare the spin state before and after running the simulation.
# HANDSON: visualize the initial spin state
# HANDSON: run the simuluation for, e.g., 100,000 steps
# HANDSON: ... and then visualize the final spin state
Animate the Ising Model#
Just for fun, let’s definte the following helper function to create animation:
# Define `animate(state, Nsub=100, N=250)` to create animation
from matplotlib.animation import FuncAnimation
from IPython.display import HTML
def animate(state, update, nsub=100, n=250, cmap='viridis'):
fig, ax = plt.subplots(1, 1, figsize=(5,5))
im = ax.imshow(state, cmap=cmap)
def func(i):
update(state, nsub)
im.set_array(state)
return [im]
an = FuncAnimation(fig, func, frames=n, interval=40, blit=True)
plt.close()
return an
# HANDSON: use the animate(state, update, nsub, n) to create an animation object `an`
# HANDSON: show the animation object as a video using HTML()
Measuring Magnetization#
Physically, magnetization measures the net magnetic moment of the system. In the 2D Ising model, each spin contributes either +1 or -1. The total magnetization is simply the sum of all spins \(M = \sum_i s_i\).
We often normalize the magnetization by dividing by the total number of spins to get a value between -1 and 1:
If all spins are aligned up, \(m = 1\).
If all spins are aligned down, \(m = -1\).
If spins are randomized, \(m \approx 0\).
This normalized magnetization tells us how “ordered” the system is.
# HANDSON: please implement the following function according to the above equation
def mag(ising):
...
With this, we can plot the magnetization as function of, e.g., every 1000 steps.
# HANDSON: plot magnetization, i.e., as a function of every 1000 steps for 1000_000 steps in total
Phase Transition#
The Ising model is important because it is the simplest model that shows phase transition.
The temperature kT
is the parameter that control the phase.
At high temperature, thermal fluctuations dominate: the system becomes disordered and \(m \approx 0\).
At low temperature, the system settles into an ordered state: \(m \approx \pm1\).
Near a critical temperature T_c, the system undergoes a phase transition from disordered to ordered.
# HANDSON: Let's implementat a function that takes the temperature
# of an Ising model and returns a list of magnetization
def getmag(kT):
...
# HANDSON: We can now use getmag(kT) to run numerical experiments.
# First, let's run getmag() and plot the history 10 times for the same temperature.
# You may use the `tqdm` package to obtain a process bar.
# HANDSON: Next, let's run getmag() and plot the history 10 times for different temperature
You should observe see the magnetization evolves very differently above critical temperature \(k T_c \approx 2.27\), which indicates a phase transition. To make it easier to see, one may plot just the absolute value of the last point of the histories:
# HANDSON: Let's only plot the absolute value of the last point of the histories
The above plot is a bit noise… Let’s average 10 realizations per temperature:
# HANDSON: Let's average 10 realizations
Phase transition is one of the most striking features of the 2D Ising model!
Hopfield Networks#
The Hopfield network, introduced by John Hopfield in 1982, is a recurrent neural network that stores and recalls patterns through energy minimization.
Pattern Storage and Recall: The network can store multiple patterns as stable configurations.
Robustness to Noise: If an input is noisy or incomplete, the network still retrieves the correct stored pattern.
Content-Addressable Memory: Unlike conventional memory, which retrieves data using addresses, Hopfield networks retrieve patterns based on similarity.
Energy Landscape#
The Hopfield network shares deep similarities with the Ising model:
Neurons replace spins, taking binary values \(+1\) or \(-1\).
Connections between neurons correspond to interactions between spins.
The system evolves toward stable, low-energy states, just like in the Ising model.
Both exhibit emergent global order despite being governed by local interactions.
The energy function of the Hopfield network is given by:
where \(W_{ij}\) represents synaptic weights and \(\sigma_i\) is the neuron state.
When a neuron updates, its energy change follows:
Why It Matters#
The Hopfield network played a foundational role in artificial intelligence by introducing energy-based models:
Inspired later developments such as Boltzmann machines and Deep Learning architectures.
Showed how memory and pattern recognition could emerge from simple neuron interactions.
Provided insights into biological neural networks, helping us understand how the brain processes and retrieves information.
Because of the similar between a Hopfield network and Ising model, similar, we subclass a numpy array to store the state of a Hopfield network.
Instead of keeping track of the model’s temperature kT
, we create an empty array W
to store the synaptic weights.m
# Subclass a numpy array to store the state of a Hopfield network
class Hopfield(np.ndarray):
def __new__(cls, shape=(64,64), seed=None):
N = shape[0] * shape[1]
if seed is not None:
np.random.seed(seed)
obj = np.random.choice([-1,1], size=shape).view(cls)
obj.W = np.zeros((N, N))
return obj
def __array_finalize__(self, obj):
if obj is None: return
self.W = getattr(obj, 'W', None)
By default, an instance of the Hopfield class is random.
# HANDSON: instanize a Hopfield network and visualize its initial random state
Hebbian Learning Rule#
The weights in a Hopfield network are learned using the Hebbian learning rule, which strengthens the connections between neurons that are activated together. It is given by: $\( W_{ij} = \frac{1}{P} \sum_p \sigma_i^{(p)} \sigma_j^{(p)}. \)$ Where:
\(N\) is the number of neurons.
\(P\) is the number of patterns.
\(\sigma_i^{(p)}\) is the state of neuron \(i\) in a pattern \(p\).
This makes Hopfield networks powerful for associative memory, allowing pattern retrieval even with noise.
Using numpy
’s function, it is straightforward to implement Hebbian’s learning rule:
# HANDSON: please implement the following function according to the above equation
def learn(hopfield, patterns):
...
To test it out, let’s download some sample images and read them into this notebook. Because downloading images does not affect our understanding of Hopfield network, let’s just run the following code cells directly:
# Download images from GitHub in case this notebook is opened in Google colab
url = "https://raw.githubusercontent.com/rndsrc/2024_nobel-phys_hs/refs/heads/main/images/"
url1 = url + "AI.jpg"
url2 = url + "cloud.jpg"
url3 = url + "computer.jpg"
! if [ ! -d images ]; then mkdir images && cd images && wget {url1} && wget {url2} && wget {url3}; fi
--2025-04-24 23:08:02-- https://raw.githubusercontent.com/rndsrc/2024_nobel-phys_hs/refs/heads/main/images/AI.jpg
Resolving raw.githubusercontent.com (raw.githubusercontent.com)... 185.199.111.133, 185.199.108.133, 185.199.109.133, ...
Connecting to raw.githubusercontent.com (raw.githubusercontent.com)|185.199.111.133|:443... connected.
HTTP request sent, awaiting response...
200 OK
Length: 14238 (14K) [image/jpeg]
Saving to: ‘AI.jpg’
AI.jpg 0%[ ] 0 --.-KB/s
AI.jpg 100%[===================>] 13.90K --.-KB/s in 0s
2025-04-24 23:08:02 (67.5 MB/s) - ‘AI.jpg’ saved [14238/14238]
--2025-04-24 23:08:02-- https://raw.githubusercontent.com/rndsrc/2024_nobel-phys_hs/refs/heads/main/images/cloud.jpg
Resolving raw.githubusercontent.com (raw.githubusercontent.com)... 185.199.111.133, 185.199.109.133, 185.199.110.133, ...
Connecting to raw.githubusercontent.com (raw.githubusercontent.com)|185.199.111.133|:443... connected.
HTTP request sent, awaiting response...
200 OK
Length: 14206 (14K) [image/jpeg]
Saving to: ‘cloud.jpg’
cloud.jpg 0%[ ] 0 --.-KB/s
cloud.jpg 100%[===================>] 13.87K --.-KB/s in 0s
2025-04-24 23:08:02 (97.1 MB/s) - ‘cloud.jpg’ saved [14206/14206]
--2025-04-24 23:08:02-- https://raw.githubusercontent.com/rndsrc/2024_nobel-phys_hs/refs/heads/main/images/computer.jpg
Resolving raw.githubusercontent.com (raw.githubusercontent.com)... 185.199.108.133, 185.199.109.133, 185.199.111.133, ...
Connecting to raw.githubusercontent.com (raw.githubusercontent.com)|185.199.108.133|:443... connected.
HTTP request sent, awaiting response...
200 OK
Length: 12092 (12K) [image/jpeg]
Saving to: ‘computer.jpg’
computer.jpg 0%[ ] 0 --.-KB/s
computer.jpg 100%[===================>] 11.81K --.-KB/s in 0s
2025-04-24 23:08:03 (95.6 MB/s) - ‘computer.jpg’ saved [12092/12092]
# Define a helper function to read images
from matplotlib import image as img
def load(filename):
im = img.imread(filename)
if im.ndim == 3:
im = np.mean(im, axis=-1) # Convert to grayscale if needed
im = np.where(im < 128, -1, 1) # Convert grayscale to binary (-1, 1)
return im
# Load images and display them
im1 = load("images/AI.jpg")
im2 = load("images/cloud.jpg")
fig, (ax0, ax1) = plt.subplots(1, 2)
ax0.imshow(im1, cmap='gray')
ax0.set_title('Pattern "AI"')
ax1.imshow(im2, cmap='gray')
ax1.set_title('Pattern "cloud"')
plt.show()
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
Cell In[28], line 6
3 im1 = load("images/AI.jpg")
4 im2 = load("images/cloud.jpg")
----> 6 fig, (ax0, ax1) = plt.subplots(1, 2)
7 ax0.imshow(im1, cmap='gray')
8 ax0.set_title('Pattern "AI"')
NameError: name 'plt' is not defined
With these sample images, let’s make our Hopfield network learn from these patterns. Note that the state of the Hopfield network has not change. It is still random.
# HANDSON: check what happen to the Hopfield state after learning the patterns
Once the patterns are learn, we can recall each pixel value by:
which can be implemented as:
# HANDSON: please implement the following function according to the above equation
def recall(hopfield, N):
...
That’s it! We’ve implemented all the functions of the Hopfield network! We can “recall” its memory and visulize the result by:
# HANDSON: recall a pattern and visualize the result
Animate the Hopfield Network#
Just for fun again, we can now combine the implemented Hopfield network with the animate(state, update, nsub, n)
helper function that we defined at the very beginning.
The recall()
function that we implemented actually has the same call signature as the update()
function.
Hence, we can use it to create the animation:
# HANDSON: animate memory recall of a Hopfield network
Discussion#
This notebook has demonstrated how concepts from physics and AI intersect through models like the Ising Model and Hopfield Networks. Here are some questions to encourage further exploration:
How do energy minimization principles in physics apply to AI and machine learning? Can you think of other real-world phenomena that could be modeled using energy-based systems? How might AI benefit from further insights from statistical mechanics and physics? What improvements could be made to the Hopfield Network to make it more efficient? For further study, students can explore more advanced models such as Boltzmann Machines, Restricted Boltzmann Networks, and Deep Learning frameworks.