How to prevent overwriting in loop

Hi, I have the following loop. Does anyone know how to avoid the output to be overwritten in each iteration, and instead get x1, x2, ..., x5 as output?

    N=5
    T=5
    for n=1:N
        x=NaN*ones(T,1)
        c=0.3
        sigma=1
        rho=0.98
        x(1)=c/(1-rho)
        for t=1:T
            x(t+1)=c+rho*x(t)+sigma*randn(1)
        end
    end

Answer 1

[Adam]

You can move a lot of that outside of your loop as a first step as follows ( never define constants in a loop ):

N = 5;
T = 5;
c = 0.3;
sigma = 1;
rho = 0.98;
x = NaN(T,1);
x(1) = c / ( 1 - rho );

for t = 1:T
   x(t+1) = c+rho*x(t)+sigma*randn(1);
end

That is not the equivalent of what you did though because it isn't obvious what you are trying to do in the outer loop once all the constants get moved outside the loop.

If you want a 2-d array then you can, as a second step, create something like:

x = NaN(T,N);
randvals = randn(T-1,N);

x(1,:) = c / ( 1 - rho );
x(2:T,:) = c + rho * x(1,1) + sigma * randvals;

instead using a vectorised approach instead of a loop.

[lily]

What I am trying to do is to simulate N sequences of a T steps AR(1) process. I am not too familiar with Monte Carlo simulations, so that is what I was trying to do in the outer loop. What I would like to get as output is N vectors of (T*1) dimension. At the moment, all I can get is one vector x with the values of the last iteration.

[Adam]

Well, my solution should give what you want (although the first row will contain all the same values because it is based on a constant for each of the n simulations).

The only difference is that my result gives you an N * T array which is far better to work with in Matlab than N individual vectors.

[Stephen23]

@lily: you should use Adam's second solution, and learn how to use vectorized code effectively.

[lily]

I do understand the need for vectorisation, as what I have now is taking ages to compute, but I cannot get my head around how to transform my code. I tried applying the examples given on the Mathworks documentation page to my own code (see below) but I cannot get it to work.

N=1000
T=100
x=NaN*ones(T,N)
sigma=1
for n=1:N
  x(1,n)=0
   for t=2:T
     x(t,n)=x(t-1,n)+sigma*randn(1)
   end
end

plot(x)

[Adam]

That looks sufficiently similar to your original that you should still be able to apply the vectorisation I gave above, with minor modifications possibly.

[lily]

OK I understand now. As you probably guessed, I am quite new to MATLAB and I'm struggling a fair deal. Anyway, it is working, thank you very much for your help.

Answer 2

[Stephen23]

Here is a fully vectorized version. First define the parameters:

N = 100;
T = 10;
sigma = 1;

Generate all random steps in an array S, then force the first row of S to be all zero, and finally use cumsum to get the "paths" of each column of S:

S = sigma*randn(T,N);
S(1,:) = 0;
x = cumsum(S,1);

Then plot it just like before:

plot(x)

to produce this figure:

This produces, accepting that the "random" steps are of course different, the same thing as what the loop does:

[lily]

Thank you very much Stephen, this really helps me understand what is going on.