[R] Random walk
Giovanni Petris
gpetris at uark.edu
Wed May 12 16:20:44 CEST 2010
Oops, I forgot the 'cumsum' stuff. Here it is again, hopefully working
this time.
> require(MASS)
> fz <- function(n, t, rho){
+ f <- array(dim = c(t, 2, n))
+ V <- matrix(c(1, rho, rho, 1), ncol = 2)
+ for(i in 1 : n){
+ f[,, i] <- apply(mvrnorm(n = t, mu = c(0,0), Sigma = V),
+ 2, cumsum)
+ }
+ return(f)
+ }
>
> set.seed(123)
> out <- fz(1000, 20, rho = 0.7)
> ts.plot(out[,,500]) # one realization of the two random walks
> apply(out, 1, function(x) cor(t(x))[1,2]) # cor at times 1,...,20
[1] 0.7306115 0.7229866 0.7103006 0.6970370 0.7016218
[6] 0.6968639 0.6933828 0.7039743 0.7006013 0.6969842
[11] 0.6935726 0.6853286 0.6809578 0.6754606 0.6804075
[16] 0.6881239 0.6896976 0.6917677 0.6929290 0.6888701
Best,
Giovanni
On Tue, 2010-05-11 at 17:29 -0400, Sergio Andrés Estay Cabrera wrote:
> Dear Giovanni,
>
> Thanks so much for your answer, but your script returns gaussian white
> noise not a random walk, at least the time series generated don't have
> the expected periodogram for a random walk. That's the reason why I use
> cumsum, the sum of a white noise is a easy way produce a random walk
> which periodogram scales at 1/f.
>
> finally I solve my problem with a brute-force algorithm. Searching in
> the web I found some pages about random walks and neural networks where
> people said that this task is very difficult and is in a exploratory
> step nowadays.
>
> Thank you very much for your help
>
> Sergio
>
>
>
>
> Giovanni Petris wrote:
> > Hi Sergio,
> >
> > Your function does not estimate what you want. In fact it does not
> > estimate anything useful. A random walk is not stationary; in
> > particular, the variance at time t is t. Therefore, estimating variances
> > based on one run, averaging over time, does not make any sense. This is
> > what you are doing with the command cor(cumsum(s[,1]),cumsum(s[,2]))
> > (the correlation is obtained by standardizing the covariance, dividing
> > by the square root of the variances).
> >
> > On the other hand, if you replicate your simulation (as you do) and
> > compute the correlation of the simulated random walks at any _fixed_
> > time, you can see that the result is what you expect. Here it is.
> >
> >
> >> require(MASS)
> >> fz <- function(n, t, rho){
> >>
> > + f <- array(dim = c(t, 2, n))
> > + V <- matrix(c(1, rho, rho, 1), ncol = 2)
> > + for(i in 1 : n){
> > + f[,, i] <- mvrnorm(n = t, mu = c(0,0), Sigma = V)
> > + }
> > + return(f)
> > + }
> >
> >>
> >> set.seed(123)
> >> out <- fz(1000, 20, rho = 0.7)
> >> apply(out, 1, function(x) cor(t(x))[1,2]) # cor at times 1,...,20
> >>
> > [1] 0.7306115 0.6810394 0.6879712 0.7095613 0.7325612 0.6922090
> > 0.7007823
> > [8] 0.6828095 0.7144861 0.6932104 0.7061211 0.6781331 0.6898726
> > 0.6926878
> > [15] 0.6949656 0.6995899 0.7524582 0.6912938 0.6931928 0.6820243
> >
> > Best,
> > Giovanni
> >
> > On Mon, 2010-05-10 at 14:55 -0400, Sergio Andrés Estay Cabrera wrote:
> >
> >> Dear R users and specially Albyn and Giovanni,
> >>
> >> thanks for your answers, but in fact I supposed the same at the
> >> beginning of my problem. However, when I generate the data seldom I
> >> obtain the expected correlation. For example using this code:
> >>
> >> fz<-function(n,t,rho){
> >> f<-NULL
> >> for(i in 1:n){
> >> s<-rmvnorm(n=t,mean=c(0,0),sigma=matrix(c(1,rho,rho,1),ncol=2),method='svd')
> >> paso<-cor(cumsum(s[,1]),cumsum(s[,2]))
> >> f<-c(f,paso)}
> >> f<-f
> >> }
> >>
> >> and then plot the histogram of the results, it is possible to observe
> >> that the distribution of the values is asymmetric with most of the
> >> simulations close to 1 when the value of rho is higher than 0.3 and
> >> looks like a uniform distribution with values below 0.3.
> >>
> >> I suspect than the only possibility is using a brute force algorithm.
> >>
> >>
> >> Any advice would be helpful
> >>
> >> Sergio A. Estay
> >> *CASEB *
> >> Departamento de Ecología
> >> Universidad Catolica de Chile
> >>
> >>
> >>
> >>
> >> Albyn Jones wrote:
> >>
> >>> Sums of correlated increments have the same correlation as the original
> >>> variables...
> >>>
> >>> library(mvtnorm)
> >>> X<- matrix(0,nrow=1000,ncol=2)
> >>> for(i in 1:1000){
> >>> Y <- rmvnorm(1000,mean=mu,sigma=S)
> >>> X[i,] <- apply(Y,2,sum)
> >>> }
> >>> cor(Y)
> >>> [,1] [,2]
> >>> [1,] 1.0000000 0.4909281
> >>> [2,] 0.4909281 1.0000000
> >>>
> >>> So, unless you meant that you want the _sample_ correlation to be
> >>> pre-specified, you are all set.
> >>>
> >>> albyn
> >>>
> >>> On Sun, May 09, 2010 at 09:20:25PM -0400, Sergio Andrés Estay Cabrera wrote:
> >>>
> >>>
> >>>> Hi everybody,
> >>>>
> >>>>
> >>>> I am trying to generate two random walks with an specific correlation, for
> >>>> example, two random walks of 200 time steps with a correlation 0.7.
> >>>>
> >>>> I built the random walks with:
> >>>>
> >>>> x<-cumsum(rnorm(200, mean=0,sd=1))
> >>>> y<-cumsum(rnorm(200, mean=0,sd=1))
> >>>>
> >>>> but I don't know how to fix the correlation between them.
> >>>>
> >>>> With white noise is easy to fix the correlation using the function rmvnorm
> >>>> in the package mvtnorm
> >>>>
> >>>> I surfed in the web in the searchable mail archives in the R web site but
> >>>> no references appears.
> >>>>
> >>>> If you have some advices to solve this problems I would be very thankful.
> >>>>
> >>>> Thanks in advance.
> >>>>
> >>>> Sergio A. Estay
> >>>> *CASEB *
> >>>> Departamento de Ecología
> >>>> Universidad Catolica de Chile
> >>>>
> >>>> --
> >>>> “La disciplina no tiene ningún mérito en circunstancias ideales. ” – Habor Mallow
> >>>>
> >>>> ______________________________________________
> >>>> R-help at r-project.org mailing list
> >>>> https://stat.ethz.ch/mailman/listinfo/r-help
> >>>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> >>>> and provide commented, minimal, self-contained, reproducible code.
> >>>>
> >>>>
> >>>>
> >>>
> >>>
> >>
> >
> >
> >
> >
> >
> >
>
>
More information about the R-help
mailing list