Enumeración exacta con sitio absorbente



In [4]:

    
function dist_primer_paso(L, t_final=1000)

    # condicion inicial
    P_vieja = zeros(L+1)
    P_vieja[2] = 1.0

    P_nueva = zeros(L+1)

    # ecuación maestra
    p = 0.4
    q = 0.4
    r = 0.2

    prob_primer_paso_a_0 = zeros(t_final)

    for t in 1:t_final

        for i in 2:L-1
            P_nueva[i] = p*P_vieja[i-1] + q*P_vieja[i+1] + r*P_vieja[i]
        end
        P_nueva[1] = q*P_vieja[2]
        P_nueva[L] = p*P_vieja[L-1] + p*P_vieja[L] + r*P_vieja[L]

        prob_primer_paso_a_0[t] = P_nueva[1]  # la prob que llegó a 0
        P_nueva[1] = 0.0  # se escapó esta probabilidad

        P_nueva, P_vieja = P_vieja, P_nueva

    end
    
    
    prob_primer_paso_a_0
end









    Out[4]:





dist_primer_paso (generic function with 2 methods)



In [5]:

    
distribucion = dist_primer_paso(10);



In [6]:

    
using PyPlot



In [7]:

    
figure(figsize=(6,4))
plot(distribucion, "o-")
xlabel(L"$t$")
ylabel(L"$P(T=t)$")









    












    Out[7]:





PyObject <matplotlib.text.Text object at 0x118155750>



In [8]:

    
figure(figsize=(6,4))
#plot(distribucion, "o-")

for L in [50, 100, 200, 1000]
    distribucion = dist_primer_paso(L, 1000000);
    loglog(distribucion, ".", label="$L")
end
xlabel(L"t")
ylabel(L"P(T=t)")
legend(loc=3)

ylim(1e-10, 1)
#yscale("log")
#xscale("log")









    












    Out[8]:





(1.0e-10,1)



In [9]:

    
figure(figsize=(6,4))
#plot(distribucion, "o-")

for L in [50, 100, 200, 1000]
    distribucion = dist_primer_paso(L, 1000000);
    semilogy(distribucion, ".", label="$L")
end
xlabel(L"t")
ylabel(L"P(T=t)")
legend(loc=0)

ylim(1e-10, 1)
xlim(0, 300000)
#yscale("log")
#xscale("log")









    












    Out[9]:





(0,300000)

Como tenemos una recta (para $L=1000$) con pendiente $-1.6$ (aproximadamente) en log--log, concluimos que

$$P(T=t) \sim t^{-1.6}$$

El resultado exacto analítico es $\sim t^{-3/2}$

Tiempo promedio (valor esperado):

$$\langle T \rangle = \sum_{t=1}^\infty t \cdot P(T=t)$$$$ = \sum_{t=1}^\infty t \cdot t^{-3/2}$$$$ = \sum_{t=1}^\infty t^{-1/2}$$

¡Diverge!

Heurísticamente, es más o menos

$$\int_{t=1}^\infty t^{-1/2} dt = [t^{1/2}]_1^\infty $$

"=$\infty$"

[Regresión lineal en Julia: linreg]



In [1]:

    
?linreg









    



INFO: Loading help data...






    



Base.linreg(x, y) -> [a; b]

   Linear Regression. Returns "a" and "b" such that "a+b*x" is
   the closest line to the given points "(x,y)". In other words,
   this function determines parameters "[a, b]" that minimize the
   squared error between "y" and "a+b*x".

   **Example**:

      using PyPlot;
      x = float([1:12])
      y = [5.5; 6.3; 7.6; 8.8; 10.9; 11.79; 13.48; 15.02; 17.77; 20.81; 22.0; 22.99]
      a, b = linreg(x,y) # Linear regression
      plot(x, y, "o") # Plot (x,y) points
      plot(x, [a+b*i for i in x]) # Plot the line determined by the linear regression

Base.linreg(x, y, w)

   Weighted least-squares linear regression.

Sin embargo, cuando hacemos la simulación Monte Carlo, no vemos que se tarde demasiado.

Podemos hablar de un tiempo "típico":

$$\langle |T|^{\alpha} \rangle^{1/\alpha}$$

$$\langle T^\alpha \rangle \sim \int t^\alpha t^{-3/2} = \int t^{\alpha - 3/2}$$

existe cuando $\alpha - 3/2 < -1$, o sea, $\alpha < 1/2$



In [ ]: