Re: [R] simulation of PowerGARCH,Threshold GARCH,and GJR GARCH

2019-03-25 Thread Eric Berger
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

2018-01-30 Thread ruipbarradas
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

2018-01-30 Thread Eric Berger
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 Barradas  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/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

2018-01-30 Thread Rui Barradas

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

2018-01-30 Thread Daniel Nordlund

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

2016-04-22 Thread David L Carlson
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-help  on 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)

2015-09-04 Thread David Winsemius
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

2014-09-15 Thread PO SU


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

2014-09-14 Thread Don McKenzie
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

2014-09-14 Thread Don McKenzie
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

2014-09-14 Thread David Winsemius


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

__
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

2014-09-14 Thread David Winsemius
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

2014-08-10 Thread William Revelle
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

2014-08-01 Thread Charles Determan Jr
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

2014-07-31 Thread Charles Determan Jr
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

2014-04-10 Thread Rui Barradas

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

2014-04-10 Thread Charles Determan Jr
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

2014-04-09 Thread Rolf Turner

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

2014-04-05 Thread Charles Determan Jr
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

2014-04-04 Thread Charles Determan Jr
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

2014-02-24 Thread Greg Snow
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

2013-08-04 Thread Rui Barradas

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

2013-08-04 Thread Jeff Newmiller
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

2013-08-04 Thread Preetam Pal
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

2013-08-04 Thread Prof Brian Ripley

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

2013-06-12 Thread Mikhail Umorin
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

2013-04-18 Thread Adams, Jean
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

2013-01-03 Thread Berend Hasselman

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

2013-01-03 Thread Berend Hasselman

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

2012-12-01 Thread Greg Snow
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

2012-11-19 Thread Uwe Ligges



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

2012-05-28 Thread Özgür Asar
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

2012-05-28 Thread R. Michael Weylandt
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

2012-04-05 Thread David Winsemius


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

2011-12-16 Thread Bert Gunter
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

2011-11-13 Thread R. Michael Weylandt
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

2011-10-26 Thread BSanders
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

2011-10-26 Thread Mehmet Suzen
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

2011-07-21 Thread Benjamin Caldwell
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

2011-07-21 Thread David Winsemius


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

2011-06-06 Thread Petr Savicky
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

2011-05-29 Thread SERDAR NESLIHANOGLU
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

2011-05-29 Thread Jorge Ivan Velez
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

2011-05-17 Thread cassie jones
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

2011-05-16 Thread Greg Snow
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

2011-05-16 Thread Greg Snow
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

2011-05-16 Thread cassie jones
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 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.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

2011-05-02 Thread Andrew Robinson
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

2011-04-08 Thread Søren Højsgaard
?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

2011-04-07 Thread Dennis Murphy
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

2011-03-01 Thread Ivan Calandra

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

2011-03-01 Thread Mike Marchywka






 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

2011-02-11 Thread Wonsang You
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

2011-02-10 Thread Kjetil Halvorsen
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

2011-01-05 Thread Mike Marchywka





 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

2011-01-05 Thread Ben Ward

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

2011-01-05 Thread Ben Ward

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

2010-11-10 Thread Giovanni Petris

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

2010-11-09 Thread Dennis Murphy
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

2010-09-10 Thread Jim Silverton
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

2010-09-10 Thread Greg Snow
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.

2010-04-14 Thread Greg Snow
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

2010-03-29 Thread Pfaff, Bernhard Dr.
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

2010-03-29 Thread Ron_M

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

2010-01-23 Thread Juliet Hannah
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

2010-01-20 Thread Rolf Turner


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

2010-01-13 Thread Tal Galili
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

2010-01-13 Thread Peter Ehlers

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

2009-08-18 Thread David Winsemius


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

2009-08-18 Thread jim holtman
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

2009-08-15 Thread Zhiliang Ma
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

2009-07-15 Thread Charles C. Berry

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

2009-07-15 Thread Roger Bivand

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

2009-05-18 Thread Kon Knafelman

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

2009-05-18 Thread Liaw, Andy
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

2009-05-18 Thread Uwe Ligges



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

2009-05-18 Thread Peter Dalgaard

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

2009-05-15 Thread Kon Knafelman

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

2009-05-15 Thread Stefan Grosse
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

2009-05-15 Thread Ben Bolker


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

2009-05-15 Thread Greg Snow
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

2009-05-15 Thread Ben Bolker
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

2009-05-15 Thread Greg Snow
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

2009-05-15 Thread Wacek Kusnierczyk
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

2009-05-14 Thread Dimitris Rizopoulos

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

2009-05-14 Thread Wacek Kusnierczyk
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

2009-05-14 Thread Wacek Kusnierczyk
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

2009-05-14 Thread Wacek Kusnierczyk
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

2009-05-14 Thread Peter Flom
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

2009-05-14 Thread Wacek Kusnierczyk
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)

2009-05-14 Thread Peter Flom
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

2009-05-14 Thread Barry Rowlingson
 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)

2009-05-14 Thread Wacek Kusnierczyk
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

2009-05-14 Thread Wacek Kusnierczyk
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

2009-05-14 Thread Ted Harding
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

2009-05-14 Thread Wacek Kusnierczyk
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

2009-05-14 Thread Debbie Zhang

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

2009-05-14 Thread Stefan Grosse
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

2009-05-14 Thread Linlin Yan
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

2009-05-14 Thread Wacek Kusnierczyk
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.


  1   2   >