extfieldsofstudy:
Physics
Computer Science
Engineering

text:
Scalable Node-Disjoint and Edge-Disjoint Multi-wavelength Routing

Probabilistic message-passing algorithms are developed for routing transmissions in multi-wavelength optical communication networks, under node and edge-disjoint routing constraints and for various objective functions. Global routing optimization is a hard computational task on its own but is made much more difficult under the node/edge-disjoint constraints and in the presence of multiple wavelengths, a problem which dominates routing efficiency in real optical communication networks that carry most of the world's Internet traffic. The scalable principled method we have developed is exact on trees but provides good approximate solutions on locally tree-like graphs. It accommodates a variety of objective functions that correspond to low latency, load balancing and consolidation of routes, and can be easily extended to include heterogeneous signal-to-noise values on edges and a restriction on the available wavelengths per edge. It can be used for routing and managing transmissions on existing topologies as well as for designing and modifying optical communication networks. Additionally, it provides the tool for settling an open and much debated question on the merit of wavelength-switching nodes and the added capabilities they provide. The methods have been tested on generated networks such as random-regular, Erd\H{o}s R\'{e}nyi and power-law graphs, as well as on the UK and US optical communication networks. They show excellent performance with respect to existing methodology on small networks and have been scaled up to network sizes that are beyond the reach of most existing algorithms.

I. INTRODUCTION
Optical communication networks underpin the global digital communications infrastructure and carry most of the Internet traffic. They comprise thousands of kilometers of optical fibers, organized in a complex web of constituent sub-networks including the Internet backbone, Metro access and Data Center farms. The exponential growth in Internet traffic and energy consumption threatens to overload the existing infrastructure and a capacity crunch is looming [1]. Not only that a matching growth in infrastructure is infeasible, it raises fundamental questions on the ultimate capacity of optical communication networks and the manner in which we could optimize their use. The next-generation digital infrastructure has to offer flexibility, low latency, high network throughput and resilience.
One of the key requirements is the routing and wavelength assignment (RWA) for all traffic demands across this complex heterogeneous network in a way that optimizes a given objective function, be it low latency, high throughput or resilience [2]. Each optical fiber carries information using light of one or many wavelengths. The latter uses, among others, dense wavelength-division multiplexing (DWDM) methods that employ as many as 80-160 channels of different laser wavelengths [3]. The main constraint in the RWA is that any complete individual route, from source to destination, uses the same single wavelength and that two separate routes using the same wavelength cannot shares the same fiber. This constraint makes the corresponding mathematical problem hard to solve in general.
Route optimization in optical communication networks can be mapped onto the hard computational problem of edge-disjoint routing on a graph, where transceivers (transmitter-receiver) are mapped to vertices (or nodes) and fibers to edges (or links). Given that routes are constrained to be contiguous and interaction between paths is non-localized, local optimization methods are insufficient and global optimization is required. Globally optimal routing of multiple messages or vehicles given a general objective function is a computationallyhard constraint satisfaction problems on its own and has been addressed in the physics literature using scalable and distributed message passing approximation techniques, inspired by statistical physics methodology [4][5][6][7]. Moreover, similar techniques have been suggested also for addressing the single-wavelength Node-Disjoint Paths (NDP) [8] and Edge-Disjoint Paths (EDP) [9] problems where multiple paths of different origin-destination pairs cannot share nodes or edges on a graph, respectively. Generally, both optimization tasks are within the class of NP-hard combinatorial problems [10][11][12][13] and the approximation offered by message passing techniques work well. However, the existing methods developed for singlewavelength routing become intractable in the presence of multiple wavelengths, making them inapplicable for realistic scenarios.
It is worthwhile noting that in some extreme cases these hard computational problems become polynomial in the single wavelength case as discussed in [9]. For instance, when the number of origin-destination pairs is low with respect to the systems size [14] and where all origindestination pairs are identical [15,16]. Nevertheless, the general problems of NDP and EDP routing on graphs are computationally hard even in the single wavelength case and a variety of methods have been used to address them in the context of optical communication networks and more general problems. Alongside established methods, such as integer/linear programming and its variants [17][18][19][20][21] bin-packing based approaches [22], Monte Carlo search [23], post-optimization [24] and greedy algorithms [25][26][27][28], a large number of heuristics have also been used to obtain approximate solutions in both problems, among them genetic algorithms [29][30][31], ant colony optimization [32] and particle swarm optimization [33]. Specifically in the area of optical communication networks, it is common to use integer/linear programming and its variants for small networks, to obtain exact results, and a variety of heuristics for larger networks. In practice, current optical networks use overprovision of capacity to compensate for sub-optimal routing, resulting in both over-engineering and underutilized capacity [34,35].
The main challenge we address here is the RWA under heavy traffic using multiple wavelengths and a very large number of origin-destination pairs under the NDP and EDP constraints and for various objective functions. Globally optimal routing is a non-localized difficult problem but adding the NDP/EDP restrictions and having a large number of different wavelengths increases the complexity considerably, making the problem intractable for existing algorithms [8,9]. We map the globally-optimal routing problem in the presence of multiple wavelengths onto a multi-layer replica of the original graph and utilize probabilistic optimization approaches. The methods developed here are based on message passing techniques, developed independently in several fields including statistical physics, computer science and information theory [36][37][38] but are closely interlinked [39,40]. The methods we develop allow for messages, in the form of conditional probability values to be passed between nodes and the replicated networks representing the different wavelengths, in a way that keeps the algorithms scalable and applicable even for a large number of wavelengths, transmissions (corresponding to source-destination pairs) and nodes.
The main result of this paper is the derivation of principled scalable algorithms, capable of obtaining approximate solutions for routing problems in large graphs, where the number of transmissions is of similar order to that of the number of free variables (vertices/edges) and a large number of wavelengths, under the NDP and EDP constraints and for various objective functions, both convex and concave. The algorithm also accommodate cases where the number of transmissions is much larger than the number of vertices (quadratic with respect to the number of vertices). The computational complexity of the NDP/EDP algorithms for sparse graphs is O(M Q(M + N + Q))/O(M Q(M/N + N + Q)), with N the number of vertices, M the number of transmissions and Q the number of wavelengths. The algorithm has been tested for a variety of sparse network topologies, both synthetic random graphs and real optical commu-nication networks, and for different objective functions, showing excellent results in obtaining high quality approximate solutions. Among the generic networks examined are random regular graphs, Erdős Rényi (ER) [41] and scale free networks [42], while the realistic networks considered include the British 22 nodes (BT22) and US 60 node (CONUS) [43]backbone optical communication networks. When tested on small networks against known results obtainable using unscalable methods like variants of integer/linear programming, it was shown to provide the optimal routing results.
The results provide the maximal number of communication pairs that could be accommodated given the network size, topology and number of wavelengths used; the minimal number of wavelengths required for a given network, topology and communication pairs; and the resulting utilization of edges. They identify the impact of using the suggested algorithm on the average path length, the utilization of wavelengths per edge and how it can be controlled using concave and convex objective functions.
In addition, our algorithms can be used for routing transmissions across networks in single instances, study the limitations of heterogeneous networks of different degree distributions, with variable edge signal-to-noise ratios and wavelength availability. Moreover, our algorithms could be employed in the design of new infrastructure, especially through the use of concave cost functions that consolidate routes, by determining the least important routes that could be removed with little effect on the network throughput or resilience. These are of both academic and practical values since the performance of optical communication networks is often directly related to their capacity limits, traffic congestion, rate of information flow and bandwidth flexibility. Moreover, we also studied a switching model, where wavelength can be converted (switched) at the vertices (transceivers) to settle an open and much debated question on the merit of wavelength-converters [44][45][46] for increasing throughput, resource utilization and resilience in optical communication networks.
While we mainly focus here on the optical communication network application and test the efficacy of the method on networks and number of wavelengths that are relevant to this application domain, one should point out that these problems are highly relevant to other domains. For instance, multi-wavelength NDP/EDP are relevant to both 5G and the future 6G wireless communication systems and wireless ad-hoc communication networks in the relay setting, where each node can act as a relay. Our algorithms can reduce path overlaps, which represent transmissions in similar wavelength, resulting in signal interference and low transmission quality, or to consolidate paths since longer paths result in signal degradation and the need for higher transmission power [26,27,47]. Another application is the design of very large system multilayer integrated circuits (VLSI), where non-overlapping wired paths to connect different components are sought to avoid cross-path interference.
In all cases, higher throughput, robustness and lower latency can be achieved for the same resource by obtaining a good approximation to the globally optimal solution. Practical algorithms for various applications often depend on the specific network topologies considered [48] and are aimed at maximizing the number of paths routed [49].
The reminder of the paper is organized as follows: in Sec. II we introduce the model used followed by the message-passing based algorithmic solutions for NDP, EDP and wavelength switching scenarios in Sec. III. Results obtained from numerical studies on a range of synthetic and real networks and a variety of objective functions are presented in Sec. IV followed by a discussion on their computational complexity. Possible extensions of the framework to accommodate real-world scenarios such as edges with different signal-to-noise ratios or wavelength availability are presented in Sec. V. Finally, we discuss the efficacy of the methods developed and point to future research directions in Sec. VI.

II. MODEL
We consider a dense wavelength-division multiplexing (DWDM) optical network G(V, E), with V ≡ {i | i ∈ G} and |V | = N , the set of nodes representing transceivers and E ≡ {(i, j) | (i, j) ∈ G} the set of edges, such that the indices (i, j) represent the optical fiber between node i and j. For a network which uses Q wavelength channels to deliver M transmissions, we introduce a variable s i,j on the link from node i to j such that s a i,j = s or −s if transmission s passes from node i to j or from j to i respectively, through link (i, j) using wavelength a; the transmission s corresponds to one of the origindestination pair {0, 1, . . . , M }. A similar variable s a i is defined for node i in the case of NDP, as explained in Sec. III. Since each transmission has to occupy an individual wavelength channel on a link, more wavelengths are generally required for more transmissions, i.e. a larger value of Q is required for a larger M . For specific network instances with Q wavelength channels, there exists a maximum number of transmissions denoted as M max which can be transmitted; alternatively, one can define the minimum number of wavelength channels, i.e. Q min , which accommodate all M transmissions on a specific instance. The relationship between M max and Q, or between M and Q min , would be highly relevant for characterizing the maximum capacity of optical networks.
This framework can accommodate a variety of objective functions; here, we consider the sum of the cost (or utility) on each link to be the objective function for optimization, given by where δ y x is the Kronecker delta such that δ y x = 1 if x = y and δ y x = 0 otherwise; the function F i,j denotes the cost on link (i, j) as a function of the argument in parenthesis; while it can take an arbitrary form we will mostly focus on simple polynomial functions. In the context of statistical physics, we introduce the inverse temperature β, and the partition function Z of the system is given by where Ω( s) is an indicator function such that Ω( s) = 1 if s satisfies all the constraints of the problem or otherwise Ω( s) = 0. The double vector notation comes to emphasize dependence on both topology and wavelength.
We now summarize the constraints of the optimization problem. Firstly, the route for each transmission must be contiguous. A loopless path is a sequence of nonrepeating nodes from origin to destination, for example, path s is constructed as O s → · · · → j → i → k → · · · → D s , where O s and D s are the origin and destination pair of transmission s. Secondly, for an intermediate node i along the path there exist only two used edges (j, i) and (i, k) for that wavelength and transmission. This constraint could be expressed as follows: if s a i,j = s = 0, then where ∂i ≡ {j | (i, j) ∈ E} is the set of the nearest neighbors of node i, and ∂i \ j ≡ ∂i − {j} is the subset of ∂i except node j. Regarding the wavelength constraints along the path, we will consider three scenarios in the subsequent analyses: 1. Node-disjoint paths (NDP) -where only a single transmission is allowed to utilize a specific wavelength channel on a node [8], but there can be multiple transmissions using different channels through the same node as shown in Fig. 1(a); this may correspond to transceivers which can only process a single transmission for each individual wavelength. The expression of the node-disjoint constraint is s a i,j = s = 0 ⇒ ∀k ∈ ∂i \ j : s a i,k ∈ {0, −s}. 2. Wavelength-switching (WS) with NDP -where all incoming or outgoing transmissions to transceivers (nodes) use different wavelength channels, but transmissions are allowed to switch between wavelength channels at the transceivers as shown in Fig. 1(b), leading to a larger routing flexibility.
3. Edge-disjoint paths (EDP) -where multiple transmissions using the same wavelength channel, are allowed to be routed through any given node [9]but cannot share an edge, as shown in Fig. 1(c). This is the typical scenario in optical communication networks.

III. MESSAGE-PASSING ALGORITHMS
To derive optimization algorithms to allocate simultaneously the optimal route and wavelength for a large number of transmissions, we solve the problem on multilayer graphs where each layer represents a different wavelength, and messages are passed within each layer for assignment of routes and between layers for allocation of wavelength. By applying the cavity approach from the study of spin glass systems [36,50], we can derive distributed message-passing algorithms for optimizing transmission routes in optical networks for NDP, WS and EDP scenarios.

A. Node-disjoint routing
First, we consider the NDP routing scenario [8] shown in Fig. 1(a). The network is represented by a factor graph with two types of variables, namely s a i defined on nodes and s a i,j defined on links. According to the Bethe-Peierls approximation [51,52], we assume only large loops exist in the network such that all neighboring nodes and edges of a node i are nearly independent in the absence of i.
To derive the multi-wavelength routing algorithm, we first define p i→j to be the message from node i to edge (i, j) and q i→j to be the message from edge (i, j) to node j. Both messages are conditional probabilities: p i→j (s) is the probability of edge (i, j) to be in state s due to the state of node i and q i→j (s) the probability of edge (i, j) to be in state s without the interaction from node i. We then write a closed set of recursion relations to express the message p i→j in terms of q k→i , as well as another set of relations representing q i→j in terms of p i→j and p j→i : where we have further defined the auxiliary quantities q b i,j given by For brevity we omit the normalization term and use the notation "∼" instead of the equality symbol. For a brief explanation of how these equations are construed: the first equation of (4), looks at the probability for no transmission on the edge i → j, as a summation of the probability of no transmission entering i (first term) and the probability of transmissions arriving at i from node m but leaving through some other node n (second term); the second equation looks at the message s arriving at node i but not leaving through any other edge but i → j; the third and fourth equations describe the probability of edge i → j to be in a given state 0 or s, given the related cost on the edge, in conjunction with all other wavelengths. In Eq. (4), node messages using wavelength a depend only on the messages from their neighbors using the same wavelength and the network is effectively mapped to a system with Q separate layers, each of which employs a different wavelength channel, as shown in Fig. 1(a); the interdependence between wavelengths is considered at the origin and destination of individual transmissions ( Fig. 2), as discussed below.
With the messages p i→j and q i→j having converged to stable values, one can express the marginal probability of node i being in state s using wavelength a as The marginal probability of edge (i, j) being in state s with wavelength a is given by where the factors after the dot symbols correspond to contributions from all neighboring links of (i, j) to the partition function, which takes into account all possible variable configurations in the trees terminated at the neighboring edges.
To simplify the algorithms, we introduce the variables φ ≡ − 1 β log q and ψ ≡ − 1 β log p, and take the zerotemperature limit for optimization β → ∞. The recursion relations of messages in Eq. (4) thus become equivalent to the min-sum belief propagation relations, given For each individual transmission, we then introduce auxiliary nodes labeled as µ = ±|µ|, which connect to the origin and destination nodes in each of the Q wavelength network layers, respectively (see in Fig. 2). These auxiliary origin-destination pairs of each transmission communicate with all wavelength replica networks to determine the path and wavelength channel allocated to each transmission, using the message passing algorithm, given by such that the transmission routes are determined independently on each wavelength channel. After introducing the auxiliary nodes, we treat the messages to and from them as to other network nodes. For example, the neighboring nodes set of node i in Fig. 2 is ∂i = {j, k, l, µ}, and the calculation of messages i → j and i → µ both follow Eq. (8). The marginal probability of edge (i, j) being in state s using wavelength a is then given by and the state s a i,j of wavelength a of edge (i, j) is deter-mined by which ultimately leads to the optimized configuration of routes for all M transmissions through the optical network with Q wavelength channels.

Linear cost
While the objective function can take many different forms we use a simple power as the cost function F i,j (x) on edges x γ [4,53], as it provides a good example of both concave and convex costs. When γ = 1, the cost function is linear and equivalent to the total length of all transmissions L ≡ (i,j);a 1 − δ 0 s a i,j , and optimizing it is equivalent to finding the shortest average path.
In this case (γ = 1), the message-passing equations Eq. (8) can be simplified to where and the variable ψ in Eq. (8) can be omitted. The marginal quantities are then given by and the messages to the auxiliary nodes at the origins and destinations are the same as Eq. (10).

