{'fraction_non_alphanumeric': 0.03654819954988747, 'fraction_numerical': 0.00897880720180045, 'mean_word_length': 4.658928097638631, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 5, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Multi-Source Data Fusion Outage Location in Distribution Systems via Probabilistic Graph Models

Efficient outage location is critical to enhancing the resilience of power distribution systems. However, accurate outage location requires combining massive evidence received from diverse data sources, including smart meter (SM) last gasp signals, customer trouble calls, social media messages, weather data, vegetation information, and physical parameters of the network. This is a computationally complex task due to the high dimensionality of data in distribution grids. In this paper, we propose a multi-source data fusion approach to locate outage events in partially observable distribution systems using Bayesian networks (BNs). A novel aspect of the proposed approach is that it takes multi-source evidence and the complex structure of distribution systems into account using a probabilistic graphical method. Our method can radically reduce the computational complexity of outage location inference in high-dimensional spaces. The graphical structure of the proposed BN is established based on the network's topology and the causal relationship between random variables, such as the states of branches/customers and evidence. Utilizing this graphical model, accurate outage locations are obtained by leveraging a Gibbs sampling (GS) method, to infer the probabilities of de-energization for all branches. Compared with commonly-used exact inference methods that have exponential complexity in the size of the BN, GS quantifies the target conditional probability distributions in a timely manner. A case study of several real-world distribution systems is presented to validate the proposed method.


I. INTRODUCTION
Frequent power outages are becoming a critical issue in the U.S. In 2018, the Department of Energy estimates that outages are costing the U.S. economy $150 billion annually [1]. 1.9 million customers in Midwest were affected by 1.4 million outages between August 10 and 13, 2020 [2]. Outage detection in distribution grids is an immediate and indispensable task after service disruptions, without which utilities cannot obtain needed situational awareness for initiating repair and restoration. This suggests an urgent need of efficient approaches to shorten the time of lateral-level outage location. Traditionally, outage location inference has been done based on manual outage mapping, which in addition to the voltage and current components measured only at the substations, has mainly depended on customers' trouble calls. However, trouble calls alone are not a reliable source for outage location inference. It is estimated that only one-third of customers report the events in the first hour of outages, which might prolong the location determination process [3]. Also, customers might contact utilities due to temporary and individual problems other than system-level outage events, which can mislead the location process and result in additional truck rolls to verify power outages.
One way of avoiding these problems is to rely on advanced metering infrastructure (AMI)-based techniques, which can send outage notifications at the grid-edge by leveraging the bidirectional communication function of smart meters (SMs) [4]. Researchers have dedicated great efforts to this topic. In [5], a hierarchical generative model is proposed that employs SM error count measurements to detect anomalies. In [4], a multi-label support vector machine model is developed that utilizes the state of customers' SMs to identify states of distribution lines. In [6], a two-stage method is presented to detect non-technical losses and outage events using realtime consumption data from SMs. In [7], a hypothesis testingbased outage location method is developed that combines the power flow measurements and SM-based load forecasts of the nodes. In [8], by using data from SMs and fault indicators, a multiple-hypothesis method with an extended protection tree is presented to detect a fault and identify the activated protective devices. A main challenge is that most AMI-based methods require full observability for distribution grids, i.e., SM installation for all customers. This assumption is not necessarily applicable to practical distribution systems, mostly due to utilities' budgetary limitations. To perform outage detection in partial observable systems, we have proposed a generative adversarial network (GAN)-based method to efficiently identify outage region [9]. Although this method is guaranteed to capture the maximum amount of information on outage location, it does not provide granular outage location estimation at the branch level due to the limitations of the single data source. This issue is further exacerbated considering that SM signal communication to the utilities' data centers can fail due to hardware/software malfunction and tampering [5].
Rather than using SM data, an alternative solution is to utilize other grid-independent data sources to identify outage events in real-time [10]. In [11], weather information data is used to detect outages in overhead distribution systems employing an ensemble learning approach. In [12], a data-driven outage identification approach is proposed that extracts textural and spatial information from social media. In [13], a mixedinteger linear program (MILP) is formulated to identify the arXiv:2012.02877v2 [eess.SP] 8 May 2021 topology under both outage and normal operating conditions using line flow measurements, forecasted load data, and ping measurements from a limited set of SMs. Nonetheless, the considerable uncertainty of these data sources can lead to erroneous outage location and additional costs for utilities. For example, only a part of SM last gasp signals can be delivered to the utility's data center due to hardware and software issues. Thus, to handle the limitations and uncertainties of individual data sources, this paper proposes a multi-source data fusion strategy to combine outage-related information from diverse sources for accurate outage location.
One fundamental challenge in multi-source outage location is the computational complexity of the problem: first, note that outage location inference is the process of computing the probabilities of topology candidates after disrupting events by leveraging available information received by utilities. Estimating these probability values requires obtaining the joint probability distribution function (PDF) of the unknown state variables and the evidence, which is a high-dimensional mathematical object. Directly quantifying this joint distribution is computationally infeasible for actual distribution systems, which requires enumerating probabilities of all possible combinations of variables. In addition, outage data sources have heterogeneous characteristics such as accuracy levels and reporting rates. Further, they may provide inconsistent and contrary information. How to integrate these data sources is a challenge.
To address these challenges and the shortcomings of the previous works in the literature, a multi-source data fusion method is presented to identify and locate the lateral-level outage events in partially observable distribution systems. To achieve this, we have adopted a probabilistic graphical modeling approach towards data fusion to reduce the computational complexity of representing high-dimensional joint PDF of the system. The basic idea of this methodology is to first build a Bayesian network (BN) for each distribution feeder and then use this probabilistic model to infer the conditional PDF of the state of network branches and the connectivity of customer switches given the observed evidence. BNs utilize a graph-based representation as the basis for analyzing statistical relationships among random variables [14]. System topology from one-line diagram and context data, such as weather data and vegetation information, from geographic information system are used to design the architecture of the BN, as shown in Fig. 1. The graph parameters are learned empirically from historical outage data. By utilizing the proposed BN-based technique, the high-dimensional joint PDF of the system is decomposed into a set of more manageable probabilistic factors obtained from conditional independencies. Consequently, the data fusion-based outage location process is efficiently transformed into online inference over the BN. This inference task is solved by leveraging a Gibbs sampling (GS) algorithm. As a Markov chain Monte Carlo (MCMC)-based algorithm, GS can provide a full characterization of the distribution of unknown variables by generating a sequence of samples. We have used several real-world feeder models from our utility partners to validate the performance of the proposed method. The main contributions of this paper can be summarized as follows: • The proposed method can seamlessly integrate heterogeneous data sources. Different data sources can complement each other to increase the amount of outage information, thus addressing low SM coverage or customer report rates in actual grids. • By exploiting the inherent conditional independencies between the evidence and state variables in distribution systems, the exponential computational complexity of the outage location task is reduced to linear complexity in the number of variables. • The proposed method is robust with respect to misreports and inconsistencies in outage evidence, as the uncertainty of each data source is explicitly modeled using graph parameters. • The proposed method introduces an intuitive way to incorporate the prior expert knowledge from utilities into the outage inference model. • The proposed method provides a principled framework to reduce overfitting under outage data scarcity by leveraging the Bayesian hypotheses.
The rest of this paper is constructed as follows: In section II, the statement of the outage location problem is described. Section III presents the proposed BN-based data fusion model, along with structure selection and parameter learning schemes. An MCMC approximate inference algorithm is given in Section IV. The numerical results are analyzed in Section V. Section VI concludes the paper with major findings.

II. OUTAGE LOCATION PROBLEM STATEMENT
Considering that outage events cause topological changes in the grid, outage location is the process of inferring the probabilities of post-event operational topology candidates. The accuracy of outage location depends on the completeness of outage information. In this word, outage-related information from diverse sources are utilized for outage location, which include SM last gasp signals, customer trouble calls, social media messages, wind speed, vegetation information, and physical parameters of the grid. The rationale behind the use of wind speed and vegetation information is that 87% of major power outages are caused by extreme weather events [4]. These data sources are commonly available to utilities. Specifically, SM will automatically send the last gasp signal to the head end system of the AMI after the power fails. Trouble calls and social media messages are reported by customer's phones and Twitter. Wind speed, vegetation information, and the physical parameters of the grid can be found from the neighboring land-based station, the geographic information system, and the grid model, respectively. Thus, our method does not require to install extra metering devices, such as SMs, line flow measurements, and micro-phasor measurement units (µPMUs), for data collection. This ensures the practicability of the proposed method for actual distribution systems. Further, since the proposed method is scalable and adaptive, utilities are able to add or remove data sources according to their situations. For example, in the absence of extreme weather events, utilities can remove vegetation information and wind speed to reduce the complexity of the model, as these two data sources may not have a significant impact on outage detection and location during normal weather. After data collection, last gasp signals, customer trouble calls, wind speed, vegetation information, and physical parameters can be directly transformed into outage evidence as input to the proposed model. For social media messages, a natural language processing tool is required to extract outage-related words, as proposed in our previous work [12]. Then, social media messages are converted into binary outage evidence, similar to customer trouble calls and last gasp signals.
Given the multi-source evidence, E E E, this inference process is mathematically formulated using the Bayes estimator [15], where the conditional PDF of network topology, Y , given the set of evidence is represented as P (Y |E E E) and calculated in terms of the joint distribution of Y and E E E, denoted by P (Y, E E E). Note that in this paper, random vectors/matrices are represented with bold letters. The most probable candidate topology, which also determines the location of the outage event, is obtained by maximizing this conditional PDF, as: where, Y * is the most likely network topology after the outage. Y is a multinomial variable which is represented in terms of the states of primary network branches (D D D) and the connection of customer switches (C C C), as Y = {D D D, C C C}.
where k is the number of branches in the feeder and D i is a binary variable representing the connectivity state for the i'th branch in the feeder: D i = 0 means that the branch is energized. In other words, there is an uninterrupted path between the branch and the substation. D i = 1 indicates that the branch is de-energized. Similarly, representing the set of connection states for all the customers that are supplied by the i'th branch.
where z i is the total number of customers that are connected to the i'th branch, and C j i is the state of the j'th customer: C j i = 0 means that the customer is energized, and C j i = 1 implies that the customer is deenergized. Note that the pre-outage topology is determined by assigning 0 to all the state variables (i.e., all branches are energized and customers are energized). Thus, P (Y |E E E) in (1) can be rewritten in terms of the joint PDF of the newly-defined variables, P (D D D, C C C, E E E), as follows [16]: . (2) Using (2), the maximization over topology candidates can be conveniently transformed into finding the best values for the individual branch/customer states belonging to {D D D, C C C} using their conditional PDFs, P (D i |e e e) and P (C j i |e e e). These conditional PDFs are obtained ∀i, j using a marginalization process over the joint PDF, as follows [17]: .
(4) Solving (3)-(4) requires quantifying the joint PDF, P (D D D, C C C, E E E). Considering the complexity of distribution grids, obtaining the explicit representation of this joint PDF is unmanageable for two reasons: (I) a complete description of P (D D D, C C C, E E E) induces an exponential complexity in the order of 2 r − 1, where r is the total cardinality of all the unknown variables, r = |D D D| + |C C C|. Hence, modeling this joint PDF using brute-force search over all possible combinations of branch/customer states is computationally infeasible for large-scale distribution systems. (II) Due to the outage data scarcity in distribution grids, it is impossible to acquire enough historical data to robustly estimate the massive number of parameters of this joint distribution. One solution is to use naive classification by assuming full independence among all evidence and unknown state variables [17]. However, this assumption is not applicable to practical distribution systems and may lead to severe misclassification due to overfitting.

III. BN-BASED DATA FUSION MODEL
To counter computational complexity and overfitting in the outage location inference, we propose a BN-based method to provide a compact representation of the high-dimensional joint PDF, P (D D D, C C C, E E E). To do this, our method exploits conditional independencies between random variables, {D D D, C C C, E E E}, to decompose the joint PDF into a set of factors with significantly smaller size. Using this computationally efficient BN-based approach, we can infer the conditional PDF of the state of each primary branch and the customer switch given outagerelated evidence from various data sources, shown in (3)-(4), to rapidly identify the location of lateral-level outage events.

A. Factorization of the Joint PDF and BN Representation
The main idea of a BN-based representation is to use conditional independencies, encoded in a graph structure, to compactly break down high-dimensional joint PDFs with a set of factors. Here, a factor refers to a low-dimensional and more manageable conditional PDF that is determined by two components: a child variable, such as D i and a number of parent variables denoted by P a(·), such as P a(D i ). Parents represent the direct causal sources of influence for a child variable. In other words, each child is a stochastic function of its parents [17]. Thus, if the values of the parents are known, then the child variable becomes conditionally independent of random variables that do not directly influence it in a causal manner. It can be shown that by using chain rule over these conditional independencies, defined by parent-child relationships, the joint PDF of a set of random variables can be simplified as the multiplication of the identified factors [17]. In the outage location problem, this factorization leads to the following data fusion representation for the joint PDF: where, u = |E E E|, and the factors are: i,j denotes the human-based evidence from the customer-side, including trouble calls and social media messages; E m i,j represents meter-based evidence from customerside, such as smart meter last gasp signals. Compared with the original model in (2) that requires 2 r − 1 independent parameters, the new representation in (5) only needs parameters. It can be observed that the number of parameters in the new representation is a function of size of parents for each variable. Considering that the number of variables' parents is typically small, the new representation achieves a radical complexity reduction in outage location inference. As a directed acyclic graph, BN offers a convenient way of representing the factorization (5). Accordingly, the random variables, {D D D, C C C, E E E}, are represented as the vertices of the BN. Using the identified factors in (5), the vertices of the BN are connected by drawing directed edges that start from parent vertices and end in child vertices. BN provides a graphical way for encoding conditional independencies defined by the factors as follows: any vertex, X, is conditionally independent of its non-descendant vertices in the graph, N d(X), if the values of its parents are known. This is symbolically denoted by (X ⊥ N d(X)|P a(X)) [18]. N d(X) is the set of the vertices of the BN, excluding parents of X, to which no directed path exists originating from X. A ⊥ B means that A and B are marginally independent.

B. BN Structure Development and Parameterization
Developing a BN requires discovering the structure of the graph and the parameters of the conditional PDFs. To do this, a knowledge discovery-based method is utilized in this paper. The basic idea of this method is to reveal the conditional independencies between variables using the grid topology information and causal relationships. Then, the conditional PDFs (i.e., factors) are parameterized based on the available statistical outage information from the literature [10], [12], [19]. For concreteness, the parent-child variables of each factor in (5) are described as follows: ( i is vegetation information, which contains vegetation species-specific constants and tree diameter; E b i represents the i'th branch's physical parameters including the length of conductors and the number of poles. Based on this parent selection scheme for branch state variables, N d(D i ) includes all the variables that are not downstream of the i'th branch in the feeder (see Fig. 2). To show the direct causal influences of these four variables on D i , two cases are described: In the first case, when the parent branch is de-energized, then D i = 1 with probability 1. Consequently, all the variables in the path from the substation to D i−1 , represented with {D 1 , ..., D i−2 }, are conditionally independent from {D i } given D i−1 = 1. The intuition behind this is that in radial networks there can only be a unique path between the substation and each branch; if this path is interrupted at any arbitrary point in {D 1 , ..., D i−2 }, we can automatically conclude D i−1 = 1 regardless of the location of outage in the path. Hence, considering the binary nature of variable D i , the conditional PDF, can be formulated as: In the second case, if the neighboring upper-stream branch is energized, then all upstream branches of the i'th branch are also energized with probability 1, and have not been impacted by outage, {D 1 = 0, ..., D i−2 = 0}. In this case, D i = 1 will only occur when this branch is damaged. As demonstrated concretely in [20], the majority of branch damage is caused by tree contacts to power lines and broken poles due to high wind speed. Thus, three context variables E w i , E v i and E b i are serve as causal evidence for the i'th branch state to estimate the probability of outage at the i'th branch. The conditional PDF, can be formulated as a Bernoulli distribution as follows: where, the probability of failure for branch i, denoted as P i l , is a function of E w i , E v i , and E b i . To formulate this function, a fragility model is leveraged [20]. Basically, the fragility model is a series model with the fragility analysis of each pole and conductor within the branch: where, L is the number of distribution poles used for supporting branch i, K is the number of conductor wires between two neighboring poles at the i'th branch, Φ is the standard normal probability integral, χ is the median of the fragility function, ξ is the logarithmic standard deviation of intensity measure, and P f (E w i , E v i ) represents the failure probability for conductor f of branch i which is modeled as follows: where, p u is the probability of conductor f being underground, F wind,f (E w i ) represents the wind force loading on the conductor and F no,f (E w i ) demonstrates the maximum perpendicular force of the conductor wire determined as shown in [21]. α describes the average tree-induced damage probability of overhead conductor, and P f (E v i ) is the fallen tree-induced failure probability of conductor f computed as in [20]. Hence, for the case D i−1 = 0, equations (8) and (9) are utilized to estimate the probability of outage for branch i given the context variables E w i , E v i , and E b i . To summarize, the conditional PDFs given in equations (6) and (7) fully determine the factors of the form P (D i |P a(D i )).
(2) Factor P (C j i |P a(C j i )) represents the conditional PDF of the status of customer j given parent variables. The parent of customer state variable is selected as P a(C j i ) = {D i } (see Fig. 2). Here, D i is the state of the immediate upper-stream branch that supplies the j'th customer. To show the casual relationship between C j i and D i , two cases are considered: D i = 1 and D i = 0.
In the first case, if the primary branch is de-energized, the probability of C j i = 1 is 1 due to the radial structure of the feeder. Utilizing this deterministic relationship, P (C j i |D i = 1) can be written as follows: In the second case, if the primary branch is energized, then the path between the substation and the i'th branch is active. In this case, customer outage, C j i = 1, can only be caused by overloading/faults at the customer-side occurring with probability π 2 . This case is represented using a Bernoulli distribution adopted from statistical outage information [22]: To account for the uncertainty of parameter π 2 , a beta distribution is defined with user-defined hyper-parameters α 2 and β 2 : π 2 ∼ Beta(α 2 , β 2 ) = γ 2 π α2−1 2 (1 − π 2 ) β2−1 (12) where, γ 2 is a normalizing constant and defined as γ 2 = Γ(α 2 + β 2 ) with Γ = ∞ 0 t x−1 e −t dt [17].
. The parents of human-based evidence, E h i,j , are selected as P a(E h i,j ) = {C j i , ∆T }, as shown in Fig. 2. ∆T refers to the time elapsed after the outage occurrence. More precisely, ∆T embodies the time period that utilities need to wait before outage reports are issued [23]. This is to avoid false alarms from customers due to temporary events. It is clear that there is a trade-off between the amount of human-based evidence and waiting time of outage location inference. For example, when the observability of the actual feeder is extremely low, the utility may increase ∆T to receive more human-based evidence for outage location inference. Within the ∆T period, the time at which the humanbased evidence is received, t, after outage occurrence at time, T 0 , is distributed according to an exponential distribution as shown in [10]: Thus, given ∆T , the probability of P (E h i,j |C j i = 1, (t−T 0 ) ≤ ∆T ) can be calculated as: Hence, the factor P (e h i,j |C j i , t ≤ ∆T ) is obtained as follows: where, π 3 denotes a small user-defined value to take into account the possibility of false positives, such as illegitimate trouble call and social media data processing errors.
. Compared to the humanbased signals E h i,j , AMI-based notification mechanism will be delivered almost instantaneously to the utilities. Thus, the parent of meter-based evidence is selected as P a(E m i,j ) = {C j i } (see Fig. 2). When the state of customer switch is known, E m i,j becomes conditionally independent of the remaining variables, as encoded by the factor: where, π 4 and π 5 represent the AMI communication reliability and the SM malfunction probability values, respectively. For concreteness, π 4 is the probability that the last gasp can be delivered to the utilities correctly for outage notification. π 5 is the probability that the SM loses power due to its own failure and sends a last gasp signal. In this work, the values of these two parameters are determined based on the historical outage report. Considering the size of the historical data is limited, beta distributions are used to model the uncertainty of these two parameters as follows:

IV. BN-BASED OUTAGE LOCATION INFERENCE USING GS
After construction and parameterization of the BN, the data fusion outage location process is efficiently transformed into a probabilistic inference over the graphical model. However, even if P (D D D, C C C, E E E) is simplified, solving (3)-(4) still requires calculating computationally expensive summation operations d d d c c c P (d d d, c c c, E E E) over all nodes of the graph simultaneously, which is not scalable for large-scale distribution grids [17]. To address this, a GS algorithm is used to perform the inference task over the BN [24].

A. GS Algorithm
GS is an MCMC-based approximate inference method, which allows one to provide a good representation of a PDF by leveraging random variable instantiations, without knowing all the distribution's mathematical properties [24].
The key advantage of this method is that it employs univariate conditional distributions for sampling, which eliminates Select one random variable X i ∈ {D D D, C C C} 5: Obtain P a(X i ) and Ch(X i ) from G 7: end for 11: end for 12: Return sample vectors: , ∀i, j 16: Select the nearest de-energized branch as the outage location the dependency on the dimension of the random variable space. Thus, compared to the commonly-used exact inference methods, such as variable elimination and clique trees, GS is insensitive to the size of BN [16]. This indicates that the GS method is especially beneficial for complex real-world applications.
When an outage occurs, the de-energization probabilities of branches/customers are inferred using the GS algorithm and the BN structure. To do this, first, all the outage evidence from the customer-side, {E h 1,1 , ..., E h z k ,k , E m 1,1 , ..., E m z k ,k }, is collected after ∆T has elapsed: if utilities receive trouble call/tweet or last gasp signal from the j'th customer at branch i, the corresponding evidence E h i,j or E m i,j is set to 1. The branch-level evidence, is obtained from the utilities' data centers and weather information systems. After collecting all evidence, arbitrary initial samples are randomly assigned to all the unknown state variables {D D D, C C C}: Then, an arbitrary state variable is selected as the sampling starting point, e.g., D i . At iteration τ + 1 of GS, following the structure of the BN, the assigned samples to the parents and children of D i are inserted into a local Bayesian estimator [16], as shown in (20), to approximate the conditional PDF of D i given the latest samples: is all the latest samples except for d i , including values of evidence variables, and: ) can be directly calculated using the determined factors, (6)- (17), in Section III-B. Note that because is a PDF over a single random variable given the samples assigned to all the others, this computation can be performed efficiently.
is drawn using the inverse transform method [17] to replace d (τ ) i . Then, the algorithm moves to a next non-evidence variable of BN to perform the local sampling process (see (20)). When all the unknown variables of the BN have been sampled once, one iteration of GS is complete. This process is able to propagate the information across the BN and combine the data from diverse sources to infer the location of outage efficiently. The sampling process is repeatedly applied until a sufficient number of random samples are generated for the unknown variables, {D D D, C C C}. It has been theoretically proved that the approximate PDFs, P Φ (·), are guaranteed to approach the target conditional PDFs, P (D i |E E E) and P (C j i |E E E), defined in (3)-(4) [17]. Thus, P (D i |E E E) and P (C j i |E E E) can be estimated by counting the samples generated by the GS algorithm. As an example, P (D i = 1|E E E) is estimated as follows: where, M is the number of iterations. After the GS process, the most likely value of each branch/customer state is determined based on the obtained approximated conditional PDFs to solve (1). To achieve this, due to the binary nature of the state variables, a 0.5 threshold is used, e.g. P (D i = 1|E E E) ≤ 0.5 indicates branch i is energized. After the connectivity states of all the branches/customers are inferred, the location of outage events are obtained by selecting the nearest de-energized branch to the substation. See Algorithm 1 for details.

B. GS Calibration Process
One challenge in GS is how to determine the number of iterations, M . In general, if the iterations have not proceeded long enough, the sampling may grossly misrepresent the target distributions, thus decreasing the inference accuracy. In contrast, if the value of M is large enough, the theory of MCMC guarantees that the stationary distribution of the samples generated using the GS algorithm [16]. However, such a strategy leads to high computational time, which increases outage duration and cost. Hence, by using GS, a trade-off exists between the accuracy and computational time of outage location. To find a reasonable maximum iteration number for a specific BN, a potential scale reduction factor, R, is utilized to diagnose the convergence of the GS at different numbers of iterations [25]. The basic idea is to measure betweenand within-sequence variances of generated sample sequences. Specifically, for each M, we start with n sample sequences produced by the GS for each unknown variable in the BN. After discarding the samples generated in the warm-up period, each sequence is divided into two halves of the same size, m, and used to complement the original sequences. All the sample sequences are concatenated into a matrix of size 2n × m, denoted as θ θ θ. Utilizing this matrix, the between-sequence and within-sequence variances are calculated as follows: where, B i is the between-sequence variance of variable i, V i is the within-sequence variance of variable i,θ θ θ .j is the within-sequence means that can be calculated usingθ θ θ .j = 1 m m i=1 θ θ θ ij .θ θ θ .. is the overall mean that can be computed usingθ θ θ .. = 1 2n 2n j=1θ θ θ .j . s 2 j denotes the j'th sample sequence variance obtained as s 2 Utilizing V i and B i , R i is defined and computed as [16]: In theory, the value of R i equals 1 as 2m → ∞. R i 1 indicates that either estimate of the variance can be further decrease by more iterations. In other words, the generated sequences have not yet made a full tour of the target PDF. Alternatively, if R i ≈ 1, the sequences are close to the target PDF. Here, following the previous work [16], a threshold R ψ = 1.1 is adopted to select the value of M . Thus, M ← 2m is set as the number of iterations that satisfy R i ≤ R ψ , ∀i for the BN.

V. NUMERICAL RESULTS
This section explores the practical effectiveness of the proposed data fusion outage location method. Three real-world distribution feeders are utilized in this case study, which are publicly available online [26]. The topological information is shown in Fig. 3. All aforementioned evidence, including trouble calls, social media messages, last gasp signal, vegetation information, and wind speed, is utilized in the simulations. Specifically, the human-based evidence is generated using an exponential PDF given ∆T . Note that the parameter of this PDF is considerably different from that of (14) to simulate the uncertainty of the BN parameterization in real-world applications. Consequently, in the outage inference task, we do not know the PDF used to generate evidence and the conditional PDF of the outage location. In this work, the value of ∆T is assigned as 10 minutes, which indicates that only a fraction of customers are active to give trouble calls or social media messages in our simulations. The number of last gasp signals is determined based on the system observability that is calculated as the ratio of customers with SMs to those without SMs. When a customer is assumed to have the SM, this indicates that he is likely to send a last gasp signal when an outage occurs. The vegetation information and the branch's physical parameters are provided by our utility partners. For some unknown parameters, such as tree diameter, we refer to previous work [27]. Further, depending on the geographical locations of the available systems, the wind speed data is obtained from national oceanic and atmospheric administration (NOVAA) [28]. For each feeder, we have tested the proposed method under three different observability levels, 25%, 50%, 75%. To validate the average performance of the proposed method, a Monte Carlo approach has been utilized to generate 1500 outage scenarios for each case. In each scenario, a portion of customers are randomly selected to install SMs. Then, the outage location is also randomly chosen. Based on the outage location and the distribution of the SMs, the human-and meter-based evidence is generated. To simulate real-world power outages, 10%, 15%, and 3% of total evidence is assumed to be wrong to simulate the illegitimate calls, natural language processing errors, and AMI communication failure.    stream of the outage location, these probabilities converge to significantly higher values compared to the branches that are not impacted by the outage event. By using the threshold, the energized branches and the de-energized branches can be easily distinguished to locate the outage. This demonstrates that the BN-based outage location inference method is able to correctly determine the state of the system.

A. Performance of the Proposed Data Fusion Model
To evaluate the performance of the proposed outage location method for 1500 generated outage cases in the test systems, several statistical metrics are applied among all branches, including accuracy, precision, recall, and F 1 score [29], [30]. These indexes are determined as follows: where, TP is the true positive (i.e., state of branch is inferred as de-energized while its actual state is also de-energized), TN is the true negative (i.e., state of branch is considered as an energized while its true state is also energized), FP is the false positive (i.e., state of branch is inferred as de-energized while its actual state is energized), FN is the false negative (i.e., state of branch is inferred as energized while its actual state is de-energized), P and N are the numbers of total positives and negatives, and β is the precision weight which is selected to be 1 in this paper. The average values of these indexes are presented in Table. I for the three different test feeders with various observability levels. In all cases, the lowest accuracy, precision, recall, and F 1 score are 98.7%, 83.47%, 98.88%, and 88.05%, respectively. For 50% and 75% observability cases, all branch-level indexes reach values over 0.9. Also, the system-level accuracy is calculated as the percentage of times that the states of all the branches/customers have been inferred correctly in outage scenarios. In other words, even though the outage location is inferred correctly, the systemlevel index may fail because of one misclassified branch. As shown in the table, when the observability is 25%, the system-level accuracy is about 70%. This could be due to the evidence scarcity. We have analyzed the failed scenarios. In more than 80% of these scenarios, the proposed method can infer the actual location of the outage but misjudged the status of one or two branches. For the cases that have 75% observability, the system-level accuracy is about 90%. This result is not surprising since we have assigned false positive and false negative alarms in each scenario. Such alarms reduce the completeness of outage information. By comparing the results of the three feeders, it can be concluded that the performance of the proposed outage location method improves as the observability increases, due to the high confidence levels of meter-based evidence. Also, the proposed algorithm shows almost the same level of performance over the different test feeders. This result demonstrates that the BN-based outage location method is nearly insensitive to the topology of the underlying network.
To further evaluate the performance of our method, coinciding multiple outage events are generated in three test systems. The locations of the outages are randomly selected. For concreteness, we have also calculated the accuracy under 25%, 50%, and 75% observability levels. Fig. 5 shows the performance indexes as a function of observability level and the number of outages. As can be seen, almost in all cases, higher observability improves the performance indexes regardless of the number of coinciding outage events. Also, the indexes have nearly similar values in cases with single and multiple outages. Hence, we can conclude that the method has a stable performance for multiple outages.

B. GS Calibration Results
Basically, the GS calibration is a trial and error process using a specific index, R. Hence, in each test feeder, we have generated 500 sample sequences for each unknown variable in the BN at different sampling iterations, M . Fig. 6 shows the values of R i in the 51-node test feeder. As can be seen, by increasing the number of M , the values of R i 's tend to converge to 1. By selecting M = 4000, all R i 's drop below the user-defined calibration threshold, R ψ = 1.1, which indicates that GS has reached a reasonable number of iterations in this BN. Note that GS calibration is a offline process; as a result, the high computational burden of the trial and error process does not impact the real-time performance of the proposed method.

C. Computational Complexity Analysis
The case study is conducted on a standard PC with an Intel(R) Xeon(R) CPU running at 4.10GHZ and with 64.0GB of RAM and an Nvidia Geforce GTX 1080ti 11.0GB GPU. To provide a comprehensive computational complexity analysis, the proposed method is conducted on two additional realworld distribution feeders: a 17-node and 164-node feeders. The detailed information of these feeders can be found in [9] and [31]. Fig. 7 shows the average computational time of outage inference for the test feeders. As described in the figure, by using our standard PC, the average computational time for outage location inference in five test feeders are {2.7s, 12.58s, 21.64s, 30.14s, 51.59s}, respectively. Further, due to significant complexity reduction by the proposed BN, it can be observed that the computational complexity and the size of the distribution feeders have almost a linear relationship rather than an exponential one. Also, the proposed model does not infer outage location in a system-wide fashion, but performs feeder-level location estimation. This strategy enables parallel computation of different feeders to further reduce the computational time. These salient features can facilitate the application of practical distribution systems.

VI. CONCLUSION
In this paper, we have presented a novel multi-source data fusion approach to detect and locate outages in partially observable distribution networks. The problem is cast as the process of inferring the probabilities of post-event operational topology candidates. Our method encodes the network's topology and the causal relationship between outage evidence and branch states into BNs by leveraging the conditional independence inherent in distribution grids. By constructing the BNs, the proposed method is able to infer the connectivity probability of individual primary branches with nearly linear complexity in the size of the network. Moreover, this method exploits data redundancy to reduce the impact of data uncertainty, and is suitable for arbitrary radial distribution systems. Based on simulation results on real-world networks, the proposed method can accurately detect and locate outage events within a short time.
Future study will seek to extend the proposed method in meshed grids with distributed energy resources. BNs alone cannot fully capture conditional independencies when there are multi-directional power flows. Hence, we plan to explore hybrid graphs that consist of both directed BNs and fully undirected Markov networks. Further, a joint Boltzmann distribution function will be investigated to embody graph parameters.