La sucesión de Fibonacci se usa para modelar crecimiento de poblaciones de bichos. Por ejemplo: conejitos. Se empieza con un granero vacío, en $t_{0}$, en $t_{1}$ una parejita de conejos se metió quién sabe cómo y usan ese tiempo para instalarse, en $t_{1}$ sigue la misma parejita, pero una nueva parejita ya se está gestando. En $t_{3}$ nace una segunda pareja, ahora son dos: la nueva usa ese tiempo para instalarse y la primera otra vez está gestando una nueva.
Uno puede imaginarse el área que ocupa la población de conejitos viendo esta figura:
También puede usarse la serie de Fibonacci para aproximar la proporción áurea:
Más información en el artículo "Sucesión de Fibonacci" en la Wikipedia.
La ecuación de recurrencia de segundo orden es así:
$F_{k+2}=F_{k+1}+F_{k}$ .......eq1
Con valores iniciales $F_{0}=1, F_{1}=1$ la secuencia es:
$0,1,1,2,3,5,8,13,...$
Hay dos formas de obtener el k-ésimo valor de esta serie. Usando la ecuación de recurrencia eq1 o usando esta fórmula cerrada:
$F_{k}=\frac{1}{\sqrt{5}}\left[\left(\frac{1+\sqrt{5}}{2}\right)^{k}-\left(\frac{1-\sqrt{5}}{2}\right)^{k}\right]$.......eq2
Sea $A$ una matriz cuadrada de dimensión n. Si $A$ tiene n valores propios distintos
$\lambda_{1},\lambda_{2},...,\lambda_{n},$
con vectores propios correspondientes
$v_{1},v_{2},...,v_{n},$
entonces $A$ es diagonalizable, o sea existen las matrices $P$ y $D$ de dimensión n tales que $D$ es una matriz diagonal y
$P^{-1}AP=D$
$D=\left[\begin{array}{cccc} \lambda_{1} & 0 & \cdots\ & 0\\ 0 & \lambda_{2} & \cdots\ & 0\\ \vdots\ & \vdots\ & \ddots\ & \vdots\\ 0 & 0 & 0 & \lambda_{n} \end{array}\right]$
$P$ es la matriz cuyas columnas son los vectores propios.
$P=\left[v_{1}|v_{2}|...|v_{n}\right]$
Se generan los núneros de la sucesión de Fibonacci en las coordenadas de los vectores en $R2$ a partir de la condición inicial
$u_{0}=\left[\begin{array}{c} F_{1}\\ F_{0} \end{array}\right]=\left[\begin{array}{c} 1\\ 0 \end{array}\right]$
Sea
$A=\left[\begin{array}{cc} 1 & 1\\ 1 & 0 \end{array}\right]$ .......eq3
entonces, si
$u_{k}=A\mathbf{u}_{k-1}$ .......eq4
para $k=1,2,3...$
$u_{k}=\left[\begin{array}{c} F_{k+1}\\ F_{k} \end{array}\right]$
$u_{k}=A^{k}\mathbf{u}_{0}$ .......eq5
Se usa la ecuación eq4 para generar algunos $u$.
$u_{0}=\left[\begin{array}{c} F_{1}\\ F_{0} \end{array}\right]=\left[\begin{array}{c} 1\\ 0 \end{array}\right]$
$u_{1}=\left[\begin{array}{cc} 1 & 1\\ 1 & 0 \end{array}\right] \left[\begin{array}{c} 1\\ 0 \end{array}\right]=\left[\begin{array}{c} 1\\ 1 \end{array}\right]$
$u_{2}=\left[\begin{array}{cc} 1 & 1\\ 1 & 0 \end{array}\right] \left[\begin{array}{c} 1\\ 1 \end{array}\right]=\left[\begin{array}{c} 2\\ 1 \end{array}\right]$
$u_{3}=\left[\begin{array}{cc} 1 & 1\\ 1 & 0 \end{array}\right] \left[\begin{array}{c} 2\\ 1 \end{array}\right]=\left[\begin{array}{c} 3\\ 2 \end{array}\right]$
$u_{4}=\left[\begin{array}{cc} 1 & 1\\ 1 & 0 \end{array}\right] \left[\begin{array}{c} 3\\ 2 \end{array}\right]=\left[\begin{array}{c} 5\\ 3 \end{array}\right]$
$u_{5}=\left[\begin{array}{cc} 1 & 1\\ 1 & 0 \end{array}\right] \left[\begin{array}{c} 5\\ 3 \end{array}\right]=\left[\begin{array}{c} 8\\ 5 \end{array}\right]$
$u_{6}=\left[\begin{array}{cc} 1 & 1\\ 1 & 0 \end{array}\right] \left[\begin{array}{c} 8\\ 5 \end{array}\right]=\left[\begin{array}{c} 13\\ 8 \end{array}\right]$
A continuación pongamos a prueba la eq5.
In [1]:
# importamos bibliotecas cómputo de matrices
import numpy as np
A = np.matrix([[1, 1],
[1, 0]])
# para k=3 A se multiplica por sí misma tres veces:
A*A*A
Out[1]:
In [2]:
u0 = [1,
0]
# u3
np.dot(A*A*A,
u0)
Out[2]:
Como $\lambda$ es un valor propio de $A$
$(\lambda I-A)\overrightarrow{X}=\overrightarrow{0}$
un eigenvector que sea su solución depende de que
$\det(\lambda I-A)=\det\left(\left[\begin{array}{cc} \lambda & 0\\ 0 & \lambda \end{array}\right]-\left[\begin{array}{cc} 1 & 1\\ 1 & 0 \end{array}\right]\right)=\overrightarrow{0}$
Luego
$\det\left(\left[\begin{array}{cc} \lambda-1 & -1\\ -1 & \lambda \end{array}\right]\right)=(\lambda-1)\lambda-(-1)(-1)=\lambda^{2}-\lambda-1$
Para obtener los valores propios usamos la fórmula general para encontrar valores de $\lambda$ que satisfagan el polinomio característico
$p(\lambda)=\lambda^{2}-\lambda-1$
Acá la fórmula general: $x=\frac{-b\pm\sqrt{b^{2}-4ac}}{2a}$
$\lambda=\frac{1\pm\sqrt{-1^{2}-4(1)(-1)}}{2(1)}=\frac{1\pm\sqrt{1+4}}{2}$
Por eso
$\begin{array}{c} \lambda_{1}=\frac{1+\sqrt{5}}{2}\\ \lambda_{2}=\frac{1-\sqrt{5}}{2} \end{array}$ .......eq6
$\left[\begin{array}{cc} \lambda-1 & -1\\ -1 & \lambda \end{array}\right]\overrightarrow{X}=\left[\begin{array}{cc} \lambda-1 & -1\\ -1 & \lambda \end{array}\right]\left[\begin{array}{c}x_{1}\\ x_{2}\end{array}\right]$
Equivale al sistema de ecuaciones
$(\lambda-1)x_{1}-x_{2}=0$....(1)
$-x_{1}+\lambda x_{2}=0$.....(2)
Despejando $x_{1}$ en (2)
$x_{1}=\lambda x_{2}$
que sustituyendo en (1)
$(\lambda-1)\lambda x_{2}-x_{2}=0$
factorizamos $x_{2}$
$x_{2}\left[(\lambda-1)\lambda-1\right]=0$
sea $x_{2}=1$ en (2)
$-x_{1}+\lambda=0$
$x_{1}=\lambda$
Por eso los vectores propios son
$v_{1}=\left[\begin{array}{c} \lambda_{1}\\ 1 \end{array}\right],\: v_{2}=\left[\begin{array}{c} \lambda_{2}\\ 1 \end{array}\right]$ .......eq7
Con $v_{1}$ y $v_{2}$ formamos la matriz $P$.
Para obtener la matriz invertida $P^{-1}$ Anton y Rorres muestran este teorema:
Si
$A=\left[\begin{array}{cc} a & b\\ c & d \end{array}\right]$
entonces
$A^{-1}=\frac{1}{ad-bc}\left[\begin{array}{cc} d & -b\\ -c & a \end{array}\right]$
Teniendo la matriz
$P=\left[\begin{array}{cc} \lambda_{1} & \lambda_{2}\\ 1 & 1 \end{array}\right]$
entonces
$P^{-1}=\frac{1}{\lambda_{1}-\lambda_{2}}\left[\begin{array}{cc} 1 & -\lambda_{2}\\ -1 & \lambda_{1} \end{array}\right]$
Para comprobar
$P^{-1}AP=\left[\begin{array}{cc} \lambda_{1} & 0\\ 0 & \lambda_{2} \end{array}\right]=D$
basta con multiplicar las matrices $P^{-1}AP$. Como está laborioso hagámoslo con numpy en la siguiente celda.
In [3]:
from math import sqrt
# eigenvalores
l1=(1+sqrt(5))/2.0
l2=(1-sqrt(5))/2.0
print("lambda1=%s, lambda2=%s" % (l1,l2))
# la matriz de transicion P
P = np.matrix([[l1, l2],
[1, 1]])
# valores como 1.66533454e-16 son como ceros! esto pasa por no hacerlo analiticamente
(P**-1)*A*P
Out[3]:
Si $A$ es diagonalizable entonces existe
$P^{-1}AP=D$
Como $PP^{-1}=I$ es la matriz identidad, y es el neutro multiplicativo.
$PP^{-1}APP^{-1}=PDP^{-1}$
luego
$A=PDP^{-1}$
y entonces
$A^{k}=\left(PDP^{-1}\right)^{k}$
ello significa que $PDP^{-1}$ se multiplica $k$ veces
$A^{k}=\left(PDP^{-1}\right)\left(PDP^{-1}\right)\left(PDP^{-1}\right)...$
cambiando la agrupación de los paréntesis se tiene
$A^{k}=PDDD...DP^{-1}$
lo cuál explica la siguiente ecuación:
$A^{k}=PD^{k}P^{-1}$ .......eq8
Sustituyendo eq8 en eq4 tenemos que
$\mathbf{u}_{k}=(PD^{k}P^{-1})\mathbf{u}_{0}$
no olvidar que
$u_{0}=\left[\begin{array}{c} F_{1}\\ F_{0} \end{array}\right]=\left[\begin{array}{c} 1\\ 0 \end{array}\right]$
Luego se trata de otra larga multiplicación de matrices.
$u_{k}=\left[\begin{array}{cc} \lambda_{1} & \lambda_{2}\\ 1 & 1 \end{array}\right]\left[\begin{array}{cc} \lambda_{1}^{k} & 0\\ 0 & \lambda_{2}^{k} \end{array}\right]\frac{1}{\lambda_{1}-\lambda_{2}}\left[\begin{array}{cc} 1 & -\lambda_{2}\\ -1 & \lambda_{1} \end{array}\right]\left[\begin{array}{c} 1\\ 0 \end{array}\right]$
$u_{k}=\left[\begin{array}{cc} \lambda_{1} & \lambda_{2}\\ 1 & 1 \end{array}\right]\left[\begin{array}{cc} \lambda_{1}^{k} & 0\\ 0 & \lambda_{2}^{k} \end{array}\right]\left[\begin{array}{c} 1\\ -1 \end{array}\right]\frac{1}{\lambda_{1}-\lambda_{2}}$
$u_{k}=\left[\begin{array}{cc} \lambda_{1} & \lambda_{2}\\ 1 & 1 \end{array}\right]\left[\begin{array}{c} \lambda_{1}^{k}\\ -\lambda_{2}^{k} \end{array}\right]\frac{1}{\lambda_{1}-\lambda_{2}}$
y finalmente
$\mathbf{u}_{k}=(PD^{k}P^{-1})\mathbf{u}_{0}=\frac{1}{\lambda_{1}-\lambda_{2}}\left[\begin{array}{c} \lambda_{1}^{k+1}-\lambda_{2}^{k+1}\\ \lambda_{1}^{k}-\lambda_{2}^{k} \end{array}\right]$
También $u_{k}=\left[\begin{array}{c} F_{k+1}\\ F_{k} \end{array}\right]$ entonces
$F_{k}=\frac{1}{\lambda_{1}-\lambda_{2}}\left(\lambda_{1}^{k}-\lambda_{2}^{k}\right)$
Por otro lado
$\frac{1}{\lambda_{1}-\lambda_{2}}=\frac{1}{\frac{1+\sqrt{5}}{2}-\frac{1-\sqrt{5}}{2}}=\frac{1}{\frac{(1+\sqrt{5})-(1-\sqrt{5})}{2}}=\frac{1}{\frac{1+\sqrt{5}-1+\sqrt{5})}{2}}=\frac{1}{\frac{2\sqrt{5}}{2}}=\frac{1}{\sqrt{5}}$
que si lo multiplicamos por $\left(\lambda_{1}^{k}-\lambda_{2}^{k}\right)$ equivale a la eq2:
$F_{k}=\frac{1}{\sqrt{5}}\left[\left(\frac{1+\sqrt{5}}{2}\right)^{k}-\left(\frac{1-\sqrt{5}}{2}\right)^{k}\right]$
In [4]:
def fibonacci(k):
if k <= 1:
return k
else:
return fibonacci(k-1) + fibonacci(k-2)
# el indice empieza en cero
for n in range(25):
print("F%s=%s" % (n,fibonacci(n)))
In [5]:
def lucas(k):
if k == 0:
return 1
elif k == 1:
return 2
else:
return lucas(k-1) + lucas(k-2)
# el indice empieza en cero
for n in range(25):
print("F%s=%s" % (n,lucas(n)))
¿Cómo se ven las series de Lucas y Fibonacci juntas?
In [6]:
# importamos bibliotecas para plotear
import matplotlib
import matplotlib.pyplot as plt
# para desplegar los plots en el notebook
%matplotlib inline
p = []
for n in range(11):
p.append([fibonacci(n), lucas(n)])
figura = plt.plot( p )
p
Out[6]:
In [7]:
def fibonacci_formula(k):
return (1.0/sqrt(5))*(((1.0+sqrt(5))/2)**k-((1.0-sqrt(5))/2)**k)
# el indice empieza en cero
for n in range(25):
print("F%i=%i" % (n,fibonacci_formula(n)))
In [8]:
import timeit
# cronometro para implementacion recursiva
tr = timeit.Timer("""def fibonacci(k):
if k <= 1:
return k
else:
return fibonacci(k-1) + fibonacci(k-2)
fibonacci(25)""")
# cronómetro para implementación con formula
tf = timeit.Timer("""from math import sqrt
def fibonacci_formula(k):
return (1.0/sqrt(5))*(((1.0+sqrt(5))/2)**k-((1.0-sqrt(5))/2)**k)
fibonacci_formula(25)""")
# ejecutarlas diez mil veces, medir su tiempo
tiempos=[tr.timeit(10000),
tf.timeit(10000)]
# armar gráfica de barras
plt.bar([0,1], tiempos, align='center', alpha=0.5)
plt.title('10,000 ejecuciones, k=25')
plt.ylabel('milisegundos')
plt.xticks([0,1],['recursivo', 'formula'])
plt.show()
tiempos
Out[8]:
¡Hay un orden de magnitud de diferencia entre los tiempos de ejecución promedio de ambas implementaciones!
Desde el punto de vista computacional la implementación recursiva es más eficiente.
Al parecer es más costoso computacionalmente hacer todas esas raices cuadradas.
Desde el punto de vista del programador, cuando ya se tiene la fórmula cerrada, el programa es mucho más fácil de escribir. ¿Qué será más costoso desde el punto de vista analítico: obtener la fórmula u obtener el algoritmo recursivo?