
I started by shortening your series to 20 realizations, so we could actually see something.
set.seed(420)
x=rnorm(20)
y=rep(NA,length(x))
y[1]=x[1]
for (i in 2:length(x)) y[i]=y[i-1]+x[i]*0.7
Then I simulated five trajectories (each of length six). First, I draw 5 normal random variables with mean 'y[15]. Then I draw another 5 normal random variables with mean
y[16]`. And so forth. Finally, I connect the first set, the second set, up to the fifth set.
pred_index <- 15:20
n_sims <- 5
sd <- 0.2
sims <- sapply(y[pred_index],FUN=function(yy)rnorm(n_sims,mean=yy,sd=sd))
This gives us five trajectories.
plot(y,type="l",lwd=2,ylim=range(c(y,sims)))
for ( jj in 1:n_sims ) lines(pred_index,sims[jj,],col="green")
lines(pred_index,y[pred_index],col="red",lwd=2)
Here are the simulations, with the last "actual" observation in the first column:
> (foo <- cbind(y[pred_index[1]-1],sims))
[,1] [,2] [,3] [,4] [,5] [,6] [,7]
[1,] -1.774497 -2.058771 -2.171860 -1.740961 -1.394201 -0.9514782 -1.6126567
[2,] -1.774497 -2.088186 -2.445334 -1.939621 -1.490676 -1.1473929 -1.5795282
[3,] -1.774497 -1.896928 -1.802883 -2.045344 -1.271107 -0.9969500 -0.9606348
[4,] -1.774497 -2.139482 -2.332128 -1.934507 -1.615830 -1.0684256 -1.1898377
[5,] -1.774497 -2.027070 -2.412320 -1.589246 -1.945217 -1.4398321 -1.1173766
Since this is supposed to be a random walk, we can look at the step-by-step increments within each simulation, by simply taking successive differences between the columns of this matrix:
> sapply(1:ncol(sims),FUN=function(jj)foo[,jj+1]-foo[,jj])
[,1] [,2] [,3] [,4] [,5] [,6]
[1,] -0.2842738 -0.11308883 0.4308988 0.3467603 0.4427225 -0.66117845
[2,] -0.3136891 -0.35714748 0.5057129 0.4489450 0.3432831 -0.43213534
[3,] -0.1224306 0.09404472 -0.2424613 0.7742371 0.2741574 0.03631521
[4,] -0.3649842 -0.19264672 0.3976211 0.3186768 0.5474049 -0.12141211
[5,] -0.2525731 -0.38524925 0.8230736 -0.3559711 0.5053851 0.32245545