This notebook shows how to set up an MRF for a vision problem.
To segment an image into foreground and background, we can try to label each pixel as '0'
(background), or '1'
(foreground).
We therefore have one random variable for each observed pixel, and one random variable for for each unobserved pixel label.
If we assume that labels are independent of pixel values and other labels given their neighbours, we get a grid of labels. Let us also assume that the observed pixel values are independent given the labels.
Since we will have many factors of the same form, let's template them.
Here is an example to show how templates can be used. Each entry in the potential table can be associated with a parameter, for example:
import numpy as np
from pyugm.factor import DiscreteFactor
from pyugm.factor import DiscreteBelief
# Specify the parameters (should be the same shape as the potential table would have been)
factor_parameters = np.array([['theta_0', 'theta_1'], ['theta_2', 'theta_3']])
variables_names_and_cardinalities = [(1, 2), (2, 2)]
# Construct the factor
factor = DiscreteFactor(variables_names_and_cardinalities, parameters=factor_parameters)
print factor
F{1, 2}
# The factor still has its default potential table (all ones)
factor.data
array([[ 1., 1.], [ 1., 1.]])
# Create a belief based on the factor
belief = DiscreteBelief(factor)
Now we can set the parameters...
# Potentials are filled with the exponent of the parameters.
belief.set_parameters({'theta_0': np.log(2), 'theta_1': np.log(0.2), 'theta_2': np.log(5), 'theta_3': np.log(1)})
belief.data
array([[ 2. , 0.2], [ 5. , 1. ]])
...and later reset them to something else:
belief.set_parameters({'theta_0': np.log(1), 'theta_1': np.log(2), 'theta_2': np.log(1), 'theta_3': np.log(1)})
belief.data
array([[ 1., 2.], [ 1., 1.]])
The same parameter name can also be specified more than once to force the table to fill different cells with the same value.
factor_parameters = np.array([['theta_0', 'theta_1'], ['theta_0', 'theta_0']])
factor = DiscreteFactor([(1, 2), (2, 2)], parameters=factor_parameters)
belief = DiscreteBelief(factor)
belief.set_parameters({'theta_0': np.log(3), 'theta_1': np.log(5)})
print belief.data
[[ 3. 5.] [ 3. 3.]]
For a baseline, let's classify pixels as background and foreground based only on their intensity.
import pickle
import matplotlib.pyplot as plt
import matplotlib
%matplotlib inline
import seaborn
seaborn.set_style("dark")
# Load pre-discretized image. Each pixel value is an integer between 0 and 32.
image = pickle.load(open('test_image.pkl'))
plt.figure(figsize=(14, 3))
plt.subplot(1, 2, 1)
_ = plt.imshow(image, cmap=matplotlib.cm.Greys_r, interpolation='nearest')
_ = plt.title('Image')
seaborn.set_style("darkgrid")
plt.subplot(1, 2, 2)
_ = plt.hist(image.flatten(), bins=32, color=(0.5, 0.5, 0.5))
_ = plt.title('Pixel intensity histogram')
Let's make a template where the probability that a pixel is foreground or background is different given that the pixels intensities are either between 13 and 17 or not.
observation_template = np.array([['obs_low'] * 32,
['obs_high'] * 32])
observation_template[0, 13:17] = 'obs_high'
observation_template[1, 13:17] = 'obs_low'
Each observed pixel variable occurs in the same factor as its label.
I, J = image.shape
factors = []
evidence = {}
# Add observation factors
for i in xrange(I):
for j in xrange(J):
label_variable_name = 'label_{}_{}'.format(i, j)
observation_variable_name = 'obs_{}_{}'.format(i,j)
factors.append(DiscreteFactor([(label_variable_name, 2), (observation_variable_name, 32)], parameters=observation_template))
evidence[observation_variable_name] = image[i, j]
from pyugm.model import Model
model = Model(factors)
Give a higher probability of seeing background pixel intensities of 13 to 17 than foreground intensities in that range.
from pyugm.infer_message import LoopyBeliefUpdateInference
from pyugm.infer_message import FloodingProtocol
order = FloodingProtocol(model, max_iterations=30)
inference = LoopyBeliefUpdateInference(model, order)
parameters = {'obs_high': 0.1, 'obs_low': -1.0}
(Doesn't really do belief updates because none of the factors are connected - but let's go through the process to see how it works)
inference.calibrate(evidence, parameters)
<pyugm.infer_message.LoopyBeliefUpdateInference at 0x12ef30c50>
labels = np.zeros(image.shape)
for i in xrange(I):
for j in xrange(J):
variable_name = 'label_{}_{}'.format(i, j)
label_factor = inference.get_marginals(variable_name)[0]
labels[i, j] = label_factor.normalized_data[0]
seaborn.set_style("dark")
plt.figure(figsize=(14, 3))
plt.subplot(1, 3, 1)
_ = plt.imshow(image, cmap=matplotlib.cm.Greys_r, interpolation='nearest')
_ = plt.title('Original image')
plt.subplot(1, 3, 2)
_ = plt.imshow(labels, cmap=matplotlib.cm.Greys, interpolation='nearest')
_ = plt.title('Inferred labels (white=foreground, black=background)')
Now we can start to think how the factors must be filled to nudge neighbouring labels to be similar and smooth out the baseline a bit.
In addition to having label and pixel-value variables in the same factor, let's also add factors where neighbouring label variables occur together.
Let's add one parameter for when these two variables' labels are the same and one for when they are different.
label_template = np.array([['same', 'different'],
['different', 'same']])
Before building a model where each label is connected to all four of its neighbours, let's first only connect the top and bottom neighbour. Pixel columns are thus connected.
evidence = {}
factors = []
# Add observation factors
for i in xrange(I):
for j in xrange(J):
label_variable_name = 'label_{}_{}'.format(i, j)
observation_variable_name = 'obs_{}_{}'.format(i, j)
factors.append(DiscreteFactor([(label_variable_name, 2), (observation_variable_name, 32)], parameters=observation_template))
evidence[observation_variable_name] = image[i, j]
# Add label factors
for i in xrange(I):
for j in xrange(J):
variable_name = 'label_{}_{}'.format(i, j)
if i + 1 < I:
neighbour_down_name = 'label_{}_{}'.format(i + 1, j)
factors.append(DiscreteFactor([(variable_name, 2), (neighbour_down_name, 2)], parameters=label_template))
Add the factors to a model.
model = Model(factors)
from pyugm.infer_message import LoopyBeliefUpdateInference
from pyugm.infer_message import FloodingProtocol
# Get some feedback on how inference is converging by listening in on some of the label beliefs.
var_values = {'label_1_1': [], 'label_10_10': [], 'label_20_20': [], 'label_30_30': [], 'label_40_40': []}
changes = []
def reporter(infe, orde):
for var in var_values.keys():
marginal = infe.get_marginals(var)[0].data[0]
var_values[var].append(marginal)
change = orde.last_iteration_delta
changes.append(change)
energy = infe.partition_approximation()
print '{:3} {:8.2f} {:5.2f} {:8.2f}'.format(orde.total_iterations, change, marginal, energy)
order = FloodingProtocol(model, max_iterations=15)
inference = LoopyBeliefUpdateInference(model, order, callback=reporter)
Lets choose parameters so that neighbouring labels are encouraged to be the same, and observations between 13 and 17 more likely to be background.
parameters = {'same': 2.0, 'different': -1.0, 'obs_high': 1.0, 'obs_low': -0.0}
inference.calibrate(evidence, parameters)
1 19715.74 0.73 8207.46 2 982.23 0.97 8668.57 3 347.46 0.99 8831.12 4 143.23 0.99 8881.90 5 68.94 1.00 8890.86 6 39.42 1.00 8893.96 7 24.94 1.00 8891.35 8 13.82 1.00 8891.82 9 8.19 1.00 8891.47 10 4.60 1.00 8892.31 11 2.94 1.00 8892.65 12 1.85 1.00 8893.02 13 1.13 1.00 8892.81 14 0.62 1.00 8892.83 15 0.34 1.00 8892.83
<pyugm.infer_message.LoopyBeliefUpdateInference at 0x131cbae90>
labels = np.zeros(image.shape)
for i in xrange(I):
for j in xrange(J):
variable_name = 'label_{}_{}'.format(i, j)
label_factor = inference.get_marginals(variable_name)[0]
labels[i, j] = label_factor.normalized_data[0]
plt.figure(figsize=(14, 3))
plt.subplot(1, 3, 1)
_ = plt.imshow(image, cmap=matplotlib.cm.Greys_r, interpolation='nearest')
_ = plt.title('Original image')
plt.subplot(1, 3, 2)
_ = plt.imshow(labels, cmap=matplotlib.cm.Greys, interpolation='nearest')
_ = plt.title('Label beliefs \n(darker=higher background belief,\n lighter=higher foreground belief')
plt.subplot(1, 3, 3)
_ = plt.imshow(labels > 0.5, cmap=matplotlib.cm.Greys, interpolation='nearest')
_ = plt.title('Thresholded beliefs')
Let's also connect up the rows, to form a grid of labels.
evidence = {}
factors = []
# Add observation factors
for i in xrange(I):
for j in xrange(J):
label_variable_name = 'label_{}_{}'.format(i, j)
observation_variable_name = 'obs_{}_{}'.format(i, j)
factors.append(DiscreteFactor([(label_variable_name, 2), (observation_variable_name, 32)], parameters=observation_template))
evidence[observation_variable_name] = image[i, j]
# Add label factors
for i in xrange(I):
for j in xrange(J):
variable_name = 'label_{}_{}'.format(i, j)
if i + 1 < I:
neighbour_down_name = 'label_{}_{}'.format(i + 1, j)
factors.append(DiscreteFactor([(variable_name, 2), (neighbour_down_name, 2)], parameters=label_template))
if j + 1 < J:
neighbour_right_name = 'label_{}_{}'.format(i, j + 1)
factors.append(DiscreteFactor([(variable_name, 2), (neighbour_right_name, 2)], parameters=label_template))
Build the model.
model = Model(factors)
from pyugm.infer_message import LoopyBeliefUpdateInference
from pyugm.infer_message import FloodingProtocol, LoopyDistributeCollectProtocol
# Get some feedback on how inference is converging by listening in on some of the label beliefs.
var_values = {'label_1_1': [], 'label_10_10': [], 'label_20_20': [], 'label_30_30': [], 'label_40_40': []}
changes = []
partitions = []
def reporter(infe, orde):
for var in var_values.keys():
marginal = infe.get_marginals(var)[0].data[0]
var_values[var].append(marginal)
change = orde.last_iteration_delta
changes.append(change)
energy = infe.partition_approximation()
partitions.append(energy)
print '{:3} {:8.2f} {:5.2f} {:8.2f}'.format(orde.total_iterations, change, marginal, energy)
order = FloodingProtocol(model, max_iterations=15)
inference = LoopyBeliefUpdateInference(model, order, callback=reporter)
parameters = {'same': 0.1, 'different': -1.0, 'obs_high': 1.0, 'obs_low': -0.0}
inference.calibrate(evidence, parameters)
1 9774.11 0.89 2175.63 2 2280.77 0.86 2816.06 3 791.25 0.83 3259.23 4 244.91 0.91 3360.01 5 93.39 0.91 3383.34 6 44.46 0.91 3390.13 7 22.15 0.91 3392.01 8 12.15 0.91 3393.14 9 8.79 0.91 3393.88 10 5.80 0.91 3394.42 11 4.05 0.91 3394.82 12 2.47 0.91 3394.97 13 1.80 0.91 3395.17 14 0.98 0.91 3395.32 15 0.70 0.91 3395.39
<pyugm.infer_message.LoopyBeliefUpdateInference at 0x133446690>
We can check that the total change in belief becomes less after each iteration, and that the free enery converges:
seaborn.set_style("darkgrid")
ax1 = plt.axes()
ax2 = ax1.twinx()
_ = ax1.plot([np.log(change) for change in changes])
seaborn.set_style('dark')
_ = ax2.plot(partitions, color='g')
ax2.set_yticks(np.linspace(ax2.get_yticks()[0], ax2.get_yticks()[-1], len(ax1.get_yticks())))
_ = ax1.set_ylabel('Log of total belief change')
_ = ax2.set_ylabel('Log of the approximation to the partition function Z')
_ = ax1.set_xlabel('Iteration number')
And here are the beliefs for a few variables:
for key, value in var_values.items():
plt.plot(value, label=key)
_ = plt.ylabel('Normalized belief that the variable equals 0')
_ = plt.xlabel('Iteration number')
plt.legend()
<matplotlib.legend.Legend at 0x13e909bd0>
labels = np.zeros(image.shape)
for i in xrange(I):
for j in xrange(J):
variable_name = 'label_{}_{}'.format(i, j)
label_factor = inference.get_marginals(variable_name)[0]
labels[i, j] = label_factor.normalized_data[0]
seaborn.set_style("dark")
plt.figure(figsize=(14, 3))
plt.subplot(1, 3, 1)
_ = plt.imshow(image, cmap=matplotlib.cm.Greys_r, interpolation='nearest')
_ = plt.title('Origin image')
plt.subplot(1, 3, 2)
_ = plt.imshow(labels, cmap=matplotlib.cm.Greys, interpolation='nearest')
_ = plt.title('Label beliefs \n(darker=higher background belief,\n lighter=higher foreground belief')
plt.subplot(1, 3, 3)
_ = plt.imshow(labels > 0.5, cmap=matplotlib.cm.Greys, interpolation='nearest')
_ = plt.title('Thresholded beliefs')
To recap: