Deterministic Blockmodeling of Two-Mode Binary Networks Using a Two-Mode KL-Median Heuristic

There are numerous examples of two-way, two-mode data in the context of social networks. In affiliation networks (Faust, 1997; Wasserman & Faust, 1994), the rows of the network matrix commonly pertain to individuals (actors), while the columns often correspond to events the actors attended (Davis, Gardner, & Gardner, 1941) or clubs of which the actors are members (Galaskiewicz, 1985). Other examples of two-mode network data include the participation of social/political groups in community projects (Mische & Pattison, 2000), the participation of individuals in electronic forums (Opsahl, 2013), the ratings viewers provided for movies (Harper & Konstan, 2015), the endorsement of political advertisements by organizations (Brusco & Doreian, 2015a), subject-by-item response data (Coombs, 1964; Hubert, 1974), U.S. Supreme Court voting data (Doreian, Batagelj, & Ferligoj, 2004, 2005), United Nations General Assembly (UNGA) voting data (Brusco, Doreian, Lloyd, & Steinley, 2013a; Doreian, Lloyd, & Mrvar, 2013), and the presence of words in documents (Xu, Liu, & Gong, 2001).

A variety of methods have been used to analyze two-mode social network data, including centrality analysis (Faust, 1997), Q-analysis (Doreian, 1979; Freeman, 1980), permutation procedures (Arabie, Hubert, & Schleutermannm 1990; Brusco & Steinley, 2006), and matrix factorization methods (Brusco, 2011; Everett & Borgatti, 2013; Pattison, 1993; Pattison & Bartlett, 1982; Pattison & Brieger, 2002). Here, we restrict our focus to deterministic two-mode blockmodeling procedures for social networks, as they have been particularly common in recent years (Brusco & Doreian, 2015a, 2015b; Brusco et al., 2013a; Brusco, Doreian, Mrvar, & Steinley, 2013b; Brusco & Steinley, 2007a, 2011; Doreian et al., 2004, 2005, 2013). More specifically, we focus on partitioning methods for deterministic blockmodeling of binary two-mode networks based on the principles of structural equivalence (Lorrain & White, 1971).

As formalized by Doreian et al. (2004, 2005), the goal of the deterministic two-mode blockmodeling problem (DTMBP) is to produce K clusters of the n row objects and L clusters of the m column objects with the goal of producing submatrices defined by the row and column clusters that consist of either all ones (or all zeros) to the greatest extent possible. Essentially, the DTMBP seeks to find partitions minimizing the inconsistency with this structure, as measured by the total number of ones in submatrices containing mostly zeros, plus the number of zeros in submatrices containing mostly ones. An exact algorithm for the DTMBP was developed by Brusco et al. (2013b), but is limited to relatively small two-mode matrices. A relocation heuristic RH for the DTMBP was described by Doreian et al. (2004). The RH is often implemented in a multistart framework, whereby the heuristic is applied to large numbers of initial randomly-generated partitions. However, the RH has also been embedded within the framework of more sophisticated heuristic procedures, such as tabu search (Brusco & Steinley, 2011) and variable neighborhood search (Brusco et al., 2013a).

One of the limitations of the RH is that its computational requirement increases significantly as the number of objects (n and m) and the number of clusters (K and L) increase. Here, we present an alternative heuristic algorithm, one dramatically more efficient than the RH. We refer to this algorithm as a two-mode KL-median heuristic (TMKLMedH). The TMKLMedH seeks to partition the row and column objects with the goal of minimizing the sum of the absolute deviations between matrix elements and the median of the elements in their submatrices. The TMKLMedH is an adaptation of a two-mode KL-means heuristic that has been popular in both social network and broader classification literature (Baier, Gaul, & Schader, 1997; Brusco & Doreian, 2015a, 2015b; Brusco & Steinley, 2007a; Gaul & Schader, 1996; Van Rosmalen, Groenen, Trejos, & Castillo, 2009; Vichi, 2001). Although TMKLMedH is not restricted to binary network matrices, our emphasis here is limited to such networks. More importantly, it generalizes naturally to valued networks, something much harder to accomplish with RH.

The next section presents the motivation for the intended contribution of our paper. This is followed by a section providing a precise formulation of the DTMBP and a section that describes two competitive algorithms for its solution: RH and TMKLMedH. A simulation-based computational comparison of RH and TMKLMedH is subsequently provided, followed by a comparison of RH and TMKLMedH for voting data from the UNGA (Brusco et al., 2013a) and a larger network corresponding to movie ratings. The paper concludes with a summary of our findings and an identification of extensions for future research.

There are two primary sources of motivation for this paper: (i) our discovery that TMKLMedH is a direct competitor of RH for solving the DTMBP; and (ii) previous findings concerning the relationships among different implementation of K-means clustering algorithms (Forgy, 1965; MacQueen, 1967; Steinhaus 1956; Steinley, 2006a). As noted previously, two-mode KL-means partitioning has received some attention in the social network literature (Brusco & Doreian, 2015a, 2015b). Pursuant to this work, we discovered that, in the case of binary networks, if the mean of submatrices is replaced with the corresponding medians, then the resulting optimization problem has the same criterion function as DTMBP, a point developed in more detail in the next section.

There are a variety of different implementations of K-means clustering (Brusco & Steinley, 2007b; Hansen & Mladenović, 2001; Späth, 1980; Steinley, 2006a, b). The most popular one begins with a partition of objects into K clusters. An iterative process consisting of the computation of cluster centroids and the reassignment of each case to its nearest centroid is then repeated until convergence. Although this procedure is commonly referred to as K-means clustering, some authors refer to it as H-means clustering (see Brusco & Steinley, 2007b; Hansen & Mladenović, 2001; Späth, 1980). Another variant that is sometimes referred to as K-means clustering is another relocation heuristic, where each object is tested in turn for reassignment from its current cluster to each of the other clusters. This algorithm, called K-means by Hansen and Mladenović (2001), converges when no object can be relocated to further improve the criterion function.

The advantage of H-means relative to K-means is that it is much faster, allowing many more restarts if the two algorithms are allotted the same amount of time. The disadvantage of H-means relative to K-means is that it does not necessarily converge to a locally-optimal solution regarding all possible object relocations. In fact, H-means solutions can even be degenerate in the form of having either empty clusters or two or more clusters with identical centroid values. To capitalize on the strengths of both methods, Hansen and Mladenović (2001) devised a procedure called HK-means, where H-means is used first to rapidly obtain a solution and K-means is then applied to refine that solution so as to improve the criterion function value. Brusco and Steinley (2007b) showed HK-means performed exceptionally well compared to a variety of different heuristic approaches.

The relationship between TMKLMedH and RH is directly analogous to the relationship between H-means and K-means. This suggested that TMKLMedH was likely to be a good alternative to RH, as it can allow for many more restarts within the same prescribed time limit. Computational studies were undertaken to investigate this claim.

To summarize, the contributions of our paper are:

We use the following notation to formulate the DTMBP:

Using this notation, a formulation of the DTMBP follows: (See also Brusco et al., 2013; Brusco & Steinley, 2007a, 2011).

where:

and

The goal of the DTMBP is finding the row-object (π) and column-object (ω) partitions minimizing the number of inconsistencies with perfect structural equivalence, whereby all blocks are either null or complete. The values of λ_kl and ρ_kl are, respectively, the number of ones and zeros in the submatrix corresponding to row cluster k and column cluster l. For each submatrix, if λ_kl ≥ ρ_kl, then the ρ_kl zeros would be the inconsistencies in the submatrix. This count is accumulated in the objective function. Likewise, if λ_kl < ρ_kl, then the λ_kl ones would be the inconsistencies in the submatrix. This count is also accumulated in the objective function, which is the sum of these inconsistency counts over all K × L blocks as expressed in equation (1).

Next, we propose an alternative formulation of the DTMBP using submatrix medians. The following additional terms are necessary:

d_kl := the median of the elements in the submatrix (or block) formed by the row objects in cluster S_k and the column objects in cluster T_l, for all 1 ≤ k ≤ K and 1 ≤ l ≤ L;

With these definitions in place, an alternative formulation of the DTMBP based on submatrix medians is as follows:

Establishing the equivalence of g(π, ω) in equation (1) and f(π, ω) in equation (5) is straightforward. For a binary network matrix, the median of the submatrix defined by row cluster k and column cluster l (∀ 1 ≤ k ≤ K and 1 ≤ l ≤ L) will be d_kl = 1 if there are more ones than zeros in the submatrix, d_kl = 0 if there are more zeros than ones in the submatrix, and d_kl = 0.5 if there are an equal number of zeros and ones in the submatrix. Furthermore, (∀ 1 ≤ k ≤ K and 1 ≤ l ≤ L), v_kl is equal to the number of zeros in the submatrix when d_kl = 1, v_kl is equal to the number of ones in the submatrix when d_kl = 0, and v_kl is equal to the number of ones (or, equivalently, the number of zeros) in the submatrix when d_kl = 0.5. It follows that, for binary data, the v_kl values reflect the number of inconsistencies in the submatrices with respect to the ideal condition under structural equivalence. That is, if d_kl = 1, then v_kl = ρ_kl, the sum of the zeros in the submatrix. Similarly, if d_kl = 0, then v_kl = λ_kl, the sum of the ones in the submatrix. Finally, if if d_kl = 0.5, then v_kl = ρ_kl = λ_kl, the sum of the zeros (or ones) in the submatrix. In summary, equivalence of the two representations stems from the fact that the value of v_kl = min{λ_kl, ρ_kl,}.