Switching wavelength channels at nodes
A generalization of the NDP scenario is to allow for transmissions to change wavelengths at nodes (WS) as shown in Fig. 1(b), which is equivalent to a singlewavelength network routing with both node and edge capacity being Q, and utilize other multiplexing techniques than wavelength division (e.g. code or time division multiplexing). A variable τ is introduced to enforce this capacity constraint on nodes. In this case, the message-passing equations are given by where τ µ i denotes the summation overψ i from a total of at most Q transmissions, excluding µ, that pass through node i, given by Specifically, τ µ i (1) could be understood as indicating that there are free wavelength channels on i that µ could take, and τ µ i (0) indicates that transmission µ does not pass through node i. The quantityψ µ i is an auxiliary variable related to the probability that node i would be a part of the path of transmission µ, which leads to The messages from the origin and the destination of transmission µ to the neighboring edges are slightly different from the previous case of NDP without switching. In Eq. (16) Q i = Q, andψ i (1) = 0 andψ(0) = ∞ for the transmissions using node i as their origin or destination. Alternatively, one can set Q i = Q − M i , where M i is the number of transmissions using node i as their origin or destination, and then exclude these transmissions when one calculates τ i . Here, we first introduce a variable µ i = ±1 to denote node i as the origin or destination of transmission µ respectively, given by The variable φ µ i,j determines the state of edge (i, j), i.e. whether it is a part of the path of transmission µ or not, and its expression is given by such that the optimized state is σ µ i,j = arg min σ φ µ i,j (σ).

B. Edge-disjoint routing
Edge-disjoint path routing (EDP) is similar to NDP routing but is less restrictive, since nodes can accommodate any number of paths with the same wavelength but edges do not [9] (Fig. 1(c)). In other words, for node i, if there exists an edge (i, j) with s a i,j = s 0 = 0, then there exists one and only one edge (i, k) with s a i,k = −s 0 continuity of transmission path, while other neighboring edges using the same wavelength channel a could be either in state 0 or take up other transmissions such that s a i,m1 = −s a i,n1 = s 1 = 0, s a i,m2 = −s a i,n2 = s 2 = 0 etc. In comparison, for NDP at most two variables s a i,j for neighboring nodes j ∈ ∂i can assume a non-zero value. In other words, single-wavelength EDP is equivalent to a generalized version of WS, where the capacity of each node is non-uniform and determined by its degree and the number of it being chosen as origins or destinations of transmissions. Consequently, the message passing equations for EDP scenarios are more complicated and their computational complexity is higher. The corresponding message passing equations are given by where the message from node i to edge (i, j) (first two equations) are different from those of NDP in Eq. (8), while Eq. (9)- (12) in the NDP scenario apply also here. It is numerically difficult to compute the message passing relations in Eq. (20), since all feasible configurations satisfying the constraints have to be considered by the first two equations, which results in a computational complexity of where K = |∂i| − 1 is the number of neighboring edges of node i, M is the total number of transmissions, x is the floor function that is equal to the greatest integer less than or equal to x, n P k = n! (n−k)! is the number of ordered permutations of k out of n elements and n C k = n! k!(n−k)! the number of unordered combinations. To simplify the computation we map the task on to the maximum weight matching problem [9,54]. As in [9], we consider all the possible pairs of considered edges {(i, k) | k ∈ ∂i \ j} and obtain the weight of each pair Then, we construct a weighted graph with |∂i| − 1 nodes where the weight of each edge is given by Eq. (22). After that, ψ a i→j in Eq. (20) could be obtained by the maximum weight matching algorithm [9,55].
In Fig. 3, we show the transformation employed for calculating the messages i → j. Links between all pairs, including auxiliary origin-destination and graph edges, are shown on the left. The right figure shows all valid links between nodes, while distinguishing between auxiliary and graph nodes: any of the ordinary graph nodes could be paired but the auxiliary nodes could only be paired with ordinary graph nodes. The maximum matching algorithm finds edge sets with the maximum sum of weights, where matched edges have no common nodes; the inverse value of the weights sum is approximately the value of minimum matched configurations in Eq. (20).

Linear cost
Similar to Sec. III A 1, if the cost function on edges is linear F i,j (x) = x, then the message-passing equations could be simplified as follows, The marginal states of nodes and edges, and messages to the auxiliary nodes at transmission origins and destinations, are the same as in the NDP case of Eq. (10) and (14).

C. Algorithmic framework
The same algorithmic procedure governs both NDP and EDP routing algorithms based on the message passing equations. The min-sum algorithm, i.e. zerotemperature optimization algorithm, follows the process outlined below: 1. initialize the messages {ψ a i→j , φ a i→j , φ µ→a }; 2. update the messages by Eq. (8)-(10) in the nodedisjoint scenarios, and Eq. (20) and (10) in the edge-disjoint scenarios until convergence or when a maximum number of iteration steps is reached; 3. calculate the marginal state of edges by Eq. (11) and the state of each wavelength channel by Eq. (12).
After obtaining the marginal states of edges and wavelengths, the allocation of transmissions to specific wavelength and path would follow. However, there are some cases when the algorithms need many iterations to provide a valid configuration, for example if the cost function F i,j (x) in Eq. (1) is concave, e.g. F (x) = √ x in these cases, a decimation procedure can be introduced to speed up convergence. As an example, here we show the NDP algorithm with decimation (fixing states at intermediate steps) incorporated as part of the process: 1. initialize the messages {ψ a i→j , φ a i→j , φ µ→a } and the state of all wavelength channels per available (undecimated) edge; 2. update the messages of the available wavelength channels and edges by Eq. (8)-(10) for a specific number of iteration steps; 3. compute the marginal states of the available wavelength channels and edges by Eq. (11); calculate the quantityφ a i,j = min s =0 φ a i,j (s) − φ a i,j (0), and fix the state of the wavelength channel a of edge (i, j) with the largest value ofφ a i,j to be s a i,j = 0.
4. iterate step 2 and 3 and determine the state of the available wavelength channels and edges by Eq. (12) until a valid solution is obtained or a maximum number of iteration steps is reached.
The algorithms with or without the decimation procedures have a similar performance in terms of the optimized cost, but for harder problem (e.g. F (x) = √ x) the decimation procedures can reduce the number of iteration steps needed for generating a valid routing solution.

A. Linear cost
We first examine the NDP, WS and EDP scenarios with linear cost F (x) = x to explore the behavior and performance of communication networks, such as capacity and average length, with impact on latency and number of wavelength channels required.

Dependence on the number of wavelength channels
We performed numerical experiments using the three algorithms on three different types of generated random networks including random regular graphs, Erdős-Rényi and scale-free networks with 100 nodes and an average node degree of 3, as well as on two real optical communication networks -CONUS60, which has 60 nodes and 79 edges [43] and BT-Core, 22 nodes and 35 edges [56] as shown in Fig. 4. Specifically, for the scale-free networks studied, the degree distribution is p(d) ∼ d −2.63 .  We define the capacity to be the maximum number of transmissions, which we denote as M max , that can be transmitted by an optical communication network with Q wavelength channels. We remark that capacity depends on network topologies and the set of origin-destination pairs, so in Fig. 5 we report the dependence of the average value of M max with standard deviation as error-bars on Q, showing the average behavior of capacities for different Q values.   5 shows the simulation results for the five types of networks -in the three different routing scenarios. For the three types of random networks and BT-Core, the average capacity per wavelength channel M max /Q keeps increasing at the beginning and become saturated as the number of wavelength channels Q increases.As for the networks CONUS and BT-Core in Fig. 5, we see that the average capacities also increase with Q but the increases are not as fast as that in the three generated networks.
The average path length of the corresponding networks are shown in Fig. 6 and 7. In Fig. 6, we see that as more wavelengths are available the shorter the average path length L/M decreases; for instance, in NDP scenarios shown in Fig. 6(a)-6(c), the average path length with Q = 1 is higher than that with Q = 4, which is higher than that of Q = 9 as the average load per wavelength channel (M/Q) increases. Moreover, with the same values of Q and M , it is clear that the wavelengthswitching model would provide us with shorter paths and a higher capacity on Erdős-Rényi and scale-free networks as shown in Fig. 6.
In comparison with NDP, the advantages of EDP with multiple transmissions using the same wavelength availability is clearer on random regular graphs as shown in Fig. 7(a), and less so on other graph types, presumably due to the graph heterogeneity and finite size effects. We observe that EDP routing increases the network capacity, i.e. average transmission load per wavelength channel, and provides shorter valid paths.

Dependence on the number of transmissions
In practice, one may wish to minimize the total number of active wavelength channels for given demand. We define the smallest number of wavelength channels which accommodate a set of transmission to be Q min ; to find Q min from for a specific instance, we gradually increase Q  We compare the values of Q min /M obtained by our algorithms and by the MGA algorithm on random regular networks, the CONUS and BT-Core networks, in both NDP and EDP scenarios. The numerical results on random regular networks presented in Fig. 8(a) and 8(c) show that our algorithms offer significant advantages over the MGA algorithm in reducing the number of wavelength channels required for specific random transmission pairs in. The improved performance is also observed on the CONUS and BT-Core networks in Fig. 9 and 10 respectively. The underlying reason for the improvement is the low success rate of random greedy assignments. The numerical results demonstrate that our MWR algorithms lead to a better use of resource. The additional flexibility provided by the wavelength-switching on transceivers leads to an even smaller value of Q min as M increases in Fig. 8(a) and 9(a), whereas the experiments on the BT-Core network do not show a significant improvement of WS over the ordinary NDP in Fig. 10(a), which may depend on the topology of the graph or our range of values tested.
Comparing the average path length found for both NDP and EDP routing on random regular graphs, we see in Fig. 8(b) that the average path lengths obtained by the MGA and MWR for NDP routing are similar but that shorter paths are found by MWR in the EDP scenario ( Fig. 8(d)). A significant reduction in average route length is not expected since typical route lengths on random graphss is O(log N ) with little variability. However, a significant reduction in the number of wavelength channels used Q min , is shown in Fig. 8(d) and a similar trend is shown for EDP scenarios in Fig. 9(a,c) and 10(a,c) for the two real networks.

B. Computational complexity
Next, we discuss the computational complexity of our multi-wavelength routing algorithms with a linear cost function F i,j (x) = x. In this case, if we use the normalization φ a i→j (s) := φ a i→j (s) − φ a i→j (0) in the min-sum equations (13) (15) (20), the computational complexities would decrease.
Under NDP routing, the complexity in computing Eq. (13) is approximately O( k 2 M ), where k is the average node degree, and there are k N Q such messages. The complexity in computing the messages from origin/destination nodes in (10)  , which linearly scales with Q as shown in Fig. 11(a).
For EDP routing, the complexity of matching a pair of incoming and outgoing transmissions with the same wavelength in Eq.(23) is O( k 3 log k ) [9]. Therefore, the total complexity of one round of update for all messages is O(2M k 4 log k × k N Q + M Q 2 ), and if the graph is sparse the complexity becomes O(M N Q + M Q 2 ). For relatively small Q, the complexity scales with Q as in NDP routing as shown in Fig. 11(b).
However, when M/N N , the effective degree becomes quite different from the average degree k ⇒ k + 2M N (see in Fig. 2, the degree of node i increases from 3 to 4 after introducing the auxiliary node µ). When M is very large, for instance M ∼ O(N 2 ) as is used in many optical communication network applications, the effect of effective degree cannot be omitted. To make the algorithm scale better we divide the messages into three types: 1. messages from auxiliary nodes to ordinary nodes, e.g. µ → i in Fig. 2 -there are 2M Q such messages; 2. messages from ordinary nodes to ordinary nodes, e.g. i → j in Fig. 2 -there are k N Q messages of this type; 3. messages from ordinary nodes to auxiliary nodes, e.g. i → µ in Fig. 2 -there are also 2M Q such entities.
The messages of type 1 obeys Eq. (10) and the resulting total complexity is O(M Q 2 ). Noticing that the messages from auxiliary nodes are sparse, only φ µ→a (0) and φ µ→a (µ) are non-trivial in Eq. (10), and only φ b→µ (−µ) is needed to generate new messages, which simplifies the calculations of Eq. (13). Under NDP routing, by separating messages of type 1 and 2, the computational complexity of generating one message of type 2 or 3 by Eq. (13) is O( k 2 M ). Then, the total complexity of the three types of messages is  (N 5 ). In Fig. 11(c), we present the one-iteration computing time on random regular graphs with degree k = 3 and N nodes, the one-iteration computing time scales approximately as N 5.14 , which fits our analysis of Under EDP routing, the computational complexity of type 2 and 3 messages depends on the complexity of the maximum weighted matching algorithm (approximately O(( k + 2M N ) α )) with α being a coefficient to be determined. Using the sparse properties of type 1 and 3 messages, the total complexity of the three types is O(M Q(N + Q + (M/N ) α )). Recalling the mapping to matching problems -the matching network has approximately k + 2M N nodes and k ( k −1 2 + 2M N ) edges, and the average node degree is approximately k when 2M N k , therefore it is a sparse network. By using the algorithm in [9], the complexity should be approxi- ). However, due to resulting topology and the special network structure our tests result in complexity scaling close to O( k + 2M N ) when tested. In other words, α = 1 and therefore the resulting total complexity is O(M Q(N + Q + M/N )). Experiments on random regular graphs presented in Fig. 11

C. Non-linear cost
In this section, we show the simulation results where a more general form of the cost function where the argument is the wavelength occupancy q i,j on edge (i, j), given by When γ > 1, the utility increases faster on heavily loaded edges, biasing solutions towards paths with uniform load, balancing the loads on edges. On the other hand, when 0 < γ < 1, the cost on edges increases slower with the load, leading to configurations which consolidate transmissions on used edges, leaving more edges and wavelength channels idle. Such configurations would be relevant for identifying less important transceiver nodes which can be switched off in hours of low usage. In Figure 12 we show the distribution p(x) of wavelength occupancy on edges of random regular graphs and the two real network CONUS and BT-Core for different values of γ. For all three networks and both NDP and EDP scenarios, we observe that when γ = 0.5 more edges were unused, i.e. a higher value at p(x = 0) in Fig. 12, and as γ increases the distribution become more evenly distributed and peaked at some values of x, which corresponds to the balancing of edge loads. We show in the insets the relative difference between p(x) obtained by γ = 1 and those obtained by γ = 0.5 or 2, i.e. δ(x; γ) ≡ p(x;γ)−p(x;γ=1) p(x;γ=1) , where the evidence for load consolidation or balancing, for γ = 0.5 or 2 respectively, become obvious.

D. Comparison with linear programming
Linear programming is a commonly used method for optimized throughput of optical networks [57]. Here, we compare the results of our proposed algorithm with those To solve the routing problem using linear programming, we find n shortest paths for each of the M transmissions, denoted by the variables { s µ k } with superscript µ = 1, . . . , M representing transmissions and subscript k = 1, . . . , n representing the k th candidate path for transmission µ. By carrying either a subscript for nodes or edges, the variables s µ k can represent a configuration of node states or edge states, and the NDP or EDP constraints can both be expressed in terms of s µ k . In this case, we introduce a variable σ a,µ k = 1 to denote that transmission µ chooses its k th path with wavelength a, σ a,µ k = 0 otherwise. The disjoint constraint can be ex-pressed as follows, In Eq. (25), if the path s µ k represents a configuration of node states, then each element of v a represents the load of a th wavelength channel on that node, and node-disjoint constraint restricts the value of the load to be no more than 1; if the path s µ k represents a configuration of edges, then the edge-disjoint constraints can be defined in a similar manner.
For each transmission, we choose one wavelength to accommodate one candidate path, which is given by The objective function for the problem with linear cost is given by The problem can be solved by linear programming with the objective function of Eq. (27) subject to the constraints in Eq. (25)- (26) where all the expressions are linear.
For the WS scenario, the variables σ µ k are introduced instead of the previous σ a,µ k , and the capacity constraint for each node is given by where s µ k is a configuration of node states. The constraint for all the transmissions is given by The objective function to be minimized is j v j . We conducted numerical experiments on four real networks [56,57] and compared the results obtained by our algorithms and linear programming in the NDP, wavelength-switching and EDP scenarios. We show the smallest numbers of wavelength channels, i.e. Q min , needed to complete all possible transmissions, and the corresponding total path length L in Tab. I. As we can see, our message-passing algorithms and linear programming yield almost identical performance in finding optimized path solutions. For relatively small network, linear programming is efficient, but it quickly becomes impractical when the size of networks increases, whereas we show in Sec. IV B that our message-passing algorithms have more practical scaling properties. Table I. The smallest number of wavelength channels Qmin required to transmit all M = |V |(|V | − 1)/2(= N (N − 1)/2) transmissions on four small real optical communication networks including NSF-Net, Google-B4, DTAG/T-systems and BT-Core, obtained by our multi-wavelength routing algorithm (MP) compared to that obtained by linear programming (LP), in edge-disjoint (ED), node-disjoint (ND) and node-disjoint wavelength-switching (WS) scenarios. The corresponding total path lengths are also shown. In real optical networks, the number of wavelength channels in different optical fibers may vary; their lengths, signal-to-noise ratios, or type of fibers used may influence the quality of communication. To make our model more general and realistic, we consider the case of optical networks with weight w i,j and number of wavelengths Q i,j defined for any individual link (i, j) (or Q i defined for node i in node-disjoint wavelength-switching scenarios). With simple modifications, our proposed algorithm can be generalized to accommodate w i,j , Q i,j or Q i .
In cases with heterogeneous Q i,j on edges, we denote the largest number of wavelength channels among all edges to be Q * , i.e. Q * = max (i,j) Q i,j , then for an edge (i, j), one can introduce Q * − Q i,j additional wavelength channels with state 0, such that all edges on the network would have virtually Q * wavelengths.
As for heterogeneous weights on edges with linear cost, the objective function becomes L = (i,j) w i,j A. Heterogeneous wavelength availability on edges In the message-passing equations (13) describing NDP scenario, the first equation φ a i→j (0) considers the condition for wavelength a of edge (i, j) to be in state 0 in the absence of node j, and the contribution to the objective function L is 0; this equation is the same in cases with heterogeneous edge weights. Nevertheless, weights w i,j should be introduced in the second equation as The marginal messages on edges, Eq. (14), should be modified as φ a i,j (s) = φ a i→j (s) + φ a j→i (−s) + w i,j (δ 0 s − 1).
The same algorithmic procedure described in Sec. III C can be applied in the present cases with heterogeneous edge weights.
For an even more general scenarios, our model can be modified to accommodate heterogeneous weights for different wavelength channels on the same edge, i.e. w a i,j = w b i,j for a = b, asymmetric directed weights such as w i,j = w j,i on directed graphs. In these cases, one only needs to replace w i,j in the above modified equations by w a i,j or directed weights.

