NETWORK Procedure

Example 2.9 Cycle Enumeration for Kidney Donor Exchange

This example looks at an application of cycle enumeration to help create a kidney donor exchange. Suppose someone needs a kidney transplant and a family member is willing to be a donor. If the donor and recipient are incompatible (because of blood type, tissue mismatch, and so on), the transplant cannot happen. Now suppose two donor-recipient pairs, i and j, are in this situation, but donor i is compatible with recipient j and donor j is compatible with recipient i. Then two transplants can take place in a two-way swap, shown in Figure 230. More generally, an n-way swap can be performed involving n donors and n recipients (CNN 2012).

Figure 230: Kidney Donor Exchange Two-Way Swap

Kidney Donor Exchange Two-Way Swap


To model this problem, define a directed graph as follows: Each node is an incompatible donor-recipient pair. Link left-parenthesis i comma j right-parenthesis exists if the donor from node i is compatible with the recipient from node j, as shown in Figure 231.

Figure 231: Kidney Donor Exchange Network

Kidney Donor Exchange Network


The link weight is a measure of the quality of the match. By introducing dummy links whose weight is 0, you can also include recipients who have no donors and altruistic donors who have no recipients. The idea is to find a maximum-weight node-disjoint union of directed cycles. You want the union to be node-disjoint so that no kidney is donated more than once, and you want cycles so that the donor from node i donates a kidney if and only if the recipient from node i receives a kidney.

Without any other constraints, the problem could be solved as a linear assignment problem, as described in SAS Optimization: The OPTNETWORK Procedure. But doing so would allow arbitrarily long cycles in the solution. For practical considerations (such as travel) and to mitigate risk, each cycle must have no more than L links. The kidney exchange problem is to find a maximum-weight node-disjoint union of short directed cycles.

One way to solve this problem is to explicitly generate all cycles whose length is at most L and then solve a set-packing problem. You can use PROC NETWORK to generate the cycles and then use PROC OPTMODEL to read the PROC NETWORK output, formulate the set-packing problem, call the mixed integer linear programming solver, and output the optimal solution. See Chapter 9, The OPTMODEL Procedure (SAS Optimization: Mathematical Optimization Procedures).

The following DATA step sets up the problem by first creating a random graph on n nodes with link probability p and Uniform(0,1) weight:

/* create random graph on n nodes with link probability p
   and uniform(0,1) weight */
%let n = 100;
%let p = 0.02;
data mylib.LinkSetIn;
   call streaminit(1);
   do from = 0 to &n - 1;
      do to = 0 to &n - 1;
         if from eq to then continue;
         else if rand('UNIFORM') < &p then do;
            weight = rand('UNIFORM');
            output;
         end;
      end;
   end;
run;

The following statements use PROC NETWORK to generate all cycles whose length is greater than or equal to 2 and less than or equal to 10:

/* generate all cycles with 2 <= length <= max_length */
%let max_length = 10;
proc network
   logLevel     = moderate
   direction    = directed
   links        = mylib.LinkSetIn;
   cycle
      minLength      = 2
      maxLength      = &max_length
      maxCycles      = all
      outCyclesLinks = mylib.CyclesLinks;
run;
%put &_NETWORK_;

PROC NETWORK finds 395 cycles of the appropriate length, as shown in Output 2.9.1.

Output 2.9.1: PROC NETWORK Log: Cycles for Kidney Donor Exchange

NOTE: ------------------------------------------------------------------------------------------
NOTE: ------------------------------------------------------------------------------------------
NOTE: Running NETWORK.                                                                          
NOTE: ------------------------------------------------------------------------------------------
NOTE: ------------------------------------------------------------------------------------------
NOTE: Reading the links data.                                                                   
NOTE: Data input used 0.00 (cpu: 0.00) seconds.                                                 
NOTE: Building the input graph storage used 0.00 (cpu: 0.00) seconds.                           
NOTE: The number of nodes in the input graph is 98.                                             
NOTE: The number of links in the input graph is 208.                                            
NOTE: Processing cycle enumeration using 16 threads across 1 machines.                          
NOTE: Processing cycle enumeration using the build algorithm.                                   
NOTE: The algorithm found 395 cycles.                                                           
NOTE: Processing cycle enumeration used 0.00 (cpu: 0.01) seconds.                               
NOTE: The Cloud Analytic Services server processed the request in 0.209434 seconds.             
NOTE: The data set MYLIB.CYCLESLINKS has 3431 observations and 5 variables.                     
STATUS=OK  PROBLEM_TYPE=CYCLE  SOLUTION_STATUS=OK  NUM_CYCLES=395  CPU_TIME=0.80  REAL_TIME=0.21


For this set of cycles, you can now formulate a mixed integer linear program (MILP) to maximize the total cycle weight. Let C define the set of cycles of appropriate length, upper N Subscript c define the set of nodes in cycle c, upper E Subscript c define the set of links in cycle c, and w Subscript e denote the link weight for link e. Define a binary decision variable x Subscript c. Set x Subscript c to 1 if cycle c is used in the solution; otherwise, set it to 0. Then, the following MILP defines the problem that you want to solve in order to maximize the quality of the kidney exchange:

The constraint (incomp_pair) ensures that each node (incompatible pair) in the graph is intersected at most once. That is, a donor can donate a kidney only once. You can use PROC OPTMODEL to solve this mixed integer linear programming problem as follows:

/* solve set-packing problem to find maximum-weight node-disjoint union
   of short directed cycles */
proc optmodel;
   /* declare index sets and parameters, and read data */
   set <num,num> LINKS;
   num weight {LINKS};
   read data mylib.LinkSetIn into LINKS=[from to] weight;
   set <num,num,num> TRIPLES;
   read data mylib.CyclesLinks into TRIPLES=[cycle from to];
   set CYCLES = setof {<c,i,j> in TRIPLES} c;
   set LINKS_c {c in CYCLES} = setof {<(c),i,j> in TRIPLES} <i,j>;
   set NODES_c {c in CYCLES} = union {<i,j> in LINKS_c[c]} {i,j};
   set NODES = union {c in CYCLES} NODES_c[c];
   num cycle_weight {c in CYCLES} = sum {<i,j> in LINKS_c[c]} weight[i,j];

   /* UseCycle[c] = 1 if cycle c is used, 0 otherwise */
   var UseCycle {CYCLES} binary;

   /* declare objective */
   max TotalWeight
      = sum {c in CYCLES} cycle_weight[c] * UseCycle[c];

   /* each node appears in at most one cycle */
   con NodePacking {i in NODES}:
      sum {c in CYCLES: i in NODES_c[c]} UseCycle[c] <= 1;

   /* call solver */
   solve;

   /* output optimal solution */
   create data Solution from
      [c]={c in CYCLES: UseCycle[c].sol > 0.5} cycle_weight;
quit;
%put &_OROPTMODEL_;

PROC OPTMODEL solves the problem by using the mixed integer linear programming solver.

Output 2.9.2: PROC OPTMODEL Log: Cycles for Kidney Donor Exchange

NOTE: There were 208 observations read from the data set MYLIB.LINKSETIN.                       
NOTE: There were 3431 observations read from the data set MYLIB.CYCLESLINKS.                    
NOTE: Problem generation will use 4 threads.                                                    
NOTE: The problem has 395 variables (0 free, 0 fixed).                                          
NOTE: The problem has 395 binary and 0 integer variables.                                       
NOTE: The problem has 64 linear constraints (64 LE, 0 EQ, 0 GE, 0 range).                       
NOTE: The problem has 3431 linear constraint coefficients.                                      
NOTE: The problem has 0 nonlinear constraints (0 LE, 0 EQ, 0 GE, 0 range).                      
NOTE: The OPTMODEL presolver is disabled for linear problems.                                   
NOTE: The initial MILP heuristics are applied.                                                  
NOTE: The MILP presolver value AUTOMATIC is applied.                                            
NOTE: The MILP presolver removed 110 variables and 31 constraints.                              
NOTE: The MILP presolver removed 1685 constraint coefficients.                                  
NOTE: The MILP presolver added 1 constraint coefficients.                                       
NOTE: The MILP presolver modified 11 constraint coefficients.                                   
NOTE: The presolved problem has 285 variables, 33 constraints, and 1746 constraint coefficients.
NOTE: The MILP solver is called.                                                                
NOTE: The parallel Branch and Cut algorithm is used.                                            
NOTE: The Branch and Cut algorithm is using up to 4 threads.                                    
          Node   Active   Sols    BestInteger      BestBound      Gap    Time                   
             0        1      3     20.7728712   1199.7787581   98.27%       0                   
             0        1      3     20.7728712     25.4194215   18.28%       0                   
             0        1      3     20.7728712     24.8914878   16.55%       0                   
             0        1      4     21.1953805     24.8914878   14.85%       0                   
             0        1      5     24.8508554     24.8508554    0.00%       0                   
             0        0      5     24.8508554     24.8508554    0.00%       0                   
NOTE: The MILP solver added 34 cuts with 5416 cut coefficients at the root.                     
NOTE: Optimal.                                                                                  
NOTE: Objective = 24.850855395.                                                                 
NOTE: The data set WORK.SOLUTION has 7 observations and 2 variables.                            
STATUS=OK ALGORITHM=BAC SOLUTION_STATUS=OPTIMAL OBJECTIVE=24.850855395 RELATIVE_GAP=0           
ABSOLUTE_GAP=0 PRIMAL_INFEASIBILITY=5.107026E-14 BOUND_INFEASIBILITY=7.338574E-14               
INTEGER_INFEASIBILITY=7.338574E-14 BEST_BOUND=24.850855395 NODES=1 SOLUTIONS_FOUND=5            
ITERATIONS=110 PRESOLVE_TIME=0.02 SOLUTION_TIME=0.05                                            


The output data table Solution, shown in Output 2.9.3, now contains the cycles that define the best exchange and their associated weight (quality).

Output 2.9.3: Maximum-Quality Solution for Kidney Donor Exchange

ccycle_weight
264.3542
814.3403
1344.9748
1695.0843
3621.9424
3852.4295
3911.7253
 24.8509


Last updated: August 07, 2026