Primary user Emulation Attack Defense in Filter Bank based Cognitive Radio

Background/Objectives: The spectrum scarcity problem arising due to reckless development in wireless communication can be addressed using cognitive radio technology. An important step in cognitive radio technology is spectrum sensing and associated with that is an attack called primary user emulation attack. Methods/Statistical Analysis: A two channel Quadrature Mirror Filter (QMF) bank is employed for a single user spectrum sensing. The energy detection method is used as reference for comparing the performance of the filter bank method in the presence of attackers. Existing approaches like Neyman-Pearson Criterion, Improved Detection Scheme with Double threshold and Location Verification Method for detecting PUE attacks were examined. A new method which combines improved detection scheme and location verification method is proposed for detecting PUE attacks. Findings: The simulation result shows that the filter bank method shows better performance than the energy detection method in the presence of attackers. The probability of detection of PUE attacks using improved detection scheme and location verification method for the filter bank method and energy detection method was analyzed and found out that the filter bank based spectrum sensing method outperforms energy detection method. By incorporating AND rule logic, OR rule logic and Alternate logic approach into the filter bank spectrum sensing technique, the probability of detection of PUE attack was investigated and it is seen that the probability of detection has increased to 0.98 in FB method at SNR = 1 dB when employing an alternate approach. Application/Improvements: The feasibility of deploying the filter bank method integrated with the alternate logic is to be tested using Universal Software Radio Peripheral (USRP) module.

Introduction
The whole world is now moving towards wireless. Almost everything we touch, every aspect of our lives nowadays, has a wireless component to it. This wireless technology has become the most stimulating areas of telecommunications and networking. As a result, there arises the problem of spectrum scarcity. In order to tackle this issue, CR network was proposed as a technology to access the spectrum in an effective and adaptable manner so as to increase spectral efficiency. A CR network consists of PU's and SU's, where SU's are allowed to access the frequency bands allotted to PU's without upsetting them in an opportunistic manner 1 . There are four main functions in CR -spectrum sensing, spectrum management, spec-trum mobility and spectrum sharing. Among the four functions, spectrum sensing is the foremost step taken to identify the white spaces. CR nodes opportunistically engage these white spaces by operating across them without interfering with the PU's. The SU's has to vacate the occupied band when the PU begins to transmit across this frequency band. All SU's has equal right to access unoccupied band when there is no active PU communication.
The SU will act as MU by modifying the air interface, thereby mimic the PU's characteristic and gets the same privilege as the PU. Therefore, the SU's has to vacate the occupied band for the MU trusting that it is a PU. Hence, the MU gets unequaled access to the PU's band. This kind of attack against CR networks is called as PUE attacks. PUE attacks are unique to CR networks in which the intimidating user takes the advantage of the integral decorum in CR networks so that the authentic SU has to evacuate the spectrum band. The presence of PUE attack may severely affect the performance of CR networks. The consequence of PUE attacks is band width waste, QoS degradation, connection unreliability, denial of service and interference with the primary user. A PUE attack can happen while the spectrum sensing is performed by using cyclostationary, energy, matched filter or FB detection method. In this work, a FB based detection method is employed as it is considered as an optimum waveform for the 5G system 2 . An Improved Energy Detection Scheme based on Channel Estimation and implementation of energy detection scheme was discussed 3,4 . Similarly Filter Bank Multicarrier (FBMC) can be used as an integrated tool for data communication as well as channel sensing. An illustration of the launching of a PUE attack on CR network is shown in Figure 1. In 5 , the author have incorporated authentication scheme to combat a PUE attack against CR users. A hybrid machine learning method for malicious activity or policy violations in a network was proposed 6 . The author proposes a new mechanism based on physical layer network coding to detect the emulators 7 . In 8 , a PUE attack on ED method using Neyman -Pearson composite hypothesis test was discussed. A survey on security aspects in software defined radio and CR Networks was elaborated 9,10 . In 11 , the author proposes an approach that estimates the attack strength and innovatively applies in a Neyman-Pearson or likelihood ratio test to improve collaborative sensing performance. In 12 , authors has focused on Cooperative Spectrum Sensing (CSS) for double threshold improved energy detector. The vulnerabilities present in CR, attack classifications and their effect on the working of cognitive radio network are discussed 13 . A threat detection technique based on localization of primary signal was proposed 14 . The inclusion of fusion logic into cooperative spectrum sensing was discussed 15 .
The paper is organized as follows. In section 2, the basics of the FB based detection method are presented. The system model for PUE attack to frame the problem is discussed in section 3. The section 4 gives the existing approaches for PUE attacks. In section 5, a new approach for defending PUE attack in CR network is proposed. In section 6, the simulation results and discussion are provided, while section 7 concludes the paper.

Filter Bank based Detection Method
The basic building blocks for Multirate signal processing are decimation and interpolation 16 . The process of reducing the sampling rate by an integer factor M is called decimation. The process of increasing the sampling rate by an integer factor L is called interpolation. The down sampling process and up sampling process is shown in Figure 2(a) and 2(b). The filter bank is set of band pass filters with either a common input or a summed output. The decomposition and reconstruction process is called the analysis filter bank and synthesis filter bank. Three basic operations are used in filter banks: linear filter, down sampling and upsampling. The function of spectrum sensing can be performed using an analysis filter bank. A M-band analysis filter bank is shown in Figure  2(c). The analysis filter bank is used to decompose the input signal into a set of subband signal, where each subband occupies a portion of the original frequency band. A two channel Quadrature Mirror Filter (QMF) bank is adequate for a single user spectrum sensing. A block diagram of spectrum sensing based on the analysis filter bank under single user scenario is shown in figure 3.The signal received from the PU is divided into two subbands and filtered using low pass and high pass filter with equal pass band respectively, and are down sampled to reduce the sample rate by an integer M. The total energy of the down sampled subband signals is added and compared with the predefined threshold to find whether PU is present or not. Vol 10 (1) | January 2017 | www.indjst.org Sabiq P.V. and D. Saraswady

System Model For PUE Attack Detection
A system model for analyzing the PUE attacks is as shown in Figure 4 (a). The malicious and SU's is confined to a circular area of radius R. The PU's and SU's are separated by a distance of atleast d p . For PU detection, the method deployed for spectrum sensing is a FB based detection method. In the FB based detection method, the total energy of the down sampled subband signals is calculated and matched with the predefined threshold to find whether PU is present or not. The coordinates of the entire MU are transformed in such a way that the SU of interest lies at the origin as shown in Figure 4  This is because the probability of PUE attack on to any user in the network is the same as there is no support between SU's. Hence, here the Probability Density Function (PDF) of the signal received from one SU only is analyzed. The primary transmitter coordinates are transformed as (d p , θ p ). The MU's are present uniformly in an annular Region (R 0 , R).The amount of power received due to the transmission from PU at the SU, P r (p) can be obtained as: The PDF of the received power is calculated as Where γ is the random variable, The amount of power received at the secondary node due to the transmissions from all m MU's is calculated as Where d j -the distance between SU and j th MU. The PDF of received power is calculated as

Neyman-Pearson Criterion
To analyze the impact of the PUE attacks on CR network, Neyman Pearson composite hypothesis testing mathematical model are used. In this method, the PDF of the received signal due to transmission by the primary and the MU at the SU is calculated. In order to calculate the decision variable, the power of the signal received from the source is measured and is given by the ratio term, ʌ Where, P (x) is the received power from the MU calculated using (5) and The ratio term is then compared with predefined threshold and from that the SU decides the following Λ ≤ λ : D 1 : Primary user transmission is taking place Λ ≥ λ : D 2 : PUE Attack is in progress There can be two possibilities • The SU may take the decision D 2 when M 1 is true.
• The SU may take the decision D 1 when M 2 is true.
The errors associated with this probability are termed as Missed Probability: P {D 2 |M 1 } = Probability of taking decision as D2 when the hypothesis M 1 is true.
False Alarm Probability: P {D 1 |M 2 } = Probability of taking decision as D 1 when the hypothesis M 2 is true.

Improved Detection Scheme with Double Threshold
By using a conventional FB method which uses single threshold value, it is difficult to differentiate between PU signal and PUE attacker. The power level of the signal received at the SU receiver is measured. Then it is compared with that from the true PUs in order for a CR network to decide whether the signal originates from an attacker or not. The test statistics in the case of FB method is calculated as However, a FB detector is not dynamic enough to challenge an advanced PUE attack. Hence, by using a double threshold value called improved detection scheme is employed. The conditions which are used to differentiate between PU and PUE attacker is given as follows.
Test statistics < Threshold1, only noise Threshold1< Test statistics < Threshold2, Primary user Test statistics > Threshold2 PUE Attacks In order to distinguish a PUE attacker from a real PU, an improved detection scheme with two energy thresholds, denoted by and are used. Here, λ1 < λ 2 and is the earlier threshold in a conventional FB method. When the received signal energy, E < λ 1 , then it is said that there is no PU or PUE attacker. If the signal energy is between the two threshold λ 1 andλ 2 , then there is presence of PU. When the signal energy is above λ 2 , then it indicates PUE attack.