B. Heterogeneous wavelength availability on nodes with wavelength-switching
For the NDP scenarios with heterogeneous wavelength availability on nodes, we do not have to introduce additional wavelengths and keep unavailable channels in state 0 if Q i are non-uniform; this is because the equation of node capacity constraint Eq. (16) has already considered the case of different Q i values. Only the second equation of Eq. (15) needs to be modified, which reads φ µ i→j (±1) ∼ w i,j + τ µ i (1)+ min

VI. CONCLUSION
Multi-wavelength NDP/EDP routing lies at the heart of the efficient running and design of optical communi-cation networks, that act as the backbone of the Internet. One of the key questions in running optical communication networks more efficiently is in the ability to carry out these routing tasks effectively for large systems. This serves for both day-to-day running of the network as well as for the design of new networks and the modification of existing infrastructure. While principled single wavelength NDP and EDP routing algorithms based on message passing have been developed already, they could not be employed in real optical networks due to the difficulty in extending the algorithms from the single wavelength to the multi-wavelength case. This essential and pivotal aspect of routing in optical communication networks makes the existing single-wavelength methods intractable and requires computational cost that grows exponentially with the system size.
To accommodate a large number of wavelengths and transmissions in large systems, we have developed algorithmic solutions that include multi-layer graphs, where each layer represents a different wavelength, and messages are passed within layer (routing assignment) and between layers (wavelength allocation). The scalable algorithm we have devised shows very good performance in manageable time scales.
We expect the algorithm to be implemented in realistic scenarios, where specific aspects of real network routing, such as heterogeneous wavelength availability and signal-to-noise ratios will have to be added in the man-ner outlined in Sec. V. We also expect the algorithm to be utilized for network design and see several possible extensions for both localized and global message passing-based implementation. Utilization of the algorithms developed here in ad-hoc network communication, multilayer VLSI design and multilayer networks will require further study of the specific requirements for the different applications.


Linguistic emergence from a networks approach: The case of modern Chinese two-character words

The models of linguistic networks and their analytical tools constitute a potential methodology for investigating the formation of structural patterns in actual language use. Research with this methodology has just started, which can hopefully shed light on the emergent nature of linguistic structure. This study attempts to employ linguistic networks to investigate the formation of modern Chinese two-character words (as structural units based on the chunking of their component characters) in the actual use of modern Chinese, which manifests itself as continuous streams of Chinese characters. Network models were constructed based on authentic Chinese language data, with Chinese characters as nodes, their co-occurrence relations as directed links, and the co-occurrence frequencies as link weights. Quantitative analysis of the network models has shown that a Chinese two-character word can highlight itself as a two-node island, i.e., a cohesive sub-network with its two component characters co-occurring more frequently than they co-occur with the other characters. This highlighting mechanism may play a vital role in the formation and acquisition of two-character words in actual language use. Moreover, this mechanism may also throw some light on the emergence of other structural phenomena (with the chunking of specific linguistic units as their basis).

Introduction
The actual use of language has attracted increasing research interest in recent decades, especially regarding its role in the formation, acquisition, evolution, etc., of linguistic structure. Viewed from a usage-based perspective, linguistic structure emerges from language experience based on exemplars of language use to which language users are exposed [1,2].
Although the structural patterns generally lack natural boundaries in actual language use, they can establish themselves and can be extracted by language users by standing out through particular relational patterns (e.g., strong associations between their component units). This has been supported, directly or indirectly, by various studies [3][4][5], which suggest that the quantitative patterns of actual language use may shed light on the emergent nature of linguistic structure. In order for the quantitative patterns of actual language use to shed light on linguistic emergence, two prerequisites need to be considered, namely, appropriate modeling and quantitative measures of actual language use.
The notion of networks often finds its way into the modeling of language experience. Hudson [6], for instance, held that language is represented as part of the conceptual network of human knowledge and employed the notion of networks in structural analysis at various language levels. Bybee [7] argued that language experience can be conceived of as networks composed of linguistic units (e.g., words) and their relations in actual language use, such as associations (e.g., phonetic and semantic similarities) and co-occurrence relations. Furthermore, Bybee [7] also sought a network-based explanation of how particular structural patterns (e.g., idioms and morphological patterns) emerge from language experience. The notion of networks constitutes a relational approach to language in that anything (be it a linguistic unit, a pattern of use, or a linguistic feature) is considered against a larger context of linguistic units and their relations. This marks an intrinsic advantage of networks thinking in usage-based inquiries into linguistic structure, for both linguistic structure and the context of actual language use manifest themselves as linguistic units and their relations.
Linguistic structure is held to emerge from the co-occurrence relations of linguistic units in actual language use [2,7]. Linguists have available to them a wide range of measures for cooccurrence strength [8], with varying computational costs. Simple measures (e.g., frequency and probability) are adopted from time to time in usage-based linguistics [7], largely due to their close bearing on language experience and thus the mechanisms of linguistic emergence. However, the mere inspection of these simple measures alone might still be insufficient for a systematic understanding of the formation of structural patterns in actual language use. For instance, a high co-occurrence frequency of two units does not necessarily mean a meaningful unit, and vice versa. However, the advantages of simple measures are worth considering. The specific roles of these measures in linguistic emergence might be further understood with appropriate models and analytical tools of actual language use. A networks approach may contribute to this understanding by taking into account the larger relational background of any co-occurrence relation and its strength.
To recap, network models plus simple measures for co-occurrence strength (such as frequency) may constitute an excellent approach to the emergence of linguistic structure. However, such networks thinking needs to be converted into substantial studies.
Rooted in a multi-disciplinary background, the models and quantitative tools of linguistic networks [9,10] constitute an operational methodology for networks thinking in linguistics. The basic form of a linguistic network model N can be generally represented as N = (V, E), with V being the set of nodes (vertices) and E the set of links (edges). Various language levels as language sub-systems can be modeled and analyzed quantitatively as linguistic networks, with the corresponding linguistic units (e.g., words) as nodes and their relations (e.g., structural and semantic relations between the words) of a particular type as links. The models can also include extra information concerning the features of linguistic units and their relations (such as co-occurrence frequencies of linguistic units) if necessary. Structural analysis of linguistic networks is supported by a wide range of quantitative tools, which rely on relational data [11] in that they always focus on the relations between linguistic units and features. These tools can hopefully shed light on structural patterns (as relational phenomena) of language at different scales of granularity.
Earlier inquiries into linguistic networks have helped to appreciate the system-level complexity of particular language sub-systems [12,13]. It has been shown later that such an understanding of the system-level complexity of language can also shed light on language development [14], special education [15], language typology [16], and so forth. In practical terms, such an understanding can contribute to various NLP tasks [17][18][19]. More recent studies have extended the scope of research. New ways of modeling have been adopted. For instance, the analysis of linguistic networks of different sub-systems of the same language has deepened the understanding of language as a multi-level system [20,21]. The use of multilayer sentence networks (distinguishing links between sentences from the same document and those connecting sentences from different documents) has facilitated NLP undertakings such as multi-document summarization [22]. Moreover, increasing interest has been devoted to the rather microscopic features of linguistic networks, which are examined with measures including those concerning node centrality [22][23][24].
Inquiries into the less macroscopic features (e.g., the various cohesive sub-networks) of linguistic networks can hopefully shed light on how structural linguistic units are formed by their component units, which in turn may contribute to a usage-based understanding of linguistic emergentism. This line of research has just started. For instance, Goh et al. [25] have shown that the motifs (cohesive sub-networks) in word co-occurrence networks based on authentic English texts function as network shortcuts, and the motif densities are drastically lowered when the texts are shuffled. More importantly, their results have revealed a link between these network motifs and linguistic constructions, which points to the mechanism of these constructions' emergence from the patterns of actual language use.
The formation of Chinese words (as structural units based on the chunking of their component characters) in actual language use constitutes an interesting case study in this line of research. This study explores the formation of Chinese two-character words (henceforth, CTCWs) in the actual use of modern Chinese, which manifests itself as continuous streams of Chinese characters and which was modeled as character co-occurrence networks. The network models consist of Chinese characters as nodes, their co-occurrence relations as directed links, and the co-occurrence frequencies as link weights. Quantitative analysis of the network models has shown that a Chinese two-character word can highlight itself as a two-node island, i.e., a cohesive sub-network with its two component characters co-occurring more frequently than they co-occur with the other characters. This highlighting mechanism may play a vital role in the formation and acquisition of CTCWs in actual language use. Moreover, this mechanism may also throw some light on the emergence of other structural phenomena (with the chunking of specific linguistic units as their basis).

