This online documentation describes the necessary steps to perform the network optimization introduced in the paper What is the Minimal Systemic Risk in Financial Exposure Networks? by Diem, Pichler & Thurner (2019). A preprint of the paper is available on arXiv:1905.05931. The goal of this paper is to show the huge systemic risk reduction potential stemming from the network structure of financial exposure networks. We showed that for ten snap shots of the Austrian interbank market the proposed optimization substantially reduces systemic risk – measured with DebtRank (Battiston et al.,2012). First, we brefiefly outline the structure of the optimization problem and its formulation as mixed integer linear program (MILP). For the details of the reformulation and interpretation of the quantities we refer to the paper.

First, we construct a randomly generated sample network as reference network. Second, we provide the R code to create the objective function and the constraints to solve optimization problem with the R package ROI (Theußl et al., 2019). Finally we visualize the results for the randomly construced sample network.

1 The Network optimization Problem

Given an interbank network, \(L\), we want to minimize its DebtRank. Because the calculation of DebtRank is an interative procedure it is difficult to optimize. Thus, we approximate DebtRank by the so called direct impacts of a bank on its neighbours (adjacent banks in the network). If bank \(i\) defaults the direct loss of bank \(j\) is \(\min \Big( \frac{L_{ij}}{e_j},1 \Big)\). In the objective function bank \(j\) has an economic weight corresponding to its share of all interbank market loans, \(\frac{a_j}{\bar{L}}\). During the optimization the total amount of loans a bank has borrowed from (row sums of \(L\)) or lent to (column sums \(L\)) other banks shall stay the same. Additionally, the optimization shall leave the sum of loans weighted by their credit riskiness \(\kappa\) unchanged. We showed that this problem can be formulated as \[\begin{align} \min_{L \in \{ M : \; M \in \mathbb{R}_+^{n \times n}, \; M_{ii}=0 \} } & \sum_{i=1}^{n}\sum_{j=1}^{n} \min \Big( \frac{L_{ij}}{e_j},1 \Big) \frac{a_j}{\bar{L}} \\ \text{s.t.} \qquad & l_i = \sum_{j=1}^{n} L_{ij} \quad \forall i \tag{1}\\ & a_i = \sum_{j=1}^{n} L_{ji} \quad \forall i \nonumber, \\ & r_i = \sum_{j=1}^{n} L_{ji}\kappa_j. \end{align}\] Details on the interpretation and definition of the quantities involved can be found in the paper.

1.1 Equivalent Mixed Integer Linear Program

We showed how the above optimization problem (Eq.1) can be rewritten into a mixed integer linear program (MILP), as

\[\begin{eqnarray} \label{eq:optim_milp} \min_{ \{z \; \in \; \mathbb{R}_+^{4N^2} \} } & & c^\top z & \\ \text{s.t.} \quad & A_1z & \leq & 0 \tag{2} \\ & A_2z & = & a \\ & A_3z & = & l \\ & A_4z & = & r. \end{eqnarray}\]

The first \(2N^2\) entrities of \(z\) correspond to the entrities of the interbank liability matrix \(L\) by \(L_{11}=z_1+z_2\), \(L_{21}=z_3+z_4\), \(\dots\),the second set of \(2N^2\) variables are binary variables, where \(z_{2N^2+k}=1\), if \(z_{k}>0\) and zero otherwise. \(A_1\) contains the constraints on the variables \(z\) to ensure the equivalence of (2) and (1). \(A_2\) and \(A_3\) represent the constraints that the row and column sums need stay constant during the optimization and \(A_4\) implents the credit risk constraint. For implementing other constraints with economic interpretations see Section 3.1.

1.2 Generate a Reference Network to Optimize

Before we can apply our optimization procedure we need an “empirically” observed financial network we can optimize. The optimization without the credit risk constraint can also be performed with the row and column sums of \(L\) and the equity vector \(e\). To illustrate the optimization, we use a randomly generated network with 25 banks. First, we sample the interbank asset vector \(a\) and use a perturbed version as the interbank liability vector \(l\). We set the equity vector \(e\) proportional to the mean of \(a\) and \(l\). The sample is constructed such that we have 2 large banks, 3 medium sized banks and 20 smaller banks. Then we use the R package systemicrisk from Gandy et al. (2016) to generate a sample network. The exact construction of the sample network, \(L\), is shwon in the code chunk below.

