Dimensión de Recuento por Cajas (Box counting)

Se define la dimensión de recuento por cajas de un conjunto $E$ como $$ dim_B E = \lim_{\delta\to0}\frac{\log N_\delta(E)}{-\log \delta} $$ (ver [1] página 45)

Cargamos un conjunto de puntos

import pandas as pd
puntos = pd.read_csv('listaPuntosSierpinski.csv',sep=';',index_col=0)
puntos.head()

	X	Y
0	0.207282	0.559435
1	0.603641	0.279718
2	0.801820	0.139859
3	0.650910	0.502942
4	0.825455	0.251471

puntos.shape

(5001, 2)

Calculamos $dim_BE$ en $[0,1]^2$

Dividimos el cuadrado unidad en cuadrados de lado mitad sucesivamente

Etapas: $k=0,1,2,\ldots,N_{etapas}$

Tamaño intervalos: $2^{-k}$

Número de cuadrados: $4^k$

Número de intervalos por eje: $2^k$

Para cada etapa $k$, dividimos los intervalos para los ejes $X$ e $Y$ en intervalos diálicos:

$$ I_{j,k}=[j\cdot2^{-k},(j+1)\cdot2^{-k}],\qquad j=0,1,\ldots,2^k-1 $$

selX = puntos.loc[((puntos.X >0.78) & (puntos.X<0.7))]
selX.shape

(0, 3)

numEtapas = 9

vecNk = []

for k in range(numEtapas):
    
    numIntervalosEje = 2**k
    sizeInt = 2**(-k)
    Nk = 0
    
    for jX in range(numIntervalosEje):
        
        limInfX = jX*sizeInt
        limSupX = (jX+1)*sizeInt
        
        selX = puntos.loc[((puntos.X >= limInfX) & (puntos.X<=limSupX))]
        if selX.shape[0]:
            
            for jY in range(numIntervalosEje):
                
                limInfY = jY*sizeInt
                limSupY = (jY+1)*sizeInt
                
                selY = selX.loc[((selX.Y >= limInfY) & (selX.Y<=limSupY))]
                
                if selY.shape[0]:
                    Nk+=1
    
    vecNk.append(Nk)

res = pd.DataFrame(columns=['Size','NumInt','NumBox','Nk'],index = range(numEtapas))
size = [2**(-k) for k in range(numEtapas)]
NumInt = [2**(k) for k in range(numEtapas)]
NumBox = [4**k for k in range(numEtapas)]
res.Size = size
res.NumInt = NumInt
res.NumBox = NumBox
res.Nk = vecNk
res

	Size	NumInt	NumBox	Nk
0	1.000000	1	1	1
1	0.500000	2	4	4
2	0.250000	4	16	13
3	0.125000	8	64	38
4	0.062500	16	256	116
5	0.031250	32	1024	346
6	0.015625	64	4096	990
7	0.007812	128	16384	2332
8	0.003906	256	65536	3756

res.to_csv('DimBoxSierpinskiTriangle.csv',sep=';')

Enlace a DimBoxSierpinskiTriangle.csv

import numpy as np
f,ax = plt.subplots(nrows=1,ncols=1,figsize =(8,8))
logX = np.log(res.Size.to_numpy())
logY = -np.log(res.Nk.to_numpy())
ax.plot(logX,logY,'.-',markersize=10)

from sklearn.linear_model import LinearRegression
reg = LinearRegression().fit(logX.reshape((-1, 1)), logY)
xmin = np.min(logX)
xmax = np.max(logX)
xx = np.linspace(xmin,xmax,100)
yy = reg.coef_*xx+reg.intercept_
ax.plot(xx,yy,'r-',label='y={:.3f}x+{:.2f}'.format(reg.coef_[0],reg.intercept_))
ax.legend()
ax.grid()

Comparación con dimensión de homotecia

El triángulo de Sierpinski es un fractal autosemejante de $N=3$ copias de tamaño $r=1/2$, por tanto la dimensión de homotecia es:

$$ dim(E) = \frac{\log N}{-\log r}=\frac{\log 3}{\log 2} $$

dimHomotecia = np.log(3)/np.log(2)
dimHomotecia

1.5849625007211563

dimCalculada = reg.coef_[0]
dimCalculada

1.512506205301025

diferencia = np.abs(dimHomotecia-dimCalculada)
diferencia

0.07245629542013132

Bibliografía

[1] Guzmán, M., Martín, M.A., Morán, M. and Reyes, M. (1993) Estructuras fractales y sus aplicaciones. Ed. Labor, Barcelona