Box counting method for fractal dimension

First I saw here that the box counting fractal dimension defined by

$D = \lim_{\epsilon \rightarrow 0}{ {\log N( \epsilon)} \over {\log { {1}\over{ \epsilon }}}}$

which makes sense for me. Then I saw here a mathematica code for calculating it. This method does a fitting of the data $\{ x= \log { {1}\over{ \epsilon }}, y= \log N( \epsilon) \}$ to a line, giving two parameter $a$ and $b$ such that

$y = a x + b$

and $a$ is interpreted as the box counting fractal dimension, which I don't understand because I don't think there is any reason to assume linear behaviour.

However, this method is "correct", because

${{y} \over {x}} = a + {{b}\over{x}} \quad$ and $\quad \lim_{\epsilon \rightarrow 0} = \lim_{x \rightarrow \infty}$

but by the same reason there are infinite other "correct" methods assuming

$y = a x + b + {{c} \over {x}} + {{d} \over {x^2}} + ...$

or other more complicated things.

Why do they choose a linear behaviour? (Instead of just plotting ${ {\log N( \epsilon)} \over {\log { {1}\over{ \epsilon }}}}$ against $\epsilon$ and taking the $y$-intercept)

migrated from mathematica.stackexchange.comNov 12 '17 at 23:09

This question came from our site for users of Wolfram Mathematica.

• The limit is really only defined for mathematically idealized objects. For real objects, we should look at a sufficiently wide range of scales, from very zoomed out up to very zoomed in, and compute the dimension at each one. A shape is considered to be fractal only when the measured dimension stays approximately constant across multiple scales. The slope computation accomplishes just that. This youtube video describes the situation nicely. – Mark McClure Nov 10 '17 at 15:39
• Precisely in that video it is noticed how slow the fractal dimension converges when calculated that way. The value $a$ in my question (the order of the polynomial in the video), only approaches the fractal dimension when $x \rightarrow 0$ (when the scaling factor goes to infinity in the video). So why not just take the computational or approximate limit? – MBolin Nov 12 '17 at 9:48

Like Euclidean objects, fractals are idealized abstractions of reality. A tree is not a fractal any more than its trunk is a line segment. As a result, we cannot compute an actual limit to find the box-counting dimension of an object. Even an attempt to approximate the limit by using a small value of $\varepsilon$ in $$\frac{\log(N({\varepsilon}))}{\log(1/\varepsilon)}$$ is hazardous because we don't know a priori what scale of $\varepsilon$ works well with the object. The process of fitting a line to the data $$\{ x= \log { {1}/ \varepsilon_k}, y= \log N( \epsilon_k) \}$$ for some terms chosen from a sequence $\varepsilon_k$ that tends down to zero allows us to estimate the dimension in a way that accounts for behavior of the set over a large scale of values of $\varepsilon$.

We can illustrate this with an example which, since this question was originally posted on mathematica.stackexchange, we present with Mathematica. Consider the following image of the Sierpinski carpet:

Import[
"https://www.marksmath.org/FractalGeometryPackages/FractalGeometry/IteratedFunctionSystems.m"
];
A = {{1, 0}, {0, 1}}/3;
IFS = {A, #} & /@ {
{-1, 1}, {0, 1},  {1, 1},
{-1, 0},          {1, 0},
{-1, -1}, {0, -1},{1, -1}
};
sierpPic = ShowIFS[IFS, 5, PlotRange -> All,
Initiator -> Polygon[{{-3, -3}, {3, -3}, {3, 3}, {-3, 3}}/2], Here's an easy way to box-count the image:

pic = Binarize[ImportString[ExportString[sierpPic, "PNG",
ImageSize -> 2^10]]];
data = ImageData[pic];
test[matrix_] := Boole[Length[Cases[matrix, 0, {2}, 1]] > 0];
count[n_] := Total[Flatten[Map[test, Partition[data, {n, n}], {2}]]];
radii = Table[2^k, {k, 0, 8}];

(* Out:
{626688, 167476, 45276, 12124, 3287, 864, 236, 60, 16}
*)


Here's the "approximate limit" computation:

Log[First[counts]]/Log[2^10] // N

(* Out: 1.92574 *)


Here's a slope fit computation:

Fit[Transpose[Log[{1/radii[[2 ;; -3]], counts[[2 ;; -3]]}]], {1, x}, x]

(* Out: 13.3465 + 1.89636 x *)


The actual dimension is $$\frac{\log(8)}{\log(3)} \approx 1.89279.$$ Thus, this choice of slope fitting does a better job.

Now, I freely admit that I fiddled a bit with the choice of points to use in my fit of the slope. The point, though, is that slope fitting can do a better job. Also, it's easy to see why the "approximate limit" approach might fail. If $\varepsilon$ is too small, the small blocks in the approximate Sierpinski carpet are solid and contribute too much to the box-count. Thus, we get too large a value. That's exactly the point when it comes to a real world object. An apparent, approximate self-similar structure is only valid down to some scale and you have no way of knowing how far down you should look.

Ultimately, though, box-counting with natural objects is necessarily approximate, error prone, and not necessarily even well defined.

If you prefer an example with Python, you find one over on Stackoverflow.