Skip to content

How to repeat a procedure

4 messages · Haiyan Chen, Thomas Lumley, Peter Wolf +1 more

Haiyan Chen
# 
Hello,

1. After I generate a 100x50 matrix by x3<-matrix(0,100,50);for (i in
1:100) {x1<-rpois(50, mu[i]);x2<-x1; x2[runif(50)<.01]<-0; x3[i,]<-x2},

2. I want to calculate means and sample variances of each row and create a
new matrix 100x2;

3. I then want to repeat above procedure 500 times so that eventually I
will have 500 100x2 matrices.

Would someone mind helping me to code 2 & 3?

Thanks ahead of time.

Heyen
Thomas Lumley
#
On Wed, 18 Feb 2004, Haiyan Chen wrote:

            

      
      
      
    

        
YOu can do this without the loop, eg

x3<-rpois(50*100, rep(mu,each=100))
x3<-ifelse(runif(50*100)<0.01, 0, x3)
x3<-matrix(x3, ncol=50)

            

      
      
      
    

        
means<-apply(x3,2,mean)
vars<-apply(x3,2,var)
cbind(means,vars)

            

      
      
      
    

        
make.a.matrix<-function(...){
	x3<-rpois(50*100, rep(mu,each=100))
	x3<-ifelse(runif(50*100)<0.01, 0, x3)
	x3<-matrix(x3, ncol=50)
	cbind(apply(x3,2,mean), apply(x3,2, var))
}

many.matrices<-lapply(1:500, make.a.matrix)

gives a list of 500 matrices.  This isn't quite the most efficient
solution, but it's not bad.


	-thomas

  
  


  


      

    

    

      
      

        


  

        

      Peter Wolf
    

    

      
      #
    

  


  

      
Haiyan Chen wrote:

        

            

      
      
      
    

        
Try:

mean.var<-function(n=10, m=5){
 mu<-(1:n)^2
 x<-matrix(rpois(n*m,mu),n,m)
 stat<-t(apply(x,1,function(x)c(mean(x),var(x))))
}

to set some elements of x to zero you can add an additional line in mean.var

            

      
      
      
    

        
n<-100; m<-50
n.comp<-500
for(i in 1:n.comp) 
eval(parse(text=paste("result",i,"<-mean.var(n,m)",sep="")))


Peter

            

      
      
      
    

        

      
      
    

  
  


  


      

    

    

      
      

        


  

        

      Patrick Burns
    

    

      
      #
    

  


  

      
I believe that Thomas got "mu" wrong.  If I understand
correctly, the line:

x3 <- rpois(50 * 100, rep(mu, each=100))

should read:

x3 <- rpois(50 * 100, rep(mu, 50))

or just

x3 <- rpois(50 * 100, mu)

Patrick Burns

Burns Statistics
patrick at burns-stat.com
+44 (0)20 8525 0696
http://www.burns-stat.com
(home of S Poetry and "A Guide for the Unwilling S User")

        
Thomas Lumley wrote: