This document describe an example of reconstructing a non-square matrix.
library(systemicrisk)
Set up vectors of row and column sums of length 10 and 4, respectively, ensuring that their sums matches up.
set.seed(1230)
l <- rexp(10,rate=5/4)
l
## [1] 0.11213094 0.24672257 0.42200570 0.59310598 0.96265285 0.90072506
## [7] 0.45717999 0.04956001 0.59984659 0.37724666
a <- rexp(4);
a <- a/sum(a)*sum(l)
a
## [1] 0.7983658 1.7702032 0.9808514 1.1717560
Construct a model with fixed matrix p and lambda.
mod <- Model.Indep.p.lambda(model.p=Model.p.constant(p=matrix(0.5,nrow=10,ncol=4)),
model.lambda=Model.lambda.constant(lambda=matrix(5,nrow=10,ncol=4)))
Run the sampler to generate 10 samples.
res <- sample_HierarchicalModel(l=l,a=a,model=mod,nsamples=10,silent=TRUE)
Here are the first two of these samples.
res$L[[1]]
## [,1] [,2] [,3] [,4]
## [1,] 0.08427217 0.0000000 0.02785877 0.0000000
## [2,] 0.24672257 0.0000000 0.00000000 0.0000000
## [3,] 0.00000000 0.0000000 0.00000000 0.4220057
## [4,] 0.37518091 0.0000000 0.02987667 0.1880484
## [5,] 0.09219012 0.0000000 0.30876088 0.5617019
## [6,] 0.00000000 0.3359300 0.56479508 0.0000000
## [7,] 0.00000000 0.4571800 0.00000000 0.0000000
## [8,] 0.00000000 0.0000000 0.04956001 0.0000000
## [9,] 0.00000000 0.5998466 0.00000000 0.0000000
## [10,] 0.00000000 0.3772467 0.00000000 0.0000000
res$L[[2]]
## [,1] [,2] [,3] [,4]
## [1,] 0.03951294 0.0000000 0.000000000 0.0726180
## [2,] 0.00000000 0.0000000 0.007797952 0.2389246
## [3,] 0.00000000 0.0000000 0.000000000 0.4220057
## [4,] 0.00000000 0.0000000 0.425118402 0.1679876
## [5,] 0.69243280 0.0000000 0.000000000 0.2702200
## [6,] 0.00000000 0.7383240 0.162401046 0.0000000
## [7,] 0.00000000 0.4571800 0.000000000 0.0000000
## [8,] 0.00000000 0.0000000 0.049560010 0.0000000
## [9,] 0.02514738 0.5746992 0.000000000 0.0000000
## [10,] 0.04127266 0.0000000 0.335974000 0.0000000