Fractals

Posted by Damian Dailisan on March 1, 2016 Tags: Python Complex Systems   3 minute read

Fractal Dimension

Fractals are geometric features in which structures are similar when viewed at small localized regions, or the entirety of the structure itself. Fractals are said to be naturally occurring in nature, which is the motivation of applying concepts of fractals in modelling physical systems.

One way to look at fractals is by their fractal dimension, given by \(D = \lim_{\epsilon \to 0} \dfrac{\log{N}}{\log{1/\epsilon}},\) where $N$ is the number of boxes that contains the structure, and $\epsilon$ is the size of the box. As such, this method is called the box counting method. This is useful for looking at the fractal dimensionality of two dimensional images.

Suppose we want to look at the fractal dimension of this image:

alt text

In [1]:

%matplotlib inline
%config InlineBackend.figure_format = 'svg'
import numpy as np
import matplotlib.pyplot as plt
from skimage import color

# img = color.rgb2gray(plt.imread('lightning-fractal.jpg'))
img = plt.imread('lightning-fractal.jpg')[:,:,0]/255.
plt.imshow(img>0.4, cmap="gray")
<matplotlib.image.AxesImage at 0x6d2ecf8>

svg

The following lines of code will count the number of blocks containing the structure $N$ for a given input of $\epsilon$.

In [2]:

def count_N(img, epsilon):
    N = 0
    for i in range(0, img.shape[0], epsilon):
        for j in range(0, img.shape[1], epsilon):
            if i==img.shape[0]:
                N += (np.sum(img[i:-1, j:j+epsilon])>0)
            elif j==img.shape[1]:
                N += (np.sum(img[i:i+epsilon, j:-1])>0)
            else:
                N += (np.sum(img[i:i+epsilon, j:j+epsilon])>0)
#             print np.sum(img[i:i+epsilon, j:j+epsilon])>0
    return N

We will then plot $\log{N}$ with $\log{1/\epsilon}$ and if the slope of the resulting plot, representing $D$ appoaches some finite value.

In [3]:

epsilon_values = np.arange(200,0,-1)
N = [count_N(img, epsilon) for epsilon in epsilon_values]

fig = plt.figure()
ax = fig.add_subplot(111)
x = np.log(1./epsilon_values)
y = np.log(N)
ax.scatter(x, y)
ax.set_ylabel(r"$\log{N}$")
ax.set_xlabel(r"$\log{1/\epsilon}$")
<matplotlib.text.Text at 0x6f50ba8>

svg

Getting the fractal dimension

From the definition of the fractal dimension, we can see that it is simply the slope of $\log{N}$ vs $\log{1/\epsilon}$ for large values of $1/\epsilon$. To obtain this, it is best to fit a line at the end of the plot. In this case, we will fit a line to the last 50 points, the slope of which is our fractal dimension.

In [4]:

from scipy import stats
import numpy as np

slope, intercept, r_value, p_value, std_err = stats.linregress(x[-50:],y[-50:])

print 'D = %s' % slope
D = 1.93467964091

Hence, the fractal dimension of the image is $D=1.935$.