[R] Matrix multiplication and random numbers

[RFish]

Hi

Sorry I don't seem to have explained what I'm trying to do very clearly.
The piece of code below multiplies the two matrices together a number of
times based on the value in the matmult(InitialPop,1) term in this case one
(year), this gives me the end population for the analysis. 

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

matmult<-function(InitialPop,N){
mat3<-matrix(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),nrow=7)

for (i in 1:N){
PVAmatrix<-matrix(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),nrow=7)
mat3<-mat3%*%PVAmatrix
}
ans<-mat3 %*% InitialPop
return(ans)
}
matmult(InitialPop,1)

The problem i have is to repeat this process say 1000 times and store this
output in a format I can export easily whilst maintaining the randomness of
the result, so that every end population is different. 

Any help would be brilliant

Tom  


 



Chris Stubben wrote:
> 
> 
> 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-tp25365184p25403899.html
Sent from the R help mailing list archive at Nabble.com.