Chinese words as structural units
Unlike other languages such as English, Chinese words lack natural boundaries in actual language use. Morphological markers such as affixes, which can serve as word boundaries, are generally non-existent in Chinese [26: 46]. Moreover, even in written Chinese, the words are not explicitly separated. In other words, they are not separated by spaces as in the case of English (see Example 1 for illustration). The actual use of Chinese language manifests itself as continuous streams of characters, which are generally understood as morphemes of the language [27], the elements of word formation. Each character is a linguistic unit corresponding to one syllable, one orthographic symbol, and usually a particular amount of meaning. To recap, spoken Chinese and written Chinese can be seen as continuous streams of syllables and orthographic symbols, respectively. Each Chinese word is a structural unit, consisting of one, two, or sometimes even more characters chunked together. CTCWs constitute a dominant portion of modern Chinese vocabulary.
Example 1 below illustrates how the actual use of modern Chinese looks like, especially the non-existence of spaces between words. It is an excerpt from the Chinese version of Bloomfield's Language [28: 21], which is followed by its English original.

Example 1
至于语言变化问题, 我们已有足够的事实可以证明, 所有的语言都同样有变化过程, 而且 都倾向于同一方向。甚至很特殊的变化类型, 在差别最大语言里也可以发生非常相同的 变化, 只不过是独立地进行而已。 As to change in language, we have enough data to show that the general processes of change are the same in all languages and tend in the same direction. Even very specific types of change occur in much the same way, but independently, in the most diverse languages.
Example 1 shows a stream of 85 Chinese characters (tokens), separated only by 7 punctuation marks. Example 2 below shows one clause in Example 1 with all the words segmented by spaces (1 st line), its phonetic transcription using the official romanization system for Standard Chinese (2 nd line), word-by-word gloss in English (3 rd line, REL = relative marker), and the English original (4 th line).

Example 2
所有 的 语言 都 同样 有 变化 过程 suáyáu de yīyán dōu tóngyàng yáu biànhuà guòchéng all REL language all same have change process [T] he general processes of change are the same in all languages. . . In actual language use, any two co-occurring (i.e., adjacent) Chinese characters may fall into either of the two cases: (1) they are within the same word (e.g., any two characters which are not separated by a space in Example 2, such as 语 (yī) and 言 (yán)); (2) they constitute the boundary between two adjacent words (e.g., any two characters separated by a space in Example 2, such as 的 (de) and 语 (yī)), i.e., the two characters belong to two different words but happen to be adjacent.
CTCWs constitute an interesting subject of inquiries into the formation of structural patterns from a usage-based perspective. Although it has proved possible to segment Chinese words largely on the basis of statistical features of the character streams [29,30], the statistics adopted by these engineering attempts generally involve complicated computation, which may be rather removed from how the patterns of actual language use affect language experience. In this study, priority is given to rather simple measures of actual language use (especially frequency), which are examined from a networks approach.

A networks approach to frequency
Frequency is held to play a fundamental role in the formation, evolution, learning, and mental representation of structural patterns at various language levels [5,31]. For instance, the frequencies of linguistic units and their combinations have a profound impact on the way how language is chunked in memory, how such chunks are connected, and how easily they are accessed [1].
Viewed from a usage-based perspective, Chinese words emerge from the co-occurrence relations of Chinese characters in actual language use. Unless otherwise specified, co-occurrence in this paper means that two linguistic units (e.g., Chinese characters) are immediate neighbors of each other and thus form a bigram as an ordered pair. Given two linguistic units u and v, the two ordered pairs they form, uv and vu, are different co-occurrence relations or bigrams. The three co-occurrence relations of characters (or character bigrams) at the beginning of Example 1, for instance, are 至于 (zhìyú), 于语 (yúyī), and 语言 (yīyán). In actual language use, any CTCW can be seen as a reusable character bigram denoting one or more meanings. This repetition of co-occurrence helps to entrench the character bigram as a structural whole. Co-occurrence frequencies of characters, therefore, may play a crucial role in the formation and acquisition of CTCWs in actual language use. For instance, the top 10 most frequent character bigrams in the Lancaster Corpus of Modern Chinese [32] are all CTCWs.
However, a more in-depth understanding of the role of frequency (e.g., the specific scope of its effects) is needed. A general law of human language is that word frequency varies drastically [33]: only few of the words being highly frequent while the majority of them rather infrequent. The same is true of CTCWs. A large number of character bigrams with low frequencies may still be CTCWs. What is worse is that many non-word character bigrams exhibit rather high frequency levels. The frequency level alone, therefore, may not be a reliable indicator of whether a character bigram forms a word or not. However, the actual language use to which language users are exposed is ever-growing and ever-changing [7]. Frequency may always play a role in the formation of structural patterns, for even the least frequently-used units or patterns according to corpus statistics can occur repeatedly in actual language use. In order for a word uv to stand out as a cohesive unit in actual language use, the frequency of uv may not have to compete with all the other character bigrams in the entire context. Instead, it may be more plausible to examine the frequency of uv against its immediate context, i.e., the bigrams which u and v respectively form with other characters.
In authentic Chinese language use, any two characters u and v forming a bigram uv also cooccur with other characters at the same time, and these co-occurrence relations constitute the immediate context of bigram uv. In Example 1, 语 (yī) and 言 (yán), in addition to co-occur with each other, also co-occur respectively with other characters, hence the immediate context of 语言. The bigram 语言 (a CTCW meaning 'language') and its immediate context can be modeled as a network in Fig 1, whereby any two characters forming a bigram in the text is connected by a directed link (an arrowed line indicating the order of characters), and the frequency of each bigram is represented by a value attached to the corresponding link and the thickness of the link.
As discussed earlier, a character bigram constitutes either a word (or part of it) or the boundary between two adjacent words. In the former case, the two characters may form a CTCW (e.g., 语言 in Example 1) and thus can be used repeatedly. A character bigram forming the boundary between two words, on the other hand, is due to the co-occurring of the two words as a result of syntactic arrangement of the utterance. Considering the flexibility of syntax, such a bigram (e.g., those formed by 语 or 言 with the other six characters as illustrated by Fig 1) may not exhibit the same degree of reusability as a CTCW. It follows that if two characters form a CTCW, they may tend to co-occur with each other more frequently than they co-occur with other characters. Given a character bigram uv found in a corpus, its own frequency can be seen as its internal tightness, and the frequencies u and v co-occur respectively with their other immediate neighbors can be seen as the external tightness of uv. A CTCW, therefore, is likely to exhibit itself as a two-character chunk with stronger internal than external tightness. In this way, a CTCW as a cohesive structural whole can stand out in its immediate context, which may contribute abundantly to its emergence from actual language use. In Example 1, for instance, 语言as a CTCW stands out as such a chunk (see Fig 1).
In sum, a networks approach to co-occurrence frequency may be more plausible in accounting for the emergence of a CTCW, which focuses on the frequency level of the corresponding character bigram relative to those of the other bigrams in its immediate context, instead of its frequency level alone. Such an approach can be operationalized using network analysis of authentic Chinese language data.

Materials and methods
In order to understand the emergence of CTCWs, this study attempts to examine (1) the twocharacter chunks with stronger internal tightness as previously defined, and (2) how such chunks relate to CTCWs. Network analysis of authentic Chinese language data can be employed to detect such chunks with stronger internal tightness.

Materials
The authentic language data adopted in this study are from two corpora, namely, the Lancaster Modern Chinese Corpus (LCMC) [32] and the Leiden Weibo Corpus (LWC) [34], which cover the formal and informal use of modern Chinese, respectively.
LCMC is a balanced corpus of modern written Chinese with about 1.6 million Chinese character tokens (roughly 1 million word tokens). It consists of 500 samples, 3,200 Chinese character tokens apiece, and covers 15 text categories. Two sub-corpora of LCMC were adopted for this study, namely, those of press reportage (LCMC_A, about 120,000 character tokens) and science (academic prose) (LCMC_J, about 230,000 character tokens). These corresponds to two representative genres of formal language use of Chinese.
LWC consists of over 5,100,000 messages posted on Weibo, China's leading social media platform. Compared with the language data in LCMC, those in LWC are rather close to daily spoken interactions. In this study, the first 8,000 messages with substantial linguistic content (i.e., consisting of at least one word) in LWC were selected as a sub-corpus (about 250,000 character tokens) concerning the informal use of Chinese.

