Re: [R] simulation of PowerGARCH,Threshold GARCH,and GJR GARCH
Doing a web search on R CRAN GJR GARCH brought up the rugarch package. The models you mentioned are discussed in the documentation to that package https://cran.r-project.org/web/packages/rugarch/vignettes/Introduction_to_the_rugarch_package.pdf On Mon, Mar 25, 2019 at 2:06 PM Amon kiregu wrote: > what is the r code for simulating PowerGARCH,Threshold GARCH,and GJR GARCH > in order to capture heteroscedasticity,volatility clustering,etc,,so that i > can have simulation of mean part and simulation on innovation part. > thanks > > [[alternative HTML version deleted]] > > __ > R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see > 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. > [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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: [R] Simulation based on runif to get mean
Hello, Right. Missed that one. Rui Barradas Enviado a partir do meu smartphone Samsung Galaxy. Mensagem original De: Eric Berger <ericjber...@gmail.com> Data: 30/01/2018 10:12 (GMT+00:00) Para: Rui Barradas <ruipbarra...@sapo.pt> Cc: Daniel Nordlund <djnordl...@gmail.com>, smart hendsome <putra_autum...@yahoo.com>, r-help@r-project.org Assunto: Re: [R] Simulation based on runif to get mean Or a shorter version of Rui's approach: set.seed(2511) # Make the results reproduciblefun <- function(n){ f <- function(){ c(mean(runif(5,1,10)),mean(runif(5,10,20))) } replicate(n, f())}fun(10) On Tue, Jan 30, 2018 at 12:03 PM, Rui Barradas <ruipbarra...@sapo.pt> wrote: Hello, Another way would be to use ?replicate and ?colMeans. set.seed(2511) # Make the results reproducible fun <- function(n){ f <- function(){ a <- runif(5, 1, 10) b <- runif(5, 10, 20) colMeans(cbind(a, b)) } replicate(n, f()) } fun(10) Hope this helps, Rui Barradas On 1/30/2018 8:58 AM, Daniel Nordlund wrote: On 1/29/2018 9:03 PM, smart hendsome via R-help wrote: Hello everyone, I have a question regarding simulating based on runif. Let say I have generated matrix A and B based on runif. Then I find mean for each matrix A and matrix B. I want this process to be done let say 10 times. Anyone can help me. Actually I want make the function that I can play around with the number of simulation process that I want. Thanks. Eg: a <- matrix(runif(5,1, 10)) b <- matrix(runif(5,10, 20)) c <- cbind(a,b); c mn <- apply(c,2,mean); mn Regards, Zuhri [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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. Here is a straight forward implementation of your code in a function with a parameter for the number simulations you want to run. sim <- function(n){ mn <- matrix(0,n, 2) for(i in 1:n) { a <- runif(5,1, 10) b <- runif(5,10, 20) c <- cbind(a,b) mn[i,] <- apply(c, 2, mean) } return(mn) } # run 10 iterations sim(10) In your case, there doesn't seem to be a need to create a and b as matrices; vectors work just as well. Also, several of the statements could be combined into one. Whether this meets your needs depends on what your real world task actually is. Hope this is helpful, Dan __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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: [R] Simulation based on runif to get mean
Or a shorter version of Rui's approach: set.seed(2511)# Make the results reproducible fun <- function(n){ f <- function(){ c(mean(runif(5,1,10)),mean(runif(5,10,20))) } replicate(n, f()) } fun(10) On Tue, Jan 30, 2018 at 12:03 PM, Rui Barradaswrote: > Hello, > > Another way would be to use ?replicate and ?colMeans. > > > set.seed(2511)# Make the results reproducible > > fun <- function(n){ > f <- function(){ > a <- runif(5, 1, 10) > b <- runif(5, 10, 20) > colMeans(cbind(a, b)) > } > replicate(n, f()) > } > > fun(10) > > Hope this helps, > > Rui Barradas > > > On 1/30/2018 8:58 AM, Daniel Nordlund wrote: > >> On 1/29/2018 9:03 PM, smart hendsome via R-help wrote: >> >>> Hello everyone, >>> I have a question regarding simulating based on runif. Let say I have >>> generated matrix A and B based on runif. Then I find mean for each matrix A >>> and matrix B. I want this process to be done let say 10 times. Anyone can >>> help me. Actually I want make the function that I can play around with the >>> number of simulation process that I want. Thanks. >>> Eg: >>> a <- matrix(runif(5,1, 10)) >>> >>> b <- matrix(runif(5,10, 20)) >>> >>> c <- cbind(a,b); c >>> >>> mn <- apply(c,2,mean); mn >>> >>> Regards, >>> Zuhri >>> >>> [[alternative HTML version deleted]] >>> >>> __ >>> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see >>> https://stat.ethz.ch/mailman/listinfo/r-help >>> PLEASE do read the posting guide http://www.R-project.org/posti >>> ng-guide.html >>> and provide commented, minimal, self-contained, reproducible code. >>> >>> >> Here is a straight forward implementation of your code in a function with >> a parameter for the number simulations you want to run. >> >> sim <- function(n){ >>mn <- matrix(0,n, 2) >>for(i in 1:n) { >> a <- runif(5,1, 10) >> b <- runif(5,10, 20) >> c <- cbind(a,b) >> mn[i,] <- apply(c, 2, mean) >> } >>return(mn) >>} >> # run 10 iterations >> sim(10) >> >> In your case, there doesn't seem to be a need to create a and b as >> matrices; vectors work just as well. Also, several of the statements could >> be combined into one. Whether this meets your needs depends on what your >> real world task actually is. >> >> >> Hope this is helpful, >> >> Dan >> >> > __ > R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posti > ng-guide.html > and provide commented, minimal, self-contained, reproducible code. > [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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: [R] Simulation based on runif to get mean
Hello, Another way would be to use ?replicate and ?colMeans. set.seed(2511)# Make the results reproducible fun <- function(n){ f <- function(){ a <- runif(5, 1, 10) b <- runif(5, 10, 20) colMeans(cbind(a, b)) } replicate(n, f()) } fun(10) Hope this helps, Rui Barradas On 1/30/2018 8:58 AM, Daniel Nordlund wrote: On 1/29/2018 9:03 PM, smart hendsome via R-help wrote: Hello everyone, I have a question regarding simulating based on runif. Let say I have generated matrix A and B based on runif. Then I find mean for each matrix A and matrix B. I want this process to be done let say 10 times. Anyone can help me. Actually I want make the function that I can play around with the number of simulation process that I want. Thanks. Eg: a <- matrix(runif(5,1, 10)) b <- matrix(runif(5,10, 20)) c <- cbind(a,b); c mn <- apply(c,2,mean); mn Regards, Zuhri [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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. Here is a straight forward implementation of your code in a function with a parameter for the number simulations you want to run. sim <- function(n){ mn <- matrix(0,n, 2) for(i in 1:n) { a <- runif(5,1, 10) b <- runif(5,10, 20) c <- cbind(a,b) mn[i,] <- apply(c, 2, mean) } return(mn) } # run 10 iterations sim(10) In your case, there doesn't seem to be a need to create a and b as matrices; vectors work just as well. Also, several of the statements could be combined into one. Whether this meets your needs depends on what your real world task actually is. Hope this is helpful, Dan __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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: [R] Simulation based on runif to get mean
On 1/29/2018 9:03 PM, smart hendsome via R-help wrote: Hello everyone, I have a question regarding simulating based on runif. Let say I have generated matrix A and B based on runif. Then I find mean for each matrix A and matrix B. I want this process to be done let say 10 times. Anyone can help me. Actually I want make the function that I can play around with the number of simulation process that I want. Thanks. Eg: a <- matrix(runif(5,1, 10)) b <- matrix(runif(5,10, 20)) c <- cbind(a,b); c mn <- apply(c,2,mean); mn Regards, Zuhri [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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. Here is a straight forward implementation of your code in a function with a parameter for the number simulations you want to run. sim <- function(n){ mn <- matrix(0,n, 2) for(i in 1:n) { a <- runif(5,1, 10) b <- runif(5,10, 20) c <- cbind(a,b) mn[i,] <- apply(c, 2, mean) } return(mn) } # run 10 iterations sim(10) In your case, there doesn't seem to be a need to create a and b as matrices; vectors work just as well. Also, several of the statements could be combined into one. Whether this meets your needs depends on what your real world task actually is. Hope this is helpful, Dan -- Daniel Nordlund Port Townsend, WA USA __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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: [R] simulation in R
I don't think we have enough information to help you with this. Do you intent the simulated values for growthrate to be selected from the values in daT$growthrate? Or do these values define a distribution of values (perhaps ranging between 0 and 1) and the simulation should use that empirical distribution? In that case any value of growthrate in the range of 0 and 1 inclusive is possible. What is the sample size of the Monte Carlo simulations? This will have a major effect on the ranges since as the sample size increases the range will be closer and closer to the range of the growthrate (0 to 1 for condition 3 and 0 to .5 for condition 1). In your example of no constraints on food (condition 2), if the growthrate is always 1, the range will always be 1 and 1. - David L Carlson Department of Anthropology Texas A University College Station, TX 77840-4352 -Original Message- From: R-help [mailto:r-help-boun...@r-project.org] On Behalf Of Kristi Glover Sent: Wednesday, April 20, 2016 1:06 AM To: R-help Subject: [R] simulation in R I realized that there was a typo error. I mean "Monte Carlo Simulation" From: R-helpon behalf of Kristi Glover Sent: April 19, 2016 11:48 PM To: R-help Subject: [R] simulation in R Hi R user, Would you mind to help me to find the range with stochastic events? For example, daT<-structure(list(sn = 1:14, growthrate = c(0.5, 0.6, 0.7, 0.99, 0.1, 0.3, 0.4, 0.5, 0.5, 0.2, 0.1, 0.4, 0.3, 0.43)), .Names = c("sn", "growthrate"), class = "data.frame", row.names = c(NA, -14L)) I want to find the ranges of growth rate of the above data using Mote corle simulation ( times) under three conditions: 1. very drought ( in that condition growth will not be more than 0.5). [what would be the range (max, min ) of the growth rate for this scenario) 2. no constraints of food (growth will be 1 or 100%) (what would be the range (max,min) of growth rate in this scenario?). 3. Control (as it is) (Range??, max.min) I tried to find whether some one had same problem but I could not find it, is it too complicated to write the code in R for this example? your help will be highly appreciated. Sincerely, KG [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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. __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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. __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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: [R] simulation in vector autoregressive model (VAR)
Before you post on Rhelp you should first read the Posting Guide. It tells you that requests for statistical advice might be answered but are not really on-topic. Furthermore requests for tutorials or extended worked examples should probably be accompanied by evidence of searching using Google or equivalent: https://www.google.com/search?q=r+vector+auto-regressive+model=utf-8=utf-8 -- David On Sep 3, 2015, at 11:53 PM, Aziz Mensah via R-help wrote: > I have a data from 4 variables ( STOCK, CPI, EXC, and CCI) from 1980 to 2012. > I want to do a forecast using VAR(12) model with a simulation of 100,000 for > 5 years. And also estimate the RMSE, MAPE, and Theil Inequality. Can anyone > help me with this problem in R? Thanks so much. > > [[alternative HTML version deleted]] > > __ > R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see > 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 Alameda, CA, USA __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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: [R] simulation data in SEM
PLZ mke sure the package installed which contains mvrnorm function. -- PO SU mail: desolato...@163.com Majored in Statistics from SJTU At 2014-09-14 09:56:34, thanoon younis thanoon.youni...@gmail.com wrote: Dear R members I want to simulate data depending on SEM and when i applied the code below i found some errors and i still cannot run it. many thanks in advance Thanoon #Do simulation for 100 replications N-1000; P-10 phi-matrix(data=c(1.0,0.3,0.3,1.0),ncol=2) #The covariance matrix of xi Ro-matrix(data=c(7.0,2.1,2.1,7.0), ncol=2) yo-matrix(data=NA,nrow=N,ncol=P) p-numeric(P); v-numeric(P) for (t in 1:100) { #Generate the data for the simulation study for (i in 1:N) { #Generate xi xi-mvrnorm(1,mu=c(0,0),phi) #Generate the fixed covariates co-rnorm(1,0,1) #Generate error term is structural equation delta-rnorm(1,0,sqrt(0.3)) #Generate eta1 according to the structural equation eta-0.8*co[i]+0.6*xi[1]+0.6*xi[2]+0.8*xi[1]*xi[2]+delta #Generate error terms in measurement equations eps-rnorm(3,0,1) #Generate theta according to measurement equations v1[1]-1.0+eta+eps[1]; v1[2]-1.0+0.7*eta+eps[2] v1[3]-1.0+0.7*eta+eps[3] v1[4]-1.0+xi[1]; v1[5]-1.0+0.8*xi[1]; v1[6]-1.0+0.8*xi[1] v1[7]-1.0+xi[2]; v1[8]-1.0+0.7*xi[2]; v1[9]-1.0+0.7*xi[2]; v1[10]-1.0+0.7*xi1[2] #transform theta to orinal variables for (j in 1:10) { if (v[j]0) yo[i,j]-1 else yo[i,j]-0 } #Input data set for WinBUGS data-list(N=200,P=10,R=Ro,z=yo) } #end [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation data in SEM
What errors? What is your output? What output did you expect? On Sep 14, 2014, at 6:56 AM, thanoon younis thanoon.youni...@gmail.com wrote: Dear R members I want to simulate data depending on SEM and when i applied the code below i found some errors and i still cannot run it. many thanks in advance Thanoon #Do simulation for 100 replications N-1000; P-10 phi-matrix(data=c(1.0,0.3,0.3,1.0),ncol=2) #The covariance matrix of xi Ro-matrix(data=c(7.0,2.1,2.1,7.0), ncol=2) yo-matrix(data=NA,nrow=N,ncol=P) p-numeric(P); v-numeric(P) for (t in 1:100) { #Generate the data for the simulation study for (i in 1:N) { #Generate xi xi-mvrnorm(1,mu=c(0,0),phi) #Generate the fixed covariates co-rnorm(1,0,1) #Generate error term is structural equation delta-rnorm(1,0,sqrt(0.3)) #Generate eta1 according to the structural equation eta-0.8*co[i]+0.6*xi[1]+0.6*xi[2]+0.8*xi[1]*xi[2]+delta #Generate error terms in measurement equations eps-rnorm(3,0,1) #Generate theta according to measurement equations v1[1]-1.0+eta+eps[1]; v1[2]-1.0+0.7*eta+eps[2] v1[3]-1.0+0.7*eta+eps[3] v1[4]-1.0+xi[1]; v1[5]-1.0+0.8*xi[1]; v1[6]-1.0+0.8*xi[1] v1[7]-1.0+xi[2]; v1[8]-1.0+0.7*xi[2]; v1[9]-1.0+0.7*xi[2]; v1[10]-1.0+0.7*xi1[2] #transform theta to orinal variables for (j in 1:10) { if (v[j]0) yo[i,j]-1 else yo[i,j]-0 } #Input data set for WinBUGS data-list(N=200,P=10,R=Ro,z=yo) } #end [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Don McKenzie Research Ecologist Pacific Wildland Fire Sciences Lab US Forest Service Affiliate Professor School of Environmental and Forest Sciences University of Washington d...@uw.edu __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation data in SEM
cc�ing to list, as requested in the posting guide, so that others may be able to help you. On Sep 14, 2014, at 8:45 AM, thanoon younis thanoon.youni...@gmail.com wrote: Thank you very much for your reply the output is #Do simulation for 100 replications N-1000; P-10 phi-matrix(data=c(1.0,0.3,0.3,1.0),ncol=2) #The covariance matrix of xi Ro-matrix(data=c(7.0,2.1,2.1,7.0), ncol=2) yo-matrix(data=NA,nrow=N,ncol=P) p-numeric(P); v-numeric(P) Error: unexpected symbol in yo-matrix(data=NA,nrow=N,ncol=P) p for (t in 1:100) { + #Generate the data for the simulation study + for (i in 1:N) { + #Generate xi + xi-mvrnorm(1,mu=c(0,0),phi) + #Generate the fixed covariates + co-rnorm(1,0,1) + #Generate error term is structural equation + delta-rnorm(1,0,sqrt(0.3)) + #Generate eta1 according to the structural equation + eta-0.8*co[i]+0.6*xi[1]+0.6*xi[2]+0.8*xi[1]*xi[2]+delta + #Generate error terms in measurement equations + eps-rnorm(3,0,1) + + #Generate theta according to measurement equations + v1[1]-1.0+eta+eps[1]; v1[2]-1.0+0.7*eta+eps[2] + v1[3]-1.0+0.7*eta+eps[3] + v1[4]-1.0+xi[1]; v1[5]-1.0+0.8*xi[1]; v1[6]-1.0+0.8*xi[1] + v1[7]-1.0+xi[2]; v1[8]-1.0+0.7*xi[2]; v1[9]-1.0+0.7*xi[2]; + v1[10]-1.0+0.7*xi1[2] + + + #transform theta to orinal variables + for (j in 1:10) { if (v[j]0) yo[i,j]-1 else yo[i,j]-0 } + + + #Input data set for WinBUGS + data-list(N=200,P=10,R=Ro,z=yo) + + } #end + also i cannot continue to get on a data. many thanks again Thanoon On 14 September 2014 18:41, Don McKenzie d...@u.washington.edu wrote: What errors? What is your output? What output did you expect? On Sep 14, 2014, at 6:56 AM, thanoon younis thanoon.youni...@gmail.com wrote: Dear R members I want to simulate data depending on SEM and when i applied the code below i found some errors and i still cannot run it. many thanks in advance Thanoon #Do simulation for 100 replications N-1000; P-10 phi-matrix(data=c(1.0,0.3,0.3,1.0),ncol=2) #The covariance matrix of xi Ro-matrix(data=c(7.0,2.1,2.1,7.0), ncol=2) yo-matrix(data=NA,nrow=N,ncol=P) p-numeric(P); v-numeric(P) for (t in 1:100) { #Generate the data for the simulation study for (i in 1:N) { #Generate xi xi-mvrnorm(1,mu=c(0,0),phi) #Generate the fixed covariates co-rnorm(1,0,1) #Generate error term is structural equation delta-rnorm(1,0,sqrt(0.3)) #Generate eta1 according to the structural equation eta-0.8*co[i]+0.6*xi[1]+0.6*xi[2]+0.8*xi[1]*xi[2]+delta #Generate error terms in measurement equations eps-rnorm(3,0,1) #Generate theta according to measurement equations v1[1]-1.0+eta+eps[1]; v1[2]-1.0+0.7*eta+eps[2] v1[3]-1.0+0.7*eta+eps[3] v1[4]-1.0+xi[1]; v1[5]-1.0+0.8*xi[1]; v1[6]-1.0+0.8*xi[1] v1[7]-1.0+xi[2]; v1[8]-1.0+0.7*xi[2]; v1[9]-1.0+0.7*xi[2]; v1[10]-1.0+0.7*xi1[2] #transform theta to orinal variables for (j in 1:10) { if (v[j]0) yo[i,j]-1 else yo[i,j]-0 } #Input data set for WinBUGS data-list(N=200,P=10,R=Ro,z=yo) } #end [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Don McKenzie Research Ecologist Pacific Wildland Fire Sciences Lab US Forest Service Affiliate Professor School of Environmental and Forest Sciences University of Washington d...@uw.edu Don McKenzie Research Ecologist Pacific Wildland Fire Sciences Lab US Forest Service Affiliate Professor School of Environmental and Forest Sciences University of Washington d...@uw.edu [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation data in SEM
On Sep 14, 2014, at 8:48 AM, Don McKenzie wrote: ccing to list, as requested in the posting guide, so that others may be able to help you. On Sep 14, 2014, at 8:45 AM, thanoon younis thanoon.youni...@gmail.com wrote: Thank you very much for your reply the output is #Do simulation for 100 replications N-1000; P-10 phi-matrix(data=c(1.0,0.3,0.3,1.0),ncol=2) #The covariance matrix of xi Ro-matrix(data=c(7.0,2.1,2.1,7.0), ncol=2) yo-matrix(data=NA,nrow=N,ncol=P) p-numeric(P); v-numeric(P) Error: unexpected symbol in yo-matrix(data=NA,nrow=N,ncol=P) p Almost any time you see an error message that says : unexpected _something_ it means you submitted a malformed expression to the interpreter that was missing a paren or a bracket or a comma or _something_. You always need to go back to the left of where the error was discovered. In this case you are missing a semicolon about here: yo-matrix(data=NA,nrow=N,ncol=P) p-numeric(P); v-numeric(P) ^ -- David. for (t in 1:100) { + #Generate the data for the simulation study + for (i in 1:N) { + #Generate xi + xi-mvrnorm(1,mu=c(0,0),phi) + #Generate the fixed covariates + co-rnorm(1,0,1) snipped -- David Winsemius, MD Alameda, CA, USA __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation data in SEM
Adding back the list address: On Sep 14, 2014, at 9:53 AM, thanoon younis wrote: thank you for your help but i still have error after putting semicolon Error: unexpected symbol in: Ro-matrix(data=c(7.0,2.1,2.1,7.0), ncol=2) yo-matrix(data=NA,nrow=N,ncol=P) p The error message shows no semicolon in the location I pointed at that was missing one. Furthermore the error message is now attaching the prior line which should not have thrown an error. Since you didn't include the actual code block that was your revision we can only guess (and I do not have a good guess why that is now occurring unless perhaps your font has two different encodings for semi-colon glyphs.) You typed (or at least that is what I see in my mail client): yo-matrix(data=NA,nrow=N,ncol=P) p-numeric(P); v-numeric(P) I suggested: yo-matrix(data=NA,nrow=N,ncol=P); p-numeric(P); v-numeric(P) PLEASE read the Posting Guide. I will not respond to offlist messages from you in the future. -- David. many thanks again On 14 September 2014 19:37, David Winsemius dwinsem...@comcast.net wrote: On Sep 14, 2014, at 8:48 AM, Don McKenzie wrote: cc’ing to list, as requested in the posting guide, so that others may be able to help you. On Sep 14, 2014, at 8:45 AM, thanoon younis thanoon.youni...@gmail.com wrote: Thank you very much for your reply the output is #Do simulation for 100 replications N-1000; P-10 phi-matrix(data=c(1.0,0.3,0.3,1.0),ncol=2) #The covariance matrix of xi Ro-matrix(data=c(7.0,2.1,2.1,7.0), ncol=2) yo-matrix(data=NA,nrow=N,ncol=P) p-numeric(P); v-numeric(P) Error: unexpected symbol in yo-matrix(data=NA,nrow=N,ncol=P) p Almost any time you see an error message that says : unexpected _something_ it means you submitted a malformed expression to the interpreter that was missing a paren or a bracket or a comma or _something_. You always need to go back to the left of where the error was discovered. In this case you are missing a semicolon about here: yo-matrix(data=NA,nrow=N,ncol=P) p-numeric(P); v-numeric(P) ^ -- David. for (t in 1:100) { + #Generate the data for the simulation study + for (i in 1:N) { + #Generate xi + xi-mvrnorm(1,mu=c(0,0),phi) + #Generate the fixed covariates + co-rnorm(1,0,1) snipped -- David Winsemius, MD Alameda, CA, USA David Winsemius Alameda, CA, USA __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation data with dichotomous varuables
Dear Thanoon, You might look at the various item simulation functions in the psych package. In particular, for your problem: R1 - sim.irt(10,1000,a=3,low = -2, high=2) R2 - sim.irt(10,1000,a=3,low = -2, high=2) R12 - data.frame(R1$items,R2$items) #this gives you 20 items, grouped with high correlations within the first 10, and the second 10, no correlation between the first and second sets. rho - tetrachoric(R12)$rho #find the tetrachoric correlation between the items lowerMat(rho) #show the correlations cor.plot(rho,numbers=TRUE) #show a heat map of the correlations Bill On Aug 4, 2014, at 8:08 PM, thanoon younis thanoon.youni...@gmail.com wrote: Dear R-users i need your help to solve my problem in the code below, i want to simulate two different samples R1 and R2 and each sample has 10 variables and 1000 observations so i want to simulate a data with high correlation between var. in R1 and also in R2 and no correlation between R1 and R2 also i have a problem with correlation coefficient between tow dichotomous var. the R- program supports just these types of correlation coefficients such as pearson, spearman,kendall. thanks alot in advance Thanoon ords - seq(0,1) p - 10 N - 1000 percent_change - 0.9 R1 - as.data.frame(replicate(p, sample(ords, N, replace = T))) R2 - as.data.frame(replicate(p, sample(ords, N, replace = T))) # pearson is more appropriate for dichotomous data cor(R1, R2, method = pearson) # subset variable to have a stronger correlation v1 - R1[,1, drop = FALSE] v1 - R2[,1, drop = FALSE] # randomly choose which rows to retain keep - sample(as.numeric(rownames(v1)), size = percent_change*nrow(v1)) change - as.numeric(rownames(v1)[-keep]) # randomly choose new values for changing new.change - sample(ords, ((1-percent_change)*N)+1, replace = T) # replace values in copy of original column v1.samp - v1 v1.samp[change,] - new.change # closer correlation cor(v1, v1.samp, method = pearson) # set correlated column as one of your other columns R1[,2] - v1.samp R2[,2] - v1.samp R1 R2 [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. William Revellehttp://personality-project.org/revelle.html Professor http://personality-project.org Department of Psychology http://www.wcas.northwestern.edu/psych/ Northwestern Universityhttp://www.northwestern.edu/ Use R for psychology http://personality-project.org/r It is 5 minutes to midnighthttp://www.thebulletin.org __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation dichotomous data
Please remember the 'reply all' for the r-help page. First Question: How can i use Pearson correlation with dichotomous data? i want to use a correlation between dichotomous variables like spearman correlation in ordered categorical variables? cor(variable1, variable2, *method = pearson*) Second Question: Would like two separate populations (1000 samples, 10 var). Variables *within* datasets highly correlated, minimal correlation *between* datasets. As I have stated in a previous response, the code you have is sufficient. You can go through as many variables as you like *for each dataset* and induce correlations. You should do this for as many variables as you require to be correlated. As the code induces these correlations randomly, there should be *minimal* correlation between datasets but still some if the datasets have the same structure (same variables correlated within). If different variables are correlated within each, then the correlation between datasets would likely be lower. It is extremely unrealistic to believe that there will be absolutely no correlation between datasets so you must decide at which point you consider it sufficiently low. One final point, in the code section # subset variable to have a stronger correlation, you can only do one at a time or you must change the name of the second object otherwise you are just overwriting the previous 'v1'. You have described what you want to me and you have the code to do it. The major hurdle here would be an implementation of some 'for loops', which is not terribly complex if you are working on your programming. However, they are not necessary if you just want to write several lines with new object names for each variable in each dataset. Give it a try, you know how to induce correlations now. Just chose which variables to correlate and do it for all of those for each dataset and compare. Regards, Dr. Charles Determan On Thu, Jul 31, 2014 at 9:10 AM, thanoon younis thanoon.youni...@gmail.com wrote: Many thanks to you firstly : how can i use Pearson correlation with dichotomous data? i want to use a correlation between dichotomous variables like spearman correlation in ordered categorical variables. secondly: i have two different population and each population has 1000 samples and 10 var. so i want to put a high correlation coefficient between variables in the first population and also put a high correlation coefficient between variables in the second population and no correlation between two populations because i want to use multiple group structural equation models. many thanks again Thanoon On 31 July 2014 16:45, Charles Determan Jr deter...@umn.edu wrote: Thanoon, You should still send the question to the R help list even when I helped you with the code you are currently using. I will not always know the best way or even how to proceed with some questions. As for to your question with the code below. Firstly, there is no 'phi' method for cor in base R. If you are using it, you must have neglected to include a package you are using. However, given that the phi coefficient is equal to the pearson coefficient for dichotomous data, you can use the 'pearson' method. Secondly, with respect to your primary concern. In this case, we have randomly chosen variables to correlate between two INDEPENDENT DATASETS (i.e. different groups of samples). The idea with this code is that R1 and R2 are datasets of 1000 samples and 10 variables. It would be miraculous if they correlated when each had variables randomly assigned as correlated. The code work correctly, the question now becomes if you want to see correlations across variables for all samples (which this does for each DATASET) or if you want two DATASETS to be correlated. ords - seq(0,1) p - 10 N - 1000 percent_change - 0.9 R1 - as.data.frame(replicate(p, sample(ords, N, replace = T))) R2 - as.data.frame(replicate(p, sample(ords, N, replace = T))) # phi is more appropriate for dichotomous data cor(R1, method = phi) cor(R2, method = phi) # subset variable to have a stronger correlation v1 - R1[,1, drop = FALSE] v1 - R2[,1, drop = FALSE] # randomly choose which rows to retain keep - sample(as.numeric(rownames(v1)), size = percent_change*nrow(v1)) change - as.numeric(rownames(v1)[-keep]) # randomly choose new values for changing new.change - sample(ords, ((1-percent_change)*N)+1, replace = T) # replace values in copy of original column v1.samp - v1 v1.samp[change,] - new.change # closer correlation cor(v1, v1.samp, method = phi) # set correlated column as one of your other columns R1[,2] - v1.samp R2[,2] - v1.samp R1 R2 On Thu, Jul 31, 2014 at 7:29 AM, thanoon younis thanoon.youni...@gmail.com wrote: dear Dr. Charles i have a problem with the following R - program in simulation data with 2 different samples and with high correlation between variables in each sample so when i applied the
Re: [R] simulation dichotomous data
Thanoon, You should still send the question to the R help list even when I helped you with the code you are currently using. I will not always know the best way or even how to proceed with some questions. As for to your question with the code below. Firstly, there is no 'phi' method for cor in base R. If you are using it, you must have neglected to include a package you are using. However, given that the phi coefficient is equal to the pearson coefficient for dichotomous data, you can use the 'pearson' method. Secondly, with respect to your primary concern. In this case, we have randomly chosen variables to correlate between two INDEPENDENT DATASETS (i.e. different groups of samples). The idea with this code is that R1 and R2 are datasets of 1000 samples and 10 variables. It would be miraculous if they correlated when each had variables randomly assigned as correlated. The code work correctly, the question now becomes if you want to see correlations across variables for all samples (which this does for each DATASET) or if you want two DATASETS to be correlated. ords - seq(0,1) p - 10 N - 1000 percent_change - 0.9 R1 - as.data.frame(replicate(p, sample(ords, N, replace = T))) R2 - as.data.frame(replicate(p, sample(ords, N, replace = T))) # phi is more appropriate for dichotomous data cor(R1, method = phi) cor(R2, method = phi) # subset variable to have a stronger correlation v1 - R1[,1, drop = FALSE] v1 - R2[,1, drop = FALSE] # randomly choose which rows to retain keep - sample(as.numeric(rownames(v1)), size = percent_change*nrow(v1)) change - as.numeric(rownames(v1)[-keep]) # randomly choose new values for changing new.change - sample(ords, ((1-percent_change)*N)+1, replace = T) # replace values in copy of original column v1.samp - v1 v1.samp[change,] - new.change # closer correlation cor(v1, v1.samp, method = phi) # set correlated column as one of your other columns R1[,2] - v1.samp R2[,2] - v1.samp R1 R2 On Thu, Jul 31, 2014 at 7:29 AM, thanoon younis thanoon.youni...@gmail.com wrote: dear Dr. Charles i have a problem with the following R - program in simulation data with 2 different samples and with high correlation between variables in each sample so when i applied the program i got on a results but without correlation between each sample. i appreciate your help and your time i did not send this code to R- help because you helped me before to write it . many thanks to you Thanoon -- Dr. Charles Determan, PhD Integrated Biosciences [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation data
Hello, At an R prompt type ?rbinom ?replicate Hope this helps, Rui Barradas Em 10-04-2014 02:28, thanoon younis escreveu: hi i want to simulate multivariate dichotomous data matrix with categories (0,1) and n=1000 and p=10. thanks alot in advance [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation data
Thanoon, My reply to your previous post should be more than enough for you to accomplish your goal. Please look over that script again: ords - seq(4) p - 10 N - 1000 percent_change - 0.9 R - as.data.frame(replicate(p, sample(ords, N, replace = T))) or alternatively as Mr. Barradas suggests with rbinom(), I leave options for you to figure out. Look at the help page, feel free to experiment with different numbers and look at the output. It is important you learn how to explore new functions you are unfamiliar with. R - as.data.frame(replicate(p, rbinom(n=#, size=#, p=#))) These lists are meant to help people with their code but not do the work for them. Given your prior questions to me as well I strongly suggest you explore some R tutorials. There are dozens online that should help you with the basics and understand the above code more clearly. Also, regarding your prior question about tetratorich correlations. You got an error previously because it is not a standard correlation within the corr() function. You can get further information about a function by checking the help pages ?corr. You will need to try and find a package that provides a function to do so or write the function yourself. This may sound daunting but if you take some time to learn how to write functions and the method for the tetratorich correlation isn't that complex you should not have too much of a problem. Summing up: 1. Find some R tutorials to get the basics down 2. Try to understand the above code for your problem 3. Find a suitable R package for your specific correlation needs 4. or Learn to write functions and find the means to calculate the tetratorich correlation. Regards, Charles Determan On Wed, Apr 9, 2014 at 8:28 PM, thanoon younis thanoon.youni...@gmail.comwrote: hi i want to simulate multivariate dichotomous data matrix with categories (0,1) and n=1000 and p=10. thanks alot in advance [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Charles Determan Integrated Biosciences PhD Candidate University of Minnesota [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation data
On 10/04/14 09:28, thanoon younis wrote: hi i want to simulate multivariate dichotomous data matrix with categories (0,1) and n=1000 and p=10. Nobody's stopping you! :-) cheers, Rolf Turner __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation data
Thanoon, Firstly, please remember to reply to the R help list as well so that other may benefit from your questions as well. Regarding your second request, I have written the following as a very naive way of inducing correlations. Hopefully this makes it perfectly clear what you change for different sample sizes. ords - seq(4) p - 10 N - 1000 percent_change - 0.9 R - as.data.frame(replicate(p, sample(ords, N, replace = T))) # spearman is more appropriate for ordinal data cor(R, method = spearman) # subset variable to have a stronger correlation v1 - R[,1, drop = FALSE] # randomly choose which rows to retain keep - sample(as.numeric(rownames(v1)), size = percent_change*nrow(v1)) change - as.numeric(rownames(v1)[-keep]) # randomly choose new values for changing new.change - sample(ords, ((1-percent_change)*N)+1, replace = T) # replace values in copy of original column v1.samp - v1 v1.samp[change,] - new.change # closer correlation cor(v1, v1.samp, method = spearman) # set correlated column as one of your other columns R[,2] - v1.samp This obviously only creates a correlation between two columns. You need to decide what you expect from this synthetic dataset. Do you want perfect correlations? Does it matter which variables are correlated? How many variables will be correlated? Are there correlations between multiple variables? Do you want negative correlations (hint: opposite values)? All of these questions would be great exercises for you to improve your R. You can also turn the above code into a function and have it randomly select two columns to be correlated if that works for you. Because of all of these possibilities I cannot provide the 'right' code but rather guide you towards something more useful. Cheers, Charles On Fri, Apr 4, 2014 at 8:37 PM, thanoon younis thanoon.youni...@gmail.comwrote: thanks alot for your help now i want two different sample size in R what should i change in previous command? and how can i get correlated simulation data (there are an interrelationships between variables) regards thanoon On 4 April 2014 18:42, Charles Determan Jr deter...@umn.edu wrote: Hi Thanoon, How about this? # replicate p=10 times random sampling n=1000 from a vector containing your ordinal categories (1,2,3,4) R - replicate(10, sample(as.vector(seq(4)), 1000, replace = T)) Cheers, Charles On Fri, Apr 4, 2014 at 7:10 AM, thanoon younis thanoon.youni...@gmail.com wrote: dear sir i want to simulate multivariate ordinal data matrix with categories (1,4) and n=1000 and p=10. thanks alot thanoon [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Charles Determan Integrated Biosciences PhD Candidate University of Minnesota -- Charles Determan Integrated Biosciences PhD Candidate University of Minnesota [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation data
Hi Thanoon, How about this? # replicate p=10 times random sampling n=1000 from a vector containing your ordinal categories (1,2,3,4) R - replicate(10, sample(as.vector(seq(4)), 1000, replace = T)) Cheers, Charles On Fri, Apr 4, 2014 at 7:10 AM, thanoon younis thanoon.youni...@gmail.comwrote: dear sir i want to simulate multivariate ordinal data matrix with categories (1,4) and n=1000 and p=10. thanks alot thanoon [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Charles Determan Integrated Biosciences PhD Candidate University of Minnesota [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation of Hierarchical design
So can the same student be associated with multiple instructors? only within the same school? or more general? by repeated student ID's do you mean that a student in school A and a student in school B can both have ID 1, but they are different students? or do you mean that instructor #1 and Instructor #2 in the same school can both have student ID 1 and that it is the same student? What do you want the effect of state, school, and instructor to be on your outcome variable? A couple of the examples for the simfun function in the TeachingDemos package do simulations for hierarchical designs (the last 3 examples), in this case simulating height for subjects nested in cities nested in states assuming that each have an effect on height. You could modify one of those to match your story (use a binomial instead of normal, or use normal and a cut-off value). On Mon, Feb 24, 2014 at 11:58 AM, farnoosh sheikhi farnoosh...@yahoo.com wrote: Hi, I want to simulate a data set with following condition: There are 6 states with 12 schools. Each two schools are coming from one states. For example school one and two from state A, school 3 and 4 from state B and etc. Each school has 10 unique instructors with random number of students meaning that there are repeated IDs( student ID). For each school, there is an indicator telling if the school has a diet program or not. There is also an indicator telling if student has been on diet or not. (possible response). There is also a variable stating the state name. state A, state B, etc. I want to create a data set in R based on this story. I really appreciate any help and direction. Regards, Farnoosh Sheikhi [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Gregory (Greg) L. Snow Ph.D. 538...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation in R
Hello, See the help page for ?sample. X - sample(0:1, 1, replace = TRUE, prob = c(0.25, 0.75)) Hope this helps, Rui Barradas Em 04-08-2013 08:51, Preetam Pal escreveu: Hi All, I want to simulate a random variable X which takes values 1 and 0 with probabilities 75% and 25% respectively and then repeat the procedure 1 times. I am sure this is trivial, I tried to look at the help pages online, but I can't quite find it. Appreciate your help. Thanks and Regards, Preetam [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation in R
So you looked at some unspecified help pages online and tried some unspecified stuff? Try being more specific next time you post. For example, try reading the footer of any R-help email. Note that it says read the Posting Guide, and provide a reproducible example (at least of what you tried that didn't work). Also note that you need to use plain text email for your reproducible examples to get through the email system undamaged. After reading the Posting Guide, the following hint may be of use: ?sample, prob argument. --- Jeff NewmillerThe . . Go Live... DCN:jdnew...@dcn.davis.ca.usBasics: ##.#. ##.#. Live Go... Live: OO#.. Dead: OO#.. Playing Research Engineer (Solar/BatteriesO.O#. #.O#. with /Software/Embedded Controllers) .OO#. .OO#. rocks...1k --- Sent from my phone. Please excuse my brevity. Preetam Pal lordpree...@gmail.com wrote: Hi All, I want to simulate a random variable X which takes values 1 and 0 with probabilities 75% and 25% respectively and then repeat the procedure 1 times. I am sure this is trivial, I tried to look at the help pages online, but I can't quite find it. Appreciate your help. Thanks and Regards, Preetam [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation in R
Thank you very much guys On Sun, Aug 4, 2013 at 3:51 PM, Prof Brian Ripley rip...@stats.ox.ac.ukwrote: On 04/08/2013 09:30, Rui Barradas wrote: Hello, See the help page for ?sample. X - sample(0:1, 1, replace = TRUE, prob = c(0.25, 0.75)) Hope this helps, ?rbinom would have been a better answer since simpler, faster, algorithms are available in that case. Or even as.integer(runif(1) 0.75) Rui Barradas Em 04-08-2013 08:51, Preetam Pal escreveu: Hi All, I want to simulate a random variable X which takes values 1 and 0 with probabilities 75% and 25% respectively and then repeat the procedure 1 times. I am sure this is trivial, I tried to look at the help pages online, but I can't quite find it. Appreciate your help. Thanks and Regards, Preetam -- Brian D. Ripley, rip...@stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~**ripley/http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 -- Preetam Pal (+91)-9432212774 M-Stat 2nd Year, Room No. N-114 Statistics Division, C.V.Raman Hall Indian Statistical Institute, B.H.O.S. Kolkata. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation in R
On 04/08/2013 09:30, Rui Barradas wrote: Hello, See the help page for ?sample. X - sample(0:1, 1, replace = TRUE, prob = c(0.25, 0.75)) Hope this helps, ?rbinom would have been a better answer since simpler, faster, algorithms are available in that case. Or even as.integer(runif(1) 0.75) Rui Barradas Em 04-08-2013 08:51, Preetam Pal escreveu: Hi All, I want to simulate a random variable X which takes values 1 and 0 with probabilities 75% and 25% respectively and then repeat the procedure 1 times. I am sure this is trivial, I tried to look at the help pages online, but I can't quite find it. Appreciate your help. Thanks and Regards, Preetam -- Brian D. Ripley, rip...@stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation from truncated skew normal
I am not aware of any such command so, I think, you may have to write one yourself: invert the CDF and use uniform random variable (runif) to sample Mikhail On Tuesday, June 11, 2013 16:18:59 cassie jones wrote: Hello R-users, I am trying to simulate from truncated skew normal distribution. I know there are ways to simulate from skewed normal distribution such as rsn(sn) or rsnorm(VGAM), but I could not find any command to simulate from a truncated skew-normal distribution. Does anyone know how to do that? Any help is appreciated. Thanks in advance. Cassie [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation\bootstrap of list factors
Tobias, I'm not sure if this is what you're after, but perhaps it will help. # create a list of 5 vectors n - 5 subsets - lapply(1:n, function(x) rnorm(5, mean=80, sd=1)) # create another list that takes 2 bootstrap samples from each of the 5 vectors and puts them in a matrix nbootstrap - 2 test - lapply(subsets, function(x) matrix(sample(x, nbootstrap*length(x), replace=TRUE), ncol=nbootstrap)) subsets test Jean On Wed, Apr 17, 2013 at 1:09 PM, Berg, Tobias van den to.vandenb...@vumc.nl wrote: Dear R experts, I am trying to simulate a list containing data matrices. Unfortunately, I don't manage to get it to work. A small example: n=5 nbootstrap=2 subsets-list() for (i in 1:n){ subsets[[i]] - rnorm(5, mean=80, sd=1) for (j in 1:nbootstrap){ test-list() test[[j]]-subsets[[i]] } } How can I get test to be 2 simulation rounds with each 5 matrices. Kind regards, Tobias [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation
On 03-01-2013, at 17:40, Simone Gogna singletontheb...@msn.com wrote: Dear R users, suppose we have a random walk such as: v_t+1 = v_t + e_t+1 where e_t is a normal IID noise pocess with mean = m and standard deviation = sd and v_t is the fundamental value of a stock. Now suppose I want a trading strategy to be: x_t+1 = c(v_t – p_t) where c is a costant. I know, from the paper where this equations come from (Farmer and Joshi, The price dynamics of common trading strategies, 2001) that the induced price dynamics is: r_t+1 = –a*r_t + a*e_t + theta_t+1 and p_t+1 = p_t +r_t+1 where r_t = p_t – p_t-1 , e_t = v_t – v_t-1 and a = c/lambda (lambda is another constant). How can I simulate the equations I have just presented? I have good confidence with R for statistical analysis, but not for simulation therefore I apologize for my ignorance. What I came up with is the following: ##general settings c-0.5 lambda-0.3 a-c/lambda n-500 ## Eq.12 (the v_t random walk) V_init_cond-0 Et-ts(rnorm(n+100,mean=0,sd=1)) Vt-Et*0 Vt[1]-V_init_cond+Et[1] for(i in 2:(n+100)) { Vt[i]-Vt[i-1]+Et[i] } Vt-ts(Vt[(length(Vt)-n+1):length(Vt)]) plot(Vt) ## Eq.13 (the strategy) Xt_init_cond-0 Xt-Xt_init_cond*0 Xt[2]-c(Vt[1]-Pt[1]) for(i in 2:(n)){ Xt[i]-c(Vt[i-1]-Pt[i-1]) } Xt-ts(Xt[(length(Xt)-n+1):length(Xt)]) plot(Xt) ## Eq. 14 (pice dynamics) P_init_cond-0 Pt-Rt*0 Pt[1]-P_init_cond+Rt[1] for(i in 2:(n+100)) { Pt[i]-Pt[i-1]+Rt[i] } Pt-ts(Pt[(length(Pt)-n+1):length(Pt)]) plot(Pt) Rt_init_cond-0 Rt-Rt_init_cond*0 Rt[2]- -a*Rt[1]+a*Et[1]+e[2] for(i in 2:(n)){ Rt[i]- -a*Rt[i-1]+a*Et[i-1]+e[i] } Rt-ts(Rt[(length(Rt)-n+1):length(Rt)]) plot(Rt) I don’t think the code above is correct, and I don’t even know if this is the approach I have to take. Any suggestion is warmly appreciated. Do not use c as a user variable. It is an R provided function. You have a formulae such as Xt[2]-c(Vt[1]-Pt[1]) for(i in 2:(n)){ Xt[i]-c(Vt[i-1]-Pt[i-1]) c is not doing here what you want. I assume you meant to multiply as in Xt[2]-c*(Vt[1]-Pt[1]) for(i in 2:(n)){ Xt[i]-c*(Vt[i-1]-Pt[i-1]) So call this constant cpar or something similar. Where has e been defined? If you reorder your equations in such a way that all initial conditions are computed first in the correct order then you simulation loops could be condensed into a single loop such as for(i in 2:(n+100)) { Vt[i] - Vt[i-1]+Et[i] Rt[i] - -a*Rt[i-1]+a*Et[i-1]+e[i] Pt[i] - Pt[i-1]+Rt[i] Xt[i] - cpar*(Vt[i-1]-Pt[i-1]) } If I am correct. Berend __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation
On 03-01-2013, at 17:40, Simone Gogna singletontheb...@msn.com wrote: Dear R users, suppose we have a random walk such as: v_t+1 = v_t + e_t+1 where e_t is a normal IID noise pocess with mean = m and standard deviation = sd and v_t is the fundamental value of a stock. Now suppose I want a trading strategy to be: x_t+1 = c(v_t – p_t) where c is a costant. I know, from the paper where this equations come from (Farmer and Joshi, The price dynamics of common trading strategies, 2001) that the induced price dynamics is: r_t+1 = –a*r_t + a*e_t + theta_t+1 and p_t+1 = p_t +r_t+1 where r_t = p_t – p_t-1 , e_t = v_t – v_t-1 and a = c/lambda (lambda is another constant). How can I simulate the equations I have just presented? I have good confidence with R for statistical analysis, but not for simulation therefore I apologize for my ignorance. What I came up with is the following: ##general settings c-0.5 lambda-0.3 a-c/lambda n-500 ## Eq.12 (the v_t random walk) V_init_cond-0 Et-ts(rnorm(n+100,mean=0,sd=1)) Vt-Et*0 Vt[1]-V_init_cond+Et[1] for(i in 2:(n+100)) { Vt[i]-Vt[i-1]+Et[i] } Vt-ts(Vt[(length(Vt)-n+1):length(Vt)]) plot(Vt) ## Eq.13 (the strategy) Xt_init_cond-0 Xt-Xt_init_cond*0 Xt[2]-c(Vt[1]-Pt[1]) for(i in 2:(n)){ Xt[i]-c(Vt[i-1]-Pt[i-1]) } Xt-ts(Xt[(length(Xt)-n+1):length(Xt)]) plot(Xt) ## Eq. 14 (pice dynamics) P_init_cond-0 Pt-Rt*0 Pt[1]-P_init_cond+Rt[1] for(i in 2:(n+100)) { Pt[i]-Pt[i-1]+Rt[i] } Pt-ts(Pt[(length(Pt)-n+1):length(Pt)]) plot(Pt) Rt_init_cond-0 Rt-Rt_init_cond*0 Rt[2]- -a*Rt[1]+a*Et[1]+e[2] for(i in 2:(n)){ Rt[i]- -a*Rt[i-1]+a*Et[i-1]+e[i] } Rt-ts(Rt[(length(Rt)-n+1):length(Rt)]) plot(Rt) I don’t think the code above is correct, and I don’t even know if this is the approach I have to take. Any suggestion is warmly appreciated. You should also have a look at package simecol which can also do discrete time models. It would certainly require some study of the manual and a bit of work on your part but I think it would be worth it. Berend __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation in R
look at functions replicate and mvrnorm functions (the later in the MASS package). On Sat, Dec 1, 2012 at 12:02 PM, mboricgs mbori...@hotmail.com wrote: Hello! How can I do 100 simulations of length 17 from bivariate bivariate normal distribution, if I know all 5 parameters? -- View this message in context: http://r.789695.n4.nabble.com/Simulation-in-R-tp4651578.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Gregory (Greg) L. Snow Ph.D. 538...@gmail.com [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation with cpm package
On 13.11.2012 15:45, Christopher Desjardins wrote: Hi, I am running the following code based on the cpm vignette's code. I believe the code is syntactically correct but it just seems to hang R. I can get this to run if I set the sims to 100 but with 2000 it just hangs. Any ideas why? No: Works for me and completes within 90 minutes. Uwe Ligges Thanks, Chris library(cpm) cpmTypes - c(Kolmogorov-Smirnov,Mann-Whitney,Cramer-von-Mises) changeMagnitudes - c(1, 2, 4, 5) changeLocations - c(50,100,300) sims - 2000 ARL0 - 500 startup - 20 results - list() for (cpmType in cpmTypes) { results[[cpmType]] - matrix(numeric(length(changeMagnitudes) * length(changeLocations)), nrow = length(changeMagnitudes)) for (cm in 1:length(changeMagnitudes)) { for (cl in 1:length(changeLocations)) { print(sprintf(cpm:%s magnitude::%s location:%s, cpmType, changeMagnitudes[cm], changeLocations[cl])) temp - numeric(sims) for (s in 1:sims) { x -c(rchisq(changeLocations[cl], df=3), rchisq(2000, df=changeMagnitudes[cm])) temp[s] -detectChangePoint(x, cpmType, ARL0=ARL0, startup=startup)$detectionTime } results[[cpmType]][cm,cl] - mean(temp[temp changeLocations[cl]]) - changeLocations[cl] } } } [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation of levene's test
Dear Dila, Try the following: library(Rcmdr) asim - 1000 pv-NULL for(i in 1:asim) { print(i) set.seed(i) g1 - rnorm(20,0,2) g2 - rnorm(20,0,2) g3 - rnorm(20,0,2) x - c(g1,g2,g3) group-as.factor(c(rep(1,20),rep(2,20),rep(3,20))) pv-c(pv,leveneTest(x,group)$Pr(F)[1]) } Best Ozgur - Ozgur ASAR Research Assistant Middle East Technical University Department of Statistics 06531, Ankara Turkey Ph: 90-312-2105309 http://www.stat.metu.edu.tr/people/assistants/ozgur/ -- View this message in context: http://r.789695.n4.nabble.com/simulation-of-levene-s-test-tp4631578p4631600.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation of levene's test
On Mon, May 28, 2012 at 12:14 PM, Özgür Asar oa...@metu.edu.tr wrote: Dear Dila, Try the following: library(Rcmdr) Or avoid the unncessary overhead of Rcmdr and use library(car) to provide levenTest instead. asim - 1000 pv-NULL It's also many orders of magnitude more efficient to preallocate pv and then simply put things into it. pv - vector(real, 1000) for(i in 1:asim) { print(i) set.seed(i) Setting the seed each loop seems excessive but I suppose it's a matter of taste. g1 - rnorm(20,0,2) g2 - rnorm(20,0,2) g3 - rnorm(20,0,2) x - c(g1,g2,g3) Is there any reason not to do this as x - rnorm(60, 0, 2) group-as.factor(c(rep(1,20),rep(2,20),rep(3,20))) and this as as.factor(rep(1:3, each = 20)) pv-c(pv,leveneTest(x,group)$Pr(F)[1]) Once you preallocate pv change this to pv[i] - leveneTest(x, group)$Pr(F)[1] But it's even better not to use the dollar sign shortcut here (defensive programming and all that -- particularly with nonstandard names which I'm pretty sure won't give a big error here but will elsewhere) pv[i] - leveneTest(x, group)[[Pr(F)]][1] And even better would be to do this all using the replicate function, but I'll leave that as an exercise to the reader. Michael } Best Ozgur - Ozgur ASAR Research Assistant Middle East Technical University Department of Statistics 06531, Ankara Turkey Ph: 90-312-2105309 http://www.stat.metu.edu.tr/people/assistants/ozgur/ -- View this message in context: http://r.789695.n4.nabble.com/simulation-of-levene-s-test-tp4631578p4631600.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation
On Apr 5, 2012, at 10:57 PM, Christopher Kelvin wrote: Hello, i need to simulate 100 times, n=40 , the distribution has 90% from X~N(0,1) + 10% from X~N(20,10) Is my loop below correct? Thank you n=40 for(i in 1:100){ x-rnorm(40,0,1) # 90% of n You are overwriting x and y and at the end of that loop you will only have two vectors of length 40 each. If you wanted a 90 10 weighting then why not lengths of 36 and 4??? To do repeated simulations you will find this help page useful: ?replicate z-rnorm(40,20,10) # 10% of n } x+z At this point you should not be using + but rather the c() function if you are trying to join those two vectors. I think you need to spend more time working through Introduction to R. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. David Winsemius, MD West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation
Suggestions? -- Yes. 1) Wrong list.. Post on R-sig-mixed-models, not here. 2) Follow the posting guide and provide the modelformula, which may well be the source of the difficulties (overfitting). -- Bert On Fri, Dec 16, 2011 at 1:56 PM, Scott Raynaud scott.rayn...@yahoo.com wrote: I'm using an R program (which I did not write) to simulate multilevel data (subjects in locations) used in power calculations. It uses lmer to fit a mixed logistic model to the simulated data based on inputs of means, variances, slopes and proportions: (fitmodel - lmer(modelformula,data,family=binomial(link=logit),nAGQ=1)) where modelformula is set up in another part of the program. Locations are treated as random and the model is random intercept only. The program is set to run 1000 simulations. I have temperature, five levels of gestational age (GA), birth wieght (BW) and four other categorical pedictors, all binary. I scaled everything so that all my slopes are in the range of -5.2 to 1.6 and variances from .01 to .08. I have a couple of categories of GA that have small probabilities (.10). I'm using a structured sampling approach looking at 20, 60, 100, and 140 locations with a total n=75k. The first looks like this: # groups n 5 800 4 2239 3 3678 3 5117 3 6557 2 7996 Total 20 75000 As the level 2 sizes increase, the cell sizes decrease. When I run this model in the simulation I get: Warning: glm.fit: algorithm did not converge every time the model is fit (I killed this long before it ran 1000 times). I tried increasing the number of iterations to no avail. I suspected linear dependencies among the predictors, so I took out GA (same result), put GA back and took out BW (same result) and then took out both GA and BW. This ran about half the time with th other half passing warnings such as: Warning: glm.fit: algorithm did not converge Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred or Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred in addition to some like the original warning. If I leave everything in but temperature, then it runs fine. I also tested the full model separately at 50 and 75 level 2 units each with total n=75k. Nothing converged. I want to include temperature, but I'm not sure what else to try. Any suggestions? __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Bert Gunter Genentech Nonclinical Biostatistics Internal Contact Info: Phone: 467-7374 Website: http://pharmadevelopment.roche.com/index/pdb/pdb-functional-groups/pdb-biostatistics/pdb-ncb-home.htm __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation over data repeatedly for four loops
Perhaps you might want to abstract your code a bit and try something like: X = rnorm(500) # Some Data replicate(1e4, mean(sample(X, 500, replace = T))) Obviously you can set up a loop over your data sets as needed. Michael On Sat, Nov 12, 2011 at 6:46 PM, Francesca francesca.panco...@gmail.com wrote: Dear Contributors, I am trying to perform a simulation over sample data, but I need to reproduce the same simulation over 4 groups of data. My ability with for loop is null, in particular related to dimensions as I always get, no matter what I try, number of items to replace is not a multiple of replacement length This is what I intend to do: replicate this operation for four times, where the index for the four groups is in the part of the code: datiPc[[1]][,2]. I have to replicate the following code 4 times, where the changing part is in the data from which I pick the sample, the data that are stored in datiPc[[1]][,2]. If I had to use data for the four samples, I would substitute the 1 with a j and replicate a loop four times, but it never worked. My desired final outcome is a matrix with 1 observations for each couple of extracted samples, i.e. 8 columns of 1 observations of means. db-c() # Estrazione dei campioni dai dati di PGG e TRUST estr1 - c(); estr2 - c(); m1-c() m2-c() tmp1- data1[[1]][,2]; tmp2- data2[[2]][,2]; for(i in 1:100){ estr1-sample(tmp1, 1000, replace = TRUE) estr2-sample(tmp2, 1000, replace = TRUE) m1[i]-mean(estr1,na.rm=TRUE) m2[i]-mean(estr2,na.rm=TRUE) } db-data.frame(cbind(m1,m2)) Thanks for any help you can provide. Best Regards -- Francesca -- [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation from discrete uniform
If you wanted a discrete uniform from 1-10 use: ceiling(10*runif(1)) if you wanted from 0-12, use: ceiling(13*runif(1))-1 -- View this message in context: http://r.789695.n4.nabble.com/Simulation-from-discrete-uniform-tp3434980p3939694.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation from discrete uniform
Why don't you use sample; sample(1:10,10,replace=TRUE) -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of BSanders Sent: 26 October 2011 08:49 To: r-help@r-project.org Subject: Re: [R] Simulation from discrete uniform If you wanted a discrete uniform from 1-10 use: ceiling(10*runif(1)) if you wanted from 0-12, use: ceiling(13*runif(1))-1 -- View this message in context: http://r.789695.n4.nabble.com/Simulation-from-discrete-uniform-tp3434980 p3939694.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. LEGAL NOTICE This message is intended for the use o...{{dropped:10}} __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation for equation with two variables
That is, run all possible combinations of the two vectors through the equation. *Ben * On Thu, Jul 21, 2011 at 10:04 AM, Benjamin Caldwell btcaldw...@berkeley.edu wrote: Hi, I'm trying to run a basic simulation and sensitivity test by running an equation with two variables and then plotting the results against each of the vectors. R is running the vectors like this : 0 with 0, 1 with 1, etc. I would like it to run them like 0 for 1:100, 1 for 1:100, and then the reverse (so 100^2*100@ iterations. How do I get it to do that? Here's what I have so far: par(mfrow=c(1,2)) t - 0:100 DBH - 0:100 M - 4000*(1-exp(-t*(1.104-(0.67*0.7)-0.163*log(DBH))^2))) plot(M~DBH) plot(M~t) *Ben* [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation for equation with two variables
On Jul 21, 2011, at 1:42 PM, Benjamin Caldwell wrote: That is, run all possible combinations of the two vectors through the equation. *Ben * For all combinations the usual route is data preparation with either expand.grid() or outer() On Thu, Jul 21, 2011 at 10:04 AM, Benjamin Caldwell btcaldw...@berkeley.edu wrote: Hi, I'm trying to run a basic simulation and sensitivity test by running an equation with two variables and then plotting the results against each of the vectors. R is running the vectors like this : 0 with 0, 1 with 1, etc. I would like it to run them like 0 for 1:100, 1 for 1:100, and then the reverse (so 100^2*100@ iterations. How do I get it to do that? Here's what I have so far: par(mfrow=c(1,2)) t - 0:100 DBH - 0:100 M - 4000*(1-exp(-t*(1.104-(0.67*0.7)-0.163*log(DBH))^2))) plot(M~DBH) plot(M~t) *Ben* David Winsemius, MD West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation
On Mon, Jun 06, 2011 at 04:50:57PM +1000, Stat Consult wrote: Dear ALL I want to simulate data from Multivariate normal distribution. GE.N-mvrnorm(25,mu,S) S -matrix(rep(0,1),nrow=100) for( i in 1:100){sigma-runif(100,0.1,10);S [i,i]=sigma[i];mu-runif(100,0,10)} for (i in 1:20){for (j in 1:20){if (i != j){S [i,j]=0.3*sigma[i]*sigma[j]}}} for (i in 21:40){for (j in 21:40){if (i != j){S [i,j]=0.3*sigma[i]*sigma[j]}}} for (i in 41:60){for (j in 41:60){if (i != j){S [i,j]=0.3*sigma[i]*sigma[j]}}} for (i in 61:80){for (j in 61:80){if (i != j){S [i,j]=0.3*sigma[i]*sigma[j]}}} for (i in 81:100){for (j in 81:100){if (i != j){S [i,j]=0.3*sigma[i]*sigma[j]}}} How should I do when S is not positive definite matrix? I saw this error: 'Sigma' is not positive definite. Hello. I am not sure, how the matrix is created. Should the command sigma-runif(100,0.1,10) be indeed inside the loop over i? I suspect that no, since otherwise, only the vector sigma used for S[100, 100] goes to the remaining part of the construction. The matrix is block diagonal. So, the corresponding distribution can be build from parts corresponding to the blocks generated independently. Let me look at the first block assuming that sigma is generated only once. The first block may be obtained also as B - 0.3*outer(sigma[1:20], sigma[1:20]) diag(B) - sigma[1:20] The result may have negative eigenvalues. For example, if all components in sigma[1:20] are 4, which is in the range used for sigma, then we have a matrix, whose diagonal elements are 4 and nondiagonal elements are 0.3*4^2 = 4.8 4. This matrix has negative eigenvalues, so it is not a covariance matrix. Is the construction of the matrix, which you sent, correct? Petr Savicky. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation from discrete uniform
Hi Also , same problem to create discrete uniform Distribution , But sample () and runif() not useful to generate discrete uniform . Ex: u-round(runif(10*10,min=1,max=10),0) table(u) u 1 2 3 4 5 6 7 8 9 10 6 10 9 10 14 6 11 14 12 8 Not useful for large number OR # for generate large number dus - sample(0:9, 100, replace = FALSE) Error in base::sample(x, size, replace = replace, prob = prob, ...) : cannot take a sample larger than the population when 'replace = FALSE' DO you have any suggestion my question ? Need to generate 1000*10 number from 1:250 with discrete uniform distribution ? Regards, Serdar [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation from discrete uniform
Hi Serdar, Take a look at the following: sample(0:9, 100, replace = FALSE) Error in sample(0:9, 100, replace = FALSE) : cannot take a sample larger than the population when 'replace = FALSE' sample(0:9, 100, replace = TRUE) [1] 5 6 5 7 3 0 8 4 8 2 2 4 7 6 0 7 0 0 0 7 5 6 3 6 0 9 6 1 2 6 9 0 0 4 7 9 8 6 4 7 0 [42] 4 6 1 8 2 5 6 3 6 5 1 7 6 0 9 5 5 3 6 3 8 7 5 1 2 3 6 6 9 3 6 5 6 2 5 9 3 6 5 0 7 [83] 8 0 8 7 3 9 9 1 4 4 1 1 0 9 8 1 9 3 HTH, Jorge On Sun, May 29, 2011 at 6:38 PM, SERDAR NESLIHANOGLU wrote: Hi Also , same problem to create discrete uniform Distribution , But sample () and runif() not useful to generate discrete uniform . Ex: u-round(runif(10*10,min=1,max=10),0) table(u) u 1 2 3 4 5 6 7 8 9 10 6 10 9 10 14 6 11 14 12 8 Not useful for large number OR # for generate large number dus - sample(0:9, 100, replace = FALSE) Error in base::sample(x, size, replace = replace, prob = prob, ...) : cannot take a sample larger than the population when 'replace = FALSE' DO you have any suggestion my question ? Need to generate 1000*10 number from 1:250 with discrete uniform distribution ? Regards, Serdar [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation from truncated poisson
It is truncated from left. On Mon, May 16, 2011 at 6:33 PM, Greg Snow greg.s...@imail.org wrote: Which direction is it truncated? (only values less than a allowed or only greater?). One simple approach is rejection sampling, just generate from a regular poisson distribution, then throw away any values in the truncated region. Another approach if the legal values are those from 0 to a, so that there is a finite number of possibilities, then you can use the sample function with replace=TRUE and using probabilities from the poisson in the legal range. -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of cassie jones Sent: Monday, May 16, 2011 5:28 PM To: r-help@r-project.org Subject: [R] simulation from truncated poisson Dear all, I need to simulate values from a Poisson distribution which is truncated at certain value 'a'. Can anyone tell me if there is in-built package in R which can simulate from a truncated Poisson? If not, what should be the steps to write a function which would do that? Thanks in advance. Cassie [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation from truncated poisson
Which direction is it truncated? (only values less than a allowed or only greater?). One simple approach is rejection sampling, just generate from a regular poisson distribution, then throw away any values in the truncated region. Another approach if the legal values are those from 0 to a, so that there is a finite number of possibilities, then you can use the sample function with replace=TRUE and using probabilities from the poisson in the legal range. -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of cassie jones Sent: Monday, May 16, 2011 5:28 PM To: r-help@r-project.org Subject: [R] simulation from truncated poisson Dear all, I need to simulate values from a Poisson distribution which is truncated at certain value 'a'. Can anyone tell me if there is in-built package in R which can simulate from a truncated Poisson? If not, what should be the steps to write a function which would do that? Thanks in advance. Cassie [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation from truncated poisson
So there is no maximum value, just a minimum (greater than 0)? Correct? In that case the sample option will not work (unless you choose some really high value and say that you won't go above that), but the rejection sampling would still work. How efficiently it works will depend on how much probability a regular poisson would put into the truncated region. Another possibility is to find the probability of being in the truncated region, then generate a uniform between that value and 1, then feed that uniform into the qpois function. From: cassie jones [mailto:cassiejone...@gmail.com] Sent: Monday, May 16, 2011 7:46 PM To: Greg Snow Cc: r-help@r-project.org Subject: Re: [R] simulation from truncated poisson It is truncated from left. On Mon, May 16, 2011 at 6:33 PM, Greg Snow greg.s...@imail.orgmailto:greg.s...@imail.org wrote: Which direction is it truncated? (only values less than a allowed or only greater?). One simple approach is rejection sampling, just generate from a regular poisson distribution, then throw away any values in the truncated region. Another approach if the legal values are those from 0 to a, so that there is a finite number of possibilities, then you can use the sample function with replace=TRUE and using probabilities from the poisson in the legal range. -Original Message- From: r-help-boun...@r-project.orgmailto:r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.orgmailto:r-help-boun...@r-project.org] On Behalf Of cassie jones Sent: Monday, May 16, 2011 5:28 PM To: r-help@r-project.orgmailto:r-help@r-project.org Subject: [R] simulation from truncated poisson Dear all, I need to simulate values from a Poisson distribution which is truncated at certain value 'a'. Can anyone tell me if there is in-built package in R which can simulate from a truncated Poisson? If not, what should be the steps to write a function which would do that? Thanks in advance. Cassie [[alternative HTML version deleted]] __ R-help@r-project.orgmailto:R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation from truncated poisson
Thanks for the help. I appreciate. On Mon, May 16, 2011 at 10:19 PM, Greg Snow greg.s...@imail.org wrote: So there is no maximum value, just a minimum (greater than 0)? Correct? In that case the sample option will not work (unless you choose some really high value and say that you wont go above that), but the rejection sampling would still work. How efficiently it works will depend on how much probability a regular poisson would put into the truncated region. Another possibility is to find the probability of being in the truncated region, then generate a uniform between that value and 1, then feed that uniform into the qpois function. *From:* cassie jones [mailto:cassiejone...@gmail.com] *Sent:* Monday, May 16, 2011 7:46 PM *To:* Greg Snow *Cc:* r-help@r-project.org *Subject:* Re: [R] simulation from truncated poisson It is truncated from left. On Mon, May 16, 2011 at 6:33 PM, Greg Snow greg.s...@imail.org wrote: Which direction is it truncated? (only values less than a allowed or only greater?). One simple approach is rejection sampling, just generate from a regular poisson distribution, then throw away any values in the truncated region. Another approach if the legal values are those from 0 to a, so that there is a finite number of possibilities, then you can use the sample function with replace=TRUE and using probabilities from the poisson in the legal range. -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of cassie jones Sent: Monday, May 16, 2011 5:28 PM To: r-help@r-project.org Subject: [R] simulation from truncated poisson Dear all, I need to simulate values from a Poisson distribution which is truncated at certain value 'a'. Can anyone tell me if there is in-built package in R which can simulate from a truncated Poisson? If not, what should be the steps to write a function which would do that? Thanks in advance. Cassie [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation Questions
Hi Shane, it sounds to me as though you have a fairly well-defined problem. You want to generate random numbers with a specific mean, variance, and correlation with another random varaible. I would reverse-enginerr the fuinctions for simple linear regression to get a result like y = beta_0 + beta_1 * x + rnorm(n, 0, sigma^2) and use that as the basis of generating random numbers. Not sure how to interpret the second question ... Cheers Andrew On Sun, May 01, 2011 at 12:33:41AM -0400, Shane Phillips wrote: I have the following script for generating a dataset. It works like a champ except for a couple of things. 1. I need the variables itbs and map to be negatively correlated with the binomial variable lunch (around -0.21 and -0.24, respectively). The binomial variable lunch needs to remain unchanged. 2. While my generated variables do come out with the desired means and correlations, the distribution is very narrow and only represents a small portion of the possible scores. Can I force it to encompass a wider range of scores, while maintaining my desired parameters and correlations? Please help... Shane Script follows... #Number the subjects subject=1:1000 #Assign a treatment condition from a binomial distribution with a probability of 0.13 treat=rbinom(1*1000,1,.13) #Assign a lunch status condition froma binomial distribution with a probability of 0.35 lunch=rbinom(1*1000,1,.35) #Generate age in months from a random normal distribution with mean of 87 and sd of 2 age=rnorm(1000,87,2) #invoke the MASS package require(MASS) #Establish the covariance matrix for MAP, ITBS and CogAT scores sigma - matrix(c(1, 0.84, 0.59, 0.84, 1, 0.56, 0.59, 0.56, 1), ncol = 3) #Establish MAP as a random normal variable with mean of 200 and sd of 9 map - rnorm(1000, 200, 9) #Establish ITBS as a random normal variable with mean of 175 and sd of 15 itbs - rnorm(1000, 175, 15) #Establish CogAT as a random normal variable with mean of 100 and sd of 16 cogat-rnorm(1000,100,16) #Create a dataframe of MAP, ITBS, and CogAT data - data.frame(map, itbs, cogat) #Draw from the multivariate distribution defined by MAP, ITBS, and CogAT means and the covariance matrix sim - mvrnorm(1000, mu=mean(data), sigma, empirical=FALSE) #Set growth at 0 growth=0 #Combine elements into a single dataset simtest=data.frame (subject=subject, treat=treat,lunch, age=round(age,0),round(sim,0),growth) #Set mean growth by treatment condition with treatd subjects having a mean growth of 1.5 and non-treated having a mean growth of 0.1 simtest-transform(simtest, growth=rnorm(1000,m=ifelse(treat==0,0.1,1.5),s=1)) simtest cor (simtest) __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Andrew Robinson Program Manager, ACERA Department of Mathematics and StatisticsTel: +61-3-8344-6410 University of Melbourne, VIC 3010 Australia (prefer email) http://www.ms.unimelb.edu.au/~andrewpr Fax: +61-3-8344-4599 http://www.acera.unimelb.edu.au/ Forest Analytics with R (Springer, 2011) http://www.ms.unimelb.edu.au/FAwR/ Introduction to Scientific Programming and Simulation using R (CRC, 2009): http://www.ms.unimelb.edu.au/spuRs/ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation from discrete uniform
?sample -Oprindelig meddelelse- Fra: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] På vegne af cassie jones Sendt: 8. april 2011 03:16 Til: r-help@r-project.org Emne: [R] Simulation from discrete uniform Dear all, I am trying to simulate from discrete uniform distribution. But I could not find any in-built code in R. Could anyone help me please? Thanks in advance for the time and help. Cassie [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation from discrete uniform
Hi: Suppose X has a discrete uniform distribution on the sample space S = {0, 1, 2, ..., 9}. Then a random sample of size 100 from this distribution, for example, would be dus - sample(0:9, 100, replace = TRUE) # Checks: table(dus) lattice::barchart( ~ table(dus), xlim = c(0, 20)) The sample space comprises the first argument of sample(), the sample size is the second argument, and the replace = TRUE argument allows an element of the sample space to be selected more than once. HTH, Dennis On Thu, Apr 7, 2011 at 6:15 PM, cassie jones cassiejone...@gmail.comwrote: Dear all, I am trying to simulate from discrete uniform distribution. But I could not find any in-built code in R. Could anyone help me please? Thanks in advance for the time and help. Cassie [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Well, knowing how your data looks like would definitely help! Say your data object is called mydata, just paste the output from dput(mydata) into the email you want to send to the list. Ivan Le 3/1/2011 04:18, bwaxxlo a écrit : I tried looking for help but I couldn't locate the exact solution. I have data that has several variables. I want to do several sample simulations using only two of the variables (eg: say you have data between people and properties owned. You only want to check how many in the samples will come up with bicycles) to estimate probabilities and that sort of thing. Now, I can only do a simulation in terms of this code: sample(1:10, size = 15, replace = TRUE). I do not know how select specific variables only. I'll appreciate the help -- Ivan CALANDRA PhD Student University of Hamburg Biozentrum Grindel und Zoologisches Museum Abt. Säugetiere Martin-Luther-King-Platz 3 D-20146 Hamburg, GERMANY +49(0)40 42838 6231 ivan.calan...@uni-hamburg.de ** http://www.for771.uni-bonn.de http://webapp5.rrz.uni-hamburg.de/mammals/eng/1525_8_1.php __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Date: Mon, 28 Feb 2011 19:18:18 -0800 From: kadodamb...@hotmail.com To: r-help@r-project.org Subject: [R] Simulation I tried looking for help but I couldn't locate the exact solution. I have data that has several variables. I want to do several sample simulations using only two of the variables (eg: say you have data between people and properties owned. You only want to check how many in the samples will come up with bicycles) to estimate probabilities and that sort of thing. Now, I can only do a simulation in terms of this code: sample(1:10, size = 15, replace = TRUE). I do not know how select specific variables only. I'll appreciate the help This is probably not the best R but you can do something like either of these. Note that this is just the easiest derivative of stuff I already had and can be fixed to your needs, I usually use runif instead of sample for example. The first example probably being much less efficient than the second, df-data.frame(a=.1*rnorm(100), b=(1:100)/100,c=(1:100)/100+.1*rnorm(100)) res=1:100; for ( i in 1:100) {res[i]=cor(df[which(runif(100).9),])[1,3] } res hist(res) res=1:100; for ( i in 1:100) {wh=which(runif(100).9); res[i]=cor(df$a[wh],df$c[wh]); } res -- View this message in context: http://r.789695.n4.nabble.com/Simulation-tp3329173p3329173.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation of Multivariate Fractional Gaussian Noise and Fractional Brownian Motion
Dear Kjetil, Thank you so much for your advice on my question. Best Regards, Wonsang 2011/2/10 Kjetil Halvorsen kjetilbrinchmannhalvor...@gmail.com What you can do to find out is to type into your R session RSiteSearch(multivariate fractional gaussian) That seems to give some usefull results. Kjetil On Tue, Feb 8, 2011 at 1:51 PM, Wonsang You y...@ifn-magdeburg.de wrote: Dear R Helpers, I have searched for any R package or code for simulating multivariate fractional Brownian motion (mFBM) or multivariate fractional Gaussian noise (mFGN) when a covariance matrix are given. Unfortunately, I could not find such a package or code. Can you suggest any solution for multivariate FBM and FGN simulation? Thank you for your help. Best Regards, Ryan - Wonsang You Leibniz Institute for Neurobiology -- View this message in context: http://r.789695.n4.nabble.com/Simulation-of-Multivariate-Fractional-Gaussian-Noise-and-Fractional-Brownian-Motion-tp3276296p3276296.html Sent from the R help mailing list archive at Nabble.com. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation of Multivariate Fractional Gaussian Noise and Fractional Brownian Motion
What you can do to find out is to type into your R session RSiteSearch(multivariate fractional gaussian) That seems to give some usefull results. Kjetil On Tue, Feb 8, 2011 at 1:51 PM, Wonsang You y...@ifn-magdeburg.de wrote: Dear R Helpers, I have searched for any R package or code for simulating multivariate fractional Brownian motion (mFBM) or multivariate fractional Gaussian noise (mFGN) when a covariance matrix are given. Unfortunately, I could not find such a package or code. Can you suggest any solution for multivariate FBM and FGN simulation? Thank you for your help. Best Regards, Ryan - Wonsang You Leibniz Institute for Neurobiology -- View this message in context: http://r.789695.n4.nabble.com/Simulation-of-Multivariate-Fractional-Gaussian-Noise-and-Fractional-Brownian-Motion-tp3276296p3276296.html Sent from the R help mailing list archive at Nabble.com. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation - Natrual Selection
Date: Wed, 5 Jan 2011 15:48:46 + From: benjamin.w...@bathspa.org To: r-help@r-project.org Subject: [R] Simulation - Natrual Selection Hi, I've been modelling some data over the past few days, of my work, repeatedly challenging microbes to a certain concentration of cleaner, until the required concentration to inhibit or kill them increaces, at which point they are challenged to a slightly higher concentration each day. I'm doing ths for two different cleaners and I'm collecting the required concentration to kill them as a percentage, the challenge number, the cleaner as a two level variable, and the lineage theyre in, because I have several different lineages. I'm expecting the values to rise for one cleaner but not the other as they aqquire resistance for one but not the other. Which has happened, but I have wide variation because one linage aqquired a very dramatic change which has made it immune to 50%, whereas the others, have exhibited a much more gradual increace, and so I have very weak p values for the cleaner variable, because it is secondary to the challenge vector, which has the most explanatory power, because without time and these challenges, the selection would no happen. I was using two bacterium species, but one was keen on giving hight erratic results, and insisted on becoming cross contaminated, BUT if I include it's data, It shoves cleaner over the p0.05 threshold, so i may just be having a problem with lack of data. So I've been asking about bootstrapping, which I plan to do to my cases, and thenfit a model to see what the confidence is like then. I assume if I bootstrap then it will re-select whole cases, and not jumble everything up, otherwise a microbe (totake the most extreme value as an example) with 50% concentration tolerance at the beginning, would make no sense at all. I'm also planning on doing models lineage by lineage, rather than putting them into one whole, just to have a look at what happens. You can't really have a p-value without a specific hypothesis to test, if you have that then all your other questions are probably easy to answer. Generally you want to sample from things that are iid or maybe you want to test the identical i. Generally you want to have done a lit search ahead of time and had some idea of likely evolution dynamics of your system given your design and things like your forcing functions etc. Most statisticians would not take seriously a posteriori designs and indeed it can be hard to avoid rationalization and selection bias ( problems that always and only effect people who disagree with me LOL) as being anything other than exploratory or hypothesis generating- you are looking for predictive value. While it is not always worthwhile doing blind tests, it may be something worth considering ( do you know which group gets what thing?) But what I really wanted to know from this email, was if there's a package or function for natrual selection simulation I could make use of, to see if I can simulate the experiment. I want to start with a http://www.google.com/#sclient=psyhl=enq=%22R+package%22+natural+selection but as implied above, R has lots of analysis stuff and maybe you would find something more useful that is not linked to the keywords you suggest. You may find, for whatever reason, you could write a differential equation to express your results but that isn't often used with natural selection. distribution of concentration tolerance values, taken from th e inhibitory concentration values from my first lot of microbes, back when term began. Draw 3000 from this. Then values in that draw that fall below the exposure concentration I did in my experiment, are removed, or have a high chance of being removed. Then, from what is left, a draw is made again - or perhaps a copy operation (rather than a random draw) until I have 3000 again, rather than have all exactly the same concentration, then a value can be added to some of them, that increaces their concentration tolerance slightly, but not by a great deal, except in a few individuals, where it may be increaced dramatically(some sort of exponential dstribution perhaps). Then when the distribution of this simulated population of microbes has reached the next concentration (possibly the mean or mode of the distribution) (I have a series of 1 in 2 dilutions, so 100% 50%, 25% and so on), then they move on to the next concentration. I know it's probably quite a heavy thing, it was just a thought that came to me, if anybody has any experience in this area of R or knows of something that allows this to be done, please let me know. Thanks, Ben. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation - Natrual Selection
On 05/01/2011 16:37, Mike Marchywka wrote: Date: Wed, 5 Jan 2011 15:48:46 + From: benjamin.w...@bathspa.org To: r-help@r-project.org Subject: [R] Simulation - Natrual Selection Hi, I've been modelling some data over the past few days, of my work, repeatedly challenging microbes to a certain concentration of cleaner, until the required concentration to inhibit or kill them increaces, at which point they are challenged to a slightly higher concentration each day. I'm doing ths for two different cleaners and I'm collecting the required concentration to kill them as a percentage, the challenge number, the cleaner as a two level variable, and the lineage theyre in, because I have several different lineages. I'm expecting the values to rise for one cleaner but not the other as they aqquire resistance for one but not the other. Which has happened, but I have wide variation because one linage aqquired a very dramatic change which has made it immune to 50%, whereas the others, have exhibited a much more gradual increace, and so I have very weak p values for the cleaner variable, because it is secondary to the challenge vector, which has the most explanatory power, because without time and these challenges, the selection would no happen. I was using two bacterium species, but one was keen on giving hight erratic results, and insisted on becoming cross contaminated, BUT if I include it's data, It shoves cleaner over the p0.05 threshold, so i may just be having a problem with lack of data. So I've been asking about bootstrapping, which I plan to do to my cases, and thenfit a model to see what the confidence is like then. I assume if I bootstrap then it will re-select whole cases, and not jumble everything up, otherwise a microbe (totake the most extreme value as an example) with 50% concentration tolerance at the beginning, would make no sense at all. I'm also planning on doing models lineage by lineage, rather than putting them into one whole, just to have a look at what happens. You can't really have a p-value without a specific hypothesis to test, if you have that then all your other questions are probably easy to answer. Generally you want to sample from things that are iid or maybe you want to test the identical i. My Hypothesis is that Cleaner A (I don't really want to go into names or brands), will exhbit a rise in concentration tolerance values, or rather, the microbial culture I keep exposed to it, will, reflecting aqquisition of antimicrobial resistance. And this has largely happened. And that in cleaner B, this will not happen, or if it does, it will not be as dramatic and take longer. So I expecting in my model, the cleaner variable to have a p below 0.05, and quite hight explanatory power, and a satisfying coefficient. The notion behind the hypothesis being that one might have a more difficult complex chemical structure, requiring more mutations to develop some resistance. I can't really do anything with genes or chemical structure at my current institution and at my level because of no equippment for that sort of thing, and that they felt it would be too far for a 3rd year project. So I'm using the concentration required to kill them - or stop them from growing, as a indication. Generally you want to have done a lit search ahead of time and had some idea of likely evolution dynamics of your system given your design and things like your forcing functions etc. Most statisticians would not take seriously a posteriori designs and indeed it can be hard to avoid rationalization and selection bias ( problems that always and only effect people who disagree with me LOL) as being anything other than exploratory or hypothesis generating- you are looking for predictive value. While it is not always worthwhile doing blind tests, it may be something worth considering ( do you know which group gets what thing?) But what I really wanted to know from this email, was if there's a package or function for natrual selection simulation I could make use of, to see if I can simulate the experiment. I want to start with a http://www.google.com/#sclient=psyhl=enq=%22R+package%22+natural+selection but as implied above, R has lots of analysis stuff and maybe you would find something more useful that is not linked to the keywords you suggest. You may find, for whatever reason, you could write a differential equation to express your results but that isn't often used with natural selection. distribution of concentration tolerance values, taken from th e inhibitory concentration values from my first lot of microbes, back when term began. Draw 3000 from this. Then values in that draw that fall below the exposure concentration I did in my experiment, are removed, or have a high chance of being removed. Then, from what is left, a draw is made again - or perhaps a copy operation (rather than a random draw) until I have 3000 again, rather than have all exactly the same concentration, then a value
Re: [R] Simulation - Natrual Selection
On 05/01/2011 17:40, Bert Gunter wrote: My hypothesis was specified before I did my experiment. Whilst far from perfect, I've tried to do the best I can to assess rise in resistance, without going into genetics as it's not possible. (Although may be at the next institution I've applied for MSc). With my hypothesis (I mentioned it below), I was of the frame of mind that a nonsignificant p-value on the cleaner variable (for now - experiment is far from over), indicated a lack of evidence for rejecting the null. And so at the minute, it looks like the type of cleaner makes no difference. I have no fundamental objection, but be careful. I would simply qualify your last sentence by saying that it means that the experimental noise is to great to precisely determine the size of the cleaner effect. Scientific reality tells us that it is never exactly 0; what your results show is that your uncertainty about the value of the difference encompasses both positive and negative values. This does NOT mean that the difference might not be scientifically large enough to be of interest -- a confidence interval for the difference (MUCH better than a P value) would help you determine that. If the interval is narrow enough that the difference, positive or negative, is too small to be of scientific interest, then you're done. If the linterval is large, then it tells you that you need more data, a better experiment (less noisy) etc. -- Bert At the moment I wouldn't call the confidence interval small, it's definately wide, and at the minute the confidence interval covers zero. My R-squared at the minite is also 0.5, this is mostly due to the few extreme cases of adaptation as I mentioned before, but I'm hesitant to remove it as papers in my literature study which also evolve bacteria show that there is often (sometimes wide) variation in the paths populations take. So whilst mathematically a bit undesirable, and makes me and the model uncertain, it does fall into place with what is known, or has been previously shown of the reality of selection. Again if I include the data from the bacteria dropped from the study, all that improves, and uncertainty is reduced. It may also be worth me mentioning, I am also taking a more traditional approach (by that I mean a more Statistics 101 approach, indeed that is all the stats tuition covered in my course as a taught element), incase what I've described above did not work or was not ideal, because we (me and my supervisor) did forsee a model report may contain a lot of uncertainty. Indeed we did expect some populations to adapt and some to not etc. So I've also been collecting data on the width of the zones of inhibition shown by putting disks of cleaner on plates of growth, and measuring the dead zone that results. I can get lots of data from this with only a few plates, and doing this at the start of the study, a few times in the middle, and at the end. Will allow me to do more traditional analysis, for example t.test on the dead zone widths at the end of the study, between cleaner a and b. Or a non parametric equivalent, maybe even a permutation test. The modelling stuff is already beyond what my supervisor expects of me, but I felt it would add value and a lot more insight to the study, allowing more variables to be accounted for, than a more short-sighted traditional test. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation from pareto distn
For a Pareto distribution, even a truncated one, the inverse CDF method should be straightforward to implement. Giovanni Petris On Tue, 2010-11-09 at 10:50 -0600, cassie jones wrote: Dear all, I am trying to simulate from truncated Pareto distribution. I know there is a package called PtProcess for Pareto distribution...but it is not for truncated one. Can anyone please help me with this? Thanks in advance. Cassie [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Giovanni Petris gpet...@uark.edu Associate Professor Department of Mathematical Sciences University of Arkansas - Fayetteville, AR 72701 Ph: (479) 575-6324, 575-8630 (fax) http://definetti.uark.edu/~gpetris/ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation from pareto distn
Hi: library(sos) findFn('truncated Pareto') On my system, it scared up 17 matches. It looks like the VGAM package would be a reasonable place to start looking. HTH, Dennis On Tue, Nov 9, 2010 at 8:50 AM, cassie jones cassiejone...@gmail.comwrote: Dear all, I am trying to simulate from truncated Pareto distribution. I know there is a package called PtProcess for Pareto distribution...but it is not for truncated one. Can anyone please help me with this? Thanks in advance. Cassie [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
I have two questions: (1) How do you 'create' an 2 x 2 table in R using say an Odd ratio of 3 or even 0.5 (2) If I have several 2 x 2 tables, how can I 'implement' dependence in the tables with say 25 of the Tables having an odds ratio of 1 and 75 of the tables having an odds ratio of 4? Jim [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Do you need a table with an odds ratio exactly equal to 3 (or other value), or a realistic sample from a population with odds ratio 3 where the sample table will have a different OR (but the various tables will cluster around the true value)? -- Gregory (Greg) L. Snow Ph.D. Statistical Data Center Intermountain Healthcare greg.s...@imail.org 801.408.8111 -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r- project.org] On Behalf Of Jim Silverton Sent: Friday, September 10, 2010 1:03 AM To: r-help@r-project.org Subject: Re: [R] Simulation I have two questions: (1) How do you 'create' an 2 x 2 table in R using say an Odd ratio of 3 or even 0.5 (2) If I have several 2 x 2 tables, how can I 'implement' dependence in the tables with say 25 of the Tables having an odds ratio of 1 and 75 of the tables having an odds ratio of 4? Jim [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting- guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation problem.
This looks like homework. If it is, you should really tell us along with what your teacher's policy is on getting help over the internet is (and note that many teachers monitor this list and can see if you are getting help). You have done the first part yourself, much better than some who have tried to get us to do the whole thing for them, so a possible hint: the new problem really has 3 groups, never sick, currently sick, and healed. Just modify your current code to allow for people to move from the currently sick to the healed group. -- Gregory (Greg) L. Snow Ph.D. Statistical Data Center Intermountain Healthcare greg.s...@imail.org 801.408.8111 -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r- project.org] On Behalf Of Piotr Arendarski Sent: Tuesday, April 13, 2010 4:27 PM To: r-help@r-project.org Subject: [R] Simulation problem. Hi, I have problem with simulating. This is my task... Suppose that there are N persons some of whom are sick with influenza. The following assumptions are made: * when a sick person meets a healthy one, the chance is á that the latter will be infected * all encounters are between two persons Write a function which simulates this model for various values of N (say, 10 000) and á (say, between 0.001 and 0.1). Monitor the history of this process, assuming that one individual is infected at the beginning. The code is: * simulation - function(number, prob){ cumulative.time - 0 current.time - 0 number.sick - 1 while(number.sicknumber){ current.time - current.time + 1 meetings - rhyper(nn=1, m=number.sick, n=number-number.sick, k=2) if(meetings==1){ one.sick - rbinom(n=1, size=1, prob) if(one.sick==1){ cumulative.time - c(cumulative.time, current.time) number.sick - number.sick + 1 }}} cumulative.time } number - 1000 prob - .05 model - simulate(number, prob)* But than add the assumption that *each infected person has a 0.01 chance of recovering at each time unit* Do you have idea how to modify the code ? Piotr Arendarski [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation of VAR
Dear Ron, have you had a look at the package dse? Here, ARMA models can be specified and simulated. The only exercise left for you, is to transform the VECM coefficients into their level-VAR values. Best, Bernhard | -Original Message- | From: r-help-boun...@r-project.org | [mailto:r-help-boun...@r-project.org] On Behalf Of Ron_M | Sent: Saturday, March 27, 2010 12:14 PM | To: r-help@r-project.org | Subject: [R] Simulation of VAR | | | Dear all, is there any package/function available which simulates a | co-integrating VAR model once the model parameters are | input over some | arbitrary horizon? Please let me know anyone aware of that. | | Thanks | -- | View this message in context: | http://n4.nabble.com/Simulation-of-VAR-tp1693295p1693295.html | Sent from the R help mailing list archive at Nabble.com. | | __ | R-help@r-project.org mailing list | https://stat.ethz.ch/mailman/listinfo/r-help | PLEASE do read the posting guide | http://www.R-project.org/posting-guide.html | and provide commented, minimal, self-contained, reproducible code. | * Confidentiality Note: The information contained in this ...{{dropped:10}} __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation of VAR
Yes I looked into dse package. Here I have implemented two approach for simulation like following : library(dse) A1 - matrix(rnorm(16),4) A2 - matrix(rnorm(16),4) mu - rnorm(4) sigma - matrix(c(0.006594712, 0.006467731, -0.000254914, 0.005939934, 0.006467731, 0.006654184, -0.000384097, 0.005658247, -0.000254914, -0.000384097, 0.000310294, 4.34141E-05, 0.005939934, 0.005658247, 4.34141E-05, 0.00574024), 4) initial.val - matrix(c(-0.2347096, -0.1803612, -0.2780356, -0.2154427 , 3.740364, 3.757908, 3.50216 , 3.57783), 2) # My approach res - matrix(NA, 4,4); res[c(1,2),] - initial.val library(mnormt); shocks - rmnorm(2, rep(0,4), sigma) for (i in 1:2) { res[i+2,] - mu + A1%*%res[i+2-1,] + A2%*%res[i+2-2,] + shocks[i,] } res # dse approach temp1 - matrix(t(cbind(diag(4), A1, A2)), ncol = 4, byrow = TRUE) model - ARMA(A=array(temp1, c(3,4,4)), B=diag(4), TREND=mu) simulate(model, y0=initial.val, sampleT=2, noise=shocks) Ideally last two rows of res and simulate() should be exactly same. However that is not what I am getting. Can anyone please tell me whether there is any mistale in any of those approaches? Am I missing somthing? Thanks -- View this message in context: http://n4.nabble.com/Simulation-of-VAR-tp1693295p1694899.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation of binary data
Check out the help page of the lrm function in the rms library. To show how lrm is used, the examples simulate data for logistic regression. This may give you some ideas. On Wed, Jan 20, 2010 at 10:41 AM, omar kairan omarkaira...@gmail.com wrote: Hi, could someone help me with dilemma on the simulation of logistic regressiondata with multicollinearity effect and high leverage point.. Thank you [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] simulation of binary data
On 21/01/2010, at 4:41 AM, omar kairan wrote: Hi, could someone help me with dilemma on the simulation of logistic regressiondata with multicollinearity effect and high leverage point.. If that is the clearest way in which you can phrase your question then I doubt that *anyone* can help you. cheers, Rolf Turner ## Attention:\ This e-mail message is privileged and confid...{{dropped:9}} __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation numbers from a probability table
If the trials are not connected then I would consider melting the table using melt() from the reshape package. And then using lapply() with the function random.function - function(my.prob, number.of.observations = 10) { sum(rbinom(number.of.observations, 1, my.prob)) } in case the trials are connected, by column, than you could use apply(the.data.table, 2, a.function) on it. Where a.function will to multinum distribution (for which I don't remember the function at the moment, but it can be searched). Best, Tal. Contact Details:--- Contact me: tal.gal...@gmail.com | 972-52-7275845 Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il (Hebrew) | www.r-statistics.com/ (English) -- On Wed, Jan 13, 2010 at 7:20 PM, Kelvin 6kelv...@gmail.com wrote: Dear friends, If I have a table like this, first row A B C D ... are different levels of the variable, first column 0 1 2 4 ... are the levels of the numbers, the numbers inside the table are the probabilities of the number occuring. A B C D... 0 0.20.30.10.05 1 0.10.10.20.2 2 0.02 0.20 0.1 4 0.30.01 0.01 0.4 ... How can I use R to do the simulation and get a table like this, first row A B C D ... are different levels of the variable, the numbers inside the table are the numbers simulated from the probailties table above? A B C D ... 0 4 2 0 2 2 0 1 0 1 4 1 2 2 0 0 ... Thanks for help! Kelvin __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation numbers from a probability table
Try this: dat - data.frame(x=11:14, pa=1:4/10, pb=4:1/10) f - function(numreps, data){ pmat - as.matrix(data[-1]) x - data[,1] result - matrix(0, nrow=numreps, ncol=ncol(pmat)) colnames(result) - c(A, B) for(i in seq_len(numreps)){ result[i,] - apply(pmat, 2, function(p) sample(x, 1, prob=p)) } result } f(5, dat) -Peter Ehlers Kelvin wrote: Dear friends, If I have a table like this, first row A B C D ... are different levels of the variable, first column 0 1 2 4 ... are the levels of the numbers, the numbers inside the table are the probabilities of the number occuring. A B C D... 0 0.20.30.10.05 1 0.10.10.20.2 2 0.02 0.20 0.1 4 0.30.01 0.01 0.4 ... How can I use R to do the simulation and get a table like this, first row A B C D ... are different levels of the variable, the numbers inside the table are the numbers simulated from the probailties table above? A B C D ... 0 4 2 0 2 2 0 1 0 1 4 1 2 2 0 0 ... Thanks for help! Kelvin __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Peter Ehlers University of Calgary 403.202.3921 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation Function - Save results
On Aug 17, 2009, at 1:40 PM, MarcioRibeiro wrote: Ok, the LIST function I understood... I didn't see how you got a random input to that function. Would seem to need one of the rdist functions as input. What I would like now is to simulate this Function A many times (S) in order to get S results for the MEAN and for the VARIANCE... ?replicate Zhiliang Ma wrote: in order to return more multiple variables, you can put them in a list and then return this list. e.g. #Function A boot-function(a,b,c){ mean_boot-(a+b)/2 var_boot-c list(mean_boot = mean_boot, var_boot = var_boot) } out - boot(1,2,3) out $mean_boot [1] 1.5 $var_boot [1] 3 On Fri, Aug 14, 2009 at 1:15 PM, MarcioRibeiromes...@pop.com.br wrote: Hi listers, I am working on a simulation... But I am having some troubles... Suppose I have a function A which produces two results (mean and variance)... Then I would like to simulate this function A with a function B many times using the results from function A For example: #Function A boot-function(a,b,c){ mean_boot-(a+b)/2 var_boot-c #list(a=a,b=b,c=c) return(a) } Then I would like to create 2 vectors with the mean and var results from S simulations #Function B simul-function(S){ teste-rep(0,S) for(i in 1:S){ teste[i]-boot(10,12,15) #ACCORDING TO FUNCTION A I AM SAVING JUST THE MEAN_BOOT, BUT I ALSO NEED THE RESULT OF VAR_BOOT } var-var(teste) mean_emp-mean(var_boot) #THIS IS NOT WORKING, BECAUSE I DONT HAVE THE VAR_BOOT AT MY VECTOR var_emp-(sum((var_boot-var)**2))/S #THIS IS NOT WORKING, BECAUSE I DONT HAVE THE VAR_BOOT AT MY VECTOR } simul(5) But my problem is that I don't know how to save my results in 2 vectors in order to use then at function B. Thanks in advance, Marcio -- View this message in context: http://www.nabble.com/Simulation-Function---Save-results-tp24977851p24977851.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- View this message in context: http://www.nabble.com/Simulation-Function---Save-results-tp24977851p25011101.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. David Winsemius, MD Heritage Laboratories West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation function
You need to store your results in a list and then access the information in the list to get the values: boot-function(a,b,c){ + media-(a+b+c)/3 + var_app-var(a) + list(media=media,var_app=var_app) + } boot(2,4,10) $media [1] 5.33 $var_app [1] NA simul-function(S){ + results-list() + for(i in 1:S){ + results[[i]]-boot(2,4,10) + } + var-var(sapply(results, '[[', 'media')) + mean_var-mean(sapply(results, '[[', 'var_app')) + var_var-var(sapply(results, '[[', 'var_app')) + list(var=var,mean_var=mean_var,var_var=var_var) + } simul(5) $var [1] 0 $mean_var [1] NA $var_var [1] NA On Tue, Aug 18, 2009 at 11:55 AM, MarcioRibeiromes...@pop.com.br wrote: Hi listers, I've been looking for a procedure, but I am not succeding... I have a function that give multiple results... Then, I would like to simulate this function n times, so I need to save/keep the n multiple results, in order to calculate my desired statistics... I have already tried with the RETURN and LIST FUNCTION, but I am not getting it right... An example of what I am looking for sounds like this... boot-function(a,b,c){ media-(a+b+c)/3 var_app-var(a) list(media=media,var_app=var_app) } boot(2,4,10) simul-function(S){ results-rep(0,S) for(i in 1:S){ results[i]-boot(2,4,10) } var-var(media) mean_var-mean(var_app) var_var-var(var_app) list(var=var,mean_var=mean_var,var_var=var_var) } simul(5) -- View this message in context: http://www.nabble.com/Simulation-function-tp25027993p25027993.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Jim Holtman Cincinnati, OH +1 513 646 9390 What is the problem that you are trying to solve? __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation Function - Save results
in order to return more multiple variables, you can put them in a list and then return this list. e.g. #Function A boot-function(a,b,c){ mean_boot-(a+b)/2 var_boot-c list(mean_boot = mean_boot, var_boot = var_boot) } out - boot(1,2,3) out $mean_boot [1] 1.5 $var_boot [1] 3 On Fri, Aug 14, 2009 at 1:15 PM, MarcioRibeiromes...@pop.com.br wrote: Hi listers, I am working on a simulation... But I am having some troubles... Suppose I have a function A which produces two results (mean and variance)... Then I would like to simulate this function A with a function B many times using the results from function A For example: #Function A boot-function(a,b,c){ mean_boot-(a+b)/2 var_boot-c #list(a=a,b=b,c=c) return(a) } Then I would like to create 2 vectors with the mean and var results from S simulations #Function B simul-function(S){ teste-rep(0,S) for(i in 1:S){ teste[i]-boot(10,12,15) #ACCORDING TO FUNCTION A I AM SAVING JUST THE MEAN_BOOT, BUT I ALSO NEED THE RESULT OF VAR_BOOT } var-var(teste) mean_emp-mean(var_boot) #THIS IS NOT WORKING, BECAUSE I DONT HAVE THE VAR_BOOT AT MY VECTOR var_emp-(sum((var_boot-var)**2))/S #THIS IS NOT WORKING, BECAUSE I DONT HAVE THE VAR_BOOT AT MY VECTOR } simul(5) But my problem is that I don't know how to save my results in 2 vectors in order to use then at function B. Thanks in advance, Marcio -- View this message in context: http://www.nabble.com/Simulation-Function---Save-results-tp24977851p24977851.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation functions for underdispered Poisson and binomial distributions
On Wed, 15 Jul 2009, Shinichi Nakagawa wrote: Dear R users I would like to simulate underdispersed Poisson and binomial distributions somehow. I know you can do this for overdispersed counterparts - using rnbinom() for Poisson and rbetabinom() for binomial. Could anyone share functions to do this? Or please share some tips for modifying existing functions to achieve this. Shinichi, You really need a model for the underdispersion. Using that model, you would calculate the probabiltities fo the binomial or Poisson counts. But you have to come up with something appropriate for your situation. For example, probs - prop.table( dbinom( 0:10, 10, .5)^3 ) or probs - prop.table( dbinom( 0:10, 10, .5) + ifelse( 0:10 == 5, 1, 0) ) will each produce probabilities for counts that are less dispersed than probs - dbinom( 0:10, 10, 0.5 ) but neither may suitably model the counts for the situation in which you are interested. --- Once you have that model in hand sample( 0:10, N, pr=probs, repl=TRUE ) will 'simulate' N such counts. HTH, Chuck Thank you very much for your help and time Shinichi Shinichi Nakagawa, PhD (Lecturer of Behavioural Ecology) Department of Zoology University of Otago 340 Great King Street P. O. Box 56 Dunedin, New Zealand Tel: +64-3-479-5046 Fax: +64-3-479-7584 http://www.otago.ac.nz/zoology/staff/academic/nakagawa.html __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Charles C. Berry(858) 534-2098 Dept of Family/Preventive Medicine E mailto:cbe...@tajo.ucsd.edu UC San Diego http://famprevmed.ucsd.edu/faculty/cberry/ La Jolla, San Diego 92093-0901 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation code error
This message is one of a number of identical copies, also cross-posted to R-sig-geo, and in fact with an obvious solution that the author gives at the end. The thread will be continued on R-sig-geo. Cross posting is advised against expressly in e.g. http://en.wikipedia.org/wiki/Crossposting, and leads to threads being left dangling if solutions do not propagate to all lists where threads were initiated. Roger Bivand Steve Hong wrote: Dear List, I have some problem with my simulation code. Here is output from R: sim.sp - function(data,CM,n,N) + { + C - matrix(rep(NA,N),ncol=1) + for(i in 1:N) + { + j - n + xx - which(colSums(CM[j,])==1) + V - names(xx) + V - paste(V, collapse=+) + V - paste(SBA~, V) + rd - round(nrow(data)*(2/3)) + d - sample(seq(1:nrow(data)),rd) + dat1 - data[d,] + dat2 - data[-d,] + crd - cbind(dat1$Longitude,dat1$Latitude) + dist80 - dnearneigh(crd,0,100,longlat=F) + dist80sw - nb2listw(dist80, style=B) + fm - errorsarlm(as.formula(V), data=dat1, listw=dist60sw) + pred - predict(fm,dat2) + C[i,1] - cor(dat2$SBA,pred) + out - cbind(C) + } + colMeans(out) + } sim.sp(df2007.5k.s2,CM,1,1000) Error in nb2listw(dist80, style = B) : Empty neighbour sets found I guess I got error message since there are some observations without neighborhoods in the process of determining neighborhood structure randomly. Is there any way to ignore neighborhood structures with observations without neighborhoods and proceed the code? Thank you in advance!! Steve Hong [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- View this message in context: http://www.nabble.com/Simulation-code-error-tp24506052p24506545.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Hi peter, Quite an insight you have there hehe. i am continuing on from the orignal problem of creating a simulation. Im now trying to find (nâ1)S2/Ï2, and fit it to a chi squared dist with 5 degrees of freedom. im having trouble with the coding for this. i think for the second part of that i need to use the fitdist function, but to get it to where i am able to do that, im not sure what to do. THis is what i have been trying to do so far, but it hasn't returned me anything good sum((x-mean(x))^2)/(length(x)-1) i am really confused, can someone please help? Cheers Date: Thu, 14 May 2009 12:05:30 +0100 From: b.rowling...@lancaster.ac.uk To: peterflomconsult...@mindspring.com CC: waclaw.marcin.kusnierc...@idi.ntnu.no; r-help@r-project.org Subject: Re: [R] Simulation As a beginner, I agree the for loop is much clearer to me. [Warning: Contains mostly philosophy] To me, the world and how I interact with it is procedural. When I want to break six eggs I do 'get six eggs, repeat break egg until all eggs broken'. I don't apply an instance of the break egg function over a range of eggs. My world is not functional (just like me, some might say...). Neither do I send a 'break yourself' message to each egg - my world is not object-oriented. That does not mean that these paradigms are not good ways of writing computer programs - they are brilliant ways of writing computer programs. But they build on procedural concepts, and we don't teach children to run before they can walk. So when someone says 'how do I do this a thousand times?' on R-help, I'll assume their knowledge level is that of a beginner, and try to map the solution to their world view. Computer scientists will write their beautiful manuscripts, but how many people who come to R because they want to do a t-test or fit a GLM will read them? That's the R-help audience now. Barry __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. _ [[elided Hotmail spam]] [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation from a multivariate normal distribution
Check out the help page for replicate(). Andy From: barbara.r...@uniroma1.it I must to create an array with dimensions 120x8x500. Better I have to make 500 simulations of 8 series of return from a multivariate normal distribution. there's the command mvrnorm but how I can do this repeating the simulation 500 times? [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Notice: This e-mail message, together with any attachme...{{dropped:12}} __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation from a multivariate normal distribution
barbara.r...@uniroma1.it wrote: I must to create an array with dimensions 120x8x500. Better I have to make 500 simulations of 8 series of return from a multivariate normal distribution. there's the command mvrnorm but how I can do this repeating the simulation 500 times? ?replicate Uwe Ligges [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation from a multivariate normal distribution
Liaw, Andy wrote: Check out the help page for replicate(). Andy Or the 'n' argument to mvrnorm (or mvtnorm::rmvnorm for that matter)... From: barbara.r...@uniroma1.it I must to create an array with dimensions 120x8x500. Better I have to make 500 simulations of 8 series of return from a multivariate normal distribution. there's the command mvrnorm but how I can do this repeating the simulation 500 times? [[alternative HTML version deleted]] -- O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - (p.dalga...@biostat.ku.dk) FAX: (+45) 35327907 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
hey guys, i've been following this discussion about the simulation, and being a beginner myself, im really unsure of the best method. I hve the same problem as the initial one, except i need 1000 samples of size 15, and my distribution is Exp(1). I've adjusted some of the loop formulas for my n=15, but im unsure how to proceed in the quickest way. Can someone please help? Much appreciated :) From: r.tur...@auckland.ac.nz Date: Thu, 14 May 2009 10:26:38 +1200 To: c...@witthoft.com CC: r-help@r-project.org Subject: Re: [R] Simulation On 14/05/2009, at 10:04 AM, Carl Witthoft wrote: So far nobody seems to have warned the OP about seeding. Presumably Debbie wants 1000 different sets of samples, but as we all know there are ways to get the same sequence (initial seed) every time. If there's a starting seed for one of the generate a single giant matrix methods proposed, the whole matrix will be the same for a given seed. If rnorm is called 1000 times (hopefully w/ different random (oops) seeds), the entire set of samples will be different. and so on. I really don't get this. The OP wanted 1000 independent samples, each of size 100. Whether she does set.seed(42) M - matrix(rnorm(100*1000),nrow=1000) # Each row is a sample. or L - list() set.seed(42) for(i in 1:1000) L[[i]] - rnorm(100) # Each list entry is a sample. she gets this, i.e. the desired result. Setting a seed serves to make the results reproducible. This works via either approach. Making results reproducible in this manner is advisable, but seed-setting is nothing that the OP needs to be *warned* about. cheers, Rolf Turner ## Attention:\ This e-mail message is privileged and confid...{{dropped:9}} __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. _ [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
On Fri, 15 May 2009 19:17:37 +1000 Kon Knafelman konk2...@hotmail.com wrote: KK I hve the same problem as the initial one, except i need 1000 KK samples of size 15, and my distribution is Exp(1). I've adjusted KK some of the loop formulas for my n=15, but im unsure how to proceed KK in the quickest way. KK Can someone please help? What exactly do you want? Please be more specific about what you did and what does not work. Stefan __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
On Fri, 15 May 2009 19:17:37 +1000 Kon Knafelman konk2...@hotmail.com wrote: KK I hve the same problem as the initial one, except i need 1000 KK samples of size 15, and my distribution is Exp(1). I've adjusted KK some of the loop formulas for my n=15, but im unsure how to proceed KK in the quickest way. KK Can someone please help? Taking a guess: matrix(rexp(15000,1),ncol=15) ? -- View this message in context: http://www.nabble.com/Simulation-tp23556274p23558953.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Another possibility (maybe more readable, gives the option of a list, probably not faster): Replicate(1000, rexp(15,1) ) -- Gregory (Greg) L. Snow Ph.D. Statistical Data Center Intermountain Healthcare greg.s...@imail.org 801.408.8111 -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r- project.org] On Behalf Of Ben Bolker Sent: Friday, May 15, 2009 6:37 AM To: r-help@r-project.org Subject: Re: [R] Simulation On Fri, 15 May 2009 19:17:37 +1000 Kon Knafelman konk2...@hotmail.com wrote: KK I hve the same problem as the initial one, except i need 1000 KK samples of size 15, and my distribution is Exp(1). I've adjusted KK some of the loop formulas for my n=15, but im unsure how to proceed KK in the quickest way. KK Can someone please help? Taking a guess: matrix(rexp(15000,1),ncol=15) ? -- View this message in context: http://www.nabble.com/Simulation- tp23556274p23558953.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting- guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Greg Snow wrote: Another possibility (maybe more readable, gives the option of a list, probably not faster): Replicate(1000, rexp(15,1) ) I think that should be replicate The matrix form is quite a bit faster, but don't know if that will matter -- times below are for doing this task (1000 x 15 replicates) 1000 times ... system.time(replicate(1000,replicate(1000,rexp(15,1 user system elapsed 12.689 0.220 12.985 system.time(replicate(1000,matrix(rexp(15000,1),ncol=15))) user system elapsed 2.512 0.452 2.976 -- Ben Bolker Associate professor, Biology Dep't, Univ. of Florida bol...@ufl.edu / www.zoology.ufl.edu/bolker GPG key: www.zoology.ufl.edu/bolker/benbolker-publickey.asc signature.asc Description: OpenPGP digital signature __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
I wrote replicate but the darn e-mail program fixed it for me. I expected replicate to be a bit slower, but not by that amount. I just wanted to include replicate as a more readable version of lapply while still improving over the loop approach. -- Gregory (Greg) L. Snow Ph.D. Statistical Data Center Intermountain Healthcare greg.s...@imail.org 801.408.8111 -Original Message- From: Ben Bolker [mailto:bol...@ufl.edu] Sent: Friday, May 15, 2009 10:19 AM To: Greg Snow Cc: r-help@r-project.org Subject: Re: [R] Simulation Greg Snow wrote: Another possibility (maybe more readable, gives the option of a list, probably not faster): Replicate(1000, rexp(15,1) ) I think that should be replicate The matrix form is quite a bit faster, but don't know if that will matter -- times below are for doing this task (1000 x 15 replicates) 1000 times ... system.time(replicate(1000,replicate(1000,rexp(15,1 user system elapsed 12.689 0.220 12.985 system.time(replicate(1000,matrix(rexp(15000,1),ncol=15))) user system elapsed 2.512 0.452 2.976 -- Ben Bolker Associate professor, Biology Dep't, Univ. of Florida bol...@ufl.edu / www.zoology.ufl.edu/bolker GPG key: www.zoology.ufl.edu/bolker/benbolker-publickey.asc __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Greg Snow wrote: Another possibility (maybe more readable, gives the option of a list, probably not faster): Replicate(1000, rexp(15,1) ) provided that simplify=FALSE: is(replicate(10, rexp(15, 1))) # matrix ... is(replicate(10, rexp(15, 1), simplify=FALSE)) # list ... vQ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Wacek Kusnierczyk wrote: Barry Rowlingson wrote: On Wed, May 13, 2009 at 5:36 PM, Wacek Kusnierczyk waclaw.marcin.kusnierc...@idi.ntnu.no wrote: Barry Rowlingson wrote: Soln - for loop: z=list() for(i in 1:1000){z[[i]]=rnorm(100,0,1)} now inspect the individual bits: hist(z[[1]]) hist(z[[545]]) If that's the problem, then I suggest she reads an introduction to R... i'd suggest reading the r inferno by pat burns [1], where he deals with this sort of for-looping lists the way it deserves ;) I don't think extending a list this way is too expensive. Not like indeed, for some, but only for some, values of m and n, it can actually be half a hair faster than the matrix and the replicate approaches, proposed earlier by others: another approach to create a matrix, a bit more efficient than using matrix() but also clean for beginners IMO, is to directly assign dimensions to a vector, e.g., library(rbenchmark) n=100; m=100 benchmark(replications=100, columns=c('test', 'elapsed'), order=NULL, list={ l=list(); for (i in 1:n) l[[i]] = rnorm(m) }, liist={ l=vector('list', n); for (i in 1:n) l[[i]] = rnorm(m) }, matrix=matrix(rnorm(n*m), n, m), matrix2 = {mat - rnorm(n*m); dim(mat) - c(n, m); mat}, replicate=replicate(m, rnorm(n)) ) #test elapsed # 1 list0.25 # 2 liist0.25 # 3matrix0.22 # 4 matrix20.17 # 5 replicate0.23 n=10; m=1000 ... #test elapsed # 1 list0.17 # 2 liist0.17 # 3matrix0.20 # 4 matrix20.15 # 5 replicate0.75 n=1000; m=10 ... #test elapsed # 1 list1.37 # 2 liist0.92 # 3matrix0.21 # 4 matrix20.17 # 5 replicate0.19 Best, Dimitris n=100; m=100 library(rbenchmark) benchmark(replications=100, columns=c('test', 'elapsed'), order=NULL, list={ l=list(); for (i in 1:n) l[[i]] = rnorm(m) }, liist={ l=vector('list', n); for (i in 1:n) l[[i]] = rnorm(m) }, matrix=matrix(rnorm(n*m), n, m), replicate=replicate(m, rnorm(n)) ) #test elapsed # 1 list 0.247 # 2 liist 0.235 # 3matrix 0.172 # 4 replicate 0.247 n=10; m=1000 ... #test elapsed # 1 list 0.169 # 2 liist 0.169 # 3matrix 0.173 # 4 replicate 0.771 n=1000; m=10 ... #test elapsed # 1 list 1.428 # 2 liist 0.902 # 3matrix 0.172 # 4 replicate 0.185 but note that when the number of replicate m-samples is sufficiently large (relative to the sample size), the list approach is orders of magnitude (here, one order in the last case above) less efficient than the other approaches. (preallocation seems only partially helpful.) which approach to choose -- if efficiency is an issue -- would be a matter of the actual parameter values. doing 1000 foo=rbind(foo,bar)s to a matrix. The overhead for extending a list should really only be adding a single new pointer to the list pointer structure. The existing list data isn't copied. i don't think it's accurate, as far as i understand lists in r. while a list is, internally, holding just pointers, and extending a list does not cause the pointed-to structures to be reallocated, the list itself is backed with a fixed-length array (sensu c arrays), not as a pairlist, so it needs to be reallocated. you don't rewrite the values, but you do rewrite the pointers. on the other hand, an appropriate implementation of pairlists would indeed allow you to extend a (pair)list in O(1) time. benchmark(columns=c('test', 'elapsed'), list={ l = list(); for (i in 1:1000) l[[i]] = 0 }, pairlist={ p = pairlist(); for (i in 1:1000) p[[i]] = 0 }) # test elapsed # 1 list 0.743 # 2 pairlist 0.523 the result is, of course, dependent on the actual implementation details (dynamic arrays, pairlists with cached end pointer and size, etc.) you may want to use the profiler (i haven't), but my guess is that extending an r list does require the pointer data (though not the pointed-to data) to be rewritten -- because it is a fixed-length c array, while extending a pairlist does require list traversal (because, presumably, there is no end pointer cache; as far as i could see in a quick look at the sources, pairlists are represented with chained SEXPs, which have previous/next pointers, but no end-of-list pointers, hence the traversal). for example: n = 10 benchmark(columns=c('test', 'elapsed'), order=NULL, short.control={ l = list() }, short={ l = list(); l[[1]] = 0 }, long.control={ l = vector('list', n) }, long={ l = vector('list', n); l[[n+1]] = 0 }) #test elapsed # 1 short.control 0.001 # 2 short 0.001 # 3 long.control 0.027 # 4 long 0.138 apparently, extending a short list by one is much more efficient
Re: [R] Simulation
Barry Rowlingson wrote: On Wed, May 13, 2009 at 9:56 PM, Wacek Kusnierczyk waclaw.marcin.kusnierc...@idi.ntnu.no wrote: Barry Rowlingson wrote: n = 1000 benchmark(columns=c('test', 'elapsed'), order=NULL, 'for'={ l = list(); for (i in 1:n) l[[i]] = rnorm(i, 0, 1) }, lapply=lapply(1:n, rnorm, 0, 1) ) # test elapsed # 1for 9.855 # 2 lapply 8.923 Yes, you can probably vectorize this with lapply or something, but I prefer clarity over concision when dealing with beginners... but where's the preferred clarity in the for loop solution? Seriously? You think: lapply(1:n, rnorm, 0, 1) is 'clearer' than: x=list() for(i in 1:n){ x[[i]]=rnorm(i,0,1) } for beginners? seriously, i do; but it does depend on who those beginners are. if they come directly from c and the like, you're probably right. Firstly, using 'lapply' introduces a function (lapply) that doesn't have an intuitive name. Also, it takes a function as an argument. The concept of having a function as a parameter to another function is something that a lot of programming beginners have trouble with - unless they were brought up on LISP of course, and few of us are. well, that's one of the first things you learn on a programming languages course that is not procedural programming-{centered,biased}. no need for prior lisp experience. if messing with closures in not involved (as here), no need for advanced discussion is needed. also, the for looping may not be as trivial stuff to explain as you might think. note, you're talking about r, not c, and the treatment of iterator variables in for loops in scripting languages differs: perl -e ' $i = 0; for $i (1..5) { # iterate with $i }; print $i\n ' # 0 ruby -e ' i = 0 for i in 1..5 # iterate with i end printf %d\n, i ' # 5 and you've gotten into explaining lexical scoping etc. I propose that the for-loop example is clearer to a larger population than the lapply version. which population have you sampled from? you may not be wrong, but give some data. Plus it's only useful in that form if the first parameter is the one you want to lapply over. If you want to work over the third parameter, say, you then need: lapply(1:n,function(i){rnorm(100,0,i)}) at which point you've introduced anonymous functions. The jump from: x[[i]] = rnorm(i,0,1) to x[[i]] = rnorm(100,0,i) is much less than the changes in the lapply version, where you have to go 'oh hang on, lapply only works on the first argument, so you have to write another function, but you can do that inline like this...'. you may be unhappy to learn that you're unaware of how the lapply solution can still be nicely adapted here: lapply(1:n, rnorm, n=100, mean=0) Okay, maybe my example is a little contrived, but I still think for a beginners context it's important not to jump too many paradigms at a time. for a complete beginner, jump into for loops may not be that trivial as you seem to think. there's still quite some stuff to be explained to clarify that i = 0 for (i in 1:n) # do stuff print(i) will print n, not 0. unless n=0, of course. vQ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Barry Rowlingson wrote: [...] Yes, you can probably vectorize this with lapply or something, but I prefer clarity over concision when dealing with beginners... but where's the preferred clarity in the for loop solution? Seriously? You think: lapply(1:n, rnorm, 0, 1) is 'clearer' than: x=list() for(i in 1:n){ x[[i]]=rnorm(i,0,1) } for beginners? Firstly, using 'lapply' introduces a function (lapply) that doesn't have an intuitive name. Also, it takes a function as an argument. The concept of having a function as a parameter to another function is something that a lot of programming beginners have trouble with - unless they were brought up on LISP of course, and few of us are. I propose that the for-loop example is clearer to a larger population than the lapply version. there's one more issue. the lapply solution effectively hides the loop variable from the user, keeping h{im,er} from doing things that should not be done, e.g., assignment to the loop variable inside the loop with the assumption that it has an impact on the looping condition. for example, this simple c code int i; for (i = 0; i 5; i++) { i = 5; printf(%d\n, i); } will print *one* 5, while this simple r code (as well as analogous code in many other scripting languages) for (i in 0:4) { i = 5 print(i) } will print *five* 5s. it's not magic, but if the 'beginner' comes from c, that's a possibly ugly surprise, which needs further ado. we have recently seen on this list an example of how a user was changing looping variables within the loop body and was surprised it didn't work. with lapply, you simply can't get at that. not so easily, at least. that said, i don't think there is anything wrong with for looping vs. lapply, but the r inferno i referred to originally *is* a source of interesting comments explaining good and bad practices. vQ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Dimitris Rizopoulos wrote: Wacek Kusnierczyk wrote: Barry Rowlingson wrote: On Wed, May 13, 2009 at 5:36 PM, Wacek Kusnierczyk waclaw.marcin.kusnierc...@idi.ntnu.no wrote: Barry Rowlingson wrote: Soln - for loop: z=list() for(i in 1:1000){z[[i]]=rnorm(100,0,1)} now inspect the individual bits: hist(z[[1]]) hist(z[[545]]) If that's the problem, then I suggest she reads an introduction to R... i'd suggest reading the r inferno by pat burns [1], where he deals with this sort of for-looping lists the way it deserves ;) I don't think extending a list this way is too expensive. Not like indeed, for some, but only for some, values of m and n, it can actually be half a hair faster than the matrix and the replicate approaches, proposed earlier by others: another approach to create a matrix, a bit more efficient than using matrix() but also clean for beginners IMO, is to directly assign dimensions to a vector, e.g., library(rbenchmark) n=100; m=100 benchmark(replications=100, columns=c('test', 'elapsed'), order=NULL, list={ l=list(); for (i in 1:n) l[[i]] = rnorm(m) }, liist={ l=vector('list', n); for (i in 1:n) l[[i]] = rnorm(m) }, matrix=matrix(rnorm(n*m), n, m), matrix2 = {mat - rnorm(n*m); dim(mat) - c(n, m); mat}, replicate=replicate(m, rnorm(n)) ) sure; you could also replace 'matrix' with 'as.matrix' in the original solution, which also gives some speedup: n=100; m=100 benchmark(replications=1000, columns=c('test', 'elapsed'), order=NULL, list={ l=list(); for (i in 1:n) l[[i]] = rnorm(m) }, liist={ l=vector('list', n); for (i in 1:n) l[[i]] = rnorm(m) }, matrix=matrix(rnorm(n*m), n, m), matrix2 = {mat - rnorm(n*m); dim(mat) - c(n, m); mat}, as.matrix=as.matrix(rnorm(n*m), n, m), replicate=replicate(m, rnorm(n)) ) # 3matrix 0.173 # 4 matrix2 0.162 # 5 as.matrix 0.169 vQ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Wacek Kusnierczyk waclaw.marcin.kusnierc...@idi.ntnu.no wrote Seriously? You think: lapply(1:n, rnorm, 0, 1) is 'clearer' than: x=list() for(i in 1:n){ x[[i]]=rnorm(i,0,1) } for beginners? Firstly, using 'lapply' introduces a function (lapply) that doesn't have an intuitive name. Also, it takes a function as an argument. The concept of having a function as a parameter to another function is something that a lot of programming beginners have trouble with - unless they were brought up on LISP of course, and few of us are. I propose that the for-loop example is clearer to a larger population than the lapply version. As a beginner, I agree the for loop is much clearer to me. Peter Peter L. Flom, PhD Statistical Consultant www DOT peterflomconsulting DOT com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Peter Flom wrote: Seriously? You think: lapply(1:n, rnorm, 0, 1) is 'clearer' than: x=list() for(i in 1:n){ x[[i]]=rnorm(i,0,1) } for beginners? Firstly, using 'lapply' introduces a function (lapply) that doesn't have an intuitive name. Also, it takes a function as an argument. The concept of having a function as a parameter to another function is something that a lot of programming beginners have trouble with - unless they were brought up on LISP of course, and few of us are. I propose that the for-loop example is clearer to a larger population than the lapply version. As a beginner, I agree the for loop is much clearer to me. well, that's quite likely. especially given that typical courses in programming, afaik, include for looping but not necessarily functional stuff -- are you an r beginner, or a programming beginner? among the perl packages i have ever downloaded from cran, it's hard to find one without a for loop, but it's easy to find one without a map. but it's not necessarily because for loops are easier; just that that's the way people are typically taught to program. the structure and interpretation of computer programs (sicp) by abelson sussman, a beautiful cs masterpiece, introduces mapping (lapplying) on p. 105, mentions a for-each control abstraction only in an exercise two pages later, and does not really discuss for looping as such. functional mapping over stateless objects is, in general, *much* easier to reason with than procedural looping over stateful objects -- an issue a beginner may not be quite aware of, and learning the basic for loop stuff without caring about, e.g., concurrent access to shared mutable state etc. may indeed make the impression that for loops are easier. anyway, once you've learned for loops, it's not a bad idea to learn lapply. and once you've learned lapply, you'll probably not go back to for loops that easily. vQ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation)
I wrote As a beginner, I agree the for loop is much clearer to me. Wacek Kusnierczyk waclaw.marcin.kusnierc...@idi.ntnu.no replied well, that's quite likely. especially given that typical courses in programming, afaik, include for looping but not necessarily functional stuff -- are you an r beginner, or a programming beginner? Both. My PhD is in psychometrics, and, both in course work and since then I've learned a good bit of statistics, but very little programming. I've picked up a little SAS programming over the years, but not much. But the loop (at least for me) translates into English more directly than the lapply statement does. among the perl packages i have ever downloaded from cran, it's hard to find one without a for loop, but it's easy to find one without a map. but it's not necessarily because for loops are easier; just that that's the way people are typically taught to program. the structure and interpretation of computer programs (sicp) by abelson sussman, a beautiful cs masterpiece, introduces mapping (lapplying) on p. 105, mentions a for-each control abstraction only in an exercise two pages later, and does not really discuss for looping as such. functional mapping over stateless objects is, in general, *much* easier to reason with than procedural looping over stateful objects -- an issue a beginner may not be quite aware of, and learning the basic for loop stuff without caring about, e.g., concurrent access to shared mutable state etc. may indeed make the impression that for loops are easier. anyway, once you've learned for loops, it's not a bad idea to learn lapply. and once you've learned lapply, you'll probably not go back to for loops that easily. Would that be a good book for a beginner? Peter Peter L. Flom, PhD Statistical Consultant www DOT peterflomconsulting DOT com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
As a beginner, I agree the for loop is much clearer to me. [Warning: Contains mostly philosophy] To me, the world and how I interact with it is procedural. When I want to break six eggs I do 'get six eggs, repeat break egg until all eggs broken'. I don't apply an instance of the break egg function over a range of eggs. My world is not functional (just like me, some might say...). Neither do I send a 'break yourself' message to each egg - my world is not object-oriented. That does not mean that these paradigms are not good ways of writing computer programs - they are brilliant ways of writing computer programs. But they build on procedural concepts, and we don't teach children to run before they can walk. So when someone says 'how do I do this a thousand times?' on R-help, I'll assume their knowledge level is that of a beginner, and try to map the solution to their world view. Computer scientists will write their beautiful manuscripts, but how many people who come to R because they want to do a t-test or fit a GLM will read them? That's the R-help audience now. Barry __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation)
Peter Flom wrote: As a beginner, I agree the for loop is much clearer to me. well, that's quite likely. especially given that typical courses in programming, afaik, include for looping but not necessarily functional stuff -- are you an r beginner, or a programming beginner? Both. My PhD is in psychometrics, and, both in course work and since then I've learned a good bit of statistics, but very little programming. I've picked up a little SAS programming over the years, but not much. don't really know sas, but i guess for looping is of essence there, while mapping is not. But the loop (at least for me) translates into English more directly than the lapply statement does. lapply easily translates to 'apply this to every item there', which is roughly an alternative version of 'for each item in there, do this with the item'. the structure and interpretation of computer programs (sicp) by abelson sussman, a beautiful cs masterpiece, introduces mapping (lapplying) on p. 105, mentions a for-each control abstraction only in an exercise two pages later, and does not really discuss for looping as such. functional mapping over stateless objects is, in general, *much* easier to reason with than procedural looping over stateful objects -- an issue a beginner may not be quite aware of, and learning the basic for loop stuff without caring about, e.g., concurrent access to shared mutable state etc. may indeed make the impression that for loops are easier. Would that be a good book for a beginner? both yes and no. this is a book that can be used by an absolute beginner in programming, but if you're focused on statistics, you're unlikely to enjoy it, at least not as a practical introduction. but it's a good read, and contains quite a lot of useful ideas anyway. vQ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Barry Rowlingson wrote: As a beginner, I agree the for loop is much clearer to me. [Warning: Contains mostly philosophy] maybe quasi ;) To me, the world and how I interact with it is procedural. When I want to break six eggs I do 'get six eggs, repeat break egg until all eggs broken'. yeah, that's the implementation level. a typical recipe would not say 'for n from 1 to 6, break the nth egg'. it'd rather say 'break the eggs', which is closer to 'apply break to the eggs'. you do of course break the eggs sequentially (or?), but that's below the abstraction useful for the recipe purpose. the example is quite good, in fact: the lapply approach is more appropriate, since what you're interested in is actually starting with six eggs and ending with six broken eggs; neither the particular order you might have chosen, nor whether you break the eggs sequentially or in parallel. in a sense, 'for n from 1 to 6, break the nth egg' is a particular operationalization of 'apply break to the eggs'. I don't apply an instance of the break egg function over a range of eggs. but you do apply the break function to the ith egg in a for loop? and actually, the procedural way fo breaking eggs can be easily written without a for loop, e.g., using a ruby-ish notation: eggs.each { |egg| break egg } which says, 'for each egg, take the egg and break it', without using a dummy index. vQ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
On 14-May-09 11:28:17, Wacek Kusnierczyk wrote: Barry Rowlingson wrote: As a beginner, I agree the for loop is much clearer to me. [Warning: Contains mostly philosophy] maybe quasi ;) To me, the world and how I interact with it is procedural. When I want to break six eggs I do 'get six eggs, repeat break egg until all eggs broken'. yeah, that's the implementation level. a typical recipe would not say 'for n from 1 to 6, break the nth egg'. it'd rather say 'break the eggs', which is closer to 'apply break to the eggs'. you do of course break the eggs sequentially (or?), but that's below the abstraction useful for the recipe purpose. But it does influence how you organise the subsequent garbage collection. :) Ted. E-Mail: (Ted Harding) ted.hard...@manchester.ac.uk Fax-to-email: +44 (0)870 094 0861 Date: 14-May-09 Time: 12:40:16 -- XFMail -- __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Barry Rowlingson wrote: Computer scientists will write their beautiful manuscripts, but how many people who come to R because they want to do a t-test or fit a GLM will read them? That's the R-help audience now. don't forget that r seems to take, maybe undeserved, the pride of being a functional programming language. it has apparently been designed by people who don't think implementing cs guys' ideas is necessarily a bad idea. you do have lapply co in r, not without a purpose. certainly not to just please computer scientists. vQ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Thanks for everyone. I think the approach below is most suitable for me, as a beginner. x=list() for(i in 1:n){ x[[i]]=rnorm(i,0,1) } Now, I am trying to obtain the sample variance (S^2) of the 1000 samples that I have generated before. I am wondering what command I should use in order to get the sample variance for all the 1000 samples. What I am capable of doing now is just typing in var(z[[1]]) var(z[[2]]). Thanks for help. Debbie Date: Thu, 14 May 2009 08:26:03 +0200 From: waclaw.marcin.kusnierc...@idi.ntnu.no To: b.rowling...@lancaster.ac.uk CC: r-help@r-project.org Subject: Re: [R] Simulation Barry Rowlingson wrote: On Wed, May 13, 2009 at 9:56 PM, Wacek Kusnierczyk waclaw.marcin.kusnierc...@idi.ntnu.no wrote: Barry Rowlingson wrote: n = 1000 benchmark(columns=c('test', 'elapsed'), order=NULL, 'for'={ l = list(); for (i in 1:n) l[[i]] = rnorm(i, 0, 1) }, lapply=lapply(1:n, rnorm, 0, 1) ) # test elapsed # 1 for 9.855 # 2 lapply 8.923 Yes, you can probably vectorize this with lapply or something, but I prefer clarity over concision when dealing with beginners... but where's the preferred clarity in the for loop solution? Seriously? You think: lapply(1:n, rnorm, 0, 1) is 'clearer' than: x=list() for(i in 1:n){ x[[i]]=rnorm(i,0,1) } for beginners? seriously, i do; but it does depend on who those beginners are. if they come directly from c and the like, you're probably right. Firstly, using 'lapply' introduces a function (lapply) that doesn't have an intuitive name. Also, it takes a function as an argument. The concept of having a function as a parameter to another function is something that a lot of programming beginners have trouble with - unless they were brought up on LISP of course, and few of us are. well, that's one of the first things you learn on a programming languages course that is not procedural programming-{centered,biased}. no need for prior lisp experience. if messing with closures in not involved (as here), no need for advanced discussion is needed. also, the for looping may not be as trivial stuff to explain as you might think. note, you're talking about r, not c, and the treatment of iterator variables in for loops in scripting languages differs: perl -e ' $i = 0; for $i (1..5) { # iterate with $i }; print $i\n ' # 0 ruby -e ' i = 0 for i in 1..5 # iterate with i end printf %d\n, i ' # 5 and you've gotten into explaining lexical scoping etc. I propose that the for-loop example is clearer to a larger population than the lapply version. which population have you sampled from? you may not be wrong, but give some data. Plus it's only useful in that form if the first parameter is the one you want to lapply over. If you want to work over the third parameter, say, you then need: lapply(1:n,function(i){rnorm(100,0,i)}) at which point you've introduced anonymous functions. The jump from: x[[i]] = rnorm(i,0,1) to x[[i]] = rnorm(100,0,i) is much less than the changes in the lapply version, where you have to go 'oh hang on, lapply only works on the first argument, so you have to write another function, but you can do that inline like this...'. you may be unhappy to learn that you're unaware of how the lapply solution can still be nicely adapted here: lapply(1:n, rnorm, n=100, mean=0) Okay, maybe my example is a little contrived, but I still think for a beginners context it's important not to jump too many paradigms at a time. for a complete beginner, jump into for loops may not be that trivial as you seem to think. there's still quite some stuff to be explained to clarify that i = 0 for (i in 1:n) # do stuff print(i) will print n, not 0. unless n=0, of course. vQ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. _ Looking to change your car this year? Find car news, reviews and more e%2Ecom%2Fcgi%2Dbin%2Fa%2Fci%5F450304%2Fet%5F2%2Fcg%5F801459%2Fpi%5F1004813%2Fai%5F859641_t=762955845_r=tig_OCT07_m=EXT [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Debbie Zhang schrieb: Now, I am trying to obtain the sample variance (S^2) of the 1000 samples that I have generated before. I am wondering what command I should use in order to get the sample variance for all the 1000 samples. What I am capable of doing now is just typing in var(z[[1]]) var(z[[2]]). Common please use the package brain 2.0. What is the difference between this and what you already had here (???): for(i in 1:n){ x[[i]]=rnorm(i,0,1) } hint: its only a small change. Please do some introductionary exercises. You find plenty of material here: http://cran.r-project.org/manuals.html and here: http://cran.r-project.org/other-docs.html Stefan __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Since you got the most suitable way to get x, why can't you get the variances in the same way? Just like: v = vector() for (i in 1:length(x)) v[i] = var(x[[i]]) BTW, it is much better to use lapply, like this: lapply(x, var) On Thu, May 14, 2009 at 8:26 PM, Debbie Zhang debbie0...@hotmail.com wrote: Thanks for everyone. I think the approach below is most suitable for me, as a beginner. x=list() for(i in 1:n){ x[[i]]=rnorm(i,0,1) } Now, I am trying to obtain the sample variance (S^2) of the 1000 samples that I have generated before. I am wondering what command I should use in order to get the sample variance for all the 1000 samples. What I am capable of doing now is just typing in var(z[[1]]) var(z[[2]]). Thanks for help. Debbie __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Simulation
Stefan Grosse wrote: Debbie Zhang schrieb: Now, I am trying to obtain the sample variance (S^2) of the 1000 samples that I have generated before. I am wondering what command I should use in order to get the sample variance for all the 1000 samples. What I am capable of doing now is just typing in var(z[[1]]) var(z[[2]]). if you have the data produced the for loop way, i.e., as a list of vectors, you can go the intuitive way: vars = list() for (i in 1:1000) vars[[i]] = z[[i]] or the unintuitive way: vars = lapply(z, var) if you have the data produced the matrix or replicate way (i.e., a matrix with columns representing samples), you can go the intuitive way: vars = c() for (i in 1:1000) vars[i] = var(z[,i]) or the unintuitive way: vars = apply(z, 2, var) consider reading an intro to r unless you like to receive responses as the one below. vQ Common please use the package brain 2.0. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.