Location Verification Method
Location verification is another defense technique which is used to distinguish PU and PUE attacker based on the distance. There are two location verification schemes -Distance Ratio Test and Distance Difference Test. To perform the analysis some assumptions are considered. The assumptions are, • All users should broadcast information regarding their location. • All the users have a predefined transmit power level and are known to each other. A ground reflection model is used for calculating the power level of the received signal. The received power level is given by equation (8).
Where, =Transmit power level, h t = Height of a transmitter h r = Height of a receiver G t = Transmitter antenna gain G r = Receiver antenna gain L = System loss factor d = Distance between transmitter and receiver Here, the distance between a SU and other users is calculated based on the coordinates location and received power level. The user is a truthful user, when the distance calculated from both these technique matches. Otherwise, it is a malicious user.

Distance Estimate based on Location Coordinates
The distance between the users can be calculated based on the location coordinates. Let us assume (x, y) is x and y coordinates of a SU and (x1,y1) is x and y coordinates of an existing primary transmitter. The distance between them d, is given as The simulation is performed based on the assumption that the location coordinates of all users are broadcasted. Hence, using eq (9) it is possible to calculate the distance between any users.

Distance Estimate based on Received Power Level
In Received Signal Strength (RSS) approach, the parameter on which the received power level is depended is the transmitting power and distance on the path between two devices. Knowing the transmit power level; it is possible to calculate the distance between users from the measured received power level. Equation (10) gives the received power level, P r , for a specific transmit power level P t . Assume h t , h r , G t , G r and L are constant and equal to one. Therefore, the received power level can be expressed as a function of transmit power level and distance.
The distance between the user can be approximated based on the received power level and known transmit power level as 4 t r P d P = (11) If d 1 , the distance estimated from location coordinates and d 2 , the distance estimated from the received signal power is approximately equal, then it is decided that the signal is from the legitimate PU. Otherwise, a malicious user is in attack. The accuracy of distance calculation will depend on the presence of noise level in the signal received. However, statistically, the distance d 1 and d 2 calculated should come close.

Proposed Approach for PUE Attack Detection
The performance of existing method which employs double threshold for FB based spectrum sensing technique to detect the presence of the PU is pitiable at low signal to noise ratio. At the same time the location verification method also gives poor performance at low SNR. In order to detect PU's in the presence of PUE attack, a new method is suggested which combines improved detection scheme and location verification method. This approach is called fusion approach. Hence the probability of detecting the presence of PU can be increased. The malicious user cannot mimic the PUs location coordinates and transmit power level. Therefore, the authenticity of PU signal can be confirmed based on the distance measurement taken out from location coordinates and the signal power level received. Three different combinations are possible by using improved detection method and location verification method. They are discussed below.

AND Rule Logic
The logic is to perform both improved detection scheme and location verification method simultaneously. If both the method says that the signal is from the attacker then the final decision is made that, it is from the attacker. The flow chart for AND rule logic is as shown in Figure 5.

OR Rule Logic
The flowchart for OR rule logic is as shown in Figure 6. Here, both improved detection scheme and location verification method are carried out simultaneously. If any one of the method say that the signal is from the attacker, the final decision is made that the signal is from the attacker.

Alternate Logic
The improved detection scheme is performed initially and local decision (d) is made. If d = 0, the final decision is made that the signal is from the PUE attacker. Otherwise, it will perform location verification and match the estimated location with the content in the database. Here if d = 0, the final decision is made that the signal is from the PUE attacker. Else the final decision is made that the signal if from real PU. The flow chart for the alternate logic is as shown in Figure 7.

Result and Discussions
The proposed method was simulated in Matlab. In order to simulate FB based spectrum sensing method, the values for the system parameters considered are listed in Table 1. The values of system parameters considered for simulating PUE attacks in CR network is listed in Table 2. It is difficult to differentiate between PU signal and PUE attacker when a conventional FB method with single threshold value is used. Therefore, an improved detection scheme which uses a double threshold value is employed in order to detect PUE attackers. Two threshold values are used instead of a single threshold to make local decision in the double threshold FB method. Double threshold FB offers a benefit over conventional FB with respect to bandwidth. The thresh-old1 value is calculated using equation (1) and threshold 2 is decided such that the difference between thresholds is minimum or no decision region shown in Figure 8 should be less. For the reference purpose, ED method is used in order to compare the performance of the FB method in the presence of attackers. Figure 9 shows the comparison between ED method and FB method in the presence of a PUE attacker using improved detection scheme. The FB method shows better performance than the ED method. The probability of correctly detecting the presence of the attacker using the location verification method is shown Figure 10. This method works only, when the distance coordinates and the transmit power of the primary transmitter is known to all the users.  The performance of the improved detection scheme and location verification method is poor in low signal to noise ratio values. So by combining the improved detec-tion scheme with the location verification method the probability of detecting the presence of the PU in the presence of the attacker can be increased. There exist three different combinations when we fuse these two methods. The probability of detection versus SNR using ED method and FB method in the presence of attackers incorporating AND rule logic is plotted in figure 11. It basically illustrates how well the AND logic rule, introduce a better defense technique towards PUE attacks. Also, it explains that the FB technique along with AND rule logic works jointly to demolish PUE attacks better than ED technique.
The final decision of PUE attack is made when any one of the method says that the signal is from PUE attacker and this is called OR rule logic. OR rule logic was performed for both ED and FB method and it was found that again FB method has a better probability of detection in the presence of attackers. Figure 12 shows the PUE attack defense in the ED and FB techniques using OR Rule Logic. Comparing Figure 11 and Figure 12, it is shown that the employment of OR Rule logic introduce better probability of detection at respective SNR. For example, at SNR = 1 dB, the probability of detection of FB method is 0.28 when incorporating AND rule logic and it is 0.5 when incorporating OR rule logic. The percentage increase in probability of detection of a PU signal in the presence of PUE attacks is 78.6. As stated earlier, an alternate approach for defending PUE attack is to perform improved detection scheme at first and decide the presence of attacker based on the out-Vol 10 (1) | January 2017 | www.indjst.org Sabiq P.V. and D. Saraswady come. If unable to decide the presence of the attacker, then go for location verification method and make the final decision on the presence of PUE attacks. Using alternate logic the performance of ED method and FB method in the presence attackers is shown in Figure 13. Again the FB method is having better performance than ED method. The probability of detection has increased to 0.98 in FB method at SNR = 1dB. Figure 13. PUEA Defense in ED and FB method using alternate logic.
A comparison of different defense techniques discussed in this paper is examined for FB method alone and a graph between probability of detection versus SNR is plotted as shown in figure 14. By employing alternate logic there is a better probability of detection in the presence of attackers. Hence, by integrating alternate logic with FB method we can shield PUE attacks in a CR network.
A numerical comparison between ED method and FB method when integrating different defense techniques is enumerated in Table 3. The table highlights the probability of detection versus SNR for different defense techniques. From the table it shows that FB method integrated with alternate logic will yield better probability of detection when compared to other combination. For example, at SNR = 1 dB, the increase in probability of detection when integrating FB method with improved detection and alternate logic is 0.75. Similarly, at SNR = 5 dB, the increase in probability of detection is 0.1. It also states that the FB method produces a better probability of detection at low SNR when compared to ED method. Hence FB method integrated with the alternate logic can be deployed in CR network for shielding PUE attacks.

Conclusions and Future Work
In this work, PUE attacks on a CR network are studied. The spectrum sensing techniques employed here is ED method and FB method. Firstly, a double threshold scheme was integrated with each spectrum sensing


s2fieldsofstudy:
Computer Science
Mathematics

text:
Graph partitioning using matrix differential equations

Given a connected undirected weighted graph, we are concerned with problems related to partitioning the graph. First of all we look for the closest disconnected graph (the minimum cut problem), here with respect to the Euclidean norm. We are interested in the case of constrained minimum cut problems, where constraints include cardinality or membership requirements, which leads to NP-hard combinatorial optimization problems. Furthermore, we are interested in ambiguity issues, that is in the robustness of clustering algorithms that are based on Fiedler spectral partitioning. The above-mentioned problems are restated as matrix nearness problems for the weight matrix of the graph. A key element in the solution of these matrix nearness problems is the use of a constrained gradient system of matrix differential equations.

1. Introduction. In this paper we present a novel approach to partitioning a connected weighted undirected graph. We consider the Frobenius-norm minimum cut problem and allow for constraints such as prescribing the minimum cardinality of connected components or assigning a priori selected vertices to a component. We use spectral graph theory as pioneered by Fiedler [6], see also the monograph by Chung [3] and the introductory articles [16,18]. We formulate and use a gradient system of matrix differential equations to drive the smallest nonzero eigenvalue of the graph Laplacian to zero. Once this eigenvalue becomes zero, the graph is disconnected and the corresponding eigenvector indicates the membership of vertices to the connected components. This approach can be extended to other partitioning problems beyond the constrained minimum cut problems considered here.
The approach of this paper takes basic ideas and techniques of recent algorithms for eigenvalue optimization via differential equations, as given for example in [9,8,11,10], to another application area. A common feature is a two-level procedure, where on the inner level a gradient flow drives perturbations to the original matrix of a fixed size into a (local) minimum of a functional that depends on eigenvalues and possibly eigenvectors, and in an outer iteration the perturbation size is determined such that the functional becomes zero. As with the previous algorithms cited above, the algorithms presented here cannot guarantee to find the global minimum of a non-smooth, non-convex optimization problem, or of an NP-hard combinatorial optimization problem. There are cases where our algorithm could get stuck in a local minimum, and we will present a contrived example where this happens. Even with this caveat, the presented algorithm performs remarkably well in the examples from the literature on which we have tested it.
As opposed to combinatorial algorithms, the algorithm presented here modifies all weights of the graph as it proceeds, and only in the end arrives at the cut and the unchanged remaining weights.
The proposed algorithm is an iterative algorithm, where in each step the second eigenvalue and the associated eigenvector of the Laplacian of a graph with perturbed weights are computed.
In the cardinality-or membership-constrained cases, additionally a linear system with an extended shifted Laplacian is solved in each step. For a large sparse connected graph (where the number of edges leaving any vertex is moderately bounded), these computations can be done in a complexity that is linear in the number of vertices. In the known (unconstrained) minimum cut algorithms, the computational complexity is at least quadratic [17]. It is thus conceivable that for large sparse connected graphs, the proposed iterative algorithm can favorably compete with the classical unconstrained minimum cut algorithms. In constrained cases, it appears that the computational complexity is even more favorable in comparison with the existing heuristic combinatorial algorithms as proposed in [2]. However, as of now no detailed comparisons of the relative merits of the conceptually and algorithmically fundamentally different approaches have been made.
In Section 2 we formulate the Frobenius-norm minimum cut problem and its cardinalityand membership-constrained variants. This is stated as a matrix nearness problem where it is asked how far, with respect to the Frobenius norm, the weight matrix of the given graph is from that of some disconnected graph which should possibly satisfy additional constraints. We give basic notation and recall Fiedler's theorem on graph connectivity. We also formulate an ambiguity problem where it is asked how far the given weight matrix is from the weight matrix of a graph for which the second and third eigenvalues of the graph Laplacian coalesce and for which therefore graph partitioning based on the Fiedler vector (the eigenvector to the second eigenvalue) becomes ambiguous.
In Section 3 we describe the two-level approach to the unconstrained Frobenius-norm minimum cut problem. This is the central section of the paper, where the basic approach is developed.
In Section 4 we extend the approach to the cardinality-and membership-constrained minimum cut problems, and in Section 5 we extend it to the ambiguity problem.
In Section 6 we describe algorithmic aspects such as the discretization of the norm-and inequality-constrained gradient flow, the choice of initial values, and stopping criteria. In particular, since it is known beforehand that the weights of the cut graph are either zero or those of the original graph, the iteration need not be carried out to full convergence. Section 7 shows numerical results of the proposed algorithm for some graphs taken from the literature.

Preparations and problem formulation.
2.1. The Frobenius-norm minimum cut problem. Consider a graph with vertex set V = {1, . . . , n} and edge set E ⊂ V × V. We assume that the graph is undirected: with (i, j) ∈ E, also (j, i) ∈ E. With the undirected graph we associate weights w ij for (i, j) ∈ E, such that The graph is connected if for all i, j ∈ V, there is a path (i 0 , i 1 ), (i 1 , i 2 ), . . . , (i −1 , i ) ∈ E of arbitrary length , such that i = i 0 and j = i and w i k−1 ,i k > 0 for all k = 1, . . . , .
The problem considered in this paper is the following: Given a connected weighted undirected graph with weights w ij , we aim to find a disconnected weighted undirected graph with the same edge set E and modified weights w ij such that The solution to this matrix nearness problem is the same as that of finding a cut C, i.e., a set of edges that yield a disconnected graph when they are removed from E, where the cut C is such that When the weights are replaced by their square roots, so that w ij instead of w 2 ij appears in the above sum, this becomes the classical minimum cut problem, for which algorithms with complexity O(|V| 2 log |V| + |V| · |E|) exist; see Stoer & Wagner [17] and references therein.

2.2.
Constrained minimum cut problems. The above problem will further be considered with additional constraints. In particular, we consider the following cases: • Membership constraint: It is required that a given set of vertices V + ⊂ V is in one connected component and another given set of vertices V − ⊂ V is in the other connected component. • Cardinality constraint: It is required that each of the connected components has a prescribed minimum number n of vertices. It is known that cardinality constraints make the problem NP-hard 2.3. Graph Laplacian and algebraic connectivity. Setting w ij = 0 for (i, j) / ∈ E, we have the symmetric weight matrix The degrees d i = n j=1 w ij are collected in the diagonal matrix The Laplacian matrix L = Lap(W ) is defined by We note that by the Gershgorin circle theorem, all eigenvalues of L are nonnegative, and L1 = 0, so that λ 1 = 0 is the smallest eigenvalue of L. Remarkably, the connectivity of the graph is characterized by the second-smallest eigenvalue of L. Theorem 2.1 (M. Fiedler [6]). Let W ∈ R n×n be the weight matrix of an undirected graph and L the corresponding Laplacian matrix. Let 0 = λ 1 ≤ λ 2 ≤ . . . ≤ λ n be the eigenvalues of L. Then, the graph is disconnected if and only if λ 2 = 0. Moreover, if 0 = λ 2 < λ 3 , then the entries of the corresponding eigenvector orthogonal to 1 assume only two different values, of different sign, which mark the membership to the two connected components.
Because of this result, the second smallest eigenvalue λ 2 of L is called algebraic connectivity of W . If λ 2 is a simple eigenvalue, then the corresponding eigenvector is known as the Fiedler vector.

2.
4. An ambiguity problem in graph partitioning. Based on Theorem 2.1, a common and computationally inexpensive strategy for partitioning a graph is to compute the Fiedler vector and to partition the graph according to the values of its entries. This becomes unreliable when a small perturbation of the weights yields a coalescence of the eigenvalues λ 2 and λ 3 . It is then interesting to know the distance of the given weight matrix from the set of weight matrices with λ 2 = λ 3 .
3. Two-level method for the Frobenius-norm minimum cut problem.
3.1. Two-level formulation. Our approach can be summarized as follows: 1. Given ε > 0, we look for a symmetric matrix E = (e ij ) ∈ R n×n with the same sparsity pattern as W (i.e., e ij = 0 if w ij = 0), of unit Frobenius norm, with W + εE ≥ 0 (with componentwise inequality) such that the second smallest eigenvalue of Lap(W + εE) is minimal. The obtained minimizer is denoted by E(ε). 2. We look for the smallest value of ε such that the second smallest eigenvalue of Lap(W + εE(ε)) equals 0.
In order to compute E(ε) for a given ε > 0, we make use of a constrained gradient system for the functional F ε (E) = λ 2 Lap(W + εE) , (3.1) under the constraints of unit Frobenius norm and W +εE ≥ 0 and the symmetry and the sparsity pattern of E.
In the outer iteration we compute the optimal ε, denoted ε , by a combined Newton-bisection method.
The algorithm computes a partition of the graph as provided by the Fiedler vector corresponding to the weight matrix W + ε E(ε ). This is not guaranteed to yield a global optimum for the Frobenius-norm minimum cut problem, since the gradient flow might converge only to a local minimum. In any case, it provides an upper bound for the distance problem (2.1).

3.2.
Constrained gradient flow for the functional F ε .
3.2.1. Eigenvalue derivatives. We will use the following standard perturbation result for eigenvalues; see, e.g., [12, Section II.1.1]. Here and in the following, we denote˙= d/dt. Lemma 3.1. Consider the differentiable symmetric n × n matrix valued function C(t) for t in a neighborhood of 0. Let λ(t) be an eigenvalue of C(t) converging to a simple eigenvalue We denote by · = · F the Frobenius norm on R n×n and by X, Y = trace(X T Y ) the corresponding inner product.
We return to the situation of the previous section. For a set of edges E, we define P E as the orthogonal projection from R n×n onto the sparsity pattern determined by E: for A = (a ij ), For a fixed given weight matrix W and for ε > 0, we call a matrix E = (e ij ) ∈ R n×n ε-feasible if the following conditions are satisfied: Consider now a regular path E(t) of ε-feasible matrices, and denote the corresponding Laplacian matrix by L(t) = Lap(W + εE(t)) and by λ 2 (t) the second smallest eigenvalue of L(t). Lemma 3.1 applied to the Laplacian matrix L(t) yields (omitting the argument t) where x(t) is a corresponding eigenvector of unit Euclidean norm. Next we rearrange (3.2) to an equationλ 2 = ε G ε (E),Ė with an appropriate matrix-valued function G ε , which is the gradient of F ε in the space of symmetric matrices with sparsity pattern E. In the following, Sym(A) = 1 2 (A + A T ) denotes the symmetric part of a quadratic matrix A, and we write x 2 = (x 2 i ) ∈ R n for the vector of squares of the entries of x = (x i ) ∈ R n . Lemma 3.2. In the above situation we havė is symmetric and has the sparsity pattern determined by the set of edges E.
Proof. We note thaṫ Combining (3.2) and (3.4), we obtaiṅ The second term is already in the desired form. We obtain for the first term SinceĖ and xx T are symmetric, this can be rewritten aṡ We then haveλ This is in the desired form: G ε (E) is symmetric and has the sparsity pattern E.