Our first example for demonstrating the TMKLMedH uses United Nations General Assembly voting data, originally compiled by Doreian et al. (2013). Two network matrices were analyzed. The first was a 141 × 276 binary network matrix pertaining to votes of 141 countries on 276 ideological resolutions. The second was a 153 × 129 matrix for the votes of 153 countries on 129 military resolutions. Brusco et al. (2013a) performed a comparison of the performance of three two-mode blockmodeling algorithms on these two networks: (i) the relocation heuristic (RH) described in Section 3; (ii) a tabu search (TS) heuristic developed by Brusco and Steinley (2011); and (iii) a variable neighborhood search (VNS) heuristic developed by Brusco et al. (2013a) for two-mode blockmodeling. Each algorithm was applied to these network matrices for all numbers of clusters in the intervals 2 ≤ K ≤ 7 and 2 ≤ L ≤ 7. For each combination of K and L, the algorithms were constrained to a time limit of 10 CPU minutes on a 2.2 GHz Pentium 4 PC with 1 GB of RAM. There was little difference in the performance of the three methods for problems where K < 4 and/or L < 4. This suggests few differences for smaller t-mode binary networks.

For comparison to the results obtained by Brusco et al. (2013a), we ran the TMKLMedH heuristic for all combinations of 4 ≤ K ≤ 7 and 4 ≤ L ≤ 7 on both networks using the same 10-minute time limit on the same 2.2 GHz Pentium 4 PC with 1 GB of RAM. These ranges for K and L resulted in 16 test problems for each of the military and ideological networks, and a total of 32 problems overall. A comparison of the results for the four methods is displayed in Table 2.

The results in Table 2 are quite surprising. Across the 32 test problems, TMKLMedH matched the best-found criterion function value for 31 of the 32 test problems. By contrast, the TS, VNS, and RH methods only matched the best-found criterion function value for 24, 24, and 17 of the test problems, respectively. It is also noteworthy that, for three test problems (military resolution network with K = 7 and L = 4, ideological network with K = 4 and L = 7, and ideological network with K = 7 and L = 5), TMKLMedH obtained a new best-found criterion function value never matched by any of the other three methods.

We also collected, for each test problem, the elapsed time for TMKLMedH to find its best-found solution. For 22 of the 32 problems, TMKLMedH found its best solution within the first minute of its allotted 10-minute time limit. Moreover, there were only three problems where the best-found solution was not realized until after five minutes of elapsed time.

To summarize, the results using the UNGA data strongly support the simulation study. The TMKLMedH yielded a better criterion function value than RH for 15 problems, and the same criterion function value for the remaining 17 problems. The RH never produced a better result than TMKLMedH within the allotted time. Although the performance of TMKLMedH relative to TS and VNS is very favorable, considerable care is warranted. The TS and VNS procedures use relocation heuristics in their implementation and, therefore, are significantly slower than TMKLMedH. If the TS and VNS algorithms were redesigned to use the same process of finding medians and reassigning cases, perhaps their performances could improve substantially.

The second comparison of RH and TMKLMedH for a real-world network was undertaken using movie-rating data from the MovieLens database (Harper & Konstan, 2015). This data was retrieved from the KONECT (the Koblenz Network Collection: http://konect.uni-koblenz.de/networks/) website. Although the edges of the available network are weighted to reflect the ratings, our binary adaptation consists of n = 943 individuals who either did (x_ij = 1) or did not (x_ij = 0) provide ratings for each of m = 1682 movies. This network is less dense (density of 6.3%) than the UNGA networks.

We ran the RH and TMKLMedH algorithms for all combinations of 2 ≤ K ≤ 7 row clusters and 2 ≤ L ≤ 7 column clusters. Again, the algorithms were implemented on the desktop computer with an Intel ® Core^TM i7-4790 CPU @ 3.6GHz with 8 GB of RAM, but with a 10-minute time limit for all combinations of K and L. The results are summarized in Table 3.

Table reveals that TMKLMedH produced a better criterion function value than RH for 30 of the 36 combinations of K and L. The two algorithms obtained the same criterion function value for the remaining six combinations. Across the 36 combinations of K and L, the average percentage improvement in the criterion function value resulting from the use of TMKLMedH instead of RH was .267%, and the maximum improvement was 2.309%. Again, the relative superiority of TMKLMedH stemmed from the fact that it accumulated far more restarts within the 10-minute time limit. Across the 36 test problems, the minimum, average, and maximum ratios of the number of TMKLMedH restarts to the number of RH restarts were 1.2, 20.7, and 35.6, respectively.

We recognized that the small percentage differences in the criterion values obtained by TMKLMedH and RH might not necessarily translate into meaningful differences with respect to blockmodel interpretation. The first step in this process was to select appropriate values for K and L to facilitate a comparison. Several criteria for selecting K and L have been proposed in the literature, such as examination of the number of local optima (Brusco et al., 2013b) and the Chull procedure (Brusco, Doreian, & Steinley, 2016; Ceulemans & Van Mechelen, 2005; Wilderjans, Ceulemans, & Meers, 2013). In this application, we use a visual aid for choosing K and L that is based on the use of heatmaps for the criterion values obtained by RH and TMKLMedH at different combinations of K and L.

Heatmaps for the criterion function values obtained by RH and TMKLMedH are displayed in Figure 1. These heatmaps were produced using the conditional formatting capabilities in Microsoft Excel using a red-yellow-green color scheme. Deep red values reflect the poorest of the criterion function values obtained, whereas deep greens signify the finest of the criterion values obtained. The progression in the heatmaps from red to reddish-yellow to yellow (and from yellow to greenish-yellow to green) visually illustrates improvement in the criterion function.

Figure 1.

An examination of the heatmaps in Figure 1 reveals some similarities and differences. Both heatmaps show deep red cells in the first row (K = 2) and first column (L = 2). This reflects the fact that very little (if any) improvement can be achieved by increasing K when L = 2, or by increasing L when K = 2. In a similar manner, both heatmaps show a deep yellow in the row for K = 3 (when L ≥ 3), as well as the column for L = 3 (when K ≥ 3). For conditions where K ≥ 4 and L ≥ 4, both heatmaps show shades of green, but to somewhat different extents. For example, a fairly sharp green begins to appear in the TMKLMedH heatmap at K = 5 and L = 5, and we used this as motivation for selecting the K = 5 and L = 5 solution for interpretation. By contrast, the green at K = 5 and L = 5 is not as sharp for the RH heatmap. Thus, if an analyst were looking at only the RH heatmap, a different K and L combination might be selected. This possibility of a different choice for K and L depending on whether RH or TMKLMedH is used is not unique to heatmaps, but could also occur for criteria such as the number of local optima or the rules used in the Chull procedure.

The criterion function values for RH and TMKLMedH at K = 5 and L = 5 were 87146 and 86498, a difference of less than one percent. However, for this same combination of K and L, the agreement between the RH and TMKLMedH partitions of the n = 943 individuals rating the movies was only ARI = .382. Similarly, the agreement between the partitions of the m = 1682 movies that were rated was ARI = .571. Both ARI measures fall to far below the 0.65 threshold for even fair agreement as recommended by Steinley (2004). Figure 2 helps to further clarify the distinctions between the RH and TMKLMedH blockmodels by providing their image matrices and cluster sizes. The image matrix in Figure 2 uses the term complete to signify blocks with greater than 50% density and null to indicate those with less than 50% density.

Figure 2.

The most striking similarity of the blockmodels obtained by RH and TMKLMedH is that each has one large cluster of individuals and one large cluster of movies that are loosely tied to the network (i.e., not associated with any complete blocks). However, once these two large clusters are accounted for, it is evident that TMKLMedH does a better job of identifying relationships between the remaining individuals and movies. The clean, interpretable structure of TMKLMedH is evidenced by the presence of all complete blocks below the main diagonal of its image matrix. This means that the five clusters of individuals associated with the TMKLMedH blockmodel were strongly tied to 0, 1, 2, 3, and 4 clusters of movies, respectively. Likewise, the five clusters of movies associated with the TMKLMedH blockmodel were strongly tied to 0, 1, 2, 3, and 4 clusters of individuals, respectively. The RH blockmodel does not have this structure. For example, TMKLMedH identified a cluster of n ₅ = 32 individuals who form complete blocks with the four latter clusters of movies. By contrast, RH did not identify such a cluster of individuals. In a similar manner, TMKLMedH obtained a cluster of m ₅ = 17 movies that form complete blocks with the four latter clusters of individuals. Again, RH did not identify such a cluster of movies.