For the illustration we choose the credit risk proxy \(\kappa\) to be interbank leverage ratio of the bank; i.e. kappa <- l/E. In practice different credit risk proxies would be desirable.

2 Creating & Solving the MILP

2.1 Create Objective Functions and Constraints for MILP

The quantities \(c, A_1, A_2, A_3, A_4\) from (Eq.2) are generated with the code chunk below. The details of the construction for the objective function can be found in Diem et al. (2019), (Eq.19), and the details for the constraint matrices are explained in Appendix A.

library(slam)
N <- 4*n^2

# Construct the objective function vector c

slopes <- rep(a/E, each=n)
c <- c(as.vector(rbind(slopes,rep(0,n*n))), rep(0, 2*n*n))
  
  
  
# Construct the four constraint matrices A1, A2, A3 & A4 
 
# A1 

  triplet_1 <- matrix(numeric(3*8*n^2), ncol=3)
  a_l <- rep(a, each = n)
  e_l <- rep(E, each = n)
  l_l <- rep(l, n)
  u_b  <- pmax( rep(0,n^2), pmin(a_l, l_l) - e_l)
  
  for(i in 1:(n^2)){                   #row     col         value
    triplet_1[(i*8-7):(i*8), ] <- rbind(c(4*i-3, 2*i-1,        1),
                                       c(4*i-3, 2*i-1+2*n^2, -e_l[i]),
                                       
                                       c(4*i-2, 2*i-1,       -1),
                                       c(4*i-2, 2*i+2*n^2,    e_l[i]),
                                       
                                       c(4*i-1, 2*i,          1),
                                       c(4*i-1, 2*i+2*n^2,   -u_b[i]),
                                       
                                       c(4*i,   2*i-1+2*n^2, -1),
                                       c(4*i,   2*i+2*n^2,    1))
  }
  
  A1 <- simple_triplet_matrix(triplet_1[,1], triplet_1[,2], triplet_1[,3])
  
  
  
  # A2

  triplet_2 <- matrix(numeric(3*(2*n^2)), ncol=3)
  
  for(i in 1:n){
    for(j in 1:n){
      v <- if(i != j){1}else{0}
      triplet_2[((2*i-2)*n+2*j-1):((2*i-2)*n+2*j), ] <- rbind( c(i, (2*j-2)*n  +2*i-1, v),
                                                               c(i, (2*j-2)*n  +2*i,   v)) 
    }
  }
  
  
  A2 <- cbind( simple_triplet_matrix(triplet_2[,1], triplet_2[,2], triplet_2[,3]),
               simple_triplet_zero_matrix(n, 2*n^2))

             
  # A3
  
  triplet_3 <- matrix(numeric(3*(2*n^2)), ncol=3)
  
  for(i in 1:n){
    for(j in 1:n){
      v <- if(i != j){1}else{0}
      triplet_3[((2*i-2)*n+2*j-1):((2*i-2)*n+2*j), ] <- rbind( c(i, (2*i-2)*n  +2*j-1, v),
                                                               c(i, (2*i-2)*n  +2*j,   v)) 
    }
  }
  
  A3 <- cbind( simple_triplet_matrix(triplet_3[,1], triplet_3[,2], triplet_3[,3]),
               simple_triplet_zero_matrix(n, 2*n^2))
         
  
  
  # A4
  
  triplet_4 <- matrix(numeric(3*(2*n^2)), ncol=3)
  
  for(i in 1:n){
    for(j in 1:n){
      v <- if(i != j){1*kappa[j]}else{0}
      triplet_4[((2*i-2)*n+2*j-1):((2*i-2)*n+2*j), ] <- rbind( c(i, (2*i-2)*n  +2*j-1, v),
                                                               c(i, (2*i-2)*n  +2*j,   v)) 
    }
  }
  
  A4 <- cbind( simple_triplet_matrix(triplet_4[,1], triplet_4[,2], triplet_4[,3]),
               simple_triplet_zero_matrix(n, 2*n^2))
       
  
  
  # Merge the constraints together
  A <- rbind(A1, A2, A3, A4)
  # Check the dimensions
  A
#> A 2782x2704 simple triplet matrix.
  
  
  # Create upper bounds for the variables
  ub <- numeric(n^2)
  
  for(i in 1:n){   # col index
    for(j in 1:n){ # row index
      ub[((i-1)*n)+j] <- max(0, min( a[i],l[j])-E[i])
    }
  }
  
  

  # Create a vector of zeros and ones to ensure that the diagonal elements L_ii=0 
  dia_constr <- matrix(rep(1,N/2), ncol = 2*n)
  
  for(j in 1:n){ # col 
    
    dia_constr[j, (2*(j-1)+1):(2*(j-1)+2)] <- c(0,0)
    
  }
  dia_constr <- as.vector(t(dia_constr))
 
  
  # Create a long vector of equity entries
  E_long <- rep(E, each = n)
  
  # Create a vector so that y_i1 is constraint by E_long and y_i2 by ub
  u_long <- as.vector(rbind(E_long, ub))
  
  # combine the upper bounds with the zero upper bounds
  u <- u_long*dia_constr
  
  # Create b, t, d for ROI OP()
  b <- c(rep(0,N), l, a, crossprod(L,kappa)) 
  t <- c(rep("C",2*n^2), rep("B", 2*n^2))
  d <- c(rep("<=", N),rep("==", 3*n))

2.2 Solving the optimization with ROI

Now we can create the ROI object with the function OP() and solve it with ROI_solve(). Then, we transform the solution vector again into a liability matrix. For the optimization we use the cplex solver, which was working best in our experience. Smaller problems could also be solved with the open source solver SYMPHONY.

Note that we can find a network with maximal direct impacts, simply by setting maximum = TRUE in the OP() function. This yields a network with a substantially higher DebtRank than the original one.

3 Comparing Original and Optimized Networks

We compare the direct impacts and DebtRanks of the three network types in Table 1. The optimized network show substantially smaller and higher DebtRank as the reference network.

Table 1, Comparison of Network DebtRanks
Minimized Empirical Maximized
DebtRank 4.48 9.69 17.91
Direct impacts 1.26 2.37 2.52

Figure 1 shows the systemic risk profile of the reference network and the minimized network. The optimization reduces the systemic risk for most banks.

Figure 1, Systemic Risk Profile

Figure 1, Systemic Risk Profile

3.1 Visualizing the Networks

The three networks are visualized in Figure 21. The node colors indicate the level of DebtRank, the node sizes correspond to bank equity and the link widths to the size of the loan.

Figure 2, Network VisualisationFigure 2, Network VisualisationFigure 2, Network Visualisation

Figure 2, Network Visualisation

4 References

Battiston, S., Puliga, M., Kaushik, R., Tasca, P., & Caldarelli, G. (2012). DebtRank: Too Central to Fail? Financial Networks, the FED and Systemic Risk. Scientific Reports, 2, 541. doi:10.1038/srep00541.

Csardi, G., & Nepusz, T. (2006). The igraph software package for complex network research. InterJournal, Complex Systems, 1695(5), 1-9.

Diem, C., Pichler, A., & Thurner, S. (2019). What is the Minimal Systemic Risk in Financial Exposure Networks?. arXiv preprint arXiv:1905.05931.

Gandy, A., & Veraart, L. A. (2016). A Bayesian Methodology for Systemic Risk Assessment in Financial Networks. Management Science, 63(12), 4428-4446. doi:10.1287/mnsc.2016.2546.

Theußl, S., Schwendinger, F., & Hornik, K. (2019) ROI: An extensible R Optimization Infrastructure. Research Report Series / Department of Statistics and Mathematics, 133. WU Vienna University of Economics and Business, Vienna.

Zeileis, A., Hornik, K., & Murrell, P. (2009). Escaping RGBland: Selecting colors for statistical graphics. Computational Statistics & Data Analysis, 53(9), 3259-3270.

  1. The networks are plotted with the R igraph package (Csardi et al., 2006) and colors are from the R colorspace package (Zeileis et al., 2009).↩