1D Random Walk

Let us consider a Bernoulli trial. In this, we will consider a bunch of particles at the origin. We will sample a random variable from [-1, +1] and then change the position of each particle. The feeling is that on an average, the particles will have a zero mean because the random variable drawn has zero mean. Let's see what the final distribution of the locations will be.

set(0, "defaultlinelinewidth", 5);
set (0, "defaulttextfontname", "TimesNewRoman")
set (0, "defaulttextfontsize", 20)
set (0, "DefaultAxesFontName", "TimesNewRoman")
set(0, 'DefaultAxesFontSize', 20)

function y=get_steps(N)
    y= randi ([0, 1], [N,1],"int8")*2-1;
end

Np = 100000; Nt = 200;
x = zeros(Np, Nt);

for i=2:1:Nt
    steps = get_steps(Np);
    x(:,i) = x(:,i-1) + steps;
    
    if (mod(i,20)==0)
        s = x(:,i);
        binsp = min(s):2:max(s);
        [counts, bins] = hist(s, binsp);
        plot(bins, counts/trapz(bins, counts));
        hold on;
    end
end
hold off;
xlabel('x'); ylabel('p(x)'); title('Probability density function of locations')

As can be seen, the distribution of the distances tends towards that of a bell-shaped curve. In fact this curve is the normal distribution. The fact that for large times, the distribution tends towards the normal distribution, is a consequence of the Central Limit Theorem. In this particular case, if the step sizes have a finite variance, then the eventual distribution tends towards the normal distribution. Let us now sample from a uniform distribution.

pkg load statistics;
function y=get_steps(N)
    y=random("uniform",-1, 1, [N,1]);
end;
Np = 100000; Nt = 200;
x = zeros(Np, Nt);

for i=2:1:Nt
    steps = get_steps(Np);
    x(:,i) = x(:,i-1) + steps;
    
      if (mod(i,20)==0)
        s = x(:,i);
        binsp = min(s):2:max(s);
        [counts, bins] = hist(s, binsp);
        plot(bins, counts/trapz(bins, counts));
        hold on;
    end
end
hold off;
xlabel('x'); ylabel('p(x)'); title('Probability density function of locations')

Let us now see how the mean and variance scales with time.

Let's also check whether we can sample the random step using some other distribution.

function y=get_steps(N)
    #y= randi ([0, 1], [N,1],"int8")*2-1;
    y= random("uniform",-(3^0.5), 3^0.5, [N,1]);
end

Np = 100000; Nt = 200;
x = zeros(Np, Nt);

for i=2:1:Nt
    steps = get_steps(Np);
    x(:,i) = x(:,i-1) + steps;
    
    if (mod(i,20)==0)
        s = x(:,i);
        binsp = min(s):2:max(s);
        [counts, bins] = hist(s, binsp);
        counts = counts/trapz(bins, counts);
        m = trapz(bins(1:end), bins(1:end).*counts);
        var = trapz(bins(1:end), (bins(1:end)-m).^2.*counts);
        plot(i, m, 'ok');
        plot(i, var, 'sr');
        hold on;
    end
end

ta =0:1:Nt;
plot(ta, ta, '--r');
hold off;
daspect([1 1 1]);

The above plot clearly elucidates a very important property of the Central Limit Theorem. The variance seems to grow linearly with time.

warning('off', 'all'); # turn of warnings.

function y=get_steps(N,a,m)
    #y= randi ([0, 1], [N,1],"int8")*2-1;
    #y= random("uniform",-(3^0.5), 3^0.5, [N,1]);
    y = (1+burrrnd(1,1,a,[N,1]))*m;
end

Np = 100000; Nt = 200;
x = zeros(Np, Nt);

for i=2:1:Nt
    steps = get_steps(Np, 1.2, 1);
    x(:,i) = x(:,i-1) + steps;
    
    if (mod(i,20)==0)
        s = x(:,i);
        binsp = min(s):2:max(s);
        [counts, bins] = hist(s, binsp);
        counts = counts/trapz(bins, counts);
        loglog(bins, counts); 
        hold on;
    end
end
hold off;

for j=1:1:11
    plot(x(j,:));
    hold on;
end
hold off;