A few comments below
On Mon, Sep 21, 2009 at 3:10 PM, drlucyasher <las...@rvc.ac.uk> wrote:
>
> Hi,
>
> I am trying to write a simulation of the movements of four animals between
> six patches. The movement between patches is based on a first-order Markov
> chain so that the next patch they visit depends on the patch they were in
> before.
>
> I have written code that allows me to simulate the movement of one animal
> but when I add more there seems to be a problem and all chains come back as
> N/A. I can fix the problem if I run the animals as separate sub-routines
> but
> whilst this is adequate for this stage of my model it won't allow me to
> progress to the next stage. The next stage will be to make the animals move
> in reference to each other using a moving average variable. Therefore I
> think I need all animals simulated within one subroutine (loop).
>
> I probably haven't explained thie very well but here is the code I have
> been
> working on....
>
> #creation of a 6*6 matrix with the transition probabilities
>
> transitions<-matrix(0,nrow=6,ncol=6); transitions[1,1]<- (0/6);
> transitions[1,2]<- (0/6); transitions[1,3]<-(0/6);
> transitions[1,4]<-(0/6);transitions[1,5]<-(0/6);
> transitions[1,6]<-(6/6);transitions[2,1]<- (1/6);transitions[2,2]<- (1/6);
> transitions[2,3]<-(1/6); transitions[2,4]<-(1/6);transitions[2,5]<-(1/6);
> transitions[2,6]<-(1/6);transitions[3,1]<- (1/6);transitions[3,2]<- (1/6);
> transitions[3,3]<-(1/6); transitions[3,4]<-(1/6);transitions[3,5]<-(1/6);
> transitions[3,6]<-(1/6);transitions[4,1]<- (1/6);transitions[4,2]<- (1/6);
> transitions[4,3]<-(1/6); transitions[4,4]<-(1/6);transitions[4,5]<-(1/6);
> transitions[4,6]<-(1/6);transitions[5,1]<- (1/6);transitions[5,2]<- (1/6);
> transitions[5,3]<-(1/6); transitions[5,4]<-(1/6);transitions[5,5]<-(1/6);
> transitions[5,6]<-(1/6);transitions[6,1]<- (1/6);transitions[6,2]<- (1/6);
> transitions[6,3]<-(1/6); transitions[6,4]<-(1/6);transitions[6,5]<-(1/6);
> transitions[6,6]<-(1/6); transitions
>
> # start position probalities
> w <- c((1/6), (1/6), (1/6), (1/6), (1/6), (1/6))
>
> #length of the chain
>
> n<-241
>
>
#matrices creation for birds 1-4
> observationb1<-matrix(nrow=n+1,ncol=6)
> observationb2<-matrix(nrow=n+1,ncol=6)
> observationb3<-matrix(nrow=n+1,ncol=6)
> observationb4<-matrix(nrow=n+1,ncol=6)
>
>
#chains starting position
> i<-1;observationb1[1,2]<-sample(1:6, 1, prob=w, repl=T)
> i<-1;observationb2[1,2]<-sample(1:6, 1, prob=w, repl=T)
> i<-1;observationb3[1,2]<-sample(1:6, 1, prob=w, repl=T)
> i<-1;observationb4[1,2]<-sample(1:6, 1, prob=w, repl=T)
>
>
the i <- 1 above is not necessary.
rather use a for loop instead of the while - easier to use.
> while(i<n+1){
>
> Here your i == 1
> stateb1<-observationb1[i,2];
>
>
> probsb1<-transitions[stateb1,1:6];
>
> nextstateb1<-rmultinom(1,1,probsb1);
>
> seqb1<-which.max(nextstateb1);
>
> j<-i+1;
>
> observationb1[j,2]<-seqb1;i<-j
> #-------------------------------------if I close the brackets here this
>
Here you set your i to j, i.e. i+1
> works but I want to inlcude this next bit of code for bird 2
>
> stateb2<-observationb2[i,2];
>
> probsb2<-transitions[stateb2,1:6];
>
> nextstateb2<-rmultinom(1,1,probsb2);
>
> seqb2<-which.max(nextstateb2);
>
> j<-i+1;
>
> observationb2[j,2]<-seqb2;i<-j
>
> #and so on for birds 3 and 4.
>
> Any help you could give me would be much appreciated.
>
>
Hope this helps,
Rainer
> Lucy
> --
> View this message in context:
> http://www.nabble.com/Four-concurrent-Markov-chains-tp25530321p25530321.html
> Sent from the R help mailing list archive at Nabble.com.
>
> ______________________________________________
> R-help@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.
>
--
Rainer M. Krug, Centre of Excellence for Invasion Biology, Stellenbosch
University, South Africa
[[alternative HTML version deleted]]
______________________________________________
R-help@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.