Admissible directions.
Since E(t) is of unit Frobenius norm by condition (i), we have Condition (iv) requires thatė ij ≥ 0 for all (i, j) ∈ E 0 , where E 0 = E 0 (εE) is the set of cut edges defined by Conditions (ii) and (iii) are satisfied if the same holds forĖ. These four conditions are in fact also sufficient for a matrix to be the time derivative of a path of ε-feasible matrices. Hence, for every ε-feasible matrix E, a matrix Z = (z ij ) ∈ R n×n is the derivative at t = 0 of some path of ε-feasible matrices starting at E if and only if the following four conditions are satisfied:

Admissible direction of steepest descent.
To determine the admissible directioṅ E of steepest descent from E, we therefore consider the following optimization problem for G = G ε (E): min Z G, Z subject to (i')-(iv') and Z, Z = 1. (3.7) The additional constraint Z = 1 just normalizes the descent direction. Problem (3.7) has a quadratic constraint. We now formulate a quadratic optimization problem with linear constraints, which is equivalent in the sense that it yields the same descent direction, provided that a strict descent direction exists, i.e., satisfying G, Z < 0 and the constraints (i')-(iv'). This is based on the fact that when G, Z < 0, there exists a scaling factor α > 0 such that G, αZ = −1. Consider the following problem: min Z Z, Z subject to (i')-(iv') and G, Z = −1.
In addition, there are conditions (i')-(iv'). It can be shown (see [5, Lemma 3.2.1]) that the symmetry and sparsity conditions (ii') and (iii') need not be imposed, but are consequences of conditions (ii) and (iii) on E and the corresponding properties of G = G ε (E).
3.2.5. Constrained gradient flow. The gradient flow of F ε under the constraints (i)-(iv) is the system of differential equationsĖ where Z(t) solves the KKT system (3.9) with G = G ε (E(t)) under the constraints (i')-(iv') with the set of edges E 0 (t) = E 0 (εE(t)). Lemma 3.3. On an interval where E 0 (t) does not change, the gradient system becomes, with P + = P E\E0 and omitting the ubiquitous argument t, Proof. The positive Lagrange multipliers µ ij > 0 just have the role to ensure thatė ij = 0. With G = G ε (E), the gradient system therefore readṡ where κ is determined from the constraint E,Ė = 0. We then have and the result follows.
In a numerical solution of the gradient system, we thus have to monitor the sets of edges where w ij + εe ij = 0 and among them further those edges where the sign of −g ij − κe ij changes. When the active set is changed, then also κ changes in a discontinuous way. Let κ − and κ + be the values of κ before and after the event of discontinuity, respectively. Then one has generically g ij + κ − e ij > 0 after the event for the critical edge (i, j), but the sign of g ij + κ + e ij may be positive or negative. In the first case, (i, j) leaves E 0 . In the latter case, only a generalized solution in the Filippov sense exists, which keeps (i, j) ∈ E 0 , i.e., w ij + εe ij = 0. This is enforced until g ij + κ + e ij changes sign. From a practical perspective, it appears reasonable just to keep (i, j) ∈ E 0 for all future time once it has entered E 0 , which means that a cut of an edge is made irreversible.
3.2.6. Monotonicity and stationary points. The following monotonicity result follows directly from the construction of the gradient system. Theorem 3.4. Let E(t) of unit Frobenius norm satisfy the differential equation (3.12) with G ε (E) of (3.3). Then, the second smallest eigenvalue λ 2 (t) of the Laplacian matrix Lap(W + εE(t)) decreases monotonically with t:λ 2 (t) ≤ 0.
Equilibrium points of (3.12) are characterized as follows. Theorem 3.5. The following statements are equivalent along solutions of (3.12): . Proof. Using Lemma 3.2 and (3.12) we obtain, with G = G ε (E), With the strong form of the Cauchy-Schwarz inequality, the result follows.
3.3. Newton-bisection outer iteration. Let E(ε) denote the minimizer of the functional F ε . In general we expect that for a given perturbation size ε < ε , the eigenvalue λ 2 (W + εE(ε)) > 0 is simple. If so, then f (ε) = F ε (E(ε)) is a piecewise smooth function of ε and we can exploit its regularity to obtain a fast iterative method to converge to ε from the left. Otherwise we can use a bisection technique to approach ε .
The following result provides an inexpensive formula for the computation of the derivative of f (ε) = F ε (E(ε)), which will be useful in the construction of the outer iteration of the method.
Assumption 3.1. We assume that the second smallest eigenvalue of Lap(W + εE(ε)) is simple. Moreover, E(ε) is assumed to be a smooth function of ε in some interval, and the set of zero-weight edges E 0 related to E(ε) is independent of ε in the interval.
We denote again P + = P E\E0 . We then have the following result. Lemma 3.6. Under Assumption 3.1, the function f (ε) = F ε (E(ε)) is differentiable and its derivative equals (with = d/dε) Proof. Under Assumption 3.1, the projection P + remains constant near ε. This means that we are effectively working on a reduced set E = E \ E 0 of edges. We set G + (ε) = P + G ε (E(ε)) and decompose εE(ε) = εP + E(ε) + R, where R = (I − P + )εE(ε) = −(I − P + )W = −P E0 W, since w ij + εe ij (ε) = 0 for all (i, j) ∈ E 0 . In particular, R is independent of ε. Differentiating f (ε) = F ε E(ε) with respect to ε we obtain f (ε) = G ε (E(ε)), P + E(ε) + εP + E (ε) = G + (ε), P + E(ε) + εP + E (ε) . (3.14) The conservation of E(ε) = 1 and of R for all ε implies Now we use the property of minimizers as stated by Theorem 3.5, which gives us From (3.14) we thus obtain the stated formula, since f (ε) ≤ 0. For ε = ε k < ε , we make use of the standard Newton iteration In a practical algorithm it is useful to couple the Newton iteration (3.15) with a bisection technique. To do this we adopt a tolerance tol which allows us to distinguish whether ε < ε , in which case we may use the derivative formula and perform the Newton step, or ε > ε , so that we have to make use of bisection. The method is formulated in Algorithm 1.

4.
The two-level method for the membership-and cardinality-constrained minimum cut problems. 4.1. Functional for the membership-constrained minimum cut problem. Our approach to the membership problem is the same two-level procedure as for the unconstrained minimum cut problem, except that the functional (3.1) is replaced by the following functional: For ε > 0 and a matrix E of unit Frobenius norm, let x = (x i ) ∈ R n be the eigenvector to the second smallest eigenvalue λ 2 of Lap(W + εE). Let V − and V + be the set of indices whose membership to different components of the cut graph is prescribed. Let 0) collect the negative and positive components of x, respectively. Let n − and n + be the numbers of negative and nonnegative components of x, respectively. We denote the averages of x − and x + by x + i .

Algorithm 1: Newton-bisection method for distance approximation
Data: Matrix W is given, k max (max number of iterations), tol (tolerance) ε 0 , ε lb and ε ub (starting values for the lower and upper bounds for ε ) Result: ε (upper bound for the distance), E(ε ) begin 1 Compute E(ε 0 ) by the inner iteration Return ε k+1 and the interval [ε lb , ε ub ] Stop else Set k = k + 1 5 Compute E(ε k ) by the inner iteration 6 Return ε = ε k Motivated by the special form of the eigenvectors as given in the Fiedler theorem (Theorem 2.1), we consider the functional where α > 0 is a weight to be chosen. The choice of the sign of the eigenvector x is such that F ε (E) takes the smaller of the two possible values. This functional is to be minimized under the inequality constraints W + εE ≥ 0, the norm constraint E = 1 and the symmetry and the sparsity pattern of E.

