Alternating between "for loops"
Locked
 
|
Dear All,
I would like to apply two different "for loops" to each set of four columns of a matrix (the loops here are simplifications of the actual loops I will be running which involve multiple if/else statements). I don't know how to "alternate" between the loops depending on which column is "running through the loop" at the time. ## Set up matrix J <- 10 N <- 10 y <- matrix(0,N,J) ## Apply this loop to the first two of every four columns ([,1:2], [,5:6],[,9:10], etc.) for (q in 1:N){ for(j in 1:J){ if(j %% 2){ y[q,j]=1 }else{y[q,j]=2} } } ## Apply this loop to the next two of every four columns ([,3:4],[,7:8],[,11:12], etc.) for (q in 1:N){ for(j in 1:J){ if(j %% 2){ y[q,j]="A" }else{y[q,j]="B"} } } I want to get something like this (preferably without the quotes): > y [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [1,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [2,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [3,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [4,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [5,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [6,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [7,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [8,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [9,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [10,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" Any help greatly appreciated! Claudia [[alternative HTML version deleted]] ______________________________________________ [hidden email] 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. |
Re: Alternating between "for loops"
|
Hello,
Maybe something along the lines of J <- 10 cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] for(i in which(cols)) { do something } for(i in which(!cols)) { do something else } Hope this helps, Rui Barradas Em 31-07-2012 00:18, Claudia Penaloza escreveu: > Dear All, > I would like to apply two different "for loops" to each set of four columns > of a matrix (the loops here are simplifications of the actual loops I will > be running which involve multiple if/else statements). > I don't know how to "alternate" between the loops depending on which column > is "running through the loop" at the time. > ## Set up matrix > J <- 10 > N <- 10 > y <- matrix(0,N,J) > ## Apply this loop to the first two of every four columns ([,1:2], > [,5:6],[,9:10], etc.) > for (q in 1:N){ > for(j in 1:J){ > if(j %% 2){ > y[q,j]=1 > }else{y[q,j]=2} > } > } > ## Apply this loop to the next two of every four columns > ([,3:4],[,7:8],[,11:12], etc.) > for (q in 1:N){ > for(j in 1:J){ > if(j %% 2){ > y[q,j]="A" > }else{y[q,j]="B"} > } > } > I want to get something like this (preferably without the quotes): > >> y > [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] > [1,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [2,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [3,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [4,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [5,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [6,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [7,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [8,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [9,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [10,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > > > > Any help greatly appreciated! > > Claudia > > [[alternative HTML version deleted]] > > ______________________________________________ > [hidden email] 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. ______________________________________________ [hidden email] 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. |
Re: Alternating between "for loops"
|
Or, assuming you only have 4 different elements :
mat<- matrix(rep(c(1,2,"A", "B"),each=10),10,10, byrow=F) mat2 <- as.data.frame(mat) mat [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [1,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [2,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [3,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [4,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [5,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [6,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [7,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [8,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [9,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" [10,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" mat2 V1 V2 V3 V4 V5 V6 V7 V8 V9 V10 1 1 2 A B 1 2 A B 1 2 2 1 2 A B 1 2 A B 1 2 3 1 2 A B 1 2 A B 1 2 4 1 2 A B 1 2 A B 1 2 5 1 2 A B 1 2 A B 1 2 6 1 2 A B 1 2 A B 1 2 7 1 2 A B 1 2 A B 1 2 8 1 2 A B 1 2 A B 1 2 9 1 2 A B 1 2 A B 1 2 10 1 2 A B 1 2 A B 1 2 Cheers, Eloi On 12-07-30 04:28 PM, Rui Barradas wrote: > Hello, > > Maybe something along the lines of > > J <- 10 > cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] > for(i in which(cols)) { do something } > for(i in which(!cols)) { do something else } > > Hope this helps, > > Rui Barradas > > Em 31-07-2012 00:18, Claudia Penaloza escreveu: >> Dear All, >> I would like to apply two different "for loops" to each set of four >> columns >> of a matrix (the loops here are simplifications of the actual loops I >> will >> be running which involve multiple if/else statements). >> I don't know how to "alternate" between the loops depending on which >> column >> is "running through the loop" at the time. >> ## Set up matrix >> J <- 10 >> N <- 10 >> y <- matrix(0,N,J) >> ## Apply this loop to the first two of every four columns ([,1:2], >> [,5:6],[,9:10], etc.) >> for (q in 1:N){ >> for(j in 1:J){ >> if(j %% 2){ >> y[q,j]=1 >> }else{y[q,j]=2} >> } >> } >> ## Apply this loop to the next two of every four columns >> ([,3:4],[,7:8],[,11:12], etc.) >> for (q in 1:N){ >> for(j in 1:J){ >> if(j %% 2){ >> y[q,j]="A" >> }else{y[q,j]="B"} >> } >> } >> I want to get something like this (preferably without the quotes): >> >>> y >> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] >> [1,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [2,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [3,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [4,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [5,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [6,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [7,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [8,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [9,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [10,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> >> >> >> Any help greatly appreciated! >> >> Claudia >> >> [[alternative HTML version deleted]] >> >> ______________________________________________ >> [hidden email] 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. > > ______________________________________________ > [hidden email] 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. > > -- Eloi Mercier Bioinformatics PhD Student, UBC Paul Pavlidis Lab 2185 East Mall University of British Columbia Vancouver BC V6T1Z4 ______________________________________________ [hidden email] 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. |
|
|
Thank you very much Rui (and Eloi for your input)... this is (the very
simplified version of) what I will end up using: J <- 10 N <- 10 y <- matrix(0,N,J) cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] for(j in which(cols)){ for (q in 1:N){ if(j %% 2){ y[q,j]=1 }else{y[q,j]=2} } } for(j in which(!cols)){ for (q in 1:N){ if(j %% 2){ y[q,j]="A" }else{y[q,j]="B"} } } Cheers, Claudia On Mon, Jul 30, 2012 at 5:38 PM, Mercier Eloi <[hidden email]> wrote: > Or, assuming you only have 4 different elements : > > mat<- matrix(rep(c(1,2,"A", "B"),each=10),10,10, byrow=F) > mat2 <- as.data.frame(mat) > > mat > > [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] > [1,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [2,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [3,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [4,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [5,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [6,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [7,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [8,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [9,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [10,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > > mat2 > V1 V2 V3 V4 V5 V6 V7 V8 V9 V10 > > 1 1 2 A B 1 2 A B 1 2 > 2 1 2 A B 1 2 A B 1 2 > 3 1 2 A B 1 2 A B 1 2 > 4 1 2 A B 1 2 A B 1 2 > 5 1 2 A B 1 2 A B 1 2 > 6 1 2 A B 1 2 A B 1 2 > 7 1 2 A B 1 2 A B 1 2 > 8 1 2 A B 1 2 A B 1 2 > 9 1 2 A B 1 2 A B 1 2 > 10 1 2 A B 1 2 A B 1 2 > > Cheers, > > Eloi > > > On 12-07-30 04:28 PM, Rui Barradas wrote: > >> Hello, >> >> Maybe something along the lines of >> >> J <- 10 >> cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] >> for(i in which(cols)) { do something } >> for(i in which(!cols)) { do something else } >> >> Hope this helps, >> >> Rui Barradas >> >> Em 31-07-2012 00:18, Claudia Penaloza escreveu: >> >>> Dear All, >>> I would like to apply two different "for loops" to each set of four >>> columns >>> of a matrix (the loops here are simplifications of the actual loops I >>> will >>> be running which involve multiple if/else statements). >>> I don't know how to "alternate" between the loops depending on which >>> column >>> is "running through the loop" at the time. >>> ## Set up matrix >>> J <- 10 >>> N <- 10 >>> y <- matrix(0,N,J) >>> ## Apply this loop to the first two of every four columns ([,1:2], >>> [,5:6],[,9:10], etc.) >>> for (q in 1:N){ >>> for(j in 1:J){ >>> if(j %% 2){ >>> y[q,j]=1 >>> }else{y[q,j]=2} >>> } >>> } >>> ## Apply this loop to the next two of every four columns >>> ([,3:4],[,7:8],[,11:12], etc.) >>> for (q in 1:N){ >>> for(j in 1:J){ >>> if(j %% 2){ >>> y[q,j]="A" >>> }else{y[q,j]="B"} >>> } >>> } >>> I want to get something like this (preferably without the quotes): >>> >>> y >>>> >>> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] >>> [1,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [2,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [3,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [4,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [5,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [6,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [7,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [8,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [9,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [10,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> >>> >>> >>> Any help greatly appreciated! >>> >>> Claudia >>> >>> [[alternative HTML version deleted]] >>> >>> ______________________________**________________ >>> [hidden email] mailing list >>> https://stat.ethz.ch/mailman/**listinfo/r-help<https://stat.ethz.ch/mailman/listinfo/r-help> >>> PLEASE do read the posting guide http://www.R-project.org/** >>> posting-guide.html <http://www.R-project.org/posting-guide.html> >>> and provide commented, minimal, self-contained, reproducible code. >>> >> >> ______________________________**________________ >> [hidden email] mailing list >> https://stat.ethz.ch/mailman/**listinfo/r-help<https://stat.ethz.ch/mailman/listinfo/r-help> >> PLEASE do read the posting guide http://www.R-project.org/** >> posting-guide.html <http://www.R-project.org/posting-guide.html> >> and provide commented, minimal, self-contained, reproducible code. >> >> >> > > -- > Eloi Mercier > Bioinformatics PhD Student, UBC > Paul Pavlidis Lab > 2185 East Mall > University of British Columbia > Vancouver BC V6T1Z4 > > [[alternative HTML version deleted]] ______________________________________________ [hidden email] 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. |
Re: Alternating between "for loops"
|
|
On 12-07-31 02:38 PM, Claudia Penaloza wrote:
> Thank you very much Rui (and Eloi for your input)... this is (the very > simplified version of) what I will end up using: Could we get the extended version ? Because right now, I don't know why you need such complicated code that can be done in 1 line. > J <- 10 > N <- 10 > y <- matrix(0,N,J) > cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] > for(j in which(cols)){ > for (q in 1:N){ > if(j %% 2){ > y[q,j]=1 > }else{y[q,j]=2} > } > } > for(j in which(!cols)){ > for (q in 1:N){ > if(j %% 2){ > y[q,j]="A" > }else{y[q,j]="B"} > } > } > J <- 10 N <- 10 y <- matrix(0,N,J) cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] for(j in which(cols)){ if(j %% 2){ y[,j]=1 }else{y[,j]=2} } for(j in which(!cols)){ if(j %% 2){ y[,j]="A" }else{y[,j]="B"} } if you really wants to use this code. Cheers, Eloi > Cheers, > Claudia > On Mon, Jul 30, 2012 at 5:38 PM, Mercier Eloi <[hidden email] > <mailto:[hidden email]>> wrote: > > Or, assuming you only have 4 different elements : > > mat<- matrix(rep(c(1,2,"A", "B"),each=10),10,10, byrow=F) > mat2 <- as.data.frame(mat) > > mat > > [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] > [1,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [2,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [3,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [4,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [5,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [6,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [7,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [8,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [9,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [10,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > > mat2 > V1 V2 V3 V4 V5 V6 V7 V8 V9 V10 > > 1 1 2 A B 1 2 A B 1 2 > 2 1 2 A B 1 2 A B 1 2 > 3 1 2 A B 1 2 A B 1 2 > 4 1 2 A B 1 2 A B 1 2 > 5 1 2 A B 1 2 A B 1 2 > 6 1 2 A B 1 2 A B 1 2 > 7 1 2 A B 1 2 A B 1 2 > 8 1 2 A B 1 2 A B 1 2 > 9 1 2 A B 1 2 A B 1 2 > 10 1 2 A B 1 2 A B 1 2 > > Cheers, > > Eloi > > > On 12-07-30 04:28 PM, Rui Barradas wrote: > > Hello, > > Maybe something along the lines of > > J <- 10 > cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] > for(i in which(cols)) { do something } > for(i in which(!cols)) { do something else } > > Hope this helps, > > Rui Barradas > > Em 31-07-2012 00 <tel:31-07-2012%2000>:18, Claudia Penaloza > escreveu: > > Dear All, > I would like to apply two different "for loops" to each > set of four columns > of a matrix (the loops here are simplifications of the > actual loops I will > be running which involve multiple if/else statements). > I don't know how to "alternate" between the loops > depending on which column > is "running through the loop" at the time. > ## Set up matrix > J <- 10 > N <- 10 > y <- matrix(0,N,J) > ## Apply this loop to the first two of every four columns > ([,1:2], > [,5:6],[,9:10], etc.) > for (q in 1:N){ > for(j in 1:J){ > if(j %% 2){ > y[q,j]=1 > }else{y[q,j]=2} > } > } > ## Apply this loop to the next two of every four columns > ([,3:4],[,7:8],[,11:12], etc.) > for (q in 1:N){ > for(j in 1:J){ > if(j %% 2){ > y[q,j]="A" > }else{y[q,j]="B"} > } > } > I want to get something like this (preferably without the > quotes): > > y > > [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] > [1,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [2,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [3,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [4,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [5,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [6,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [7,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [8,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [9,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [10,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > > > > Any help greatly appreciated! > > Claudia > > [[alternative HTML version deleted]] > > ______________________________________________ > [hidden email] <mailto:[hidden email]> 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. > > > ______________________________________________ > [hidden email] <mailto:[hidden email]> 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. > > > > > -- > Eloi Mercier > Bioinformatics PhD Student, UBC > Paul Pavlidis Lab > 2185 East Mall > University of British Columbia > Vancouver BC V6T1Z4 > > -- Eloi Mercier Bioinformatics PhD Student, UBC Paul Pavlidis Lab 2185 East Mall University of British Columbia Vancouver BC V6T1Z4 [[alternative HTML version deleted]] ______________________________________________ [hidden email] 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. |
Re: Alternating between "for loops"
|
|
Hello,
You're right, and there's a way without loops. Summarizing: J <- 10 N <- 10 cols <- rep(c(TRUE, TRUE, FALSE, FALSE), ceiling(J / 4))[seq_len(J)] #--- 1st way: nested loops y <- matrix(0, N, J) for(j in which(cols)){ for (q in 1:N) y[q, j] <- if(j %% 2) 1 else 2 } for(j in which(!cols)){ for (q in 1:N) y[q, j] <- if(j %% 2) "A" else "B" } #--- 2nd way: one loop only z <- matrix(0, N, J) for(j in which(cols)) z[, j] <- if(j %% 2) 1 else 2 for(j in which(!cols)) z[, j] <- if(j %% 2) "A" else "B" #--- 3rd way: no loops w <- matrix(0, J, N) # reversed order w[cols, ] <- ifelse(which(cols) %% 2, 1, 2) w[!cols, ] <- ifelse(which(!cols) %% 2, "A", "B") w <- t(w) identical(y, z) # compare original to each identical(y, w) # of the other results Em 31-07-2012 23:54, Mercier Eloi escreveu: > On 12-07-31 02:38 PM, Claudia Penaloza wrote: >> Thank you very much Rui (and Eloi for your input)... this is (the very >> simplified version of) what I will end up using: > Could we get the extended version ? Because right now, I don't know why > you need such complicated code that can be done in 1 line. >> J <- 10 >> N <- 10 >> y <- matrix(0,N,J) >> cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] >> for(j in which(cols)){ >> for (q in 1:N){ >> if(j %% 2){ >> y[q,j]=1 >> }else{y[q,j]=2} >> } >> } >> for(j in which(!cols)){ >> for (q in 1:N){ >> if(j %% 2){ >> y[q,j]="A" >> }else{y[q,j]="B"} >> } >> } >> > There is no need for a double loop (on N) : > > J <- 10 > N <- 10 > y <- matrix(0,N,J) > cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] > for(j in which(cols)){ > if(j %% 2){ > y[,j]=1 > }else{y[,j]=2} > } > for(j in which(!cols)){ > if(j %% 2){ > y[,j]="A" > }else{y[,j]="B"} > } > > if you really wants to use this code. > > Cheers, > > Eloi >> Cheers, >> Claudia >> On Mon, Jul 30, 2012 at 5:38 PM, Mercier Eloi <[hidden email] >> <mailto:[hidden email]>> wrote: >> >> Or, assuming you only have 4 different elements : >> >> mat<- matrix(rep(c(1,2,"A", "B"),each=10),10,10, byrow=F) >> mat2 <- as.data.frame(mat) >> >> mat >> >> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] >> [1,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [2,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [3,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [4,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [5,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [6,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [7,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [8,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [9,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [10,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> >> mat2 >> V1 V2 V3 V4 V5 V6 V7 V8 V9 V10 >> >> 1 1 2 A B 1 2 A B 1 2 >> 2 1 2 A B 1 2 A B 1 2 >> 3 1 2 A B 1 2 A B 1 2 >> 4 1 2 A B 1 2 A B 1 2 >> 5 1 2 A B 1 2 A B 1 2 >> 6 1 2 A B 1 2 A B 1 2 >> 7 1 2 A B 1 2 A B 1 2 >> 8 1 2 A B 1 2 A B 1 2 >> 9 1 2 A B 1 2 A B 1 2 >> 10 1 2 A B 1 2 A B 1 2 >> >> Cheers, >> >> Eloi >> >> >> On 12-07-30 04:28 PM, Rui Barradas wrote: >> >> Hello, >> >> Maybe something along the lines of >> >> J <- 10 >> cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] >> for(i in which(cols)) { do something } >> for(i in which(!cols)) { do something else } >> >> Hope this helps, >> >> Rui Barradas >> >> Em 31-07-2012 00 <tel:31-07-2012%2000>:18, Claudia Penaloza >> escreveu: >> >> Dear All, >> I would like to apply two different "for loops" to each >> set of four columns >> of a matrix (the loops here are simplifications of the >> actual loops I will >> be running which involve multiple if/else statements). >> I don't know how to "alternate" between the loops >> depending on which column >> is "running through the loop" at the time. >> ## Set up matrix >> J <- 10 >> N <- 10 >> y <- matrix(0,N,J) >> ## Apply this loop to the first two of every four columns >> ([,1:2], >> [,5:6],[,9:10], etc.) >> for (q in 1:N){ >> for(j in 1:J){ >> if(j %% 2){ >> y[q,j]=1 >> }else{y[q,j]=2} >> } >> } >> ## Apply this loop to the next two of every four columns >> ([,3:4],[,7:8],[,11:12], etc.) >> for (q in 1:N){ >> for(j in 1:J){ >> if(j %% 2){ >> y[q,j]="A" >> }else{y[q,j]="B"} >> } >> } >> I want to get something like this (preferably without the >> quotes): >> >> y >> >> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] >> [1,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [2,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [3,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [4,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [5,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [6,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [7,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [8,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [9,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [10,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> >> >> >> Any help greatly appreciated! >> >> Claudia >> >> [[alternative HTML version deleted]] >> >> ______________________________________________ >> [hidden email] <mailto:[hidden email]> 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. >> >> >> ______________________________________________ >> [hidden email] <mailto:[hidden email]> 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. >> >> >> >> >> -- >> Eloi Mercier >> Bioinformatics PhD Student, UBC >> Paul Pavlidis Lab >> 2185 East Mall >> University of British Columbia >> Vancouver BC V6T1Z4 >> >> > ______________________________________________ [hidden email] 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. |
Re: Alternating between "for loops"
|
|
In reply to this post by Mercier Eloi
Your code almost gave me a headache. :-/
There are a lot of unnecessary tests and conditions. I tried to break down the code and write the tests that have been done when assigning a variable. I simplified your the first part but cannot guarantee that all criteria are met. #####COMMENTS##### for(j in which(dtp)){ for (q in 1:N){ (observable=1) #prefer TRUE/FALSE for boolean operators if(j %% 2) { survive=runif(1,0,1) if(survive<=S) { stay=runif(1,0,1) ####Is observable ?### if (observable==1) #if (observable); but can observable!=1 here ???? It is defined as observable=1 above { if(stay<=PsiAA) { observable=1 #stay<=PsiAA && observable && survive <= S && j%%2 #not needed; observable is already equals to 1 }else{ observable=0 #stay>PsiAA && observable && survive <= S && j%%2 #ie. observable = !observable; ie. was 1, now is 0 } }else{ return=runif(1,0,1) if(return<=PsiaA) { observable=1 #return<=PsiAA && !observable && survive <= S && j%%2 #ie. observable = !observable; ie. was 0, now is 1 }else{ observable=0 #return>PsiAA && !observable && survive <= S && j%%2 #not needed; observable is already equals to 0 } } ####### if(observable==1) #you could move this block above, where you assign the value to observable { capture=runif(1,0,1) if(capture<=p) {y[q,j]="A"} #capture<=p && observable && survive <= S && j%%2 } ####### }else{ deado=runif(1,0,1) if(deado<=PsiAd) { y[q,j]="d" #deado<=PsiAd && survive > S && j%%2 }else{ (y[q,j]="D") #deado<PsiAd && survive > S && j%%2 } #use ifelse instead ####### if(y[q,j]%in%c("D","d")) # test not needed; y[q,j] will always be either D or d according the assignment above { break #survive > S && j%%2 #the use of break is not recommended, especially in this case with so many loops - it's hard to tell where the break will exit } ####### } }else{ if(observable==1) { recap=runif(1,0,1) if(recap<=c) { y[q,j]="A"#recap<=c && observable && !(j%%2) } } } } } ####SUGGESTED CODE###### dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] for(j in which(dtp)){ for (q in 1:N){ observable <-TRUE if(j %% 2) { survive=runif(1,0,1) if(survive<=S) { stay=runif(1,0,1) return=runif(1,0,1) ####Is observable ?### if (return<=PsiAA || stay>PsiAA) { observable = !observable } #otherwise, observable stays identical ####### if(observable) { capture=runif(1,0,1) if(capture<=p) {y[q,j]="A"} #capture<=p && observable && survive <= S && j%%2 } ####### }else{ deado=runif(1,0,1) y[q,j]=ifelse(deado<=PsiAd, "d", "D") #deado<=PsiAd && survive > S && j%%2 break ####### } }else{ recap=runif(1,0,1) if(recap<=c && observable) {y[q,j]="A"}#recap<=c && observable && !(j%%2) } } } Regards, Eloi On 12-08-01 09:47 AM, Claudia Penaloza wrote: > I will accept any help you can give me, especially to free myself of > > SAS-brain inefficiencies... > My supervisor knows SAS not R, which may explain my code. > What I'm actually trying to do is simulate mark-recapture data with a > peculiar structure. > It is a multistate robust design model with dummy time periods... it > will basically be a matrix with a succession of letters (for different > states/age classes) and zeros that are generated following a certain > set of conditions (probability statements; drawing from a random > uniform distribution if an event happens). > Normally, survival and transition to another state occur between > primary sampling periods (in my R example, every two columns.. between > [,2] and [,3]) but in the "dummy time period" model survival occurs > first and then transition to another state, which is why I need my > "conditions" to alternate. Additionally, the robust design has > secondary sampling sessions that are within the same year, i.e., > survival is assumed = 1, which is why my columns are paired. Each > state can also be in an unobservable state at any given time... all of > these details complicate the coding. > Below I've pasted what I have thus far... I have not debugged it yet > (the second loop isn't working yet and the first loop still has some > glitches). > Further below is properly working code for a robust design without > dummy time periods (it also lacks the dead states the dummy time > period model has, which happen to be part of the glitchy-ness). > Is there a better (more R-ish) way of doing this? > Thank you for all your help, > Claudia > ################################################################ > # Robust design with Markovian emigration and dummy time periods > ################################################################ > J = 24 > N = 10 > S = .9 > PsiAD = 0.06 > PsiAd = 0.04 > PsiAA = .4 > PsiaA = .3 > p = .5 > c = p > y <- matrix(0,N,J) > y[,1]="A" > dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] > for(j in which(dtp)){ > for (q in 1:N){ > (observable=1) > if(j %% 2){survive=runif(1,0,1) > if(survive<=S){ > stay=runif(1,0,1) > if (observable==1){ > if(stay<=PsiAA){ > observable=1 > }else{observable=0} > }else{ > return=runif(1,0,1) > if(return<=PsiaA){ > observable=1 > }else{observable=0}} > if(observable==1){ > capture=runif(1,0,1) > if(capture<=p){ > y[q,j]="A"} > }}else{ > deado=runif(1,0,1) > if(deado<=PsiAd){ > y[q,j]="d" > }else(y[q,j]="D") > if(y[q,j]%in%c("D","d")){ > break > } > } > }else{ > if(observable==1){ > recap=runif(1,0,1) > if(recap<=c){ > y[q,j]="A" > } > } > } > } > } > > for(j in which(!dtp)){ > for (q in 1:N){ > if(j %% 2){ > stay=runif(1,0,1) > if(observable==1){ > if(stay<=PsiAA){ > observable=1 > }else{observable=0} > }else{ > return=runif(1,0,1) > if(return<=PsiaA){ > observable=1 > }else{observable=0} > } > if(observable==1){ > capture=runif(1,0,1) > if(capture<=p){ > y[q,j]="A"} > }}else { > if(observable==1){ > recap=runif(1,0,1) > if(recap<=c){ > y[q,j]="A" > } > } > } > } > } > ########################################### > ### Robust design with markovian emigration > ########################################### > J = 24 > N = 1000 > S = .9 > PsiOO = .4 > PsiUO = .3 > p = .5 > c = p > y <- matrix(0,N,J) > y[,1]=1 > for (q in 1:N){ > (observable=1) > for(j in 2:J){ > if(j %% 2){surviv=runif(1,0,1) > if(surviv<=S){ > stay=runif(1,0,1) > if(observable==1){ > if(stay<=PsiOO){ > observable=1 > }else{observable=0} > }else{ > return=runif(1,0,1) > if(return<=PsiUO){ > observable=1 > }else{observable=0}} > if(observable==1){ > cap=runif(1,0,1) > if(cap<=p){ > y[q,j]=1} > }else y[q,j]=0 > }else{break} > }else{ > if (observable==1){ > recap=runif(1,0,1) > if (recap<=c){ > y[q,j]=1} > else{y[q,j]=0} > } > } > } > } > On Tue, Jul 31, 2012 at 7:16 PM, David Winsemius > <[hidden email] <mailto:[hidden email]>> wrote: > > Claudia; > > The loop constructs will keep you mired in SAS-brain > inefficiencies forever. Please, listen more carefully to Mercier. > > -- > David > > > > On Jul 31, 2012, at 3:54 PM, Mercier Eloi wrote: > > On 12-07-31 02:38 PM, Claudia Penaloza wrote: > > Thank you very much Rui (and Eloi for your input)... this > is (the very > simplified version of) what I will end up using: > > Could we get the extended version ? Because right now, I don't > know why > you need such complicated code that can be done in 1 line. > > J <- 10 > N <- 10 > y <- matrix(0,N,J) > cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] > for(j in which(cols)){ > for (q in 1:N){ > if(j %% 2){ > y[q,j]=1 > }else{y[q,j]=2} > } > } > for(j in which(!cols)){ > for (q in 1:N){ > if(j %% 2){ > y[q,j]="A" > }else{y[q,j]="B"} > } > } > > There is no need for a double loop (on N) : > > J <- 10 > N <- 10 > y <- matrix(0,N,J) > cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] > for(j in which(cols)){ > if(j %% 2){ > y[,j]=1 > }else{y[,j]=2} > } > for(j in which(!cols)){ > if(j %% 2){ > y[,j]="A" > }else{y[,j]="B"} > } > > if you really wants to use this code. > > Cheers, > > Eloi > > Cheers, > Claudia > On Mon, Jul 30, 2012 at 5:38 PM, Mercier Eloi > <[hidden email] <mailto:[hidden email]> > <mailto:[hidden email] > <mailto:[hidden email]>>> wrote: > > Or, assuming you only have 4 different elements : > > mat<- matrix(rep(c(1,2,"A", "B"),each=10),10,10, byrow=F) > mat2 <- as.data.frame(mat) > > mat > > [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] > [1,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [2,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [3,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [4,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [5,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [6,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [7,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [8,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [9,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > [10,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" > > mat2 > V1 V2 V3 V4 V5 V6 V7 V8 V9 V10 > > 1 1 2 A B 1 2 A B 1 2 > 2 1 2 A B 1 2 A B 1 2 > 3 1 2 A B 1 2 A B 1 2 > 4 1 2 A B 1 2 A B 1 2 > 5 1 2 A B 1 2 A B 1 2 > 6 1 2 A B 1 2 A B 1 2 > 7 1 2 A B 1 2 A B 1 2 > 8 1 2 A B 1 2 A B 1 2 > 9 1 2 A B 1 2 A B 1 2 > 10 1 2 A B 1 2 A B 1 2 > > Cheers, > > Eloi > > > On 12-07-30 04:28 PM, Rui Barradas wrote: > > Hello, > > Maybe something along the lines of > > J <- 10 > cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] > for(i in which(cols)) { do something } > for(i in which(!cols)) { do something else } > > Hope this helps, > > Rui Barradas > > Em 31-07-2012 00 <tel:31-07-2012%2000> > <tel:31-07-2012%2000>:18, Claudia Penaloza > > escreveu: > > Dear All, > I would like to apply two different "for loops" > to each > set of four columns > of a matrix (the loops here are simplifications > of the > actual loops I will > be running which involve multiple if/else > statements). > I don't know how to "alternate" between the loops > depending on which column > is "running through the loop" at the time. > ## Set up matrix > J <- 10 > N <- 10 > y <- matrix(0,N,J) > ## Apply this loop to the first two of every > four columns > ([,1:2], > [,5:6],[,9:10], etc.) > for (q in 1:N){ > for(j in 1:J){ > if(j %% 2){ > y[q,j]=1 > }else{y[q,j]=2} > } > } > ## Apply this loop to the next two of every > four columns > ([,3:4],[,7:8],[,11:12], etc.) > for (q in 1:N){ > for(j in 1:J){ > if(j %% 2){ > y[q,j]="A" > }else{y[q,j]="B"} > } > } > I want to get something like this (preferably > without the > quotes): > > y > > [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] > [1,] "1" "2" "A" "B" "1" "2" "A" "B" > "1" "2" > [2,] "1" "2" "A" "B" "1" "2" "A" "B" > "1" "2" > [3,] "1" "2" "A" "B" "1" "2" "A" "B" > "1" "2" > [4,] "1" "2" "A" "B" "1" "2" "A" "B" > "1" "2" > [5,] "1" "2" "A" "B" "1" "2" "A" "B" > "1" "2" > [6,] "1" "2" "A" "B" "1" "2" "A" "B" > "1" "2" > [7,] "1" "2" "A" "B" "1" "2" "A" "B" > "1" "2" > [8,] "1" "2" "A" "B" "1" "2" "A" "B" > "1" "2" > [9,] "1" "2" "A" "B" "1" "2" "A" "B" > "1" "2" > [10,] "1" "2" "A" "B" "1" "2" "A" "B" > "1" "2" > > Any help greatly appreciated! > > Claudia > > > -- > Eloi Mercier > Bioinformatics PhD Student, UBC > Paul Pavlidis Lab > 2185 East Mall > University of British Columbia > Vancouver BC V6T1Z4 > > > > > David Winsemius, MD > Alameda, CA, USA > > -- Eloi Mercier Bioinformatics PhD Student, UBC Paul Pavlidis Lab 2185 East Mall University of British Columbia Vancouver BC V6T1Z4 [[alternative HTML version deleted]] ______________________________________________ [hidden email] 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. |
Re: Alternating between "for loops"
|
|
Hello,
Em 01-08-2012 20:02, Mercier Eloi escreveu: > Your code almost gave me a headache. :-/ Agree. There's a good way of avoiding headaches, package formatR, function tidy.source. My simplification is different in many places. A common one is to treat 'observable' as a logical variable, not as a numeric one. I've simplified some tests, the syntax in some places (namely, the condition 'if' is a function). And eliminated some calls to runif(), whenever they could be done only once. I think my code is the equivalent of Claudia's and I've worked it all. here it goes but without comments. # # # Robust design # # # with Markovian # # # emigration and # # # dummy time # # # periods J <- 24 N <- 10 S <- 0.9 PsiAD <- 0.06 PsiAd <- 0.04 PsiAA <- 0.4 PsiaA <- 0.3 p <- 0.5 c <- p y <- matrix(0, N, J) y[, 1] <- "A" dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] for (j in which(dtp)) { for (q in 1:N) { (observable <- TRUE) if (j %% 2) { survive <- runif(1, 0, 1) decide <- runif(1, 0, 1) if (survive <= S) { if (observable) { observable <- (decide <= PsiAA) } else { observable <- (decide <= PsiaA) } if (observable) { if (runif(1, 0, 1) <= p) y[q, j] <- "A" } } else { y[q, j] <- if (decide <= PsiAd) "d" else "D" break } } else { if (observable) { if (runif(1, 0, 1) <= c) y[q, j] <- "A" } } } } for (j in which(!dtp)) { for (q in 1:N) { if (j %% 2) { decide <- runif(1, 0, 1) if (observable) { observable <- decide <= PsiAA } else { observable <- decide <= PsiaA } if (observable) { if (runif(1, 0, 1) <= p) y[q, j] <- "A" } } else { if (observable) { if (runif(1, 0, 1) <= c) y[q, j] <- "A" } } } } # # # Robust design with markovian # # # emigration J <- 24 N <- 1000 S <- 0.9 PsiOO <- 0.4 PsiUO <- 0.3 p <- 0.5 c <- p y <- matrix(0, N, J) y[, 1] <- 1 for (q in 1:N) { (observable <- TRUE) for (j in 2:J) { if (j %% 2) { surviv <- runif(1, 0, 1) if (surviv <= S) { decide <- runif(1, 0, 1) if (observable) { observable <- (decide <= PsiOO) } else { observable <- (decide <= PsiUO) } if (observable) { y[q, j] <- if (runif(1, 0, 1) <= p) 1 else 0 } } else { break } } else { if (observable) { y[q, j] <- if (runif(1, 0, 1) <= c) 1 else 0 } } } } Hope this helps, Rui Barradas > There are a lot of unnecessary tests and conditions. I tried to break > down the code and write the tests that have been done when assigning a > variable. I simplified your the first part but cannot guarantee that all > criteria are met. > > #####COMMENTS##### > for(j in which(dtp)){ > for (q in 1:N){ > (observable=1) #prefer TRUE/FALSE for boolean operators > if(j %% 2) > { > survive=runif(1,0,1) > if(survive<=S) > { > stay=runif(1,0,1) > ####Is observable ?### > if (observable==1) #if (observable); but can > observable!=1 here ???? It is defined as observable=1 above > { > if(stay<=PsiAA) > { > observable=1 #stay<=PsiAA && observable && > survive <= S && j%%2 > #not needed; observable is already equals to 1 > }else{ > observable=0 #stay>PsiAA && observable && > survive <= S && j%%2 > #ie. observable = !observable; ie. was 1, now is 0 > } > }else{ > return=runif(1,0,1) > if(return<=PsiaA) > { > observable=1 #return<=PsiAA && !observable && > survive <= S && j%%2 > #ie. observable = !observable; ie. was 0, now is 1 > }else{ > observable=0 #return>PsiAA && !observable && > survive <= S && j%%2 > #not needed; observable is already equals to 0 > } > } > ####### > if(observable==1) #you could move this block above, > where you assign the value to observable > { > capture=runif(1,0,1) > if(capture<=p) > {y[q,j]="A"} #capture<=p && observable && > survive <= S && j%%2 > } > ####### > }else{ > deado=runif(1,0,1) > if(deado<=PsiAd) > { > y[q,j]="d" #deado<=PsiAd && survive > S && j%%2 > }else{ > (y[q,j]="D") #deado<PsiAd && survive > S && j%%2 > } #use ifelse instead > ####### > if(y[q,j]%in%c("D","d")) # test not needed; y[q,j] will > always be either D or d according the assignment above > { > break #survive > S && j%%2 > #the use of break is not recommended, especially in > this case with so many loops - it's hard to tell where the break will exit > } > ####### > } > }else{ > if(observable==1) > { > recap=runif(1,0,1) > if(recap<=c) > { > y[q,j]="A"#recap<=c && observable && !(j%%2) > } > } > } > } > } > > > ####SUGGESTED CODE###### > > dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] > for(j in which(dtp)){ > for (q in 1:N){ > observable <-TRUE > if(j %% 2) > { > survive=runif(1,0,1) > if(survive<=S) > { > stay=runif(1,0,1) > return=runif(1,0,1) > ####Is observable ?### > if (return<=PsiAA || stay>PsiAA) > { > observable = !observable > } #otherwise, observable stays identical > ####### > if(observable) > { > capture=runif(1,0,1) > if(capture<=p) > {y[q,j]="A"} #capture<=p && observable && > survive <= S && j%%2 > } > ####### > }else{ > deado=runif(1,0,1) > y[q,j]=ifelse(deado<=PsiAd, "d", "D") #deado<=PsiAd && > survive > S && j%%2 > break > ####### > } > }else{ > recap=runif(1,0,1) > if(recap<=c && observable) > {y[q,j]="A"}#recap<=c && observable && !(j%%2) > } > } > } > > Regards, > > Eloi > > On 12-08-01 09:47 AM, Claudia Penaloza wrote: >> I will accept any help you can give me, especially to free myself of >>> SAS-brain inefficiencies... >> My supervisor knows SAS not R, which may explain my code. >> What I'm actually trying to do is simulate mark-recapture data with a >> peculiar structure. >> It is a multistate robust design model with dummy time periods... it >> will basically be a matrix with a succession of letters (for different >> states/age classes) and zeros that are generated following a certain >> set of conditions (probability statements; drawing from a random >> uniform distribution if an event happens). >> Normally, survival and transition to another state occur between >> primary sampling periods (in my R example, every two columns.. between >> [,2] and [,3]) but in the "dummy time period" model survival occurs >> first and then transition to another state, which is why I need my >> "conditions" to alternate. Additionally, the robust design has >> secondary sampling sessions that are within the same year, i.e., >> survival is assumed = 1, which is why my columns are paired. Each >> state can also be in an unobservable state at any given time... all of >> these details complicate the coding. >> Below I've pasted what I have thus far... I have not debugged it yet >> (the second loop isn't working yet and the first loop still has some >> glitches). >> Further below is properly working code for a robust design without >> dummy time periods (it also lacks the dead states the dummy time >> period model has, which happen to be part of the glitchy-ness). >> Is there a better (more R-ish) way of doing this? >> Thank you for all your help, >> Claudia >> ################################################################ >> # Robust design with Markovian emigration and dummy time periods >> ################################################################ >> J = 24 >> N = 10 >> S = .9 >> PsiAD = 0.06 >> PsiAd = 0.04 >> PsiAA = .4 >> PsiaA = .3 >> p = .5 >> c = p >> y <- matrix(0,N,J) >> y[,1]="A" >> dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] >> for(j in which(dtp)){ >> for (q in 1:N){ >> (observable=1) >> if(j %% 2){survive=runif(1,0,1) >> if(survive<=S){ >> stay=runif(1,0,1) >> if (observable==1){ >> if(stay<=PsiAA){ >> observable=1 >> }else{observable=0} >> }else{ >> return=runif(1,0,1) >> if(return<=PsiaA){ >> observable=1 >> }else{observable=0}} >> if(observable==1){ >> capture=runif(1,0,1) >> if(capture<=p){ >> y[q,j]="A"} >> }}else{ >> deado=runif(1,0,1) >> if(deado<=PsiAd){ >> y[q,j]="d" >> }else(y[q,j]="D") >> if(y[q,j]%in%c("D","d")){ >> break >> } >> } >> }else{ >> if(observable==1){ >> recap=runif(1,0,1) >> if(recap<=c){ >> y[q,j]="A" >> } >> } >> } >> } >> } >> >> for(j in which(!dtp)){ >> for (q in 1:N){ >> if(j %% 2){ >> stay=runif(1,0,1) >> if(observable==1){ >> if(stay<=PsiAA){ >> observable=1 >> }else{observable=0} >> }else{ >> return=runif(1,0,1) >> if(return<=PsiaA){ >> observable=1 >> }else{observable=0} >> } >> if(observable==1){ >> capture=runif(1,0,1) >> if(capture<=p){ >> y[q,j]="A"} >> }}else { >> if(observable==1){ >> recap=runif(1,0,1) >> if(recap<=c){ >> y[q,j]="A" >> } >> } >> } >> } >> } >> ########################################### >> ### Robust design with markovian emigration >> ########################################### >> J = 24 >> N = 1000 >> S = .9 >> PsiOO = .4 >> PsiUO = .3 >> p = .5 >> c = p >> y <- matrix(0,N,J) >> y[,1]=1 >> for (q in 1:N){ >> (observable=1) >> for(j in 2:J){ >> if(j %% 2){surviv=runif(1,0,1) >> if(surviv<=S){ >> stay=runif(1,0,1) >> if(observable==1){ >> if(stay<=PsiOO){ >> observable=1 >> }else{observable=0} >> }else{ >> return=runif(1,0,1) >> if(return<=PsiUO){ >> observable=1 >> }else{observable=0}} >> if(observable==1){ >> cap=runif(1,0,1) >> if(cap<=p){ >> y[q,j]=1} >> }else y[q,j]=0 >> }else{break} >> }else{ >> if (observable==1){ >> recap=runif(1,0,1) >> if (recap<=c){ >> y[q,j]=1} >> else{y[q,j]=0} >> } >> } >> } >> } >> On Tue, Jul 31, 2012 at 7:16 PM, David Winsemius >> <[hidden email] <mailto:[hidden email]>> wrote: >> >> Claudia; >> >> The loop constructs will keep you mired in SAS-brain >> inefficiencies forever. Please, listen more carefully to Mercier. >> >> -- >> David >> >> >> >> On Jul 31, 2012, at 3:54 PM, Mercier Eloi wrote: >> >> On 12-07-31 02:38 PM, Claudia Penaloza wrote: >> >> Thank you very much Rui (and Eloi for your input)... this >> is (the very >> simplified version of) what I will end up using: >> >> Could we get the extended version ? Because right now, I don't >> know why >> you need such complicated code that can be done in 1 line. >> >> J <- 10 >> N <- 10 >> y <- matrix(0,N,J) >> cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] >> for(j in which(cols)){ >> for (q in 1:N){ >> if(j %% 2){ >> y[q,j]=1 >> }else{y[q,j]=2} >> } >> } >> for(j in which(!cols)){ >> for (q in 1:N){ >> if(j %% 2){ >> y[q,j]="A" >> }else{y[q,j]="B"} >> } >> } >> >> There is no need for a double loop (on N) : >> >> J <- 10 >> N <- 10 >> y <- matrix(0,N,J) >> cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] >> for(j in which(cols)){ >> if(j %% 2){ >> y[,j]=1 >> }else{y[,j]=2} >> } >> for(j in which(!cols)){ >> if(j %% 2){ >> y[,j]="A" >> }else{y[,j]="B"} >> } >> >> if you really wants to use this code. >> >> Cheers, >> >> Eloi >> >> Cheers, >> Claudia >> On Mon, Jul 30, 2012 at 5:38 PM, Mercier Eloi >> <[hidden email] <mailto:[hidden email]> >> <mailto:[hidden email] >> <mailto:[hidden email]>>> wrote: >> >> Or, assuming you only have 4 different elements : >> >> mat<- matrix(rep(c(1,2,"A", "B"),each=10),10,10, byrow=F) >> mat2 <- as.data.frame(mat) >> >> mat >> >> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] >> [1,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [2,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [3,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [4,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [5,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [6,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [7,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [8,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [9,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [10,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> >> mat2 >> V1 V2 V3 V4 V5 V6 V7 V8 V9 V10 >> >> 1 1 2 A B 1 2 A B 1 2 >> 2 1 2 A B 1 2 A B 1 2 >> 3 1 2 A B 1 2 A B 1 2 >> 4 1 2 A B 1 2 A B 1 2 >> 5 1 2 A B 1 2 A B 1 2 >> 6 1 2 A B 1 2 A B 1 2 >> 7 1 2 A B 1 2 A B 1 2 >> 8 1 2 A B 1 2 A B 1 2 >> 9 1 2 A B 1 2 A B 1 2 >> 10 1 2 A B 1 2 A B 1 2 >> >> Cheers, >> >> Eloi >> >> >> On 12-07-30 04:28 PM, Rui Barradas wrote: >> >> Hello, >> >> Maybe something along the lines of >> >> J <- 10 >> cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] >> for(i in which(cols)) { do something } >> for(i in which(!cols)) { do something else } >> >> Hope this helps, >> >> Rui Barradas >> >> Em 31-07-2012 00 <tel:31-07-2012%2000> >> <tel:31-07-2012%2000>:18, Claudia Penaloza >> >> escreveu: >> >> Dear All, >> I would like to apply two different "for loops" >> to each >> set of four columns >> of a matrix (the loops here are simplifications >> of the >> actual loops I will >> be running which involve multiple if/else >> statements). >> I don't know how to "alternate" between the loops >> depending on which column >> is "running through the loop" at the time. >> ## Set up matrix >> J <- 10 >> N <- 10 >> y <- matrix(0,N,J) >> ## Apply this loop to the first two of every >> four columns >> ([,1:2], >> [,5:6],[,9:10], etc.) >> for (q in 1:N){ >> for(j in 1:J){ >> if(j %% 2){ >> y[q,j]=1 >> }else{y[q,j]=2} >> } >> } >> ## Apply this loop to the next two of every >> four columns >> ([,3:4],[,7:8],[,11:12], etc.) >> for (q in 1:N){ >> for(j in 1:J){ >> if(j %% 2){ >> y[q,j]="A" >> }else{y[q,j]="B"} >> } >> } >> I want to get something like this (preferably >> without the >> quotes): >> >> y >> >> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] >> [1,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [2,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [3,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [4,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [5,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [6,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [7,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [8,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [9,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [10,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> >> Any help greatly appreciated! >> >> Claudia >> >> >> -- >> Eloi Mercier >> Bioinformatics PhD Student, UBC >> Paul Pavlidis Lab >> 2185 East Mall >> University of British Columbia >> Vancouver BC V6T1Z4 >> >> >> >> >> David Winsemius, MD >> Alameda, CA, USA >> >> > ______________________________________________ [hidden email] 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. |
Re: Alternating between "for loops"
|
In reply to this post by Mercier Eloi
On Aug 1, 2012, at 12:02 PM, Mercier Eloi wrote: > Your code almost gave me a headache. :-/ I had a similar reaction. However, my approach might have been to request a more complete verbal description of the data structures being operated on and the methods and assumptions being used. Generally it is going to be easier to go from a procedural description to good R code than it is to go from a SAS Data Proc ... even if it were tested and debugged... and yours was not even debugged. > ################################################################ > # Robust design with Markovian emigration and dummy time periods > ################################################################ > J = 24 > N = 10 > S = .9 > PsiAD = 0.06 > PsiAd = 0.04 > PsiAA = .4 > PsiaA = .3 > p = .5 > c = p > y <- matrix(0,N,J) # So is 'y' a matrix of states where the N rows are levels and the J columns are discrete times? # Actually I decided that context suggested you were using letters to represent states. > y[,1]="A" # So we start with N <something>'s in state "A"? # It seems as though it might be the case that every row is independent of the others. # .. and you are performing ( replicate()-ing in R terminology) this test N times # states: #("A" "D") #("a" "d") #transitions: matrix( c( runif(1, 0, 1) <= 0.4, # A -> A runif(1,0,1) <= 0.3, # a -> A runif(1,0,1) <= 0.06. # A -> D runif(1,0,1) <= 0.04 # A -> d) ) # What is state "a"? # How do you get from A to a? # can you get from D or d to A or a? # Not sure I have gotten the model assumptions down. # Is D" or "d" an absorbing state such as "dead or "Dead"? > dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] #Presumably the "dummy time periods" > for(j in which(dtp)){ # j will be c(1,2, 5,6, 9,10 , ..., 21,22) > for (q in 1:N){ # This can almost surely become vectorized > (observable=1) > if(j %% 2){survive=runif(1,0,1) > if(survive<=S){ > stay=runif(1,0,1) > if (observable==1){ # It would help if the concept of "observable" were explained # Is this something to do with the capture-recapture observables? > if{ > observable=1 > }else{observable=0} > }else{ # After allowing for Mercier's astute observation that observable will always be 1, I'm wondering if it can be replaced with just this code (and not need the loop surrounding it.) observable <- (stay<=PsiAA) #-------------- > return=runif(1,0,1) # better NOT to use the name "return" since it is a crucial R function name. > if(return<=PsiaA){ > observable=1 > }else{observable=0}} #------- perhaps: randret <- return=runif(N,0,1) observable <- randret <= PsiaA # ------------------- > if(observable==1){ > capture=runif(1,0,1) > if(capture<=p){ #---------- perhaps: randcap <- runif(N, 0,1) y[ randcap <= p, j] <- "A" #That would replace with "A" (within the j-column) ... # only the rows in 'y' for which randcap were less than the randcap threshold. #------------ > y[q,j]="A"} > }}else{ > deado=runif(1,0,1) > if(deado<=PsiAd){ > y[q,j]="d" > }else(y[q,j]="D") # -------Perhaps deado <- runif(N, 0,1) y[ , j] <- ifelse( deado<=PsiAd, "d", "D") #------------------------ > if(y[q,j]%in%c("D","d")){ > break # ----------Really? I thought that condition was assured at this point? > } > } > }else{ > if(observable==1){ > recap=runif(1,0,1) > if(recap<=c){ > y[q,j]="A" > } > } > } > } > } > There are a lot of unnecessary tests and conditions. I tried to break > down the code and write the tests that have been done when assigning a > variable. I simplified your the first part but cannot guarantee that > all > criteria are met. > > Regards, > > Eloi > > On 12-08-01 09:47 AM, Claudia Penaloza wrote: >> I will accept any help you can give me, especially to free myself of >>> SAS-brain inefficiencies... >> My supervisor knows SAS not R, which may explain my code. >> What I'm actually trying to do is simulate mark-recapture data with a >> peculiar structure. >> It is a multistate robust design model with dummy time periods... it >> will basically be a matrix with a succession of letters (for >> different >> states/age classes) and zeros that are generated following a certain >> set of conditions (probability statements; drawing from a random >> uniform distribution if an event happens). >> Normally, survival and transition to another state occur between >> primary sampling periods (in my R example, every two columns.. >> between >> [,2] and [,3]) but in the "dummy time period" model survival occurs >> first and then transition to another state, which is why I need my >> "conditions" to alternate. Additionally, the robust design has >> secondary sampling sessions that are within the same year, i.e., >> survival is assumed = 1, which is why my columns are paired. Each >> state can also be in an unobservable state at any given time... all >> of >> these details complicate the coding. >> Below I've pasted what I have thus far... I have not debugged it yet >> (the second loop isn't working yet and the first loop still has some >> glitches). >> Further below is properly working code for a robust design without >> dummy time periods (it also lacks the dead states the dummy time >> period model has, which happen to be part of the glitchy-ness). >> Is there a better (more R-ish) way of doing this? >> Thank you for all your help, >> Claudia >> ################################################################ >> # Robust design with Markovian emigration and dummy time periods >> ################################################################ >> J = 24 >> N = 10 >> S = .9 >> PsiAD = 0.06 >> PsiAd = 0.04 >> PsiAA = .4 >> PsiaA = .3 >> p = .5 >> c = p >> y <- matrix(0,N,J) >> y[,1]="A" >> dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] >> for(j in which(dtp)){ >> for (q in 1:N){ >> (observable=1) >> if(j %% 2){survive=runif(1,0,1) >> if(survive<=S){ >> stay=runif(1,0,1) >> if (observable==1){ >> if(stay<=PsiAA){ >> observable=1 >> }else{observable=0} >> }else{ >> return=runif(1,0,1) >> if(return<=PsiaA){ >> observable=1 >> }else{observable=0}} >> if(observable==1){ >> capture=runif(1,0,1) >> if(capture<=p){ >> y[q,j]="A"} >> }}else{ >> deado=runif(1,0,1) >> if(deado<=PsiAd){ >> y[q,j]="d" >> }else(y[q,j]="D") >> if(y[q,j]%in%c("D","d")){ >> break >> } >> } >> }else{ >> if(observable==1){ >> recap=runif(1,0,1) >> if(recap<=c){ >> y[q,j]="A" >> } >> } >> } >> } >> } >> >> for(j in which(!dtp)){ >> for (q in 1:N){ >> if(j %% 2){ >> stay=runif(1,0,1) >> if(observable==1){ >> if(stay<=PsiAA){ >> observable=1 >> }else{observable=0} >> }else{ >> return=runif(1,0,1) >> if(return<=PsiaA){ >> observable=1 >> }else{observable=0} >> } >> if(observable==1){ >> capture=runif(1,0,1) >> if(capture<=p){ >> y[q,j]="A"} >> }}else { >> if(observable==1){ >> recap=runif(1,0,1) >> if(recap<=c){ >> y[q,j]="A" >> } >> } >> } >> } >> } >> ########################################### >> ### Robust design with markovian emigration >> ########################################### >> J = 24 >> N = 1000 >> S = .9 >> PsiOO = .4 >> PsiUO = .3 >> p = .5 >> c = p >> y <- matrix(0,N,J) >> y[,1]=1 >> for (q in 1:N){ >> (observable=1) >> for(j in 2:J){ >> if(j %% 2){surviv=runif(1,0,1) >> if(surviv<=S){ >> stay=runif(1,0,1) >> if(observable==1){ >> if(stay<=PsiOO){ >> observable=1 >> }else{observable=0} >> }else{ >> return=runif(1,0,1) >> if(return<=PsiUO){ >> observable=1 >> }else{observable=0}} >> if(observable==1){ >> cap=runif(1,0,1) >> if(cap<=p){ >> y[q,j]=1} >> }else y[q,j]=0 >> }else{break} >> }else{ >> if (observable==1){ >> recap=runif(1,0,1) >> if (recap<=c){ >> y[q,j]=1} >> else{y[q,j]=0} >> } >> } >> } >> } >> On Tue, Jul 31, 2012 at 7:16 PM, David Winsemius >> <[hidden email] <mailto:[hidden email]>> wrote: >> >> Claudia; >> >> The loop constructs will keep you mired in SAS-brain >> inefficiencies forever. Please, listen more carefully to Mercier. >> >> -- >> David >> >> >> >> On Jul 31, 2012, at 3:54 PM, Mercier Eloi wrote: >> >> On 12-07-31 02:38 PM, Claudia Penaloza wrote: >> >> Thank you very much Rui (and Eloi for your input)... this >> is (the very >> simplified version of) what I will end up using: >> >> Could we get the extended version ? Because right now, I don't >> know why >> you need such complicated code that can be done in 1 line. >> >> J <- 10 >> N <- 10 >> y <- matrix(0,N,J) >> cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] >> for(j in which(cols)){ >> for (q in 1:N){ >> if(j %% 2){ >> y[q,j]=1 >> }else{y[q,j]=2} >> } >> } >> for(j in which(!cols)){ >> for (q in 1:N){ >> if(j %% 2){ >> y[q,j]="A" >> }else{y[q,j]="B"} >> } >> } >> >> There is no need for a double loop (on N) : >> >> J <- 10 >> N <- 10 >> y <- matrix(0,N,J) >> cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] >> for(j in which(cols)){ >> if(j %% 2){ >> y[,j]=1 >> }else{y[,j]=2} >> } >> for(j in which(!cols)){ >> if(j %% 2){ >> y[,j]="A" >> }else{y[,j]="B"} >> } >> >> if you really wants to use this code. >> >> Cheers, >> >> Eloi >> >> Cheers, >> Claudia >> On Mon, Jul 30, 2012 at 5:38 PM, Mercier Eloi >> <[hidden email] <mailto:[hidden email]> >> <mailto:[hidden email] >> <mailto:[hidden email]>>> wrote: >> >> Or, assuming you only have 4 different elements : >> >> mat<- matrix(rep(c(1,2,"A", "B"),each=10),10,10, >> byrow=F) >> mat2 <- as.data.frame(mat) >> >> mat >> >> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [, >> 10] >> [1,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [2,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [3,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [4,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [5,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [6,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [7,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [8,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [9,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> [10,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >> >> mat2 >> V1 V2 V3 V4 V5 V6 V7 V8 V9 V10 >> >> 1 1 2 A B 1 2 A B 1 2 >> 2 1 2 A B 1 2 A B 1 2 >> 3 1 2 A B 1 2 A B 1 2 >> 4 1 2 A B 1 2 A B 1 2 >> 5 1 2 A B 1 2 A B 1 2 >> 6 1 2 A B 1 2 A B 1 2 >> 7 1 2 A B 1 2 A B 1 2 >> 8 1 2 A B 1 2 A B 1 2 >> 9 1 2 A B 1 2 A B 1 2 >> 10 1 2 A B 1 2 A B 1 2 >> >> Cheers, >> >> Eloi >> >> >> On 12-07-30 04:28 PM, Rui Barradas wrote: >> >> Hello, >> >> Maybe something along the lines of >> >> J <- 10 >> cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3) >> [seq_len(J)] >> for(i in which(cols)) { do something } >> for(i in which(!cols)) { do something else } >> >> Hope this helps, >> >> Rui Barradas >> >> Em 31-07-2012 00 <tel:31-07-2012%2000> >> <tel:31-07-2012%2000>:18, Claudia Penaloza >> >> escreveu: >> >> Dear All, >> I would like to apply two different "for loops" >> to each >> set of four columns >> of a matrix (the loops here are simplifications >> of the >> actual loops I will >> be running which involve multiple if/else >> statements). >> I don't know how to "alternate" between the >> loops >> depending on which column >> is "running through the loop" at the time. >> ## Set up matrix >> J <- 10 >> N <- 10 >> y <- matrix(0,N,J) >> ## Apply this loop to the first two of every >> four columns >> ([,1:2], >> [,5:6],[,9:10], etc.) >> for (q in 1:N){ >> for(j in 1:J){ >> if(j %% 2){ >> y[q,j]=1 >> }else{y[q,j]=2} >> } >> } >> ## Apply this loop to the next two of every >> four columns >> ([,3:4],[,7:8],[,11:12], etc.) >> for (q in 1:N){ >> for(j in 1:J){ >> if(j %% 2){ >> y[q,j]="A" >> }else{y[q,j]="B"} >> } >> } >> I want to get something like this (preferably >> without the >> quotes): >> >> y >> >> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] >> [,10] >> [1,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [2,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [3,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [4,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [5,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [6,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [7,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [8,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [9,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> [10,] "1" "2" "A" "B" "1" "2" "A" "B" >> "1" "2" >> >> Any help greatly appreciated! >> >> Claudia >> >> >> -- >> Eloi Mercier >> Bioinformatics PhD Student, UBC >> Paul Pavlidis Lab >> 2185 East Mall >> University of British Columbia >> Vancouver BC V6T1Z4 >> >> >> >> >> David Winsemius, MD >> Alameda, CA, USA >> >> > > > -- > Eloi Mercier > Bioinformatics PhD Student, UBC > Paul Pavlidis Lab > 2185 East Mall > University of British Columbia > Vancouver BC V6T1Z4 > > > [[alternative HTML version deleted]] > > ______________________________________________ > [hidden email] 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. David Winsemius, MD Alameda, CA, USA ______________________________________________ [hidden email] 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. |
|
|
Thank you very much for all your suggestions. I am very sorry my code is so
crude (it gives me a headache too!), I'm very new to R programing. I will have to look at your suggestions/questions very carefully (I'm barely holding on at this point) and get back to you (Dr. Winsemius) with some answers. Thank you! Claudia On Wed, Aug 1, 2012 at 2:11 PM, David Winsemius <[hidden email]>wrote: > > On Aug 1, 2012, at 12:02 PM, Mercier Eloi wrote: > > Your code almost gave me a headache. :-/ >> > > I had a similar reaction. However, my approach might have been to request > a more complete verbal description of the data structures being operated on > and the methods and assumptions being used. Generally it is going to be > easier to go from a procedural description to good R code than it is to go > from a SAS Data Proc ... even if it were tested and debugged... and yours > was not even debugged. > > > ##############################**##############################**#### >> # Robust design with Markovian emigration and dummy time periods >> ##############################**##############################**#### >> J = 24 >> N = 10 >> S = .9 >> PsiAD = 0.06 >> PsiAd = 0.04 >> PsiAA = .4 >> PsiaA = .3 >> p = .5 >> c = p >> y <- matrix(0,N,J) >> > > # So is 'y' a matrix of states where the N rows are levels and the J > columns are discrete times? > # Actually I decided that context suggested you were using letters to > represent states. > > y[,1]="A" >> > > # So we start with N <something>'s in state "A"? > # It seems as though it might be the case that every row is independent of > the others. > # .. and you are performing ( replicate()-ing in R terminology) this test > N times > # states: > #("A" "D") > #("a" "d") > #transitions: > matrix( c( runif(1, 0, 1) <= 0.4, # A -> A > runif(1,0,1) <= 0.3, # a -> A > runif(1,0,1) <= 0.06. # A -> D > runif(1,0,1) <= 0.04 # A -> d) ) > > # What is state "a"? > # How do you get from A to a? > # can you get from D or d to A or a? > > # Not sure I have gotten the model assumptions down. > # Is D" or "d" an absorbing state such as "dead or "Dead"? > > > > dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] >> > > #Presumably the "dummy time periods" > > > for(j in which(dtp)){ # j will be c(1,2, 5,6, 9,10 , ..., 21,22) >> > > for (q in 1:N){ # This can almost surely become vectorized >> >> (observable=1) >> if(j %% 2){survive=runif(1,0,1) >> if(survive<=S){ >> stay=runif(1,0,1) >> if (observable==1){ >> > > # It would help if the concept of "observable" were explained > # Is this something to do with the capture-recapture observables? > > > if{ >> observable=1 >> }else{observable=0} >> }else{ >> > > # After allowing for Mercier's astute observation that observable will > always be 1, I'm wondering if it can be replaced with just this code (and > not need the loop surrounding it.) > > observable <- (stay<=PsiAA) > > #-------------- > >> return=runif(1,0,1) >> > > # better NOT to use the name "return" since it is a crucial R function > name. > > > if(return<=PsiaA){ >> observable=1 >> }else{observable=0}} >> > > #------- perhaps: > > randret <- return=runif(N,0,1) > observable <- randret <= PsiaA > # ------------------- > > if(observable==1){ >> capture=runif(1,0,1) >> if(capture<=p){ >> > > #---------- perhaps: > randcap <- runif(N, 0,1) > y[ randcap <= p, j] <- "A" > #That would replace with "A" (within the j-column) ... > # only the rows in 'y' for which randcap were less than the randcap > threshold. > > #------------ > > y[q,j]="A"} >> }}else{ >> > > > deado=runif(1,0,1) >> if(deado<=PsiAd){ >> y[q,j]="d" >> }else(y[q,j]="D") >> > > # -------Perhaps > deado <- runif(N, 0,1) > y[ , j] <- ifelse( deado<=PsiAd, "d", "D") > #------------------------ > > > if(y[q,j]%in%c("D","d")){ >> break >> > > # ----------Really? I thought that condition was assured at this point? > > > > } >> } >> }else{ >> if(observable==1){ >> recap=runif(1,0,1) >> if(recap<=c){ >> y[q,j]="A" >> } >> } >> } >> } >> } >> > > > > There are a lot of unnecessary tests and conditions. I tried to break >> down the code and write the tests that have been done when assigning a >> variable. I simplified your the first part but cannot guarantee that all >> criteria are met. >> >> > Regards, >> >> Eloi >> >> On 12-08-01 09:47 AM, Claudia Penaloza wrote: >> >>> I will accept any help you can give me, especially to free myself of >>> >>>> SAS-brain inefficiencies... >>>> >>> My supervisor knows SAS not R, which may explain my code. >>> What I'm actually trying to do is simulate mark-recapture data with a >>> peculiar structure. >>> It is a multistate robust design model with dummy time periods... it >>> will basically be a matrix with a succession of letters (for different >>> states/age classes) and zeros that are generated following a certain >>> set of conditions (probability statements; drawing from a random >>> uniform distribution if an event happens). >>> Normally, survival and transition to another state occur between >>> primary sampling periods (in my R example, every two columns.. between >>> [,2] and [,3]) but in the "dummy time period" model survival occurs >>> first and then transition to another state, which is why I need my >>> "conditions" to alternate. Additionally, the robust design has >>> secondary sampling sessions that are within the same year, i.e., >>> survival is assumed = 1, which is why my columns are paired. Each >>> state can also be in an unobservable state at any given time... all of >>> these details complicate the coding. >>> Below I've pasted what I have thus far... I have not debugged it yet >>> (the second loop isn't working yet and the first loop still has some >>> glitches). >>> Further below is properly working code for a robust design without >>> dummy time periods (it also lacks the dead states the dummy time >>> period model has, which happen to be part of the glitchy-ness). >>> Is there a better (more R-ish) way of doing this? >>> Thank you for all your help, >>> Claudia >>> ##############################**##############################**#### >>> # Robust design with Markovian emigration and dummy time periods >>> ##############################**##############################**#### >>> J = 24 >>> N = 10 >>> S = .9 >>> PsiAD = 0.06 >>> PsiAd = 0.04 >>> PsiAA = .4 >>> PsiaA = .3 >>> p = .5 >>> c = p >>> y <- matrix(0,N,J) >>> y[,1]="A" >>> dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] >>> for(j in which(dtp)){ >>> for (q in 1:N){ >>> (observable=1) >>> if(j %% 2){survive=runif(1,0,1) >>> if(survive<=S){ >>> stay=runif(1,0,1) >>> if (observable==1){ >>> if(stay<=PsiAA){ >>> observable=1 >>> }else{observable=0} >>> }else{ >>> return=runif(1,0,1) >>> if(return<=PsiaA){ >>> observable=1 >>> }else{observable=0}} >>> if(observable==1){ >>> capture=runif(1,0,1) >>> if(capture<=p){ >>> y[q,j]="A"} >>> }}else{ >>> deado=runif(1,0,1) >>> if(deado<=PsiAd){ >>> y[q,j]="d" >>> }else(y[q,j]="D") >>> if(y[q,j]%in%c("D","d")){ >>> break >>> } >>> } >>> }else{ >>> if(observable==1){ >>> recap=runif(1,0,1) >>> if(recap<=c){ >>> y[q,j]="A" >>> } >>> } >>> } >>> } >>> } >>> >>> for(j in which(!dtp)){ >>> for (q in 1:N){ >>> if(j %% 2){ >>> stay=runif(1,0,1) >>> if(observable==1){ >>> if(stay<=PsiAA){ >>> observable=1 >>> }else{observable=0} >>> }else{ >>> return=runif(1,0,1) >>> if(return<=PsiaA){ >>> observable=1 >>> }else{observable=0} >>> } >>> if(observable==1){ >>> capture=runif(1,0,1) >>> if(capture<=p){ >>> y[q,j]="A"} >>> }}else { >>> if(observable==1){ >>> recap=runif(1,0,1) >>> if(recap<=c){ >>> y[q,j]="A" >>> } >>> } >>> } >>> } >>> } >>> ##############################**############# >>> ### Robust design with markovian emigration >>> ##############################**############# >>> J = 24 >>> N = 1000 >>> S = .9 >>> PsiOO = .4 >>> PsiUO = .3 >>> p = .5 >>> c = p >>> y <- matrix(0,N,J) >>> y[,1]=1 >>> for (q in 1:N){ >>> (observable=1) >>> for(j in 2:J){ >>> if(j %% 2){surviv=runif(1,0,1) >>> if(surviv<=S){ >>> stay=runif(1,0,1) >>> if(observable==1){ >>> if(stay<=PsiOO){ >>> observable=1 >>> }else{observable=0} >>> }else{ >>> return=runif(1,0,1) >>> if(return<=PsiUO){ >>> observable=1 >>> }else{observable=0}} >>> if(observable==1){ >>> cap=runif(1,0,1) >>> if(cap<=p){ >>> y[q,j]=1} >>> }else y[q,j]=0 >>> }else{break} >>> }else{ >>> if (observable==1){ >>> recap=runif(1,0,1) >>> if (recap<=c){ >>> y[q,j]=1} >>> else{y[q,j]=0} >>> } >>> } >>> } >>> } >>> On Tue, Jul 31, 2012 at 7:16 PM, David Winsemius >>> <[hidden email] <mailto:[hidden email]**>> wrote: >>> >>> Claudia; >>> >>> The loop constructs will keep you mired in SAS-brain >>> inefficiencies forever. Please, listen more carefully to Mercier. >>> >>> -- >>> David >>> >>> >>> >>> On Jul 31, 2012, at 3:54 PM, Mercier Eloi wrote: >>> >>> On 12-07-31 02:38 PM, Claudia Penaloza wrote: >>> >>> Thank you very much Rui (and Eloi for your input)... this >>> is (the very >>> simplified version of) what I will end up using: >>> >>> Could we get the extended version ? Because right now, I don't >>> know why >>> you need such complicated code that can be done in 1 line. >>> >>> J <- 10 >>> N <- 10 >>> y <- matrix(0,N,J) >>> cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] >>> for(j in which(cols)){ >>> for (q in 1:N){ >>> if(j %% 2){ >>> y[q,j]=1 >>> }else{y[q,j]=2} >>> } >>> } >>> for(j in which(!cols)){ >>> for (q in 1:N){ >>> if(j %% 2){ >>> y[q,j]="A" >>> }else{y[q,j]="B"} >>> } >>> } >>> >>> There is no need for a double loop (on N) : >>> >>> J <- 10 >>> N <- 10 >>> y <- matrix(0,N,J) >>> cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] >>> for(j in which(cols)){ >>> if(j %% 2){ >>> y[,j]=1 >>> }else{y[,j]=2} >>> } >>> for(j in which(!cols)){ >>> if(j %% 2){ >>> y[,j]="A" >>> }else{y[,j]="B"} >>> } >>> >>> if you really wants to use this code. >>> >>> Cheers, >>> >>> Eloi >>> >>> Cheers, >>> Claudia >>> On Mon, Jul 30, 2012 at 5:38 PM, Mercier Eloi >>> <[hidden email] <mailto:[hidden email]> >>> <mailto:[hidden email] >>> >>> <mailto:[hidden email]>**>> wrote: >>> >>> Or, assuming you only have 4 different elements : >>> >>> mat<- matrix(rep(c(1,2,"A", "B"),each=10),10,10, byrow=F) >>> mat2 <- as.data.frame(mat) >>> >>> mat >>> >>> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] >>> [1,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [2,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [3,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [4,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [5,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [6,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [7,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [8,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [9,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> [10,] "1" "2" "A" "B" "1" "2" "A" "B" "1" "2" >>> >>> mat2 >>> V1 V2 V3 V4 V5 V6 V7 V8 V9 V10 >>> >>> 1 1 2 A B 1 2 A B 1 2 >>> 2 1 2 A B 1 2 A B 1 2 >>> 3 1 2 A B 1 2 A B 1 2 >>> 4 1 2 A B 1 2 A B 1 2 >>> 5 1 2 A B 1 2 A B 1 2 >>> 6 1 2 A B 1 2 A B 1 2 >>> 7 1 2 A B 1 2 A B 1 2 >>> 8 1 2 A B 1 2 A B 1 2 >>> 9 1 2 A B 1 2 A B 1 2 >>> 10 1 2 A B 1 2 A B 1 2 >>> >>> Cheers, >>> >>> Eloi >>> >>> >>> On 12-07-30 04:28 PM, Rui Barradas wrote: >>> >>> Hello, >>> >>> Maybe something along the lines of >>> >>> J <- 10 >>> cols <- rep(c(TRUE, TRUE, FALSE, FALSE), 3)[seq_len(J)] >>> for(i in which(cols)) { do something } >>> for(i in which(!cols)) { do something else } >>> >>> Hope this helps, >>> >>> Rui Barradas >>> >>> Em 31-07-2012 00 <tel:31-07-2012%2000> >>> <tel:31-07-2012%2000>:18, Claudia Penaloza >>> >>> escreveu: >>> >>> Dear All, >>> I would like to apply two different "for loops" >>> to each >>> set of four columns >>> of a matrix (the loops here are simplifications >>> of the >>> actual loops I will >>> be running which involve multiple if/else >>> statements). >>> I don't know how to "alternate" between the loops >>> depending on which column >>> is "running through the loop" at the time. >>> ## Set up matrix >>> J <- 10 >>> N <- 10 >>> y <- matrix(0,N,J) >>> ## Apply this loop to the first two of every >>> four columns >>> ([,1:2], >>> [,5:6],[,9:10], etc.) >>> for (q in 1:N){ >>> for(j in 1:J){ >>> if(j %% 2){ >>> y[q,j]=1 >>> }else{y[q,j]=2} >>> } >>> } >>> ## Apply this loop to the next two of every >>> four columns >>> ([,3:4],[,7:8],[,11:12], etc.) >>> for (q in 1:N){ >>> for(j in 1:J){ >>> if(j %% 2){ >>> y[q,j]="A" >>> }else{y[q,j]="B"} >>> } >>> } >>> I want to get something like this (preferably >>> without the >>> quotes): >>> >>> y >>> >>> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] >>> [1,] "1" "2" "A" "B" "1" "2" "A" "B" >>> "1" "2" >>> [2,] "1" "2" "A" "B" "1" "2" "A" "B" >>> "1" "2" >>> [3,] "1" "2" "A" "B" "1" "2" "A" "B" >>> "1" "2" >>> [4,] "1" "2" "A" "B" "1" "2" "A" "B" >>> "1" "2" >>> [5,] "1" "2" "A" "B" "1" "2" "A" "B" >>> "1" "2" >>> [6,] "1" "2" "A" "B" "1" "2" "A" "B" >>> "1" "2" >>> [7,] "1" "2" "A" "B" "1" "2" "A" "B" >>> "1" "2" >>> [8,] "1" "2" "A" "B" "1" "2" "A" "B" >>> "1" "2" >>> [9,] "1" "2" "A" "B" "1" "2" "A" "B" >>> "1" "2" >>> [10,] "1" "2" "A" "B" "1" "2" "A" "B" >>> "1" "2" >>> >>> Any help greatly appreciated! >>> >>> Claudia >>> >>> >>> -- >>> Eloi Mercier >>> Bioinformatics PhD Student, UBC >>> Paul Pavlidis Lab >>> 2185 East Mall >>> University of British Columbia >>> Vancouver BC V6T1Z4 >>> >>> >>> >>> >>> David Winsemius, MD >>> Alameda, CA, USA >>> >>> >>> >> >> -- >> Eloi Mercier >> Bioinformatics PhD Student, UBC >> Paul Pavlidis Lab >> 2185 East Mall >> University of British Columbia >> Vancouver BC V6T1Z4 >> >> >> [[alternative HTML version deleted]] >> >> ______________________________**________________ >> [hidden email] mailing list >> https://stat.ethz.ch/mailman/**listinfo/r-help<https://stat.ethz.ch/mailman/listinfo/r-help> >> PLEASE do read the posting guide http://www.R-project.org/** >> posting-guide.html <http://www.R-project.org/posting-guide.html> >> and provide commented, minimal, self-contained, reproducible code. >> > > David Winsemius, MD > Alameda, CA, USA > > [[alternative HTML version deleted]] ______________________________________________ [hidden email] 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. |
|
|
Answers inserted in BLUE below
On Thu, Aug 2, 2012 at 5:34 PM, Claudia Penaloza <[hidden email]>wrote: > Thank you very much for all your suggestions. I am very sorry my code is > so crude (it gives me a headache too!), I'm very new to R programing. I > will have to look at your suggestions/questions very carefully (I'm barely > holding on at this point) and get back to you (Dr. Winsemius) with some > answers. > > Thank you! > Claudia > > On Wed, Aug 1, 2012 at 2:11 PM, David Winsemius <[hidden email]>wrote: > >> >> On Aug 1, 2012, at 12:02 PM, Mercier Eloi wrote: >> >> Your code almost gave me a headache. :-/ >>> >> >> I had a similar reaction. However, my approach might have been to request >> a more complete verbal description of the data structures being operated on >> and the methods and assumptions being used. Generally it is going to be >> easier to go from a procedural description to good R code than it is to go >> from a SAS Data Proc ... even if it were tested and debugged... and yours >> was not even debugged. >> >> >> ##############################**##############################**#### >>> # Robust design with Markovian emigration and dummy time periods >>> ##############################**##############################**#### >>> J = 24 >>> N = 10 >>> S = .9 >>> PsiAD = 0.06 >>> PsiAd = 0.04 >>> PsiAA = .4 >>> PsiaA = .3 >>> p = .5 >>> c = p >>> y <- matrix(0,N,J) >>> >> >> # So is 'y' a matrix of states where the N rows are levels and the J >> columns are discrete times? >> # Actually I decided that context suggested you were using letters to >> represent states. >> > Yes, 'y' is a matrix of state, N rows are levels (independent of each Yes, I am using letters to represent states. > >> y[,1]="A" >>> >> >> # So we start with N <something>'s in state "A"? >> # It seems as though it might be the case that every row is independent >> of the others. >> # .. and you are performing ( replicate()-ing in R terminology) this test >> N times >> # states: >> #("A" "D") >> #("a" "d") >> #transitions: >> matrix( c( runif(1, 0, 1) <= 0.4, # A -> A >> runif(1,0,1) <= 0.3, # a -> A >> runif(1,0,1) <= 0.06. # A -> D >> runif(1,0,1) <= 0.04 # A -> d) ) >> >> # What is state "a"? >> # How do you get from A to a? >> # can you get from D or d to A or a? >> > Yes, we start with N "individuals" in state "A", each individual is State "a" is an unobserved (!observable) version of state "A" (biologically, an individual that has temporarily left the sampling area and is not available for capture) An individual in state "A" transitions to state "a" if it is observable and doesn't stay (stay >= AA) "D" and "d" are dead (absorbing) states. Once in either "D" or "d", the individual can no longer transition and is no longer "captured" (the row should only have zeros after a "D" or a "d") > > > >> # Not sure I have gotten the model assumptions down. >> # Is D" or "d" an absorbing state such as "dead or "Dead"? >> > Yes. > dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] >> >> >> #Presumably the "dummy time periods" >> > Yes. > >> >> for(j in which(dtp)){ # j will be c(1,2, 5,6, 9,10 , ..., 21,22) >> >> >> >> for (q in 1:N){ # This can almost surely become vectorized >> >> (observable=1) >> if(j %% 2){survive=runif(1,0,1) >> if(survive<=S){ >> stay=runif(1,0,1) >> if (observable==1){ >> >> >> # It would help if the concept of "observable" were explained >> # Is this something to do with the capture-recapture observables? >> > first time). Individuals can then stay in their current state or transition to another one, if they are observable they can continue to be observable, or they can become unobservable and viceversa. Transition depends on what state you are in at any given time (A->A != A->a != a->A != a->a). > >> >> if{ >> observable=1 >> }else{observable=0} >> }else{ >> >> >> # After allowing for Mercier's astute observation that observable will >> always be 1, I'm wondering if it can be replaced with just this code (and >> not need the loop surrounding it.) >> >> >> > The individual will always "start" as observable... but as the loop and unobservable (at least that is what I wanted to code). > >> >> observable <- (stay<=PsiAA) >> >> #-------------- >> >> return=runif(1,0,1) >> >> >> # better NOT to use the name "return" since it is a crucial R function >> name. >> >> >> > Will change. Thank you. > >> >> >> >> if(return<=PsiaA){ >> observable=1 >> }else{observable=0}} >> >> >> #------- perhaps: >> >> randret <- return=runif(N,0,1) >> observable <- randret <= PsiaA >> # ------------------- >> >> >> if(observable==1){ >> capture=runif(1,0,1) >> if(capture<=p){ >> >> >> #---------- perhaps: >> randcap <- runif(N, 0,1) >> y[ randcap <= p, j] <- "A" >> #That would replace with "A" (within the j-column) ... >> # only the rows in 'y' for which randcap were less than the randcap >> threshold. >> >> #------------ >> >> >> y[q,j]="A"} >> }}else{ >> >> >> >> >> deado=runif(1,0,1) >> if(deado<=PsiAd){ >> y[q,j]="d" >> }else(y[q,j]="D") >> >> >> # -------Perhaps >> deado <- runif(N, 0,1) >> y[ , j] <- ifelse( deado<=PsiAd, "d", "D") >> #------------------------ >> >> >> >> if(y[q,j]%in%c("D","d")){ >> break >> >> >> # ----------Really? I thought that condition was assured at this >> point? >> >> >> > I thought so too but non-zero elements appeared toward the right in rows > >> } >> } >> }else{ >> if(observable==1){ >> recap=runif(1,0,1) >> if(recap<=c){ >> y[q,j]="A" >> } >> } >> } >> } >> } >> >> >> >> >> >> >> There are a lot of unnecessary tests and conditions. I tried to break >> down the code and write the tests that have been done when assigning a >> variable. I simplified your the first part but cannot guarantee that all >> criteria are met. >> >> >> > I will run the three modified versions and see how things go. I hope my [[alternative HTML version deleted]] ______________________________________________ [hidden email] 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. |
|
|
This is what my code looks like now. However, there is one thing I
can't/don't know how to fix. I can't get it to be "once dead always dead", i.e., in any given row, after a "D" or a "d" there should be only zeros. I've tried applying a flag to break the loop if dead but I can't get it to work. Could you please help me with this? Thank you for your time, Claudia J = 24 N = 10 S = .9 PsiADd = 0.4 PsiAA = .4 PsiaA = .3 p = .5 c = p y <- matrix(0,N,J) y[,1]="A" dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] for(j in which(dtp)){ for (q in 1:N){ (observable <- TRUE) if(j %% 2){ if (runif(1,0,1)<=S) { if (observable) { observable <- (runif(1,0,1)<=PsiAA) } else { observable <- (runif(1,0,1)<=PsiaA) } if (observable) { if (runif(1,0,1) <= p) y[q,j] <- "A" } } else { y[q,j] <- ifelse(runif(1,0,1) <= PsiADd,"d","D") break } } else { if(observable){ if(runif(1,0,1) <= c) y[q,j] <- "A" } } } } for (j in which(!dtp)) { for (q in 1:N) { if (j %% 2) { if (observable) { observable <- runif(1, 0, 1) <= PsiAA } else { observable <- runif(1, 0, 1) <= PsiaA } if (observable) { if (runif(1, 0, 1) <= p) y[q, j] <- "A" } } else { if (observable) { if (runif(1, 0, 1) <= c) y[q, j] <- "A" } } } } On Wed, Aug 8, 2012 at 2:04 PM, Claudia Penaloza <[hidden email]>wrote: > Answers inserted in BLUE below > > On Thu, Aug 2, 2012 at 5:34 PM, Claudia Penaloza < > [hidden email]> wrote: > >> Thank you very much for all your suggestions. I am very sorry my code is >> so crude (it gives me a headache too!), I'm very new to R programing. I >> will have to look at your suggestions/questions very carefully (I'm barely >> holding on at this point) and get back to you (Dr. Winsemius) with some >> answers. >> >> Thank you! >> Claudia >> >> On Wed, Aug 1, 2012 at 2:11 PM, David Winsemius <[hidden email]>wrote: >> >>> >>> On Aug 1, 2012, at 12:02 PM, Mercier Eloi wrote: >>> >>> Your code almost gave me a headache. :-/ >>>> >>> >>> I had a similar reaction. However, my approach might have been to >>> request a more complete verbal description of the data structures being >>> operated on and the methods and assumptions being used. Generally it is >>> going to be easier to go from a procedural description to good R code than >>> it is to go from a SAS Data Proc ... even if it were tested and debugged... >>> and yours was not even debugged. >>> >>> >>> ##############################**##############################**#### >>>> # Robust design with Markovian emigration and dummy time periods >>>> ##############################**##############################**#### >>>> J = 24 >>>> N = 10 >>>> S = .9 >>>> PsiAD = 0.06 >>>> PsiAd = 0.04 >>>> PsiAA = .4 >>>> PsiaA = .3 >>>> p = .5 >>>> c = p >>>> y <- matrix(0,N,J) >>>> >>> >>> # So is 'y' a matrix of states where the N rows are levels and the J >>> columns are discrete times? >>> # Actually I decided that context suggested you were using letters to >>> represent states. >>> >> Yes, 'y' is a matrix of state, N rows are levels (independent of each > other) and J columns are discrete times. > Yes, I am using letters to represent states. > >> >>> y[,1]="A" >>>> >>> >>> # So we start with N <something>'s in state "A"? >>> # It seems as though it might be the case that every row is independent >>> of the others. >>> # .. and you are performing ( replicate()-ing in R terminology) this >>> test N times >>> # states: >>> #("A" "D") >>> #("a" "d") >>> #transitions: >>> matrix( c( runif(1, 0, 1) <= 0.4, # A -> A >>> runif(1,0,1) <= 0.3, # a -> A >>> runif(1,0,1) <= 0.06. # A -> D >>> runif(1,0,1) <= 0.04 # A -> d) ) >>> >>> # What is state "a"? >>> # How do you get from A to a? >>> # can you get from D or d to A or a? >>> >> Yes, we start with N "individuals" in state "A", each individual is > independent of each other. > State "a" is an unobserved (!observable) version of state "A" > (biologically, an individual that has temporarily left the sampling area > and is not available for capture) > An individual in state "A" transitions to state "a" if it is observable > and doesn't stay (stay >= AA) > "D" and "d" are dead (absorbing) states. Once in either "D" or "d", the > individual can no longer transition and is no longer "captured" (the row > should only have zeros after a "D" or a "d") > > >> >> >> >>> # Not sure I have gotten the model assumptions down. >>> # Is D" or "d" an absorbing state such as "dead or "Dead"? >>> >> Yes. > >> dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] >>> >>> >>> #Presumably the "dummy time periods" >>> >> Yes. > > >> >>> >>> for(j in which(dtp)){ # j will be c(1,2, 5,6, 9,10 , ..., 21,22) >>> >>> >>> >>> for (q in 1:N){ # This can almost surely become vectorized >>> >>> (observable=1) >>> if(j %% 2){survive=runif(1,0,1) >>> if(survive<=S){ >>> stay=runif(1,0,1) >>> if (observable==1){ >>> >>> >>> # It would help if the concept of "observable" were explained >>> # Is this something to do with the capture-recapture observables? >>> >> > Yes. All individuals start out "observable" (we need to capture them a > first time). Individuals can then stay in their current state or transition > to another one, if they are observable they can continue to be observable, > or they can become unobservable and viceversa. Transition depends on what > state you are in at any given time (A->A != A->a != a->A != a->a). > >> >>> >>> if{ >>> observable=1 >>> }else{observable=0} >>> }else{ >>> >>> >>> # After allowing for Mercier's astute observation that observable will >>> always be 1, I'm wondering if it can be replaced with just this code (and >>> not need the loop surrounding it.) >>> >>> >>> >> The individual will always "start" as observable... but as the loop > progresses through the columns in a row, it can go between being observable > and unobservable (at least that is what I wanted to code). > >> >>> >>> observable <- (stay<=PsiAA) >>> >>> #-------------- >>> >>> return=runif(1,0,1) >>> >>> >>> # better NOT to use the name "return" since it is a crucial R function >>> name. >>> >>> >>> >> Will change. Thank you. > >> >>> >>> >>> >>> if(return<=PsiaA){ >>> observable=1 >>> }else{observable=0}} >>> >>> >>> #------- perhaps: >>> >>> randret <- return=runif(N,0,1) >>> observable <- randret <= PsiaA >>> # ------------------- >>> >>> >>> if(observable==1){ >>> capture=runif(1,0,1) >>> if(capture<=p){ >>> >>> >>> #---------- perhaps: >>> randcap <- runif(N, 0,1) >>> y[ randcap <= p, j] <- "A" >>> #That would replace with "A" (within the j-column) ... >>> # only the rows in 'y' for which randcap were less than the randcap >>> threshold. >>> >>> #------------ >>> >>> >>> y[q,j]="A"} >>> }}else{ >>> >>> >>> >>> >>> deado=runif(1,0,1) >>> if(deado<=PsiAd){ >>> y[q,j]="d" >>> }else(y[q,j]="D") >>> >>> >>> # -------Perhaps >>> deado <- runif(N, 0,1) >>> y[ , j] <- ifelse( deado<=PsiAd, "d", "D") >>> #------------------------ >>> >>> >>> >>> if(y[q,j]%in%c("D","d")){ >>> break >>> >>> >>> # ----------Really? I thought that condition was assured at this >>> point? >>> >>> >>> >> I thought so too but non-zero elements appeared toward the right in rows > that had already "died". > >> >>> } >>> } >>> }else{ >>> if(observable==1){ >>> recap=runif(1,0,1) >>> if(recap<=c){ >>> y[q,j]="A" >>> } >>> } >>> } >>> } >>> } >>> >>> >>> >>> >>> >>> >>> There are a lot of unnecessary tests and conditions. I tried to break >>> down the code and write the tests that have been done when assigning a >>> variable. I simplified your the first part but cannot guarantee that all >>> criteria are met. >>> >>> >>> >> I will run the three modified versions and see how things go. I hope my > answers are helpful. > > > [[alternative HTML version deleted]] ______________________________________________ [hidden email] 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. |
Re: Alternating between "for loops"
|
|
Hello,
I'm not sure it works, but try the following. for(j in which(dtp)){ for (q in 1:N){ if(y[q, j] %in% c("d", "D")) break [...etc...] and in the other loop the same, for (j in which(!dtp)) { for (q in 1:N) { if(y[q, j] %in% c("d", "D")) break [...etc...] Em 10-08-2012 20:42, Claudia Penaloza escreveu: > This is what my code looks like now. However, there is one thing I > can't/don't know how to fix. > I can't get it to be "once dead always dead", i.e., in any given row, after > a "D" or a "d" there should be only zeros. > I've tried applying a flag to break the loop if dead but I can't get it to > work. > Could you please help me with this? > > Thank you for your time, > Claudia > > J = 24 > N = 10 > S = .9 > PsiADd = 0.4 > PsiAA = .4 > PsiaA = .3 > p = .5 > c = p > > y <- matrix(0,N,J) > y[,1]="A" > dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] > for(j in which(dtp)){ > for (q in 1:N){ > (observable <- TRUE) > if(j %% 2){ > if (runif(1,0,1)<=S) { > if (observable) { > observable <- (runif(1,0,1)<=PsiAA) > } else { > observable <- (runif(1,0,1)<=PsiaA) > } > if (observable) { > if (runif(1,0,1) <= p) y[q,j] <- "A" > } > } else { > y[q,j] <- ifelse(runif(1,0,1) <= PsiADd,"d","D") > break > } > } else { > if(observable){ > if(runif(1,0,1) <= c) y[q,j] <- "A" > } > } > } > } > > for (j in which(!dtp)) { > for (q in 1:N) { > if (j %% 2) { > if (observable) { > observable <- runif(1, 0, 1) <= PsiAA > } else { > observable <- runif(1, 0, 1) <= PsiaA > } > if (observable) { > if (runif(1, 0, 1) <= p) y[q, j] <- "A" > } > } else { > if (observable) { > if (runif(1, 0, 1) <= c) y[q, j] <- "A" > } > } > } > } > On Wed, Aug 8, 2012 at 2:04 PM, Claudia Penaloza > <[hidden email]>wrote: > >> Answers inserted in BLUE below >> >> On Thu, Aug 2, 2012 at 5:34 PM, Claudia Penaloza < >> [hidden email]> wrote: >> >>> Thank you very much for all your suggestions. I am very sorry my code is >>> so crude (it gives me a headache too!), I'm very new to R programing. I >>> will have to look at your suggestions/questions very carefully (I'm barely >>> holding on at this point) and get back to you (Dr. Winsemius) with some >>> answers. >>> >>> Thank you! >>> Claudia >>> >>> On Wed, Aug 1, 2012 at 2:11 PM, David Winsemius <[hidden email]>wrote: >>> >>>> On Aug 1, 2012, at 12:02 PM, Mercier Eloi wrote: >>>> >>>> Your code almost gave me a headache. :-/ >>>> I had a similar reaction. However, my approach might have been to >>>> request a more complete verbal description of the data structures being >>>> operated on and the methods and assumptions being used. Generally it is >>>> going to be easier to go from a procedural description to good R code than >>>> it is to go from a SAS Data Proc ... even if it were tested and debugged... >>>> and yours was not even debugged. >>>> >>>> >>>> ##############################**##############################**#### >>>>> # Robust design with Markovian emigration and dummy time periods >>>>> ##############################**##############################**#### >>>>> J = 24 >>>>> N = 10 >>>>> S = .9 >>>>> PsiAD = 0.06 >>>>> PsiAd = 0.04 >>>>> PsiAA = .4 >>>>> PsiaA = .3 >>>>> p = .5 >>>>> c = p >>>>> y <- matrix(0,N,J) >>>>> >>>> # So is 'y' a matrix of states where the N rows are levels and the J >>>> columns are discrete times? >>>> # Actually I decided that context suggested you were using letters to >>>> represent states. >>>> >>> Yes, 'y' is a matrix of state, N rows are levels (independent of each >> other) and J columns are discrete times. >> Yes, I am using letters to represent states. >> >>>> y[,1]="A" >>>> # So we start with N <something>'s in state "A"? >>>> # It seems as though it might be the case that every row is independent >>>> of the others. >>>> # .. and you are performing ( replicate()-ing in R terminology) this >>>> test N times >>>> # states: >>>> #("A" "D") >>>> #("a" "d") >>>> #transitions: >>>> matrix( c( runif(1, 0, 1) <= 0.4, # A -> A >>>> runif(1,0,1) <= 0.3, # a -> A >>>> runif(1,0,1) <= 0.06. # A -> D >>>> runif(1,0,1) <= 0.04 # A -> d) ) >>>> >>>> # What is state "a"? >>>> # How do you get from A to a? >>>> # can you get from D or d to A or a? >>>> >>> Yes, we start with N "individuals" in state "A", each individual is >> independent of each other. >> State "a" is an unobserved (!observable) version of state "A" >> (biologically, an individual that has temporarily left the sampling area >> and is not available for capture) >> An individual in state "A" transitions to state "a" if it is observable >> and doesn't stay (stay >= AA) >> "D" and "d" are dead (absorbing) states. Once in either "D" or "d", the >> individual can no longer transition and is no longer "captured" (the row >> should only have zeros after a "D" or a "d") >> >> >>> >>> >>>> # Not sure I have gotten the model assumptions down. >>>> # Is D" or "d" an absorbing state such as "dead or "Dead"? >>>> >>> Yes. >>> dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] >>>> >>>> #Presumably the "dummy time periods" >>>> >>> Yes. >> >>>> for(j in which(dtp)){ # j will be c(1,2, 5,6, 9,10 , ..., 21,22) >>>> >>>> >>>> >>>> for (q in 1:N){ # This can almost surely become vectorized >>>> >>>> (observable=1) >>>> if(j %% 2){survive=runif(1,0,1) >>>> if(survive<=S){ >>>> stay=runif(1,0,1) >>>> if (observable==1){ >>>> >>>> >>>> # It would help if the concept of "observable" were explained >>>> # Is this something to do with the capture-recapture observables? >>>> >> Yes. All individuals start out "observable" (we need to capture them a >> first time). Individuals can then stay in their current state or transition >> to another one, if they are observable they can continue to be observable, >> or they can become unobservable and viceversa. Transition depends on what >> state you are in at any given time (A->A != A->a != a->A != a->a). >> >>>> if{ >>>> observable=1 >>>> }else{observable=0} >>>> }else{ >>>> >>>> >>>> # After allowing for Mercier's astute observation that observable will >>>> always be 1, I'm wondering if it can be replaced with just this code (and >>>> not need the loop surrounding it.) >>>> >>>> >>>> >>> The individual will always "start" as observable... but as the loop >> progresses through the columns in a row, it can go between being observable >> and unobservable (at least that is what I wanted to code). >> >>>> observable <- (stay<=PsiAA) >>>> >>>> #-------------- >>>> >>>> return=runif(1,0,1) >>>> >>>> >>>> # better NOT to use the name "return" since it is a crucial R function >>>> name. >>>> >>>> >>>> >>> Will change. Thank you. >>>> >>>> >>>> if(return<=PsiaA){ >>>> observable=1 >>>> }else{observable=0}} >>>> >>>> >>>> #------- perhaps: >>>> >>>> randret <- return=runif(N,0,1) >>>> observable <- randret <= PsiaA >>>> # ------------------- >>>> >>>> >>>> if(observable==1){ >>>> capture=runif(1,0,1) >>>> if(capture<=p){ >>>> >>>> >>>> #---------- perhaps: >>>> randcap <- runif(N, 0,1) >>>> y[ randcap <= p, j] <- "A" >>>> #That would replace with "A" (within the j-column) ... >>>> # only the rows in 'y' for which randcap were less than the randcap >>>> threshold. >>>> >>>> #------------ >>>> >>>> >>>> y[q,j]="A"} >>>> }}else{ >>>> >>>> >>>> >>>> >>>> deado=runif(1,0,1) >>>> if(deado<=PsiAd){ >>>> y[q,j]="d" >>>> }else(y[q,j]="D") >>>> >>>> >>>> # -------Perhaps >>>> deado <- runif(N, 0,1) >>>> y[ , j] <- ifelse( deado<=PsiAd, "d", "D") >>>> #------------------------ >>>> >>>> >>>> >>>> if(y[q,j]%in%c("D","d")){ >>>> break >>>> >>>> >>>> # ----------Really? I thought that condition was assured at this >>>> point? >>>> >>>> >>>> >>> I thought so too but non-zero elements appeared toward the right in rows >> that had already "died". >> >>>> } >>>> } >>>> }else{ >>>> if(observable==1){ >>>> recap=runif(1,0,1) >>>> if(recap<=c){ >>>> y[q,j]="A" >>>> } >>>> } >>>> } >>>> } >>>> } >>>> >>>> >>>> >>>> >>>> >>>> >>>> There are a lot of unnecessary tests and conditions. I tried to break >>>> down the code and write the tests that have been done when assigning a >>>> variable. I simplified your the first part but cannot guarantee that all >>>> criteria are met. >>>> >>>> >>>> >>> I will run the three modified versions and see how things go. I hope my >> answers are helpful. >> >> >> ______________________________________________ [hidden email] 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. |
|
|
Thank you all for your help... below is the code that worked as I intended:
J <- 48 N <- 1000 S <- 0.9 PsiADD <- 0.6 PsiAA <- 0.4 PsiaA <- 0.3 p <- 0.5 c <- p y <- matrix(0,N,J) y[,1]="A" dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 12)[seq_len(J)] for (q in 1:N){ (observable <- TRUE) for(j in which(!dtp)){ if(j>1) if(y[q, j-1] %in% c("d", "D")) break { if(j %% 2){ if (runif(1,0,1)<=S) { if (observable) { observable <- (runif(1,0,1)<=PsiAA) } else { observable <- (runif(1,0,1)<=PsiaA) } if (observable) { if (runif(1,0,1) <= p) y[q,j] <- "A" } } else { y[q,j] <- ifelse(runif(1,0,1) <= PsiADD,"D","d") break } } else { if(observable){ if(runif(1,0,1) <= c) y[q,j] <- "A" } } } } } for (q in 1:N) { for (j in which(dtp)) { if(j>2) if(y[q, j-2] %in% c("d", "D")) break { if (j %% 2) { if (observable) { observable <- runif(1, 0, 1) <= PsiAA } else { observable <- runif(1, 0, 1) <= PsiaA } if (observable) { if (runif(1, 0, 1) <= p) y[q, j] <- "A" } } else { if (observable) { if (runif(1, 0, 1) <= c) y[q, j] <- "A" } } } } } On Sat, Aug 11, 2012 at 12:15 PM, Rui Barradas <[hidden email]> wrote: > Hello, > > I'm not sure it works, but try the following. > > > for(j in which(dtp)){ > for (q in 1:N){ > if(y[q, j] %in% c("d", "D")) break > [...etc...] > > and in the other loop the same, > > > for (j in which(!dtp)) { > for (q in 1:N) { > if(y[q, j] %in% c("d", "D")) break > [...etc...] > > > Em 10-08-2012 20:42, Claudia Penaloza escreveu: > >> This is what my code looks like now. However, there is one thing I >> can't/don't know how to fix. >> I can't get it to be "once dead always dead", i.e., in any given row, >> after >> a "D" or a "d" there should be only zeros. >> I've tried applying a flag to break the loop if dead but I can't get it to >> work. >> Could you please help me with this? >> >> Thank you for your time, >> Claudia >> >> J = 24 >> N = 10 >> S = .9 >> PsiADd = 0.4 >> PsiAA = .4 >> PsiaA = .3 >> p = .5 >> c = p >> >> y <- matrix(0,N,J) >> y[,1]="A" >> dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] >> for(j in which(dtp)){ >> for (q in 1:N){ >> (observable <- TRUE) >> if(j %% 2){ >> if (runif(1,0,1)<=S) { >> if (observable) { >> observable <- (runif(1,0,1)<=PsiAA) >> } else { >> observable <- (runif(1,0,1)<=PsiaA) >> } >> if (observable) { >> if (runif(1,0,1) <= p) y[q,j] <- "A" >> } >> } else { >> y[q,j] <- ifelse(runif(1,0,1) <= PsiADd,"d","D") >> break >> } >> } else { >> if(observable){ >> if(runif(1,0,1) <= c) y[q,j] <- "A" >> } >> } >> } >> } >> >> for (j in which(!dtp)) { >> for (q in 1:N) { >> if (j %% 2) { >> if (observable) { >> observable <- runif(1, 0, 1) <= PsiAA >> } else { >> observable <- runif(1, 0, 1) <= PsiaA >> } >> if (observable) { >> if (runif(1, 0, 1) <= p) y[q, j] <- "A" >> } >> } else { >> if (observable) { >> if (runif(1, 0, 1) <= c) y[q, j] <- "A" >> } >> } >> } >> } >> On Wed, Aug 8, 2012 at 2:04 PM, Claudia Penaloza >> <[hidden email]>**wrote: >> >> Answers inserted in BLUE below >>> >>> On Thu, Aug 2, 2012 at 5:34 PM, Claudia Penaloza < >>> [hidden email]> wrote: >>> >>> Thank you very much for all your suggestions. I am very sorry my code is >>>> so crude (it gives me a headache too!), I'm very new to R programing. I >>>> will have to look at your suggestions/questions very carefully (I'm >>>> barely >>>> holding on at this point) and get back to you (Dr. Winsemius) with some >>>> answers. >>>> >>>> Thank you! >>>> Claudia >>>> >>>> On Wed, Aug 1, 2012 at 2:11 PM, David Winsemius <[hidden email] >>>> >wrote: >>>> >>>> On Aug 1, 2012, at 12:02 PM, Mercier Eloi wrote: >>>>> >>>>> Your code almost gave me a headache. :-/ >>>>> I had a similar reaction. However, my approach might have been to >>>>> request a more complete verbal description of the data structures being >>>>> operated on and the methods and assumptions being used. Generally it is >>>>> going to be easier to go from a procedural description to good R code >>>>> than >>>>> it is to go from a SAS Data Proc ... even if it were tested and >>>>> debugged... >>>>> and yours was not even debugged. >>>>> >>>>> >>>>> ##############################****############################** >>>>> ##**#### >>>>> >>>>>> # Robust design with Markovian emigration and dummy time periods >>>>>> ##############################****############################** >>>>>> ##**#### >>>>>> >>>>>> J = 24 >>>>>> N = 10 >>>>>> S = .9 >>>>>> PsiAD = 0.06 >>>>>> PsiAd = 0.04 >>>>>> PsiAA = .4 >>>>>> PsiaA = .3 >>>>>> p = .5 >>>>>> c = p >>>>>> y <- matrix(0,N,J) >>>>>> >>>>>> # So is 'y' a matrix of states where the N rows are levels and the J >>>>> columns are discrete times? >>>>> # Actually I decided that context suggested you were using letters to >>>>> represent states. >>>>> >>>>> Yes, 'y' is a matrix of state, N rows are levels (independent of each >>>> >>> other) and J columns are discrete times. >>> Yes, I am using letters to represent states. >>> >>> y[,1]="A" >>>>> # So we start with N <something>'s in state "A"? >>>>> # It seems as though it might be the case that every row is independent >>>>> of the others. >>>>> # .. and you are performing ( replicate()-ing in R terminology) this >>>>> test N times >>>>> # states: >>>>> #("A" "D") >>>>> #("a" "d") >>>>> #transitions: >>>>> matrix( c( runif(1, 0, 1) <= 0.4, # A -> A >>>>> runif(1,0,1) <= 0.3, # a -> A >>>>> runif(1,0,1) <= 0.06. # A -> D >>>>> runif(1,0,1) <= 0.04 # A -> d) ) >>>>> >>>>> # What is state "a"? >>>>> # How do you get from A to a? >>>>> # can you get from D or d to A or a? >>>>> >>>>> Yes, we start with N "individuals" in state "A", each individual is >>>> >>> independent of each other. >>> State "a" is an unobserved (!observable) version of state "A" >>> (biologically, an individual that has temporarily left the sampling area >>> and is not available for capture) >>> An individual in state "A" transitions to state "a" if it is observable >>> and doesn't stay (stay >= AA) >>> "D" and "d" are dead (absorbing) states. Once in either "D" or "d", the >>> individual can no longer transition and is no longer "captured" (the row >>> should only have zeros after a "D" or a "d") >>> >>> >>> >>>> >>>> # Not sure I have gotten the model assumptions down. >>>>> # Is D" or "d" an absorbing state such as "dead or "Dead"? >>>>> >>>>> Yes. >>>> dtp <- rep(c(TRUE, TRUE, FALSE, FALSE), 6)[seq_len(J)] >>>> >>>>> >>>>> #Presumably the "dummy time periods" >>>>> >>>>> Yes. >>>> >>> >>> for(j in which(dtp)){ # j will be c(1,2, 5,6, 9,10 , ..., 21,22) >>>>> >>>>> >>>>> >>>>> for (q in 1:N){ # This can almost surely become vectorized >>>>> >>>>> (observable=1) >>>>> if(j %% 2){survive=runif(1,0,1) >>>>> if(survive<=S){ >>>>> stay=runif(1,0,1) >>>>> if (observable==1){ >>>>> >>>>> >>>>> # It would help if the concept of "observable" were explained >>>>> # Is this something to do with the capture-recapture observables? >>>>> >>>>> Yes. All individuals start out "observable" (we need to capture them a >>> first time). Individuals can then stay in their current state or >>> transition >>> to another one, if they are observable they can continue to be >>> observable, >>> or they can become unobservable and viceversa. Transition depends on what >>> state you are in at any given time (A->A != A->a != a->A != a->a). >>> >>> if{ >>>>> observable=1 >>>>> }else{observable=0} >>>>> }else{ >>>>> >>>>> >>>>> # After allowing for Mercier's astute observation that observable will >>>>> always be 1, I'm wondering if it can be replaced with just this code >>>>> (and >>>>> not need the loop surrounding it.) >>>>> >>>>> >>>>> >>>>> The individual will always "start" as observable... but as the loop >>>> >>> progresses through the columns in a row, it can go between being >>> observable >>> and unobservable (at least that is what I wanted to code). >>> >>> observable <- (stay<=PsiAA) >>>>> >>>>> #-------------- >>>>> >>>>> return=runif(1,0,1) >>>>> >>>>> >>>>> # better NOT to use the name "return" since it is a crucial R function >>>>> name. >>>>> >>>>> >>>>> >>>>> Will change. Thank you. >>>> >>>>> >>>>> >>>>> if(return<=PsiaA){ >>>>> observable=1 >>>>> }else{observable=0}} >>>>> >>>>> >>>>> #------- perhaps: >>>>> >>>>> randret <- return=runif(N,0,1) >>>>> observable <- randret <= PsiaA >>>>> # ------------------- >>>>> >>>>> >>>>> if(observable==1){ >>>>> capture=runif(1,0,1) >>>>> if(capture<=p){ >>>>> >>>>> >>>>> #---------- perhaps: >>>>> randcap <- runif(N, 0,1) >>>>> y[ randcap <= p, j] <- "A" >>>>> #That would replace with "A" (within the j-column) ... >>>>> # only the rows in 'y' for which randcap were less than the randcap >>>>> threshold. >>>>> >>>>> #------------ >>>>> >>>>> >>>>> y[q,j]="A"} >>>>> }}else{ >>>>> >>>>> >>>>> >>>>> >>>>> deado=runif(1,0,1) >>>>> if(deado<=PsiAd){ >>>>> y[q,j]="d" >>>>> }else(y[q,j]="D") >>>>> >>>>> >>>>> # -------Perhaps >>>>> deado <- runif(N, 0,1) >>>>> y[ , j] <- ifelse( deado<=PsiAd, "d", "D") >>>>> #------------------------ >>>>> >>>>> >>>>> >>>>> if(y[q,j]%in%c("D","d")){ >>>>> break >>>>> >>>>> >>>>> # ----------Really? I thought that condition was assured at this >>>>> point? >>>>> >>>>> >>>>> >>>>> I thought so too but non-zero elements appeared toward the right in >>>> rows >>>> >>> that had already "died". >>> >>> } >>>>> } >>>>> }else{ >>>>> if(observable==1){ >>>>> recap=runif(1,0,1) >>>>> if(recap<=c){ >>>>> y[q,j]="A" >>>>> } >>>>> } >>>>> } >>>>> } >>>>> } >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> There are a lot of unnecessary tests and conditions. I tried to break >>>>> down the code and write the tests that have been done when assigning a >>>>> variable. I simplified your the first part but cannot guarantee that >>>>> all >>>>> criteria are met. >>>>> >>>>> >>>>> >>>>> I will run the three modified versions and see how things go. I hope >>>> my >>>> >>> answers are helpful. >>> >>> >>> >>> > [[alternative HTML version deleted]] ______________________________________________ [hidden email] 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. |
«
Return to R help
|
1 view|%1 views
Loading...
| Powered by Nabble | Edit this page |