Network modeling
Considering the advantages of linguistic networks in dealing with relational phenomena (such as linguistic structure), it is necessary to construct and analyze linguistic network models based on authentic language data. These networks have linguistic units as nodes and their relations (e.g., co-occurrence, syntactic dependency, or semantic dependency relations) observed in actual language use as links. Network models constructed in this way, therefore, are potentially capable of modeling language users' exemplar-based language experience. Quantitative analysis of such networks can help to detect various cohesive clusters formed by linguistic units in actual language use, which can hopefully shed empirical light on the emergent mechanisms of particular structural patterns, such as CTCWs.
As previously noted, the actual use of Chinese language can be seen as continuous streams of characters (as the immediately-observable linguistic units). In network analysis, these can be modeled by directed and weighted character co-occurrence networks. Such a network can be represented as N = (V, E, w), with V being the set of all the characters (types) in the corpus, E the set of co-occurrence relations of characters, and w the set of link weights representing frequency values of the corresponding co-occurrence relations. Only Chinese characters were counted as network nodes. Others symbols, such as Arabic numerals and punctuation marks, were all ruled out. Any two co-occurring characters u and v forming a bigram uv were connected by a directed link pointing from u to v in the corresponding network. In rare cases, a character u might co-occur with itself to form a bigram uu. For instance, some CTCWs are formed by reduplicating the same character, such as 叔叔 (shūshu, 'uncle'), 狒狒 (fèifèi, 'baboon'), and 天天 (tiāntiān, 'every day'). A co-occurrence relation of this type was represented by a link pointing from one character to itself, which is termed as a self-loop. Note that two characters separated by a punctuation mark representing a pause (such as comma, and period 。) did not count as co-occurring. The frequency of any co-occurrence relation of characters is the weight of the corresponding link. In Fig 2, the number attached to each directed link represents its weight, that is, the frequency of the corresponding bigram. For instance, as the character bigram (in fact a CTCW) 语言 (yīyán, 'language') occurred 3 times in the text, the weight of the corresponding link is 3 (see also Fig 1).
With the scheme introduced above, the three sub-corpora were converted into three network models (henceforth, Networks LCMC_A, LCMC_J, and LWC), respectively.

Data analysis
As previously discussed, in a body of authentic language data, if two characters co-occur more frequently than they co-occur with the other characters, they form a chunk with stronger

