[R] Matrix multiplication and random numbers

Chris Stubben stubben at lanl.gov
Wed Sep 9 19:17:32 CEST 2009



RFish wrote:
> 
> I new to using R and am struggling with some matrix multiplication. 
> 

I'm not sure what you're trying to print, but you could place this vector in
an expression

mat3<-expression(c(0,rnorm(1,0.6021,0.0987),0,0,0,0,0,0,0,rnorm(1,0.6021,0.0987),0,0,0,0,1.9,0,0,rnorm(1,0.6021,0.0987),0,0,0,4.8,0,0,0,rnorm(1,0.6021,0.0987),0,0,9.7,0,0,0,0,rnorm(1,0.6021,0.0987),0,18,0,0,0,0,0,rnorm(1,0.6021,0.0987),32.6,0,0,0,0,0,0))

# and then evaluate to get a new matrix each time
matrix(eval(mat3), nrow=7)

#I think this may be easier to follow. First create a matrix of zeros, stick
in fertilities and then add random survival probabilities each time


mat3<-diag(0,7)
#fertilities
mat3[1,3:7]<-c(1.9, 4.8, 9.7, 18, 32.6)
# random survival on sub-diagonal
mat3[row(mat3)==col(mat3)+1]<-rnorm(6,0.6021,0.0987)


# and if you want to project the population over 10 time steps in a loop ?

n<-matrix(c(500,0,0,0,0,0,0))

popsize <- matrix(numeric(7 * 10), nrow = 7)
 for (i in 1:10) {
        popsize[, i] <- n
         mat3[row(mat3)==col(mat3)+1]<-rnorm(6,0.6021,0.0987)
         n <- mat3 %*% n
    }


     [,1]     [,2]     [,3]     [,4]      [,5]      [,6]       [,7]     
[,8]
[1,]  500   0.0000   0.0000 531.6256 709.89940 940.19337 1697.52862
3403.6610
[2,]    0 352.5116   0.0000   0.0000 298.97874 424.71160  561.32525
1027.1605
[3,]    0   0.0000 279.8029   0.0000   0.00000 231.45988  316.83352 
424.8883
[4,]    0   0.0000   0.0000 147.8957   0.00000   0.00000  136.36804 
220.7370
[5,]    0   0.0000   0.0000   0.0000  96.92715   0.00000    0.00000 
108.6551
[6,]    0   0.0000   0.0000   0.0000   0.00000  69.87527    0.00000   
0.0000
[7,]    0   0.0000   0.0000   0.0000   0.00000   0.00000   65.86229   
0.0000


Chris Stubben

-- 
View this message in context: http://www.nabble.com/Matrix-multiplication-and-random-numbers-tp25365184p25369495.html
Sent from the R help mailing list archive at Nabble.com.




More information about the R-help mailing list