Functional for the cardinality-constrained minimum cut problem.
For the cardinality-constrained problem we use the same functional F ε , except that the sets V − and V + are not given a priori, but are chosen depending on E in the following way: V − and V + collect the indices of the smallest and largest n components of the eigenvector x, respectively, augmented by those indices for which the components of x do not differ by more than a threshold δ from the average of the smallest and largest n components, respectively.  x(t) be the associated eigenvector normalized such that x(t) 2 = 1. Moreover, let M (t) = C(t) − λ(t)I and let M (t) † be the Moore-Penrose pseudoinverse of M (t). Then, the derivative of the eigenvector is given byẋ We remark that in [15] this is formulated with the group inverse, which in the symmetric case is the same as the Moore-Penrose pseudoinverse.
4.3.2. Gradient of F ε . Consider a differentiable path E(t) of ε-feasible matrices, and denote the corresponding Laplacian matrix by L(t) = Lap(W + εE(t)), by λ 2 (t) the second smallest eigenvalue of L(t) , and by x(t) the associated eigenvector. We set and, with e i denoting the ith standard unit vector, We define which is computed as the solution of the linear system We denote by x • y = (x i y i ) the vector obtained by componentwise multiplication of the entries of x and y. We then have the following result. Lemma 4.2. In the above situation we have is symmetric and has the sparsity pattern determined by the set of edges E. Proof. We have d dt and similarly for the sum over V + . With K = (L − λ 2 I) † we obtain from Lemma 4.1 thaṫ Using the expression forL given in (3.4) and proceeding as in the proof of Lemma 3.2 gives the result.
The computational cost of computing G ε (E) lies in computing the second eigenvalue and its eigenvector and in solving the linear system (4.3). For a sparse weight matrix, these computations have a complexity that is linear in the number of vertices. With this gradient G ε (E), the further procedure is now exactly the same as in Section 3.

5.
The two-level method for the ambiguity problem.

Two-level formulation.
For the ambiguity problem we proceed similarly as in Section 3.
1. Given ε > 0, we look for a symmetric matrix E ∈ R n×n with the same sparsity pattern as W (i.e., e ij = 0 if w ij = 0), of unit Frobenius norm, with W + εE ≥ 0 (with componentwise inequality) such that the difference between the third and second smallest eigenvalues of Lap(W + εE) is minimized. The obtained minimizer is denoted by E(ε). 2. We look for the smallest value of ε such that the second and third eigenvalues of Lap(W + εE(ε)) coalesce. In order to compute E(ε) for a given ε > 0, we make use of a constrained gradient system for the functional under the inequality constraints W + εE ≥ 0, the norm constraint E = 1 and the symmetry and the sparsity pattern of E.
In the outer iteration we compute the optimal ε, denoted ε , by a combined Newton-bisection method as in Section 3.

5.2.
Gradient of F ε . Consider a regular path E(t) of ε-feasible matrices, and denote the corresponding Laplacian matrix by L(t) = Lap(W +εE(t)) and by λ 2 (t) and λ 3 (t) the second and third smallest eigenvalues of L(t), respectively. We denote by x(t) and y(t) the corresponding eigenvectors of unit Euclidean norm.
Lemma 5.1. In the above situation we have is symmetric and has the sparsity pattern determined by the set of edges E.
With this gradient we then proceed further as in Section 3.
6.1. Discretizing the constrained gradient flow. We use a modified explicit Euler method for the approximate integration of the differential equation (3.10). For a given ε > 0, a step-size h > 0 and from the ε-feasible perturbation matrix E n of the nth time step, we compute E n+1 as follows. We compute G ε (E n ) = g n ij and define E n+1 with θ ∈ [0, 1) such that w ij + ε e n ij − θhg n ij = 0, that is, with θ = (w ij + εe n ij )/(εhg n ij ). We would ideally take the new perturbation matrix E n+1 such that E n+1 − E n+1 → min subject to E n+1 = 1 and W + εE n+1 ≥ 0.
We approximate this optimization problem by treating the two constraints one after the other in an alternating way. Let E 0 be the set of edges (i, j) ∈ E for which w ij + εẽ n+1 ij = 0 (cut edges), and let P 0 = P E0 and P + = P E\E0 be the complementary projections as defined in Section 3.2.2. We first normalize by choosing ρ > 0 such that has unit Frobenius norm, i.e., (In case that P 0 E n+1 is larger than 1 or very close to 1, we replace P 0 E n+1 by P 0 E n in the two lines above.) We denote by E − the set of inadmissible edges (i, j) for which w ij + εê n+1 ij < 0. We then resetê n+1 augment E 0 := E 0 ∪ E − and consider the updated projection P 0 = P E0 . We then normalize the so obtained matrix E n+1 in the same way as above by leaving the entries of P 0 E n+1 unchanged, reset the entries for inadmissible edges, normalize, and so on. As there are only finitely many edges, this iteration terminates after finitely many steps (typically after the first step). Finally, we have obtained an ε-feasible perturbation matrix E n+1 .
6.2. Choice of step-size. The step-size h can, for example, be selected by the following adaptive algorithm. Here the objective is to reduce the function F ε , not to follow accurately a trajectory of the constrained gradient differential equation.
6.3. Stopping criterion. Let F n = F ε (E n ). In order to stop the integration when F n has approximately reached a stationary value, we use a criterion of the following type: where tol is a tolerance parameter (e.g., tol = 10 −6 ) and β and δ are further parameters. We had good experience with the choice β = 10 · tol and δ = tol/100; see further [5] where also the sensitivity of the algorithm to the chosen parameters is discussed.
6.4. Initial value of the constrained gradient flow for a new ε. When we change to a new value of ε in the outer iteration, we need an initial value for the constrained gradient flow. A first idea might be to take the terminal perturbation matrix E 0 = E(ε old ) as the initial value, but usually this does not satisfy the nonnegativity constraints W + εE 0 ≥ 0 if ε ≥ ε old . We therefore modify E 0 to E 0 by solving approximately E 0 − E 0 → min subject to E 0 = 1 and W + εE 0 ≥ 0, alternating between normalization and enforcing the nonnegativity constraints as in Section 6.1, but this time beginning with the empty set E 0 = ∅. Compute E n+1 (h) and its function value F ε (E n+1 (h)) with the step-size h if F ε (E n+1 (h)) ≥ F ε (E n ) then halve the step-size, h := h/2 and repeat from 2 else if h = h n−1 then compute E n+1 (2h) and its function value Return E n+1 and h n 6.5. Choice of the inital perturbation size ε 0 and the initial perturbation matrix. While one might just start with a random perturbation, a more educated guess starts from the normalized free gradient E 0 = −G ε (0)/ G ε (0) and determines ε 0 as the largest number ε such that W + εE 0 ≥ 0. 6.6. Stopping the outer iteration before convergence. In the exact solution W = W + ε E to the constrained minimum cut problem, the entries of W are either zero or those of W . To decide about the cut, it is therefore not necessary to iterate towards W with very high accuracy, but instead the cut can be inferred earlier from a moderately accurate approximation W + εE(ε), for example using the following criterion, with a small threshold parameter ϑ > 0: Stop if for every edge (i, j) ∈ E, either w ij + εe ij ≤ ϑw ij or |εe ij | ≤ ϑw ij .
In the first case one would then cut to w ij = 0 and in the second case one would leave the weight unchanged: w ij = w ij . It can finally be checked if the so obtained cut graph is indeed disconnected, by computing λ 2 (Lap(W )). Instead, in the ambiguity problem such a shortcut is not feasible.

Numerical examples.
We consider a few illustrative examples for both the cardinality and the membership constraints. At the end we will also consider graphs to illustrate the ambiguity problem.
The standard Fiedler spectral partitioning algorithm, to which we refer below, is simply based on the sign of the components of the eigenvector of Lap(W ) associated to λ 2 .
Example 1 (Zachary's karate club). This weighted graph consisting of 34 vertices describes the relationship between 34 members of a karate club (for a detailed description see [19]). Using the Fiedler spectral partitioning we obtain two connected components of 16 and 18 vertices as shown in Figure 7.1a.
According tho the standard Fiedler partitioning the first component is led by the vertex labeled as 1, while the second one is led by the vertex labeled as 34.
We next consider the following constraints: (i) Cardinality constraint with threshold equal ton = 17 vertices; (ii) Membership constraint; (iii) Both constraints. In more detail: (i) By asking for a cardinality constraint withn = 17 vertices in each component, the approximate computed distance is ε ≈ 10.03631. The results in Table 7.1 are obtained by setting a tolerance tol = 10 −5 and the weight α = 3 in (4.1). With respect to the standard partition obtained by the Fiedler eigenvector of Lap(W ), the vertex that changes partition is the vertex labeled as 9.  (iii) Finally we consider both constraints, that is, we ask for a cardinality constraint with thresholdn = 17 vertices, and we require that vertex 9 is in the same connected component as the vertex 34. Table 7.3 shows values of ε k , f (ε k ) = F ε k (E(ε k )) computed with the tolerance tol = 10 −5 and the weights α c = 3, α m = 10 in the functional F ε that combines the cardinality and membership functionals. We obtain that the vertex 10 further changes the connected component in order to satisfy the cardinality constraint. Example 2 (A misbehavior of the algorithm). We present an example where the algorithm fails.
Consider an unweighted graph with N vertices, such that each vertex 2, . . . , N is connected to the following two vertices, i. e.
The first vertex is connected to the second one, i.e., w 1,2 = 1, but not to the third one. It is clear that the minimum cut is obtained by removing the edge (1, 2), since this is the only possibility to obtain a disconnected graph when only one edge is removed. We tried to solve this example for different values of N . The algorithm works correctly when N = 8 but fails when N ≥ 12.
In the latter case we obtain a disconnected graph, but with the wrong edges removed. For  [14]. Two books are linked if they were purchased by the same person, and the books are colored red, blue or green, based on book buying data (see Fig.(7.3b)). The links determine the grouping and coloring of the vertices. The vertices are colored by Fiedler spectral partitioning in Figure 7.3a. By this partitioning we obtain two connected components of 52 and 53 vertices. We ask for the   Table 7.4) for each one of the vertex in the green group of Figure 7.3b to belong to the blue group (ε b ) or to belong to the red group (ε r ). The results reported in the tables are obtained with the weight α = 1. Example 4 (Les Miserables). Figure 7.4 shows the graph of character co-occurence in Les Miserables [13]. This graph consists of 77 vertices (representing characters). According to the Fiedler partitioning, 22 of these belong to one part and the remaining 55 belong to the other part. Asking for the partitioning of the graph with the cardinality constraint with threshold n = 35, we obtain the result shown in Figure 7  Example 5 (Planted Partition Model -ambiguity problem). We consider here a class of graphs for which we investigate the distance to ambiguity (that is, a coalescence of the second and third eigenvalues in the associated Laplacian matrix). We consider Planted Partition Models [4], which constitute a special case of Stochastic Block models, a commonly used generative model for social and biological networks. The probability matrix consists of a constant value p in on the diagonal and a different constant value p out off the diagonal; in addition, the number of vertices is n = 40, while the communities are 4, each one made of 10 vertices. In Figure 7.3 we can see two Planted Partition Models; on the left when p in = 0.8 and p out = 0.2, on the right when p in = 0.9 and p out = 0.1. In Table 7.5 a comparison among various values of p in and p out is shown. These results are obtained by setting the tolerance to tol = 10 −5 .  The results show that relatively small perturbations may yield an ambiguity in the Fiedler partitioning for these graphs.