PLOS ONE
internal than external tightness in its immediate context. A chunk of this type is in line with the notion of islands [35] in network analysis. An island is defined as a sub-network, whose inner-links, those connecting its member nodes, all have greater weight values than its outerlinks, those connecting its member nodes with nodes outside it. A two-character chunk with stronger internal tightness due to its repeated use (e.g., 语言 (yīyán, 'language') as illustrated by Fig 1), therefore, corresponds to a two-node island (henceforth, TNI) in a linguistic network constructed as previously described.
The classical definition of islands is based on networks without self-loops. In the actual use of Chinese, the co-occurring characters do form self-loops, which are rare yet worth considering. Some of these self-loops are CTCWs, e.g., 叔叔 (shūshu, 'uncle'). For the purpose of this study, an adjusted definition of TNIs is formulated, so that the self-loops formed by co-occurring characters are considered. Let the directed link between two nodes u and v be e (allowing u = v, and hence a self-loop), and any other directed link involving u or v e'. If the weight of e is greater than the maximal weight of the latter, i.e., w(e) > max(w(e')), then <u, v> forms a TNI.
Based on the above definition of TNIs, data analysis of this study took two major steps: (1) extraction of TNIs from the linguistic networks, and (2) determination of the relationship between the TNIs and CTCWs.
Chinese words are rather fluid units in that it is sometimes difficult to draw the boundary between a word and phrase. Dictionaries, lexicons, and lexical databases published officially or academically constitute an operational (though not perfect) tool for wordhood judgment. The sources adopted by this study are: (1) Contemporary Chinese Dictionary (7 th Edition) [36] (henceforth, CCD), the most authoritative reference work on modern Chinese words; (2) Lexicon of Common Words in Contemporary Chinese [37] (henceforth, LCWCC), released by the State Language Commission of China; (3) the Chinese Lexical Database [38] (henceforth, CLD), a large-scale lexical database for simplified modern Chinese which well reflects the linguistic experience of modern Chinese users; (4) Wiktionary (https://www.wiktionary.org); and (5) Baidu Encyclopedia (https://baike.baidu.com/). Both (4) and (5) are powerful online dictionaries with a good coverage of Chinese neologisms, slang words, and dialect words. Technical terms which are not listed in the above sources were judged on the basis of (1) the terms published by China National Committee for Terms in Sciences and Technologies (partially searched on http://www.termonline.cn/index.htm), and (2) Dacihai [39], a large-scale dictionary and encyclopedia of modern Chinese.
In order to ensure the accuracy of wordhood judgment, the TNIs extracted were examined manually based on the above sources and their contexts. A TNI has to meet one major criterion and two minor criteria to be counted as a CTCW. The major criterion is that the TNI should be listed in at least one of the above sources as an entry. The two minor criteria are that (1) not all occurrences of the TNI are boundaries between two adjacent words (i.e., co-occurring but belong to two adjacent words), and (2) not all occurrences of the TNI are parts of longer words. The major criterion is generally sufficient for wordhood judgment, while the two minor criteria help to rule out two very rare cases, which further improves the accuracy of judgment.

Results
The TNIs extracted from each network were classified into three types: (1) CTCWs, (2) wordlike chunks, and (3) non-word chunks. Table 1 displays the number and percentage of TNIs of each type (see Appendices A-C in S1 File). The few TNIs with weight 1 (altogether 3, 1, and 3 in LCMC_A, LCMC_J, and LWC, respectively) have been ruled out, for they were not formed through repeated use of the corresponding character bigrams. Those (altogether 2) in the three sub-networks of Network LCMC_J were treated in the same way.
The Type-1 TNIs are those identified as words according to the abovementioned criteria for wordhood judgement. As can be seen from Table 1, CTCWs constitute an overwhelming majority of the TNIs in each network. The CTCWs extracted as TNIs from Network LWC cover some slang words widely used in daily conversations and on social media platforms, such as 粉丝 (fěnsī, 'fans'), 童鞋 (tóngxié, 'schoolmate' or 'classmate'), and 坑爹 (kēngdiē, 'cheating' or 'deceiving'). These words reflect the informal use of modern Chinese.
The Type-2 TNIs do not count as words according to the criteria of this study. However, they constitute either a two-character phrase (e.g., 万元 (wànyuán, 'ten thousand Yuan RMB') or part a word/phrase of three or more characters (e.g., 铃虫 (língchóng) as part of 棉铃虫 (miánlíngchóng, 'bollworm') and 红铃虫 (hónglíngchóng, 'red bollworm')). As previously noted, the boundary between a two-character phrase and a CTCW is not always clear-cut and even LCWCC and CLD list a number of CTCWs which are sometimes treated as two-character phrases. A two-character phrase behaves very much like a CTCW in that they are both used as a structural and meaningful whole. For instance, 万元 (with 万 (wàn) meaning 'ten thousand' and 元 (yuan) 'Yuan RMB') can also be seen as a CTCW for a unit of currency.
A Type-3 TNI consists of two characters which respectively belong to two adjacent words and cannot be treated as a phrase. As can be seen from Table 1, the proportion of non-word islands in each network can be almost ignored.
In sum, the above results have preliminarily shown that the two-character chunks with stronger internal than external tightness due to their repeated use are generally CTCWs.
In network analysis, the minimal weight value of the links between the component nodes of an island is termed as the height of the island [35]. In this study, the height of a TNI is the frequency of the corresponding chunk. In Fig 1, for instance, 语言 (yīyán, 'language') as a TNI exhibits a height of 3. The distribution of word frequency in natural language can be captured by power law [33], which is known as Zipf's Law. In other words, when the words in a language are ranked in descending order of frequency, word frequency (f) and rank (r) generally follows a power-law relationship f~r -γ . The rank-frequency distribution obtained is extremely uneven: a rather limited number of words having extremely high levels of frequency while the majority of words having rather low levels of frequency. Most CTCWs have very low levels of frequency, which however does not deny their status as words. Fig 3 displays the rank-frequency distributions of Type-1 TNIs, i.e., the CTCWs, in the three networks (see Appendices A-C in S1 File). All the plots in Fig 3 exhibit Zipf-like distributions, with word frequency dropping abruptly to a low level and then decreasing rather slowly as word rank increases. It might be unnecessary to test statistically how well the distributions fit power law. The most important message from the distributions, however, is that the TNIs extracted covered CTCWs with vastly different frequency levels. In other words, there is no necessary connection between the formation of a TNI and the frequency level of the corresponding CTCW.
The definition of islands constitutes a rather rigorous criterion for detection of two-character chunks with stronger internal than external tightness, for an island is a local sub-network with the greatest link weight values. The same character (e.g., 语 (yī)) may occur in more than one word (e.g., 语言 (yīyán, 'language') and语法 (yīfǎ, 'grammar')). However, a TNI uv in a given network will prevent the other character bigrams containing u or v in the same network from forming islands, even though some of them might also be CTCWs. In addition, a CTCW which occurs in particular immediate contexts (e.g., one that happens to occur always between the same neighboring characters in a particular text) may also fail to form a TNI. It follows that island extraction as conducted in this study will inevitably miss a number of CTCWs in a network. The CTCWs that are missed by island extraction are worthy of investigation, especially whether they can form TNIs with the change of their immediate contexts. Actual language use is rich and flexible. The same is true of the immediate context of a CTCW. A body of authentic language use, regardless of its size, only constitutes a fragment of language experience. Given a particular fragment of language use (e.g., the sub-corpora LCMC_A and LCMC_J), a CTCW may fail to form a TNI in its immediate context. However, other fragments of language use may provide appropriate immediate contexts where the CTCW may well form a TNI. This is due to the nature of CTCWs as reusable character bigrams. This reusability gives rise to a frequency effect, which makes a CTCW a structural whole with stronger internal than external tightness in appropriate immediate contexts. In order to illustrate this point, the TNIs extracted from three sub-networks of Network LCMC_J were examined. These three sub-networks were network models constructed respectively on the basis of three text samples selected at random from sub-corpus LCMC_J (with 583, 1,495, and 3,353 character tokens, respectively). The three sub-networks are labeled as Networks LCMC_J_1, LCMC_J_2, and LCMC_J_3 according to their sizes (with Network LCMC_J_1 being the smallest). Table 2 displays statistics of the TNIs in the three sub-networks (see Appendix D in S1 File). The TNIs are generally CTCWs in the three networks, especially Networks LCMC_J_2 and LCMC_J_3. The relatively low proportion of CTCWs in Network LCMC_J_1 is probably

PLOS ONE
due to the small size of the text sample on the basis of which it was constructed. In each of the three sub-networks, there were CTCWs which failed to form TNIs in Network LCMC_J (see Table 2 for their numbers). For instance, 数量 (shùliàng, 'quantity') failed to form a TNI in Network LCMC_J due to the formation of 数据 (shùjù, 'data') as a TNI. However, 数量 (shùliàng, 'quantity') formed a TNI due to the change of its immediate context in Network LCMC_J_1. If text samples of various sizes can be selected repeatedly from a large Chinese corpus, then the TNIs extracted from these samples will cover the overwhelming majority, if not all, of the CTCWs (with a frequency of at least 2) in the corpus. Therefore, even with the rather rigorous criterion for island extraction, a CTCW still has chances to form a TNI thanks to its repeated use in appropriate immediate contexts.
To recap, CTCWs have been examined, in their immediate contexts, through network analysis based on the notion of islands. There are three major findings: (1) the TNIs in the linguistic networks based on Chinese characters and their co-occurrence relations are generally CTCWs; (2) a CTCW of any frequency level (usually at least 2) may form a TNI; and (3) any CTCW (usually with a frequency of at least 2) has the potential to form a TNI in appropriate immediate contexts.

Discussion
In a wide range of research fields, the notion of islands has helped researchers to identify, in systems of various types, cohesive groups which are meaningful in one way or another [40][41][42]. In this study, linguistic networks as models of actual language use of Chinese have been constructed and analyzed, with Chinese characters as nodes, their co-occurrence relations as links, and the co-occurrence frequencies as link weights. It has been found that the TNIs extracted from these networks are generally CTCWs. Furthermore, it has been shown that a CTCW always has the potential to form a TNI in appropriate immediate contexts. The findings of this study, therefore, suggest a highlighting mechanism of CTCWs in their immediate contexts. While this highlighting mechanism is by itself a statistical (frequency-based) effect, it may have a bearing on the emergence of CTCWs, especially when viewed from the perspective of attention/activation. An island as a cohesive whole is 'a local summit in the network', which is 'raised above its immediate surroundings' [35: 129]. A CTCW, by forming an island in appropriate immediate contexts, can highlight itself as a structural whole with stronger internal tightness (as can be seen from the case of 语言 (yīyán, 'language') in Fig 1). In this way, it can form a local peak in terms of attention/activation level in language experience and can be readily perceived and represented as a whole. The ever-changing active areas of language experience [6: 74] give rise to great flexibility of immediate contexts of the CTCWs, so that any CTCW can have the opportunities to form a local peak of attention/activation level in language experience. The repeated highlighting of a CTCW in appropriate immediate contexts can entrench a CTCW as a local peak in language experience, so that it can it can establish and maintain its status as a structural whole. Similarly, this highlighting mechanism may also contribute to the holistic representation of a two-character phrase (e.g., those extracted from the network models).

PLOS ONE
Considering that two-character phrases constitute an important source of new CTCWs in Chinese [43], this mechanism may also play a role in the future lexicalization of some of these phrases. Moreover, the highlighting mechanism points to a contextual/networks approach to cooccurrence frequency. The findings indicate that for a CTCW to be highlighted as a cohesive whole in its immediate context, its component characters' co-occurrence frequency does not need to reach a particular threshold. Such a non-threshold effect of frequency is also found in other studies [44][45][46]. What matters most to the highlighting of a CTCW has been found to be the advantage of its component characters' co-occurrence frequency in the immediate context. Such an advantage suggests that a CTCW can be used repeatedly in diversified contexts, highlighting it as cohesive whole in its immediate surroundings. In other words, co-occurrence frequency studied in isolation is insufficient in accounting for the chunking of Chinese characters into CTCWs; rather, it needs to be examined in its immediate contexts. Such importance of contextual information has been appreciated by a growing body of usage-based research. For instance, empirical inquiries [47,48] have shown that contextual diversity facilitates the processing and acquisition of lexical and sub-lexical units better than their frequencies alone.
This contextual/networks approach may also apply to other structural patterns (based on the chunking of lower-level units), with other quantitative measures of co-occurrence strength possibly involved. The above discussion is generally speculative. In future research, it is necessary to investigate the cognitive reality of the highlighting mechanism. For instance, will the two-character sequences that form TNIs differ in the shape of P300 from those that do not, considering that P300 is related to attention and activation of working memory [49,50]?

Conclusions
With appropriate models and analytical techniques of linguistic networks, this study focuses on the emergence of CTCWs (as structural units based on the chunking of their component characters) in the actual use of modern Chinese, which manifests itself as continuous streams of Chinese characters. Based on a relational view of co-occurrence frequency, two-character chunks with their component characters co-occurring more frequently than they do with the other characters (i.e., two-node islands in terms of network analysis) were extracted from linguistic network models of authentic language use of modern Chinese, and the relationship between the TNIs extracted and CTCWs was determined. It has been shown that the TNIs extracted are generally CTCWs and a CTCW, regardless of its frequency level (usually above 2), always has the potential to form a TNI in appropriate immediate contexts.
The findings of this study have shed some light on the emergent mechanism of CTCWs in language use. This mechanism helps to highlight a CTCW as a structural whole in its immediate context, and thus may play a vital role in the formation and acquisition of CTCWs in actual language use. This mechanism may also throw some light on the emergence of other structural patterns (with the chunking of specific linguistic units as their basis). Moreover, the findings may help to further understand the role of frequency (and probably other types of co-occurrence strength) in linguistic emergence, especially the specific scope of the role.
Methodologically speaking, this study has further shown the value of linguistic networks in usage-based research into the structural patterns of language. To be specific, linguistic network models based on authentic language data may provide an effective means to researching the emergence of structural patterns of human language. The findings of this study have also confirmed the advantages of linguistic networks in handling relational linguistic phenomena in authentic language use. This study is still preliminary and further researches are wanted. For instance, the highlighting mechanism is still at the computational level [51] in that it is based on statistics of actual language use and its cognitive realism needs to be determined by researches concerning the psychological and neural levels. In addition, considering the advantages of linguistic networks in dealing with relational linguistic phenomena, network-based investigations can be conducted to unveil the formation of other structural patterns of human language.