Search code examples
image-processingpython-3.xbayesianpymc

Examples/tutorials of Adaptive Metropolis for images using PyMC


I am looking for examples or tutorials of the AdaptiveMetropolis step method used for image processing.

The only vaguely image-related resource that I have found until now is this astronomy dissertation and the related GitHub repo.

This wider question does not seem to provide PyMC example code.

What about finding the peak on this simulated array?

import numpy as np
from matplotlib import pyplot as plt

sz = (12,18)
data_input = np.random.normal( loc=5.0, size=sz )
data_input[7:10, 2:6] = np.random.normal( loc=100.0, size=(3,4) )

fig = plt.figure()
ax = fig.add_subplot(1,1,1)
im = ax.imshow( data_input )
ax.set_title("input")

enter image description here


Solution

  • The closest I know of is here: http://nbviewer.ipython.org/github/CamDavidsonPilon/Probabilistic-Programming-and-Bayesian-Methods-for-Hackers/blob/master/Chapter5_LossFunctions/LossFunctions.ipynb (see Example: Kaggle contest on Observing Dark World)

    On this thread you asked a specific question: https://github.com/pymc-devs/pymc/issues/653 about finding an array in an image. Here is a first attempt at a model:

    In that case it seems like you are trying to estimate a 2-D uniform distribution with gaussian noise. You'll have to translate this into an actual model but this would be one idea:

    lower_x ~ DiscreteUniform(0, 20)

    upper_x ~ DiscreteUniform(0,20)

    lower_y ~ DiscreteUniform(0, 20)

    upper_y ~ DiscreteUniform(0, 20)

    height ~ Normal(100, 1)

    noise ~ InvGamma(1, 1)

    means = zeros((20, 20))

    means[[lower_x:upper_x,lower_y: upper_y]] = height # this needs to be a deterministic

    data ~ Normal(mu=means, sd=noise)

    It might be better to code upper_x as an offset and then do lower_x:lower_x+offset_x, otherwise you need a potential to enforce lower_x < upper_x.