extfieldsofstudy:
Mathematics
Computer Science

text:
Learning Semantic Maps with Topological Spatial Relations Using Graph-Structured Sum-Product Networks

We introduce Graph-Structured Sum-Product Networks (GraphSPNs), a probabilistic approach to structured prediction for problems where dependencies between latent variables are expressed in terms of arbitrary, dynamic graphs. While many approaches to structured prediction place strict constraints on the interactions between inferred variables, many real-world problems can be only characterized using complex graph structures of varying size, often contaminated with noise when obtained from real data. Here, we focus on one such problem in the domain of robotics. We demonstrate how GraphSPNs can be used to bolster inference about semantic, conceptual place descriptions using noisy topological relations discovered by a robot exploring large-scale office spaces. Through experiments, we show that GraphSPNs consistently outperform the traditional approach based on undirected graphical models, successfully disambiguating information in global semantic maps built from uncertain, noisy local evidence. We further exploit the probabilistic nature of the model to infer marginal distributions over semantic descriptions of as yet unexplored places and detect spatial environment configurations that are novel and incongruent with the known evidence.

Introduction
It is essential for a mobile robot to maintain a representation of spatial knowledge, a framework that organizes the understanding of the environment. Mobile robots have access to information at both local and global scale. Therefore, it is desirable for a representation to enable integration of knowledge across spatial scales and levels of abstraction with the help of discovered spatial relations. Topological maps are an established framework for representing spatial relations between local places that enables anchoring high-level conceptual information and easy access for a planning algorithm. As a result, several semantic mapping approaches rely on topological graphs as part of their representation and associate topological nodes with semantic place attributes (Friedman, Pasula, and Fox 2007;Zender et al. 2008;Pronobis and Jensfelt 2012).
In order to integrate the collected spatial knowledge, resolve ambiguities, and make predictions about unobserved places, such frameworks often employ structured prediction algorithms. Unfortunately, the relations discovered by a robot exploring a real-world environment tend to be complex and noisy, resulting in difficult inference problems. Topological maps are dynamic structures, growing as the robot explores its environment, and containing a different number of nodes and relations for every environment. At the same time, many approaches to structured prediction place strict constraints on the interactions between inferred variables to achieve tractability (Bach and Jordan 2002) and require that the number of output latent variables be constant and related through a similar global structure (Belanger and McCallum 2016). This makes them either inapplicable or impractical in robotics settings and require compromising on the structure complexity (Mozos et al. 2007), introducing prior structural knowledge (Friedman, Pasula, and Fox 2007) or making hard commitments about values of semantic attributes (Pronobis and Jensfelt 2012). These problems are not unique to robotics and often present themselves in other domains, such as computer vision .
In this paper, we present Graph-Structured Sum-Product Networks (GraphSPNs), a general probabilistic framework for modeling graph-structured data with complex, noisy dependencies between a varying number of latent variables. Our framework builds on Sum-Product Networks (SPNs) (Poon and Domingos 2011), a probabilistic deep architecture with solid theoretical foundations (Peharz et al. 2017). SPNs can learn probabilistic models capable of representing context-specific independence directly from high-dimensional data and perform fast, tractable inference on high-treewidth models. GraphSPNs learn template SPN models representing distributions over attributes of subgraphs of arbitrary complexity. Then, to perform inference for a specific, potentially expanding graph, they assemble a mixture model over multiple decompositions of the graph into sub-graphs, with each sub-graph modeled by an instantiation of an appropriate template.
We apply GraphSPNs to the problem of modeling largescale, global semantic maps with noisy topological spatial relations built by robots exploring multiple office environ-ments. We make no assumptions about the structure of the topological map that would simplify the inference over semantic attributes. Our approach is capable of disambiguating uncertain and noisy local information about semantic attributes of places as well as inferring distributions over semantic attributes for yet unexplored places, for which local evidence is not available. Furthermore, we illustrate the benefits of the probabilistic representation by relying on likelihood of global semantic maps to detect novel and incongruent spatial configurations. We compare the performance of our model with the traditional approach based on Probabilistic Graphical Models assembled according to the structure of the topological graph. We show that GraphSPNs significantly outperforms Markov Random Fields built from pairwise and higher-order potentials relying on the same uncertain evidence.

Related Work
Probabilistic graphical models (PGMs) (Koller and Friedman 2009) provide a flexible framework for structured prediction. PGMs express the conditional dependence between random variables in terms of a graph. This results in a factorization of the joint distribution into a normalized products of factors, with factors defined over subsets of variables. The factors can be considered templates from which a distribution over an arbitrary number of variables can be assembled. As a result, PGMs were used in the past for modeling graph-structured data (Friedman, Pasula, and Fox 2007;Pronobis and Jensfelt 2012). Unfortunately, inference in PGMs is generally intractable, with the exception of low treewidth models (Bach and Jordan 2002). Even models with pairwise potentials require approximate inference techniques with no guarantee of convergence, such as Loopy Belief Propagation (BP) (Murphy, Weiss, and Jordan 1999), when the graph structure involves loops (as in this work). Higher-order models can pose a challenge even to the approximate methods. Additionally, many distributions which can be represented compactly, cannot be represented in terms of the factorization performed by PGMs (Poon and Domingos 2011).
An increasing number of structure prediction approaches utilize deep architectures. Several approaches build on graphical models for representing dependencies between output variables, while relying on deep models to realize the factors Chen et al. 2015). Other types of models focus specifically on sequential data by combining local convolutions with classification (Collobert et al. 2011) or by building sequence to sequence predictors based on multiple LSTM layers (Sutskever, Vinyals, and Le 2014). Finally, fully convolutional models exist for such problems as pixel-wise segmentation (Shelhamer, Long, and Darrell 2017). Other deep models rely on deep architectures to define an energy function over candidate labels and then predictions are produced by using backpropagation to iteratively optimize the energy with respect to the labels. Structured Prediction Energy Networks (SPENs) (Belanger and McCallum 2016) take that approach and can capture dependencies that would lead to intractable graphical models. Unfortunately, many of the deep approaches are not probabilis-tic and are mostly applicable to data of the same global structure and number of output labels as the training examples.
Sum-Product Networks (SPNs) are unique in being a deep and probabilistic architecture capable of representing hightreewidth models that would result in unmanageable PGMs. Yet, they provide tractable exact inference by exploiting context-specific independence and determinism (Gens and Domingos 2012). In contrast to SPENs and approximate inference over PGMs, marginal and MPE inference in SPNs does not require iterations and can be achieved using a single pass through the network (Poon and Domingos 2011). SPNs are quite expressive and have been used to solve difficult problems in several domains (Poon and Domingos 2011;Gens and Domingos 2012;Amer and Todorovic 2015), including semantic place classification in robotics . These benefits resulted in several approaches to structured prediction using SPNs. In (Ratajczak, Tschiatschek, and Pernkopf 2014), Linear-chain Conditional Random Fields (LC-CRFs) extended with local factors modeled using SPNs were used to model sequences. In (Nath and Domingos 2015), Relational-SPNs were proposed that model graph-based relational data based on first-order logic. This method models graphs with potentially varying sizes by summarizing multiple variables with an aggregate statistic. In contrast, we directly model each output variable associated with nodes of the graph, and construct an SPN structure specific to each graph instance.
There have been numerous attempts to employ structured prediction to modeling semantic maps with topological spatial relations. In (Mozos et al. 2007), Hidden Markov Models were used to smooth sequences of AdaBoost classifications of place observations into semantic categories. (Friedman, Pasula, and Fox 2007) proposed Voronoi Random Fields (VRFs) which are CRFs constructed according to a Voronoi graph extracted from an occupancy grid map. VRFs utilize pairwise potentials to model dependency between neighboring graph nodes and 4-variable potentials to model junctions. In (Pronobis and Jensfelt 2012), Markov Random Fields were used to model pairwise dependencies between semantic categories of rooms according to a topological map. The categorical variables were connected to Bayesian Networks that reasoned about local environment features, forming a chain graph. This approach relied on a door detector to segment the environment into a topological graph with only one node per room. While these approaches are probabilistic, they rely on approximate inference using Loopy BP leading to problems with convergence (Friedman, Pasula, and Fox 2007). Moreover, in both cases, additional prior knowledge or hard commitments about the semantics of some places were employed in order to obtain a clean and manageable topological graph structure. In contrast, in this work, we rely on a graph built primarily to support navigation and execution of actions by the robot. Such graph provides a better coverage, but results in more noisy structure. Furthermore, we make no hard commitments about the semantics of the places at time of structure creation and defer such inference to the final model. Our experiments show, that under such conditions, graphical models with pairwise or higher-order potentials deteriorate quickly. Figure 1: An simple SPN for a naive Bayes mixture model P (X 1 , X 2 ), with three components over two binary variables. The bottom layer consists of indicators for X 1 and X 2 . Weighted sum nodes, with weights attached to inputs, are marked with +, while product nodes are marked with ×.

