systemicrisk (version 0.4.2)

calibrate_FitnessEmp: Calibrate empirical fitness model to a given density

Description

The model is an empirical fitness based model for the existence of links (more details below) which contains one fixed parameter and the weight of each existing link follows an exponential distribution with a fixed rate parameter. This function chooses the two parameters such that the density of the network (the average proportion of existing links) with these given row and column sums is a certain desired value.

Usage

calibrate_FitnessEmp(l, a, targetdensity, L_fixed = NA,
  nsamples_calib = 100, thin_calib = 100)

Arguments

l

row sums of matrix to be reconstructed

a

column sum of matrix to be reconstructed

targetdensity

desired proportion of reconstructed entries to be positive

L_fixed

Matrix containing known values of L, where NA signifies that an element is not known. If L_fixed equates to NA (the default) then no values are assumed to be known.

nsamples_calib

number of matrices to generate during calibration.

thin_calib

amount of thinning to use during calibration

Value

Model that can be used to generate the desired samples using sample_HierarchicalModel.

Details

The empirical fitness model assumes that every node \(f[i]=log(l[i]+a[i])\) has a fitness given by the observered row and column sum and that the existence probability of a link between node i and j is then given by \(p[i,j]=1/(1+exp(-(alpha+f[i]+f[j])))\), where alpha is an additional parameter. The resulting model uses observed quantities (the row and column sums of the matrix) as input to the model and is thus an empirical Bayes approach.

Examples

# NOT RUN {
## first generate a true network
n <- 10 # size of network
ftrue <- rnorm(n) # vector of underlying fitnesses
p <- outer(ftrue,ftrue,FUN=function(x,y) 1/(1+exp(-(x+y))))
lambda <- 0.1
L <- matrix(nrow=n,rbinom(n*n,prob=p,size=1)*rexp(n*n,rate=lambda))

# then reconstruct with a target density of 0.7
model <- calibrate_FitnessEmp(l=rowSums(L),a=colSums(L),
                              targetdensity=0.7,nsamples_calib=10,thin_calib=50)
Lsamp <- sample_HierarchicalModel(l=rowSums(L),a=colSums(L),model=model,
                                    nsamples=10,thin=1e2)
# check row sums
rowSums(L)
rowSums(Lsamp$L[[10]])
# check calibration
mean(Lsamp$L[[10]]>0)

# now an example with some fixed entries
L_fixed <- L
L_fixed[1:(n/2),] <- NA
# then reconstruct with a target density of 0.9
model <- calibrate_FitnessEmp(l=rowSums(L),a=colSums(L),L_fixed=L_fixed,
                              targetdensity=0.9,nsamples_calib=10,thin_calib=50)
Lsamp <- sample_HierarchicalModel(l=rowSums(L),a=colSums(L),L_fixed=L_fixed,
                                  model=model,nsamples=10,thin=1e2)
mean(Lsamp$L[[10]][-(1:(n/2)),]>0) # known entries
mean(Lsamp$L[[10]][(1:(n/2)),]>0)  #reconstructed entries
# }