Preliminaries
We begin by giving a brief introduction to SPNs. For details, the reader is referred to (Peharz et al. 2017;Poon and Domingos 2011). Then, we describe the topological mapping framework we used to build the topological graphs.

Sum-Product Networks
One of the primary limitations of probabilistic graphical models is the complexity of their partition function, often requiring complex approximate inference in the presence of non-convex likelihood functions. In contrast, SPNs represent joint or conditional distributions with partition functions that are guaranteed to be tractable and involve a polynomial number of sum and product operations, permitting exact inference.
As shown in Fig. 1, an SPN is a directed acyclic graph composed of weighted sum and product operations. The sums can be seen as mixture models over subsets of variables, with weights representing mixture priors. Products can be viewed as combinations of features. SPNs can be defined for both continuous and discrete variables, with evidence for categorical variables often specified in terms of binary indicators.
Formally, following (Poon and Domingos 2011), we can define an SPN as follows: Definition 1 An SPN over variables X 1 , . . . , X V is a rooted directed acyclic graph whose leaves are the indicators (X 1 1 , . . . , X I 1 ), . . . , (X 1 V , . . . , X I V ) and whose internal nodes are sums and products. Each edge (i, j) emanating from a sum node i has a non-negative weight w ij . The value of a product node is the product of the values of its children. The value of a sum node is j∈Ch(i) w ij v j , where Ch(i) are the children of i and v j is the value of node j. The value of an SPN S[X 1 , . . . , X V ] is the value of its root.
Not all architectures consisting of sums and products result in a valid probability distribution. While a less constraining condition on validity has been derived in (Poon and Domingos 2011), a simpler condition, which does not limit the power of the model is to guarantee that the SPN is complete and decomposable (Peharz et al. 2017).
Definition 2 A sum-product network is complete iff all children of the same sum node have the same scope.
Definition 3 A sum-product network is decomposable iff no variable appears in more than one child of a product node.
The scope of a node is defined as the set of variables that have indicators among the descendants of the node.
A valid SPN will compute unnormalized probability of evidence expressed in terms of indicators. However, the weights of each sum can be normalized, in which case the value of the SPN S[X 1 1 , . . . , X I V ] is equal to the normalized probability P (X 1 , . . . , X V ) of the distribution modeled by the network.
Partial or missing evidence can be expressed by setting the appropriate indicators to 1. Inference is then accomplished by an upwards pass which calculates the probability of the evidence and a downwards pass which obtains gradients for calculating marginals or MPE state of the missing evidence. The latter can be obtained by replacing sum operations with weighted max operations (Peharz et al. 2017).
Parameters of an SPN can be learned generatively (Poon and Domingos 2011) or discriminatively (Gens and Domingos 2012) using Expectation Maximization (EM) or gradient descent. Additionally, several algorithms were proposed for simultaneous learning of network parameters and structure (Hsu, Kalra, and Poupart 2017;Gens and Domingos 2013). In this work, we use a simple structure learning technique  to build template SPNs representing each sub-graph. We begin by initializing the SPN with dense structure by recursively generating nodes based on multiple random decompositions of the set of variables into multiple subsets until each subset is a singleton. The resulting structure consists of products combining the subsets in each decomposition and sums mixing different decompositions at each level. Then, we employ hard EM to learn the model parameters, which was shown to work well for generative learning  and overcomes the diminishing gradient problem. After parameter learning, the generated structure can be pruned by removing edges associated with weights close to zero.

Topological Graphs
GraphSPNs are applicable to arbitrary graphs. However, in this work, we apply them specifically to topological graphs built by a mobile robot exploring a large-scale environment (Pronobis, Riccio, and Rao 2017). The primarily purpose of our topological graph is to support the behavior of the robot. As a result, nodes in the graph represent places the robot can visit and the edges represent navigability. The graph nodes are associated with latent variables representing semantics and the edges can be seen as spatial relations forming a global semantic map. Local evidence about the semantics of a place might be available and we assume that such evidence is inherently uncertain and noisy. Additional nodes in the graph are created to represent exploration frontiers, possible places the robot has not yet visited, but can navigate to. We call such nodes placeholders, and assume that the robot has not yet obtained any evidence about their semantics.
The topological graph is assembled incrementally based on dynamically expanding 2D occupancy map. The 2D map is built from laser range data captured by the robot using a grid mapping approach based on Rao-Blackwellized particle filters (Grisetti, Stachniss, and Burgard 2007). Placeholders are added at neighboring, reachable, but unexplored locations and connected to existing places. Then, once the robot performs an exploration action, a placeholder is converted into a place and local evidence captured by the robot about the semantic place category is anchored to the graph node.
We formulate the problem of finding placeholder locations as sampling from a distribution that models location relevance and suitability. Specifically, the distribution is specified as: represents the existence of a new place at location i in the occupancy map, G is the occupancy grid, and E is the set of locations of all existing places. The potential φ S ensures that placeholders are located in areas that are safe and preferred for navigation (are within safe distance from obstacles, with the preference towards centrally located places). The potential φ N models the neighborhood of a place, and guarantees sufficient coverage of the environment by promoting positions at a certain distance d n from existing places. Final location of a new placeholder is chosen through MPE inference in P (E|G, E). An edge is then created to represent navigability. It connects the placeholder to an existing place in the graph based on A* search directly over the potential φ S . An example of such semantictopological map is shown in Fig. 4.

GraphSPNs
GraphSPNs learn a template model over arbitrary graphstructured data, with local evidence X i and latent variables Y i = {Y i1 , · · · , Y iM } for each graph node or edge i, with dependencies between the latent variables expressed in terms of the graph structure. Then, an instance Graph-SPN distribution P (X 1 , Y 1 , . . . , X N , Y N ) is assembled for a specific graph to perform inference.
To this end, we define a set S of template sub-graphs, and associate each template sub-graph S ∈ S with a separate template SPN modeling the distribution over variables X i and Y i corresponding to the nodes and edges of the template sub-graph. The structure and parameters of each template SPN can be learned directly from data obtained by decomposing training graphs into sub-graphs corresponding to S.
Given a set of trained template SPNs, and a specific graph to be modeled, an instance GraphSPN is assembled as illustrated in Fig. 2. First, the graph is decomposed multiple times, each time differently, into sub-graphs using subgraph templates S in descending order of the template size (i.e. more complex templates have priority). The sub-graphs should not overlap in each decomposition and the corresponding template SPNs should cover all variables X i , Y i in the model. This condition guarantees completeness and decomposability resulting in a valid instance GraphSPN. For each decomposition and each sub-graph, we instantiate the corresponding template SPN resulting in multiple SPNs sharing weights and structure. The instantiations for a single graph decomposition are combined with a product node and the product nodes for all decompositions become children of a root sum node realizing the complete mixture model. Figure 2: An instance GraphSPN modeling a simple 5-node graph (red) with variables X i and Y i associated with graph nodes. Solid lines illustrate one decomposition of the graph based on two template sub-graphs and SPNs (green and blue), while dashed lines illustrate another decomposition.
In order to incorporate the latent variables Y ij , we include an intermediate layer of product nodes into the template SPNs. As shown in Fig. 2, each such product node combines an arbitrary distribution D k ij (X i ) with an indicator λ Yij =c k j for a specific value c k j of Y ij . The template SPN built on top of the product nodes can be learned from data and the distributions D k ij (X i ) can be arbitrary, potentially also realized with an SPN with data-driven structure.
In our experiments, we assumed only one latent variable (semantic place category) Y i per graph node i, with V al(Y i ) = {c 1 , . . . , c L }, and we defined D k i (X i ) for a single hypothetical binary observation x i , which we assumed to be observed: Such simplification allows us to thoroughly evaluate Graph-SPNs for the problem of learning topological semantic maps by directly simulating hypothetical evidence about the semantic category of varying uncertainty and under various noise conditions. Furthermore, it allows us to compare GraphSPNs with Markov Random Fields using the same α k i as the value of local potentials, i.e. φ i (Y i = c k ) = α k i . The proposed approach naturally extends to the case where a more complex distribution is used to model semantic place categories based on robot observations, such as the SPNbased approach presented in . Note, that we still learn the structure of the template SPNs built on top of distributions D k i (X i ).

Experimental Procedure Dataset
The semantic maps with topological relations used in our experiments were obtained by deploying the topological mapping on sequences of laser range data and odometry captured by a mobile robot exploring multiple large-scale environments (Pronobis and Caputo 2009). The dataset contains 99 sequences, and as a result 99 topological graphs, captured on 11 floors of 3 buildings in different cities. We identified 10 semantic place classes that are common for all buildings (e.g. a corridor, a doorway, a 1-person office, see Fig. 4 for a complete list) and annotated each topological graph node with its groundtruth class. Our goal in this work is to evaluate the ability of Graph-SPNs to disambiguate semantic place classes despite noisy and uncertain local evidence by exploiting spatial relations captured in a noisy topological graph. Probabilistic place classification algorithms, such as the SPN-based approach in  associate decisions based on local observations with probability estimates. However, the certainty of a decision can be low or the decision can be incorrect. In order to measure how sensitive the evaluated approaches are to such noise, we simulate local evidence attached to topological graph nodes by adding increasing noise to groundtruth information.
To this end, for each node i in each topological graph, we generated a local evidence distribution with values P (Y i = c k , X i = x i ) = α k i . For each graph, we first randomly selected 20% of all nodes for which the most likely local result should be incorrect. For those, we selected a random incorrect class to be associated with the highest probability value. Then, we randomized the value D incorrect , which is a difference between the highest probability and the probability of the groundtruth class, from a uniform distribution in a range depending on the noise level. For the remaining 80% of nodes, we ensured that the groundtruth class is associated with the highest probability. However, we simulated uncertainty by randomizing the value D groundtruth , which is a difference between the probability of the groundtruth class and the second highest probability. With these constraints, we used random values for the remaining likelihoods and made sure that each distribution is normalized. Intuitively, lower D groundtruth indicates higher uncertainty and higher D incorrect indicates stronger noise. The statistics of the values of D incorrect and D groundtruth for the final evidence at different noise levels are shown in Tab. 1.

Learning GraphSPNs
We learned GraphSPNs from a simple set of sub-graph templates shown in Fig. 3, matching from 1 to 5 nodes and simple edge configurations. We assumed that each node is associated with a single latent variable Y i representing the semantic place class. For each template with at least 2 nodes, we learned a template SPN of specific structure and parameters from sub-graph examples in the training set in a supervised fashion (Y i was set to the groundtruth). Each training graph was partitioned in 10 different ways to obtain sub-graphs. For a single-node template we simply assumed a uniform SPN. During testing, we built the instance GraphSPN based on 5 different graph decompositions (see Fig. 4). In our experiments, we always learned GraphSPN from all graphs from two buildings in the dataset and tested on graphs with different evidence noise levels from the remaining building. GraphSPNs are implemented using Lib-SPN (Pronobis, Ranganath, and Rao 2017).

Constructing Markov Random Fields
We compared GraphSPNs to a traditional approach based on MRFs structured according to the represented graph. The MRF was constructed from two types of potentials: potential φ i (Y i = c k ) = α k i used to provide local evidence, and potentials modeling latent variable dependencies. For the latter, we tried two models: using pairwise potentials for each pair of variables associated with connected nodes or defined over three variables for three connected nodes in any configuration. In each case, the potentials were obtained by generating co-occurrence statistics of variable values in the training graphs used for learning GraphSPNs. Inference in the MRF was performed using Loopy BP implemented in the libDAI library (Mooij 2010).

Experimental Results
We performed several experiments comparing the learned GraphSPN model to the MRF models with pairwise potentials (marked as MRF-2) and three-variable potentials (marked as MRF-3). First, we tasked all models with disambiguating noisy local evidence about semantic place class for places visited by the robot. For each topological graph in the test set, we performed marginal inference 1 , based on which we selected the final classification result argmax k P (Y i = c k |X = x). The percentage of correctly classified nodes in the graph averaged over all test graphs is shown in Tab. 2 and a visualization of results for a single graph together with the decompositions used to build the instance GraphSPN is shown in Fig. 4. We evaluated all assignments of the three buildings into training and test sets as well as different noise levels listed in Tab. 1. Each test set consisted of 5 topological maps, each with 3 random sets of noisy local evidence resulting in 15 different test graphs. Since the local evidence for 20% of nodes in each graph indicates an incorrect class as the most likely one, only accuracy greater than 80% demonstrates that the model was able to recover from the noise using learned spatial relations. Lower accuracy suggests that the incorrect evidence was too strong or that the correct evidence was too uncertain to influence the semantic class of a place.
Analyzing the performance reported in Tab. 2, we see that pairwise MRF performs well when there is little noise in the evidence, however it deteriorates quickly with increasing noise levels. At the same time, GraphSPN generating robust results (accuracy greater than 80%) despite substantial noise and uncertainty. With substantial noise, approximate Loopy BP inference for MRF converges to a solution consisting primarily of the dominant class (the corridor). At the same time, we see that using higher-order potentials with MRF actually hurts performance.
In the second experiment, we tasked the models with inferring marginal distributions over the semantic classes of places not yet visited by the robot (placeholders) for which local evidence is unavailable. We used the same setup as in the previous experiment, with local evidence for explored places including noise and uncertainty. Examples of such marginal distributions are shown in Fig. 5, while the classification accuracy when considering the most likely class is reported in Tab. 3. Again, GraphSPN significantly outperformed the MRF for this inference task. If we analyze the marginal distributions for three representative placeholders shown in Fig. 5, we see that GraphSPN is confident about the correct class for the two placeholders for which nearby nodes provide correct (albeit uncertain) evidence. For the placeholder connected to nodes for which evidence indicates an incorrect class the marginal distribution is almost uniform. This indicates the ability of GraphSPN to generate useful confidence signal in presence of noisy inputs.
Our final experiment was designed to evaluate the quality of likelihood produced by GraphSPN instantiated over a graph of arbitrary size. In this experiment, we tasked the models with detecting novel spatial environment configurations that are incongruent with the evidence available during training. To this end, we produced novel graphs by swapping groundtruth semantic class labels between corridors   Table 3: Accuracy of semantic class inference for placeholders without local evidence, for all models and test buildings, and at two representative noise levels.
and doorways as well as cooridors and 1-person offices. We contrasted these graphs with the original groundtruth as well as swaps that are consistent with the training data (swapping 1-person offices with 2-person offices). To detect novel configurations, we measured the likelihood of the complete graph with certain label assignments (without any noise) and thresholded the obtained likelihood to decide whether a graph is likely to be generated from the distribution. Since the likelihood depends on the size of the graph, we normalized it by the number of variables in the test graph before performing thresholding. We used the models learned in the previous experiments and produced false and true positive rates for various threshold values over all trained models and test sets.
The results for such detection task are shown as ROC curves in Fig. 6. All three models perform well on this Figure 5: Visualization of results for the experiment involving placeholders without local evidence, at noise level 2. Left, bottom: the semantic map with groundtruth semantic place classes (including for placeholders). Left, top: the 20% of nodes for which the most likely evidence indicates an incorrect class (black) and placeholders with no evidence (gray). Right: inferred marginal distributions over semantic classes of placeholders (pie charts). Figure 6: ROC curve for the novelty detection task with certain semantic label assignments. task (as measured by the area under the curve). We see that GraphSPN performs comparably to the MRF models demonstrating its ability to produce useful likelihood values. At the same time, MRFs perform slightly better on this task. This is a result of lack of noise and uncertainty in the provided evidence. It shows that MRFs can capture relevant spatial relations and suggests that their performance deteriorates due to the limitations of approximate inference.

Conclusions
We presented GraphSPNs, a probabilistic deep model for graph-structured data that learns a template distribution allowing for making inferences about graphs of different global structure. While existing works applied SPNs to data organized as fixed-size grids or sequences, this paper presents a novel attempt at deploying SPNs on arbitrary graphs of varying size. Based on GraphSPNs, we proposed a method for learning topological spatial relations in semantic maps constructed by a mobile robot. Our method is robust to noise and uncertainty inherent in real-world problems where information about the environment is captured with robot sensors. Our framework is universal and compatible with any distributions defined over local evidence. However, it is particularly well suited for integration with other