Thu Truong

Ta Phuong Huong

Bac

Ngan Kieu

Ha

Viet Nguyen

Hoang

Xuan Nguyen

Nguyen Tai Hoang

Kim Phuc Hung

Tran

Huong

Ta Phuong Bac

Ngan Ha

Nguyen Tai Hung
email: hung.nguyentai@hust.edu.vn

Kim Phuc Tran

Federated Learning-Based Explainable Anomaly Detection for Industrial Control Systems

Keywords: Anomaly Detection, ICS, Federated Learning, XAI, VAE

des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

SVDD

Introduction

AI and Bigdata present excellent potential in migrating the manufacturing paradigm to smart manufacturing, as it enables AI-driven IIoT systems to operate in real-time and be more precise and efficient [START_REF] Lu | Smart manufacturing process and system automation a critical review of the standards and envisioned scenarios[END_REF].

Within the context of smart manufacturing, Industrial Control Systems (ICS) are an essential component of industrial systems, and their safety and security are becoming increasingly important in the Industrial Internet of Things (IIoT) landscape. However, the exponential rise of IIoT brings not only enormous benefits but also significant obstacles in terms of developing and deploying secured ICSs [START_REF] Haddadpajouh | A survey on internet of things security: Requirements, challenges, and solutions[END_REF], [START_REF] Tuptuk | Security of smart manufacturing systems[END_REF]. In reality, a contemporary ICS is no longer a stand-alone system but rather linked to the Internet. As a result, if hackers were to acquire control of a network and steal security-critical data, or viruses and infections would infiltrate and damage a production line's operating system, the effects would be severe and costly. Industrial Control Systems based on the IIoT is currently one of the top industries attacked by various threats. As threats are becoming more complex, an anomaly detection method that can identify attacks quickly and correctly while being lightweight enough to be used in Internet of Things (IoT) devices with limited processing capacity in industrial environments is required.

From another perspective, Federated Learning (FL) -a distributed machine learning mechanism [START_REF] Abdulrahman | A survey on federated learning: The journey from centralized to distributed on-site learning and beyond[END_REF] is a promising candidate for communication costs in a distributed environment. Therefore, in this paper, we develop a FL-based anomaly detection for ICSs right at edge sites. This way helps aggregating distributed local learning models for a global model update, thereby giving each single local model knowledge of other data patterns at other edge zones. Therefore, Federated Learning can achieve efficiency similar to the centralized learning manner while distributing the detection task to different local edge zones, thereby not offloading the central cloud. While in another side, the benefit of deploying anomaly detection tasks at the edge with FL is that it definitely improves the system response time upon attack arrivals since the detection model is executed right at the edge which is near to attack/anomaly sources. The integration of the FL technique with an MLbased detection scheme helps to achieve the advantages produced by both techniques in an efficient and lightweight way.

Besides, although deploying FL enables distributed deep learning algorithms to work efficiently for anomaly detection in IIoT-based ICSs, anomaly detection techniques can only help detect abnormalities. The output of the Machine Learning-based detection model is difficult to explain or interpret, especially in ICSs where information is often abstract. Interpretability is the degree to which a human can understand the cause of a decision [START_REF] Miller | Explanation in artificial intelligence: Insights from the social sciences[END_REF]. An explanation denotes the subset of elements in a sample that has the highest impact on predicting a label output of an ML-based detection model. Note that this is a very machine-learning-centric definition. In the domain of cybersecurity analysts, a satisfying explanation would also need a description of why those attributes are critical. Because of this limitation, persuading experts to accept and use anomaly detection technologies is difficult. Such ML-based model' outputs may contain abnormal cases that the systems analyst was previously unaware of, and an explanation of why an instance is abnormal might boost the analyst's confidence in the algorithm. The higher the interpretability of an ML model is, the more easily administrators can comprehend why certain predictions have been made. Furthermore, explanations might be contradictory, which is valuable and important for explaining anomalies. To overcome this drawback, the concept of eXplainable Artificial Intelligence (XAI) came into play for ICSs. XAI has been developed to explain predictions from anomaly detection algorithms.

Motivated by these potentials, in this paper, we propose a Federated learning-based Explainable Anomaly Detection for Industrial Control Systems -called FedeX as the whole architecture to detect and analyze anomalies in ICSs and to enable detection in a distributed environment with FL.

In FedeX, we propose to:

• Distribute the anomaly detection task to the edge. Distributed edge computing approach brings a few benefits: detection response can be faster since the detection task is carried out right at the edge which is near attack sources. In case of attacks, closing an edge zone to seal attacks within does not affect the other zones' operation. Practitioners also can more easily allocate the source of attacks within each small monitored area with a limited number of connecting devices.

• Solve the detection problem by using the so-called FedVAE-SVDD model to guarantee high detection performance and real-time operation (i.e. minute-time scale). FedVAE-SVDD combines the advantage of Variational Autoencoder (VAE) [START_REF] Kingma | Auto-encoding variational bayes[END_REF] and Support Vector Data Description (SVDD) [START_REF] Tax | Support vector data description applied to machine vibration analysis[END_REF] implemented at the Edge. While most current approaches determine outlier thresholds by either using heuristic methods or normally-distributed data assumptions that are unrealistic and increase false-positive rates. Therefore, FedVAE-SVDD, by using SVDD to seek for an optimal threshold, can deliver an automatic anomaly detection solution with high performance that is proved to outperform some other learning models such as Logistic Regression (LR), Linear Discriminant Analysis (LDA), K-nearest Neighbours (KNN), Nave Bayes (NB), Support Vector Machine (SVM), and Classification and Regression Tree (CART) in the same ICS context [START_REF] Selim | Anomaly events classification and detection system in critical industrial internet of things infrastructure using machine learning algorithms[END_REF].

• Leverage FL to reduce huge work offloading in the cloud. Also, exchanging only model information from the edge devices to the cloud also solves the problem of data shortage in each device, especially with high-dimensional datasets. In addition, by using FL, each single detection model based on each local training dataset at each zone will be updated globally.

• Integrate XAI (i.e SHAP [START_REF] Lundberg | A unified approach to interpreting model predictions[END_REF]) to deal with the current "black-box" detection approaches. SHAP can thoroughly explain anomalies detected by the FedVAE-SVDD model. This type of XAI aims to provide a comprehensive explanation and remove drawbacks in understanding the output of the FedVAE-SVDD anomaly detection. This XAI method can evaluate the model's faithfulness as well as provide the ability to detect a specific component causing the system problem. SHAP is the most exact agnostic XAI method available today, as a cause-ofinterpretation tool. SHAP plots and outputs recorded values and accommodate practitioners to analyze and interpret the results, further evaluating the reliability of the proposed ML-based training model.

In summary of FedeX's benefits, FedeX is the framework that applies XAI to explain anomaly for ICSs in a liquid-storage infrastructure. FedeX is also proved to have high detection performance while being able to be deployed on top of such weak hardware of edge nodes in a distributed computing architecture. The FL-based architecture of Fedex enables faster system response capability upon attacks. FL also assures that the detection model in each single distributed zone can have global knowledge by updating the global model aggregated from all distributed learning models. Moreover, with real-time training capability (i.e in minute times-scale), FedeX is able to retrain its learning model constantly in order to cope with any drift in the normal/abnormal behavior of data coming from devices (for example, drift caused by device aging inside a smart factory).

The rest of our paper is structured as follows. Section 2 discovers related and cutting-edge researches in the field of anomaly detection for ICSs and XAI for Anomaly detection in ICSs. The FedeX anomaly detection architecture will be elaborated in detail in Section 3. The evaluation of the FedeX performance in terms of detection capability, system response time, edge computing capability, and anomalies explanation is presented in Section 4. Finally, the conclusion of our findings is presented in Section 5.

Related Work

For the non time series data type, we can find anomaly detection approaches for ICSs in [START_REF] Das | Anomaly detection in industrial control systems using logical analysis of data[END_REF] , [START_REF] Ghaeini | State-aware anomaly detection for industrial control systems[END_REF], [START_REF] Perales Gmez | Madics: A methodology for anomaly detection in industrial control systems[END_REF]. In [START_REF] Das | Anomaly detection in industrial control systems using logical analysis of data[END_REF], the author proposes a Logical Analysis of Data -LAD-ADS solution using a rule-based method to detect anomalous behaviors in ICS systems over the SWat dataset. LAD-ADS performs detection by extracting rules from a huge of data in the past. However, rule-based systems are often complex and challenging to manage and determine the cause of detected anomalies. In addition, it requires experience from experts and operational engineers in deploying and operating the system. In the same scenario of ICSs, [START_REF] Ghaeini | State-aware anomaly detection for industrial control systems[END_REF] proposes a state-aware anomaly detection method that uses the CUSUM (Cumulative Sum) control chart to the statedependent detection threshold. The training process is done centrally on a large amount of data. In fact, the data in ICSs are often distributed; using a centralized solution can cause disadvantages such as latency for sending all raw data to the central cloud and huge computer resource consumption for training. Moreover, in CUSUM [START_REF] Ghaeini | State-aware anomaly detection for industrial control systems[END_REF], no detection performance metric such as Accuracy, Precision, Recall, or F1-score is revealed except the false alarm rate. In another aspect, training in a distributed environment and the privacy of data is addressed in [START_REF] Perales Gmez | Madics: A methodology for anomaly detection in industrial control systems[END_REF]. In [START_REF] Perales Gmez | Madics: A methodology for anomaly detection in industrial control systems[END_REF], the authors present a methodology called MADICS for Anomaly Detection in Industrial Control Systems using a semi-supervised anomaly detection paradigm with five main steps. The performance of MADICS in terms of Recall is slightly low over its testing dataset. In addition, this mechanism requires a large amount of data for semi-supervised learning, which faces data privacy issues when transmitting a large amount of raw data for training, resource capacity, and computational resources of the system. In terms of detection performance, Fedex is also proved to outperform the previously proposed solutions MADICs [START_REF] Perales Gmez | Madics: A methodology for anomaly detection in industrial control systems[END_REF], LAD-ADS [START_REF] Das | Anomaly detection in industrial control systems using logical analysis of data[END_REF] in the same factory contexts. Additionally, in [START_REF] Kravchik | Detecting cyber attacks in industrial control systems using convolutional neural networks[END_REF], a statistical window-based anomaly detection method is adopted by using various deep-neural network architectures, showing effectiveness in detecting the attacks in a Secure Water Treatment (SWaT) infrastructure which is, in our opinion, not a purely time-series data scenario as well. However, the authors also indicated that their work needs to be improved with the interpretability of the outcomes and the behavior detection of fault ICS components.

For the time-series data type, we have observed various proposed AD solutions. Training in a distributed environment and the privacy of data is also an issue that needs to be addressed in [START_REF] Nguyen | A gradient-based explainable variational autoencoder for network anomaly detection[END_REF]. From the aspect of using Federated Learning to implement an anomaly detection solution in a distributed ICS system, ensuring high accuracy while protecting data privacy, we can find some researches such as [START_REF] Liu | Deep anomaly detection for time-series data in industrial iot: A communication-efficient on-device federated learning approach[END_REF], [START_REF] Mothukuri | Federated learning-based anomaly detection for iot security attacks[END_REF]. In [START_REF] Liu | Deep anomaly detection for time-series data in industrial iot: A communication-efficient on-device federated learning approach[END_REF], the author proposes an FL framework that allows decentralized edge devices to cooperate in training an anomaly detector with an attention mechanism-based convolutional neural network long short-term memory (AMCNN-LSTM) model. Although designing an FL-based approach, but the experiments lack insight analysis in the performance of deploying such a learning model in an edge environment (i.e in weak hardware of an edge node). Due to the complexity of AMCNN-LSTM caused by using multi-layer CNN and LSTM, according to our experience, it is hard to feasibly deploy such a learning model on edge devices, much less for an expectation of achieving low computing complexity for running the learning model in minute-time scale and low power consumption. With a similar lack of performance testing on the edge hardware, work [START_REF] Mothukuri | Federated learning-based anomaly detection for iot security attacks[END_REF] proposes an FL-based anomaly detection approach for IoT networks based on the combination between Gated Recurrent Units (GRUs) and Long short term memory toward detecting anomalies with decentralized on-device data. However, the performance of the proposed method is not good enough in the distributed scenario; accuracy in each FL client is just around 90% on average.

In several other studies [START_REF] Chang | Anomaly detection for industrial control systems using k-means and convolutional autoencoder[END_REF], [START_REF] Al-Abassi | An ensemble deep learning-based cyber-attack detection in industrial control system[END_REF], the authors develop and investigate attack detection solutions in ICS cyberspaces. In [START_REF] Al-Abassi | An ensemble deep learning-based cyber-attack detection in industrial control system[END_REF], the authors propose an attack detection model that uses a Deep Neural Network and a Decision Tree classifier to identify cyber-attacks in the ICS context with an F1-Score of 93.83% with the ICS gas pipeline dataset, which is higher than other algorithms such as Support Vector Machine (SVM), Long Short-term Memory (LSTM), Nave Bayes (NB), Decision Tree (DT), Deep Neural Network (DNN), Random Forest (RF). Study [START_REF] Chang | Anomaly detection for industrial control systems using k-means and convolutional autoencoder[END_REF] uses the semi-supervised techniques by leveraging K-means and Convolutional Autoencoder to protect the ICS system from cyberattack. Like in [START_REF] Al-Abassi | An ensemble deep learning-based cyber-attack detection in industrial control system[END_REF] , the experiments of the proposed methods were performed with the gas pipeline dataset and the water storage tank dataset. However, the anomaly detection performance of the proposed method needs to be improved. In contrast to our study, these studies only focused on evaluating the performance of detection algorithms and did not consider other important metrics when being implemented in the edge environments of an ICS such as detection time and power consumption. Thus, in this paper, we present a comprehensive evaluation of both detection performance and system performance.

Although the studies described above solve challenges surrounding cyberattack detection in ICSs, all of them have not concerned the interpretability of the models detected results up to now. As stated in [START_REF] Amarasinghe | Toward explainable deep neural network based anomaly detection[END_REF], the interpretability of an anomaly detection model is almost as crucial as the prediction accuracy of the model. In the field of explaining the detection outcomes (XAI -Explainable AI), Kasun et al. in [START_REF] Amarasinghe | Toward explainable deep neural network based anomaly detection[END_REF] used a method named Layer-wise Relevance Propagation (LRP) to calculate the input features relevance to explain the trained Deep Neural Network model with DoS attacks detection task. The evaluation is conducted with a subset of NSL-KDD Dataset an old network intrusion detection dataset released in 1999. Even though the combination of solutions to solve the black-box problem of DNN helps domain experts intuitively access the insight of the DNN algorithms, classification accuracy improvement is required when producing predictions in the test set. Very recently, the authors in [START_REF] Hwang | E-sfd: Explainable sensor fault detection in the ics anomaly detection system[END_REF] have proposed to use XAI to interprete anomaly detection outcomes of the multiple Bi-LSTM learning model in an ICS ecosystem. The scope of the ICS is the smart factory of steam-turbine power generation and pumped-storage hydropower generation. This paper can be considered as the forefront of interpreting anomaly detection in the ICS ecosystem.

As a result, in this paper, we design and investigate the FedeX solution for not just ensuring high detection performance and lightweight implementation at the edge devices, but also providing a detailed explanation for the detection model deployed in a liquid storage infrastructure.

Federated learning-based Explainable Anomaly Detection for

ICSs -FedeX

FedeX Overview

In this paper, an architecture using FL for anomaly detection is proposed for ICSs. As Fig. 1 shows, ICSs in smart factories can be organized in various zones (i.e. Zone 1, Zone 2, Zone 3...), and each of which is monitored by a local unit (i.e. Edge 1, Edge 2...) to detect anomalies. Those local monitoring units serve as edge computing stations that run an anomaly detection function based on their own incoming local data. Computing can be carried out at the distributed edges as long as detection algorithms running on top of it require a reasonable computing capacity. If this requirement is fulfilled, then this architecture becomes effective since the detection module is implemented near attack sources which makes the whole detection process respond faster. Moreover, this solution reduces the workload offloading up on the central cloud server as the traditional centralized computing architecture does.

As illustrated in Fig. 1, the FedeX workflow consists of 6 steps, as follows:

Step 1. Anomaly Detection Model on Edges: VAE [START_REF] Kingma | Auto-encoding variational bayes[END_REF] is a tuned Autoencoder architecture to run on top of the edge device for effective anomaly detection. The benefit of VAE is the ability to minimize over-fitting by ensuring that features from its latent space are good enough for data generation. A basic idea here is, if a model was only trained on normal data, then when encountered with anomalous data, the inability to reconstruct data or, more precisely, the range of the reconstruction error that it entails, can signal the presence of anomalous data.

2. Distributed Learning Mechanism for Efficient Resource Sharing: Moreover, this architecture expands to an FL-based VAE (i.e. FedVAE) anomaly detection solution that can solve the problem of missing training data at each edge device when deep learning models often need large amounts of data to train. With the FL technique, the central cloud can federate information with different characteristics from various zones to improve the detection performance of the overall network without the need of having knowledge of original raw data.

Accurate, Fast and Automatic Determination Mechanism of Threshold for Anomaly Detection:

In addition, the SVDD method [START_REF] Tax | Pump failure detection using support vector data descriptions[END_REF] is proposed to automatically determine the threshold for efficient anomalies detection (FedVAE-SVDD). Every sample at the output of the VAEs block with the loss larger than the chosen threshold will be considered anomalous.

Explainable Artificial Intelligent to Interpret the Outputs of Black-Box Learning Models:

Finally, a XAI-based explanation method called SHAP is integrated into our architecture. SHAP enables us to explain the reason why an instance is predicted as an anomaly by showing the contribution of the features to the prediction, thereby enhancing the reliability of the black-box model. As a result, FedeX effectively assists domain engineers in finding the physical cause of anomalies quickly and making responses timely.

In the following sections, we will present a detailed design and deployment of the VAEs, FL, and SVDD as the hybrid-anomaly-detection model at the edge. And finally we will elaborate how to deploy Explainable-Artificial-Intelligent (XAI) for better understanding of anomalies occurred in the ICS of a liquid storage infrastructure.

Design and Development of Variational AutoEncoder (VAE) as Local Training Model on Edge

Since the detection module is implemented on Edge hardware, the overall design is supposed to be lightweight, whilst still ensuring the detection accuracy requirement. Therefore, in our proposed detection method, we try our best to reduce the computing complexity of the algorithms. As elaborated in Figure . 1, the AI-based detection learning model each is deployed locally at each single Edge. Then the model update of each local edge will be done through a global updating, facilitated by Federated Learning. cau nay la de giai thich cho edge computing tai comment 1.3, can bo sung them neu can.. In this research, we propose to utilize VAE for anomaly detection purposes. In fact, many different VAE architectures have been proposed, with different types of layers such as Dense, LSTM, and CNN. We design the VAE encoder and decoder with only two fully connected layers each because this approach aims to achieve the model's simplicity and lightweight, allowing it to be trained on top of edge devices with limited hardware resources, as part of an IoT-based ICS, while also lowering communication costs and providing sufficient detection performance. We also emphasize the real-time training guarantee for this model which will be illustrated in Section 4.

For the background, an autoencoder (AE) is a symmetrical unsupervised neural network, slightly different from other network architectures in that: The network uses the input itself as the ground truth. It consists of 3 main parts: encoder, latent representation, and decoder. Usually, the centre hidden layer has fewer nodes than the input and output layer (a "bottleneck"). Thus, the network learns to compress the input to the bottleneck layer and then from which subsequently restore the input. This middle layer thus becomes the "latent representation" of the input, retaining most information about the input using fewer features. The part of the network before this layer becomes the encoder, and the part after becomes the decoder.

A variational autoencoder (VAE) [START_REF] Kingma | Auto-encoding variational bayes[END_REF] is a combination of the AE with the Variational Bayesian method. But instead of generating a representation in the hidden space for a data point in the original space, the underlying principle behind VAE is to find a probability distribution for that data point. Let X be data generated by inputting a latent variable Z through a random process with network parameters θ. The goal is to model the data as a distribution, Q θ (X). Since the computation cost to calculate Q θ (X) is expensive, we translate the problem into the autoencoder domain by defining the "probabilistic encoder" P φ (Z|X) to approximate the posterior distribution Q θ (Z|X).

The VAE module is designed to combine with FL as follows:

• In each iteration of the training, one option is to let the participating nodes synchronize their VAE models over the parameter server. This, however, necessitates many rounds of communication between the federated nodes and the parameter server, possibly resulting in network communication congestion.

• Instead, we let the participating nodes run a number of local modifications and occasionally synchronize with the parameter server. Specifically, after the nodes receive an updated model from the server, they update it locally by running ρ iterations of the clients' optimizer algorithm, and then send the appropriate data to the server for updating the aggregated model. Finding the optimal choice of ρ for lowering the overall communication cost of the process is a critical trade-off problem that we have to solve.

To update the network weights by backpropagation, as with any deep learning problem, a differentiable loss function must be defined. For VAEs, the objective is to minimize the generative model parameters θ to reduce the reconstruction error between the network's input and output, and also φ, to have P φ (Z|X) as close to Q θ (Z|X) as possible.

For calculating the distance between two distributions, the commonly used function is the Kullback-Leibler divergence function. The loss function called evidence lower bound (ELBO), is obtained:

-L (θ,φ) = log(Q θ (X)) -D KL (P φ (Z|X) || Q θ (Z|X)) ≤ log(Q θ (X)) (1)
The parameters of the model can thus be expressed as:

(θ * , φ * ) = argmax θ,φ L (θ,φ) (x) (2)
As the idea mentioned above, we will sample in the hidden space and feed to the encoder to generate new data. There should be a certain rule for this sampling. A simple yet efficient way is suggested to use the normal distribution N (0, 1). That is, we add a constraint that each data point will be represented by a normal distribution that approximates the N (0, 1) distribution. This is the idea of reparameterization [START_REF] Kingma | Auto-encoding variational bayes[END_REF]. To approximate two distributions, again the KL-Divergence loss of two Gaussian distributions with univariate function is shown below:

D KL [N (µ(X), Σ(X)) N (0, 1)] = 1 2 k exp(Σ(X)) + µ 2 (X) -1 -Σ(X)
(3) with k is the dimension of our Gaussian. Although the cloud has the storage and computing power to manage the volume of data generated in manufacturing, the computationally intensive operations and vast data storage that are hosted on cloud servers may cause a delay. Because this delay is caused by the time required to send, transfer, and process massive amounts of data from IoT devices at production sites. This is a significant issue in a smart factory that must undertake huge monitoring and detection in real-time. Within this context, the concept of Edge-Cloud Computing combined with FL shall arise to circumvent this constraint.

• Firstly, the initial model is created by the Cloud Server as a Weight "Federator".

• The VAE model was then applied to solve anomaly detection. It then subscribes to numerous MQTT topics to which the zones will send the weights of their models.

• After the first model's weights are published to the aggregated model topics, the Cloud Server awaits requests from the VAE model configuration from each zone.

• Local models are trained at each edge based on their own dataset.

• In each communication round, the weights of the trained models ω

(c) r,ρ
are sent to the Cloud Server for FL.

• The Cloud then uses the formula (4) to calculate the weight of the federated global model:

ω r+1 = 1 C C c=1 ω (c) r,ρ (4)
Where:

C: the number of zones ω • Finally, the weight from the federated global model is sent downward to update the local model of each zone.

For every communication round, the server, at first, randomly picks up C zones, then sends the model ω

(c)
r,ρ to each zone. The model, which is the VAE model, will be run ρ iterations locally with the aim to minimize loss function L θ,φ (X). f (ω) stands for the neural network with parameter ω. Moreover, since Adam [START_REF] Kingma | A method for stochastic optimization[END_REF] is the best among the adaptive optimizers in most of the cases, we set it for client optimizer Opt by default. The FedVAE model is presented in Algorithm 1.

Design and Development of FedVAE-SVDD to Determine Thresholds

As aforementioned, in Phase 2, to monitor and automatically determine the threshold for the anomaly detection model -the FedVAE model, we propose to use the Support Vector Data Description (SVDD) method to go over this request accurately while keeping the real-time assurance.

Usually, experts in the industry establish the threshold after attempting a range of values, then select the one that best balances the requirements (performance, true positive, or false negative ,.etc). SVDD also works well as an outlier detection algorithm, especially with high-dimensional datasets, but just like all SVMs, it does not scale to large datasets. Therefore, we suggest a combination of FedVAE and SVDD, as a moderate addition: FedVAE serves as the main anomaly detection model for the distributed system, while SVDD, trained with a small set of error vectors from the output of the FedVAE model, can correspond to finding a small region that encompasses all instances.

Support Vector Data Description or SVDD ( [START_REF] Tax | Support vector data description applied to machine vibration analysis[END_REF]) is a type of support vector method used for single-class classification and outlier detection. The primary idea behind SVDD is to wrap samples in a high-dimensional space with the smallest volume. For the anomaly detection task in which most of the collected data are normal, the hypersphere is usually taken as the boundary around normal samples, separating them from outliers.

Primal Form: Ojective Function:

minR 2 + C n i=1 ς i (5)
Subject to:

x i -I 2 ≤ R 2 + ς i , ∀i = 1, ..., nς i ≥ 0, ∀i = 1, ..., n (6)
Where: The dual formulation is obtained using the Lagrange multipliers.

x i ∈ R m , i =
Objective Function:

M ax n i=1 α i (x i • x j ) - n i,j=1 α i α j (x i • x j ). ( 7
)
s.t 0 ≤ α i ≤ C and n i=1 α i = 1 (8)
Where: α i ∈ R, i = 1, . . . , n are the Lagrange coefficients Duality Information:

The following results are valid depending on the observation position:

• Position of Centre I:

n i=1 α i x i (9)
• Inside Position:

x i -I < R → α i = 0 (10)
• Boundary Position:

x i -I = R → 0 < α i < C (11)
• Outside Position:

x i -I > R → α i = C (12)
The circular data boundary can include amount of very sparse distribution of training observations space that can increase the probability of false positives. Hence, the support-vector-based boundary is usually used rather, since it is more flexible to cover the data with the volume as small as possible.

The SVDD becomes more flexible by replacing the inner product (x i • x j ) with an appropriate kernel function K(x i , x j ). Results 9 through 12 hold true when the kernel function is used in the mathematical formulation. Then the threshold R is calculated using the Kernel function as follows:

r = K(x l , x l ) -2 i (x i , x l ) - i,j α i α j (x i , x j ) ( 13
)
using any x k ∈ S where S is the set of support vectors that have α k < C. Any function that meets the Mercer condition can be used as a kernel function. Some commonly used kernels are Gaussian kernel

K (x, y) = exp -x -y 2 2σ 2 (14)
σ is the Gaussian kernel width Exponential kernel

K (x, y) = exp -x -y 2σ 2 (15)
The Gaussian kernel and the exponential kernel are very similar, with only the square of the norm left out.

Laplacian kernel

K(x, y) = exp(- x -y σ ) (16)
In fact, the Laplace Kernel is completely equivalent to the exponential kernel, except for being less sensitive for changes in the σ parameter. Although both radial basis function kernels above give excellent performance, in our case study, we decide to use the Laplace kernel as default since it allows the SVDD model to run faster.

Scoring:

For each new observation Z, the distance d(Z, I) is calculated as

d(Z, I) = (Z • Z) -2 1 i=1 α i K (Z • x i ) + l i,j=1 α i α j K (x i • x j ) (17)
Thus, any dataset point with d 2 (Z, I) > r 2 are indicated as an outlier.

More specifically in the process of Phase 2, another portion of normal data will be passed through the trained FedVAE model (at the end of Phase 1, which will return a set of loss value vectors. We use these vectors to train the SVDD model in Phase 2: The hypersphere's radius as calculated by Equation ( 24) will be considered as the threshold of the FedVAE-SVDD model .

Each sample belonging to the hypersphere needs to satisfy the condition:

(x i -I) T (x i -I) ≤ r 2 (18)
Where:

x i ∈ R m , i = 1, . . . , n represents the training data r: the radius that represents the decision variable I: the center, a decision variable The dual formulation is obtained using the Lagrange multipliers. Also, the support-vector-based boundary is usually used rather than the hypersphere, since it is more flexible to cover the data with the volume as small as possible. Thus it can reduce the amount of very sparse distribution of training observations space that can increase the probability of false positives. Any function that meets the Mercer condition ( [START_REF] Scholkopf | Input space versus feature space in kernel-based methods[END_REF]) can be used as a kernel function. Therefore, the Objective Function is obtained as follows:

M ax n i=1 α i K(x i • x i ) - n i,j=1 α i α j K(x i • x j ). ( 19
)
s.t 0 ≤ α i ≤ C and n i=1 α i = 1 (20)
Where: α i ∈ R, i = 1, . . . , n are the Lagrange coefficients C: the penalty constant that controls the trade-off between the volume and the errors. K(.) stands for the kernel function. Some commonly used kernels are Gaussian kernel

K (x, y) = exp -x -y 2 2σ 2 (21)
σ is the Gaussian kernel width Exponential kernel

K (x, y) = exp -x -y 2σ 2 (22)
The Gaussian kernel and the exponential kernel are very similar, with only the square of the norm left out.

Laplacian kernel

K(x, y) = exp(- x -y σ ) (23)
In fact, the Laplace Kernel is completely equivalent to the exponential kernel, except for being less sensitive for changes in the σ parameter. Although both radial basis function kernels above give excellent performance, in our case study, we decide to use the Laplace kernel as default since it allows the SVDD model to run faster. The threshold r is calculated as:

r = K(x l , x l ) -2 i (x i , x l ) - i,j α i α j (x i , x j ) (24)
using any x k ∈ S where S is the set of support vectors that have α k < C.

For each new observation Z, the distance d(Z, I) is calculated as

d(Z, I) = (Z • Z) -2 1 i=1 α i K (Z • x i ) + l i,j=1 α i α j K (x i • x j ) (25)
Thus, any dataset point with d 2 (Z, I) > r 2 are indicated as an outlier.

More specifically in the process of Phase 2, another portion of normal data will be passed through the trained FedVAE model (at the end of Phase 1, which will return a set of loss value vectors. We use these vectors to train the SVDD model in Phase 2: The hypersphere's radius as calculated by Equation (24) will be considered as the threshold of the FedVAE-SVDD model .

Algorithm 2: Phase-2: FedVAE-SVDD Input: X -NormalTrainData; C -penalty constant; K(.) -Kernel function; Output: Threshold;

1 Recontruction Data X = F edV AE(X); 2 Get set of error vectors S X = {∆X|∆X = X -X};

3 for x i , x j ∈ S X do 4 if (α i < C and α j < C) then calculate x = K(x l , x l ) -2 i (x i , x l ) -i,j α i α j (x i , x j ); 5 end if 6 end for 7 T hreshold = √
x ; 8 return Threshold ;

Integration of Explainable Artificial Intelligence (XAI SHAP) -FedeX

In anomaly detection, although algorithms related to neural network models tend to benefit rather than signature-based methods and techniques, its drawback is insufficient interpretability. Therefore, the reason why an instance is predicted to be abnormal can not be easily discovered in such cases, which could render researchers vague in analyzing anomaly or outlier outcomes. To overcome this limitation, a widely-used approach called Explainable Artificial Intelligence (XAI) can be adopted.

The aim of XAI is to assist humans to understand the results of solutions using black-box models by the assessment of feature attributions, thereby demonstrating how much each feature participated in making a decision for each data point of the model. With simple machine learning models, such as logistic regression and linear regression, the importance of features can be assessed via the coefficient of each feature in the data set. Meanwhile, as aforementioned, for several complicated models related to neural networks such as VAE, it is difficult to measure or compute the influence level of each feature on an output decision. This is simply because there are a large number of parameters engaging in the model. In fact, with the advent of some state-of-the-art XAI frameworks, this problem has been handled.

There are several effective XAI frameworks such as Local Interpretable Model-agnostic Explanations (LIME) [START_REF] Ribeiro | why should I trust you?[END_REF] and Deep Learning Important Fea-tures (DeepLIFT) [START_REF] Shrikumar | Learning important features through propagating activation differences[END_REF]. However, with the scope of this study, the main XAI approach is based on SHapley Additive exPlanations (SHAP) [START_REF] Lundberg | A unified approach to interpreting model predictions[END_REF] method using the Shapley values, which comes from the theory of the cooperative game [START_REF] Petrosyan | Game theory[END_REF]. By computing an individual player's contribution for each coalition P (a possible subset from feature set) and then averaging over all of these contributions, the Shapley value tells us the payout the one is assigned fairly in received payouts. Similarly, in terms of XAI, for an individual decision, each feature value can be considered as a player, and the payouts can be treated as the decision. Mathematically, the Shapley value of a feature in a particular prediction model can be defined as:

ξ n (f, a) = P⊆a |P|!(m -|P| -1)! m! [f a (P) -f a (P \ n)] (26)
Where, ξ n is the Shapley value for feature n f is a "black-box" model that needs explaining a is an input datapoint a is the simplified data input P is one of all possible subsets of feature set, considered as a coalition m is the number of features in the dataset Due to the fixed input size, commonly, the features of a model that are omitted in the Eq.( 26) are substituted with random input values from the background dataset. Looking at this formula, it can be seen that the total possible subsets of an m-feature set used for interpretation is 2 m , which leads to the complexity of computing Shapley values being massive if m increase more. Therefore, to deal with this issue without calculating all combinations, Kernel SHAP [START_REF] Lundberg | A unified approach to interpreting model predictions[END_REF] can be employed by sampling feature subsets and then fitting them into a linear regression model:

K(a) = γ 0 + γ 1 a 1 + γ 2 a 2 + γ 3 a 3 + ...γ m a m (27)
where, a i is a encoded feature and γ i is the corresponding coefficient representing the contribution of the feature to the model, i = 1, 2, ..., m. In this linear regression model, the variables a i are encoded according to the presence or absence of ones. Thanks to this, Sharply values can be approximated as the output values of the trained linear regressing model. In this work, our FedeX architecture is applied with SHAP to account for the impact level of features on the anomalies that are predicted from the FedVAE-SVDD model, via their SHAP values. In this case, as depicted in Fig. 3, the Kernel SHAP is fed in with the FedVAE-SVDD model and the test data to construct a local linear regression explanation model and compute the SHAP values. Subsequently, the explanatory model computes the SHAP value of classified anomalies and displays them visually. As feature values are measured by sensors, by using this explanation, operators or domain engineers can easily determine the sensors likely causing the abnormality and make a faster detection response.

Performance Evaluation

In this section, we evaluate the FedeX architecture in various aspects. From the detection performance perspective, the FedVAE-SVDD learning model is evaluated in comparison with different cutting-edge solutions for Anomaly detection in ICSs. The resource requirement of FedVAE-SVDD over the embedded edge device is also taken into account. Moreover, the results of using the XAI-SHAP technique to interpret the prediction results of FedVAE-SVDD are also described with the main case study in a SCADA liquid storage infrastructure dataset [START_REF] Laso | Dataset of anomalies and malicious acts in a cyber-physical subsystem[END_REF].

Experiment Setup

For the evaluation, we implement FedeX in a small-scale IoT testbed including:

• 4 Raspberry-Pi-4-Model-B kits acting as edge devices; Raspberry-Pi-4 equipped with quad-core 1.5 GHz ARM Cortex-A72 processor and 4

GB RAM with 32-bit Raspbian OS

• 01 Dell Precision 3640 Tower Workstation serves as Cloud Server; the workstation with Intel Core i710700K 3.8 GHz (up to 5.1 GHz), 16 GB RAM, working on Linux operating system.

• All edge devices and the Cloud Server are connected by a router through a WIFI interface At the edge devices (i.e. Raspberry-Pi-4), we implement our FedeX framework in Python 3 with the Tensorflow 2 platform, which is built with the support of the FL framework -FedML [START_REF] Li | Detecting cyberattacks in industrial control systems using online learning algorithms[END_REF]. In the FedeX architecture, the edge devices and cloud server exchange the weights and bias matrix of the VAE model using the standardized MQTT protocol for an IoT environment [START_REF] Light | Mosquitto: server and client implementation of the mqtt protocol[END_REF]. EMQ X Broker (2021) is hosted on the cloud server as an MQTT broker for better long-term performance. We discover EMQ X Broker as the most scalable open-source broker that could accept more advantageous devices linked to the server.

A Case Study for ICSs

In our main case study, we consider the SCADA liquid storage infrastructure dataset [START_REF] Laso | Dataset of anomalies and malicious acts in a cyber-physical subsystem[END_REF], which simulates a fuel storage system supplying an automated production line monitored by an ICS system. The high-level overview of the testbed system is shown in Fig. 4. As depicted in Fig 4 , the system is composed of the main tank and secondary tank with a capacity of 9 and 7 liters, respectively. Data is collected by connecting the sensors to a PLC; four discrete sensors in the main tank (i.e., IN0, IN1, IN2, IN3) and one in the secondary tank (i.e., R4) are used to measure the level of fuel in tanks. Besides, pump1 and pump2 control the flow of fuel between two tanks, it is also connected with two sensors (i.e., PP, PG). PLC registers 2 through 4 provided output data defining the system's state used to analyze the data obtained and register 2 provides the bits that indicate the discrete sensors' binary status. To extract the state of each sensor separately, a population count can be performed on the register. Register 3 holds the pump's active or inactive state, whereas Register 4 holds the ultrasound sensors' step value from 0 to 10,000. (e.g. Step 3,000 represents 2.1 liters of liquid in the tank).

As described in [START_REF] Laso | Dataset of anomalies and malicious acts in a cyber-physical subsystem[END_REF], the data set consists of 14 distinct scenarios. Each scenario includes one of 5 operational situations (such as sabotage, breakdown, accident, or cyber-attack) as well as six affected components. Affected components are those parts of the system that are directly impacted by the abnormality.

The data set has 10 features, and we use the normal data to perform the training model for all scenarios. In order to simulate the 4 distributed zones, we split the original data set into 4 independent subsets, each of which is used for local training corresponding to each of the 4 edge devices. After that, a test set containing both normal and abnormal data is utilized to evaluate the model performance. Accordingly, the proportion of the training set and test set in the whole data set is 0.7 and 0.3, respectively.

FedeX Performance Evaluation 4.3.1. Detection Capability

As aforementioned, in our FedeX testbed, we set up 4 distributed edge zones. To evaluate the detection performance of the FedVAE-SVDD model, the common detection metrics such as F1 score, Accuracy, Recall, Precision are measured in various testing scenarios:

• Scenario 1: FedVAE-SVDD versus its centralized counterpart in our main case study of ICS (i.e. the SCADA liquid storage infrastructure dataset [START_REF] Laso | Dataset of anomalies and malicious acts in a cyber-physical subsystem[END_REF])

• Scenario 2: FedVAE-SVDD and its centralized counterpart versus other previously proposed AD solutions in our main case study of ICS and different SCADA datasets.

As elaborated above, the SVDD model is used to determine the optimal detection threshold for the FedVAE model depending on each different training dataset, in order to achieve good detection performance. Here, we use the Laplacian kernel as default to calculate the kernel distance. Any loss value greater than this threshold is considered an outlier. In our experiments, the optimal thresholds found for 4 distributed zones: Zone 1, Zone 2, Zone 3, and Zone 4 are 0.11, 0.09, 0.09, and 0.09 respectively. Scenario 1: FedVAE-SVDD vs. its centralized VAE-SVDD over our main case study of a non-time series SCADA dataset

In this scenario, we measure the detection performance of our FedVAE-SVDD solution with the Centralized VAE-SVDD in which the training process is supposed to be carried out at the Central Cloud. Since the learning model converges after 3 communication rounds, the results retrieved after the rounds are shown in Table . 1.

In fact, in the context of seeking a detection scheme in a distributed manner that can provide a fast system response upon attacks or anomalies, and can work lightweight to cope with the limited computing capacity of edge devices in an IoT environment, we always have to think about the trade-off due to its distributed nature. This can be considered as an advantage of Decentralized Learning, so the results when comparing our model with the Centralized learning method are slightly higher.

Considering the performance of FedVAE-SVDD only, we can see that all detection metrics are very good. Only Precision in Zone 4 gets a bit low at 0.864. However, in a smart factory, even the smallest abnormal incident can adversely affect the entire factory. So in general, we need to avoid discarding anomalies (i.e. Recall is important) and accept that sometimes the model can miss detecting a normal sample to be abnormal (i.e. Precision). Because engineers can easily test it and then operate the factory properly. Therefore, the Recall results of our model prove that this model can be a very good candidate to be deployed in a smart factory. First of all, as the main case study of our detection design, we compare the performance of FedVAE-SVDD with the results of other AD solutions for ICSs found by a recent research work [START_REF] Selim | Anomaly events classification and detection system in critical industrial internet of things infrastructure using machine learning algorithms[END_REF] who work in the same ICS context (i.e the same SCADA liquid storage infrastructure data set [START_REF] Laso | Dataset of anomalies and malicious acts in a cyber-physical subsystem[END_REF]). The results can be seen in Table . 2 in which our FedVAE-SVDD detection solution is shown to outperform other machine-learning algorithms LR, LDA, KNN, CART, NB, SVM in all metrics (i.e. Accuracy, Precision, Recall, and F1-Score). Even its centralized counterpart (i.e the centralized VAE-SVDD) performs better than those ones, which shows that this learning model is suitable for such a case study. Note that the detection results of our proposed solution are retrieved after 3 communication rounds as the convergence point of the learning model. In the extended experiments, we run our FedVAE-SVDD model and the centralized VAE-SVDD model over the SWaT dataset [START_REF] Goh | A dataset to support research in the design of secure water treatment systems[END_REF] which is the case study of several other researches such as LAD-ADS [START_REF] Das | Anomaly detection in industrial control systems using logical analysis of data[END_REF], MADICS [START_REF] Perales Gmez | Madics: A methodology for anomaly detection in industrial control systems[END_REF], and 1D CNN ensembled attacks [START_REF] Kravchik | Detecting cyber attacks in industrial control systems using convolutional neural networks[END_REF]. As it can be seen in Table . 3, since SWaT represents a data type that is not purely time-series, therefore similar to the case of SCADA liquid storage dataset, FedVAE-SVDD and its centralized counterparts (VAE-SVDD) perform quite well in comparison with LAD-ADS [START_REF] Das | Anomaly detection in industrial control systems using logical analysis of data[END_REF], MADICS [START_REF] Perales Gmez | Madics: A methodology for anomaly detection in industrial control systems[END_REF], and 1D CNN ensembled attacks [START_REF] Kravchik | Detecting cyber attacks in industrial control systems using convolutional neural networks[END_REF] in terms of Recall and F1-score. In the FL-manner, the detection model FedVAE-SVDD provides higher detection performance which is varied from one zone to another one depending on the local data of each zone.

Again, let us note that Recall and F1-score are the 2 important metrics for ICSs. In the case of imbalanced datasets like these considered datasets in which the number of abnormal samples is so much different from the number of normal samples, a good F1-score figure is necessary. Although the above results demonstrate that FedeX achieves good anomaly detection performance, we want to investigate the reasons why they are predicted so. Since our case study is based on the data set gathered in a liquid storage infrastructure [START_REF] Laso | Dataset of anomalies and malicious acts in a cyber-physical subsystem[END_REF] as described in Section 4.2, we expect that FedeX could support domain engineers quickly and visually in finding and checking abnormal behavior of those sensors or actuators . Therefore, SHAP is employed to identify how features contribute to the anomalies predicted by the FedVAE-SVDD model. Thanks to this, decisions and priorities in checking and maintaining systems can be made effectively, allowing operators to save more time.

From a practical perspective, anomalies can arise from various threats such as accidents, sabotage, breakdown, and cyber-attack. This promotes us to perform two explanation scenarios, where SHAP is employed to explain two sets, corresponding to two different intervals, drawn randomly from anomalous samples predicted in the test set at Zone 1. The results of both scenarios are visualized in Fig. 5, a summary plot for the distribution of SHAP values over whole computed data points, pointing out the impor-tance of features through their impact. In the visualizations, the dots in each feature correspond to the SHAP values of each data point, piling up to depict density. The position on the x-axis is denoted by Shapley values and on the y-axis by features ordered as per importance. Besides, the value of features from low to high is displayed by color gradation. As characterized in Fig. 5(a), we can observe that overwhelmingly, Slope and R4 are critical features, while the other features do not contribute to the anomaly. Consequentially, it can be inferred that the ultrasound sensor which measures the physical values of the R4 and Slope feature may break down. This incident can come from some weather factors like humidity. Therefore, by checking the ultrasound sensor quickly, domain engineers can make reasonable solutions, without verifying other physical components in the system. On the other hand, if the sensor still works properly, i.e false alarm occurs, the operator can consider retraining the model for higher anomaly detection accuracy. For the remaining scenario, Fig. 5(b) shows that IN3 is the most crucial feature, while Slope, IN2, R4, and IN1 have a remarkable influence. Based on these signs, as an engineer, we could determine that the anomaly is most likely to arise from sabotage impacting physical components such as the discrete sensors in the main tank and the ultrasound sensor, thereby prioritizing checking them.

In both of these scenarios, SHAP suggests that the R4 feature has a significant impact on predicted anomalies, similar to the analyses mentioned in a SCADA dataset research [START_REF] Hindy | Improving siem for critical scada water infrastructures using machine learning[END_REF]. The authors confirm that the ultrasound sensor badly affects most of the abnormal scenarios in the dataset and the value of R4 measured by this sensor is most significant. Accordingly, it can be seen that our XAI-based explanatory solution is capable of precisely identifying the primary cause related to the anomaly in reality.

In conclusion, based on these positive findings, we would like to make some comments and recommendations. Firstly, our scheme can make a comprehensive explanation for detected anomalies, boosting the reliability of FedeX. Besides, if there are the occurrence of unknown threats, FedeX will still support operators to determine affected physical components and come up with timely responses rather than inspecting the entire system. This issue may not be solved by other multi-class classification-based anomaly detection solutions. Furthermore, based on data records, we recommend that domain engineers should run SHAP periodically, for example, once per week, to check and schedule system maintenance depending on attack types, or to retrain the model for higher detection performance.

Edge Computing Capacity

Deploying a learning model at the edge is challenging due to limited capacity of those embedded devices. Therefore, to get insight into the efficiency and feasibility of the FedeX architecture, we conduct a few experiments for the FedVAE-SVDD training phase to evaluate the edge performance during the training, based on some metrics such as: power consumption, CPU usage, memory usage, and model running time. load the power consumption of the range from 5000 to 6000 mW on average to afford the worst case. These real-world metrics give us a better idea of how deploying distributed machine learning models on edge devices will consume more energy for that computation. It is conspicuous that the running time in the whole process is very fast, but in the worst case, the SVDD phase accounts for 100% of the CPU usage while this ratio of the FedVAE phase is just over 20%. Based on these findings, we would like to make a few recommendations. Firstly, in reality, with a runtime of only 70 seconds, the threshold update process in the SVDD phase can be retrained during system maintenance time or the night on schedule, rather than implemented on a real-time scale (i.e minutes or seconds scale). Thanks to this, other services would not be interrupted on the edge device every update time. From another perspective, these findings seem to be an acceptable tradeoff between the running time for high detection performance and the hardware resource. Furthermore, it is possible to consider upgrading to edge hardware devices with higher processing capacity than Raspberry-Pi-4. With more powerful edge hardware, the FedVAE-SVDD model will be the effective detection model for such a factory. • Memory usage: In the same experimental setup with the power measurement, the percentage of memory usage in the VAE-SVDD phase at an edge device is demonstrated in Fig. 8. Obviously, throughout the period of 150 seconds of the VAE-SVDD training process with one communication round, the memory usage of the VAE and SVDD phase is quite steady, with just over 9% and 14%, respectively. With these ratios, it can be inferred that in the training process, the memory resource is still available for other tasks.

• Model running time: Using deep learning to detect anomalies inside ICS is common worldwide, but we always have to face the fact that training models take quite a long time to train with the whole data set, even with smaller data sets from distributed zones. We have usually observed the computing running time on the scale of an hour or a few hours; and those figures mean that the system should be only Basically, this shows that FedeX not only responds to anomalies in real-time but also handles the running time problem in a trade-off for a high performance described in a previous work [START_REF] Huong | Detecting cyberattacks 960 using anomaly detection in industrial control systems: A federated learning approach[END_REF], especially for the smart factory scenario as presented in this paper.

Conclusion

In this paper, we have elaborated our proposed hybrid model which combines an effective and fast detection scheme based on VAE and SVDD with the Federated-Learning technique that enables the hybrid model to perform efficiently on the weak computing hardware of distributed edge devices installed in the IoT-based system of a Smart Factory. With the FL architecture design, the detection task is distributed to smaller local zones located in the last premise of traffic senders. Therefore, anomalies or attacks can be quickly identified and quarantined in each separate zones. This FL architecture also helps to deal with Big Data created from a variety of devices inside a huge smart Factory 4.0 of the future.

In addition, this detection model stands out with a very fast training time in the minute-time scale (i.e 7.5 minutes). In case IoT devices in a Factory are aging, leading to changes of data patterns over time (i.e. concept drift), FedeX still works well since it can be retrained quickly every 7.5 minutes.

Moreover, the scheme has been proved to achieve high anomaly detection metrics such as Accuracy, Precision, Recall, F1-score, especially in the context of distributed training environment where the different edge of numerous zones trains their own model from different datasets.

In this paper, we also introduce XAI to help enhancing anomaly detection by interpreting how features contribute to the anomalies predicted by a "black-box" ML-based learning model. This way paves the way for engineers to have a deeper outlook on checking the systems more effectively.

1 ○:

 1 The edge device uses the sensing data collected from nodes within a zone as a local dataset. Step 2 ○: The edge device performs the local model (i.e., VAE model) and the mechanism for determining thresholds at the last communication round (i.e., FedVAE-SVDD model) training on the local dataset. Step 3 ○: The edge device uploads the weight matrix to the cloud aggregator. Step 4 ○: The cloud aggregator obtains a new global model by aggregating the weights uploaded by the edge device. Step 5 ○: The cloud aggregator sends the new global model to each edge device. The steps above are repeated until the global model achieves optimal convergence. This ideal global model can be used by decentralized devices to conduct anomaly detection tasks. Step 6○: Periodically, the XAI-SHAP model will be run to interpret and verify the anomaly detection model; and identify the anomalycausing elements in ICSs.

Figure 1 :

 1 Figure 1: Our FedeX Model in ICS

Figure 2 :

 2 Figure 2: VAE Structure: the data are represented via a normal distribution, and the data dimensions are independent random variables.

Algorithm 1 : 17 V

 117 Phase-1: FedVAE Input: Initial model w 0 , Client optimizer Opt Output: V AEcomplete -Trained VAEs model in each client 1 N -the number of zones; 2 Rounds -number of communication rounds; 3 for r = 1 to Rounds -1 do 4 Server randomly picks C zones; 5 Server sends ω r to C zones ; 6 for node c ∈ C do 7 ω r,0 ← ω r ; 8 for t = 0 to ρ -1 do 9 Compute stochastic gradient 10 ∇f (X, ω) = ∇L θ,φ (X); 11 set ω (c) r,t+1 ← Opt ( ∇f c (ω (c) r,t ), ω (c) r,t , α, t); AEcomplete = f (ω Rounds ); 18 return V AEcomplete 3.3. Design and Development of Federated-Learning based VAE -FedVAE In our Federated-Learning based VAE model (or called FedVAE), each edge device performs the training and detection process with local data from each manufacturing area, and the Edge device only sends information of the weight matrix of the trained model to the cloud server, rather than the entire raw data, as a traditional cloud-based training system would.

  : the weight of the local model of zone c at round r ω r+1 : the federated global model's weight at round r + 1.

1 ,

 1 ..., n indicates the training data R: radius represents the threshold ς i : the slack of each variable I: the centre C = 1 ne : the penalty constant that controls the trade-off between the volume and the errors e: the expected outlier fraction Dual Form:

Figure 3 :

 3 Figure 3: FedVAE-SVDD phase

Figure 4 :

 4 Figure 4: High level architecture of the SCADA liquid storage infrastructure system.

Scenario 2 :

 2 FedVAE-SVDD vs. other Anomaly detection solutions over different SCADA data sets

2 Figure 5 :

 25 Figure 5: Summary plot of SHAP values

Figure 6 :

 6 Figure 6: Power consumption of an edge device in one communication round

Figure 7 :

 7 Figure 7: CPU usage of an edge device in one communication round

Figure 8 :

 8 Figure 8: Memory usage of an edge device in one communication round

  retrained periodically on the scale of an hour, day, or week, since it can not capture any sudden change of traffic patterns in real-time. For example, if IoT devices inside a factory are aging over time, resulting in a change in the data set characteristic, we will need to wait until the next period of the training model operation to retrain it again. In contrast to this setback, our outcomes in the testbed indicate that in terms of time, Raspberry-Pi-4 takes relatively little total time of 150 seconds to run the FedVAE-SVDD model (with just about 70 and 80 seconds in the FedVAE phase and SVDD phase respectively) in each communication round, whereas still ensuring the high detection performance depicted in Table.1. In our case study, it just needs to run 3 communication rounds for the training to converge and the detection results shown in the previous section are retrieved after 3 rounds. Therefore, the FedVAE-SVDD model takes only roughly 450 seconds (i.e 7.5 minutes) overall to produce such high detection performance.

Table 1 :

 1 FedVAE-SVDD performance measured in 4 zones vs. Centralized VAE-SVDD over the SCADA liquid storage infrastructure dataset[START_REF] Laso | Dataset of anomalies and malicious acts in a cyber-physical subsystem[END_REF] which is supposed to be slightly lower than the detection performance of the centralized monitoring and training manner.However, as we can see, in our ICS main case study, the hybrid FedVAE-SVDD solution even outperforms the Centralized learning manner. It can be explained that Federated learning offers the improvement of generalizability of the VAE-SVDD model through the collaboration of multiple edge devices by taking advantage of separate data sources when compared to a single global model under data heterogeneity. FL eliminates a single point of failure

		Zone 1 Zone 2 Zone 3 Zone 4 Centralized
	Threshold 0.11	0.09	0.09	0.09	0.26
	Accuracy 1	0.9587 0.9992 0.9210 0.9017
	Precision 1	0.9237 0.9985 0.864	0.9059
	Recall	1	1	1	0.999	0.9806
	F1	1	0.96	0.9992 0.9269 0.9418
	AUC	1	1	1	0.92	0.9
	with detection performance

Table 2 :

 2 Our proposal vs. other anomaly detection solutions over the SCADA liquid storage dataset

	Learning model	Accuracy Precision Recall F1-Score
	LR	0.87	0.78	0.51	0.49
	LDA	0.88	0.87	0.53	0.53
	KNN	0.91	0.85	0.7	0.75
	CART	0.94	0.86	0.86	0.86
	NB	0.67	0.63	0.80	0.60
	SVM	0.91	0.90	0.68	0.74
	Centralized VAE-SVDD	0.9017	0.9059	0.9806 0.9418
	FedVAE-SVDD @ Zone 1 1	1	1	1
	FedVAE-SVDD @ Zone 2 0.9587	0.9237	1	0.96
	FedVAE-SVDD @ Zone 3 0.9992	0.9985	1	0.9992
	FedVAE-SVDD @ Zone 4 0.9210	0.864	0.999 0.9269

Table 3 :

 3 Our proposal vs. other anomaly detection solutions over the SWaT dataset

	Learning model	Accuracy Precision Recall F1-Score
	LAD-ADS[10]	N.A	0.939	0.891 0.914
	MADICS[12]	0.9659	0.984	0.75	0.851
	1D CNN ensembled attacks[13] N.A	1.0	0.853 0.920
	1 FedVAE-SVDD @Zone 1	0.942	0.942	0.9999 0.97
	FedVAE-SVDD @Zone 2	0.9727	0.9718	1	0.9857
	FedVAE-SVDD @Zone 3	0.9427	0.9427	1	0.9705
	FedVAE-SVDD @Zone 4	0.9433	0.9433	1	0.9708
	Centralized VAE-SVDD	0.9725	0.9751	0.9962 0.9855
	4.3.2. Explainable AI -SHAP

Acknowledgment

This work was supported by Hanoi University of Science and Technology (HUST) under Project T2021-PC-010.


Fernando F Dos Santos

Paolo Rech

Performance-Reliability Trade-Off in Graphics Processing Units

We show that most performance improvements in GPUs increase the number of executions correctly completed before experiencing a failure. We consider four different performance improvements: architectural solutions, software implementations, compiler optimizations, and degree of parallelism.

I. INTRODUCTION

Graphics Processing Units (GPUs) have evolved from being devices dedicated to gaming, graphics, and video rendering to flexible accelerators for a variety of High Performance Computing (HPC) and safety-critical applications, including autonomous vehicles. This market shift led to a burst in GPUs' computing capabilities and efficiency, significant improvements in the programming frameworks, and a sudden concern about their hardware reliability.

A common belief is that reliability and performance are conflicting proprieties with opposite requirements that should be traded-off, trying to achieve a sufficient reliability level without losing too much on performance or ensuring that the recovery from faults does not significantly impact the performance. In this paper, we demystify this false myth on GPUs. Our goal is to show how and why better performances, in most configurations, despite increasing the error rate, can still lead to a higher amount of correctly processed data.

We consider four common and effective ways to improve GPUs performance that are available to the programmer and do not require extra resources, specific hardening solutions, nor extra effort to be implemented: (1) software optimizations, i.e., algorithm efficient implementations, (2) degree of parallelism improvement, i.e., increase the number of parallel threads reducing the operations in each thread, (3) architectural solutions, such as tensor core and mixed-precision, and (4) compiler optimizations. Overall, we consider more than 50 different configurations.

Intuitively, the best performance is achieved on a GPU when its parallel capabilities are fully exploited. Nevertheless, when more operations are executed in parallel, the device error rate might increase. Thus a higher occupation of the GPU hardware delivers higher performance at the cost of a higher error rate. However, in the vast majority of the cases, for GPUs the performance gain grows faster than the error rate. In other words, the benefit of better using additional resources (shorter execution time) is, most of the time, higher than the possible drawback (higher FIT rate).

II. RELIABILITY-PERFORMANCE EVALUATION METRICS

The primary metric to measure the reliability of a device is the Failure In Time (FIT) rate. The FIT rate depends on the amount of resources used for computation, i.e., the cross section, and their corruption probability of affecting the output, i.e., the Architectural or Program Vulnerability Factors (AVF or PVF) [START_REF] Mukherjee | Proceedings of the 36th Annual IEEE/ACM International Symposium on Microarchitecture[END_REF], [START_REF] Sridharan | IEEE 15th International Symposium on High Performance Computer Architecture[END_REF]. The FIT rate, being by definition a rate, does not depend on the execution time. If the same amount of memory is exposed for a time interval t or 2t, its FIT rate will not change. In fact, under a constant flux (neutrons/cm 2 /sec), in 2t, we expect twice the errors and a double neutrons fluence (n/cm 2 ) to hit the device, giving the same error rate. Similarly, executing x or 2x sequential (independent, for simplicity) instructions does not change the code FIT rate. On the contrary, if the additional x instructions are executed in parallel with the original sequence, the FIT rate is expected to double (same execution time, same fluence, but doubled the error rate). We use these premises to comment on how performance improvement can impact the GPU FIT rate and reliability. The interesting aspect, only apparently playing against GPUs, is that a slower execution does not increase the FIT rate while using more parallel resources or bigger hardware cores, with a potential benefit on GPUs performance, increases the FIT rate. We will demonstrate in Section V that on GPUs most of the common optimizations increase the performance faster than the error rate.

To combine error rate and performance, the Mean Instructions, Executions, or Work Between Failure metrics (MIBF, MEBF. MWBF) were introduced [START_REF] Weaver | Proceedings. 31st Annual International Symposium on Computer Architecture[END_REF], [START_REF] Reis | 32nd International Symposium on Computer Architecture (ISCA'05)[END_REF]. The idea is to consider how many instructions, executions, or workloads can be correctly completed before the error occurrence. Considering a constant FIT rate, then the faster configuration will have a higher MIBF, MEBF or MWBF.

III. GPUS PERFORMANCE IMPROVEMENT STRATEGIES

In this Section, we discuss the reasons why, in GPUs, better performances can provide a higher number of corrected executions before the error occurrence even if they increase the error rate. The optimizations we consider, as shown in the cartoon of Figure 1, are:

Algorithm Implementation: The performance of an algorithm is improved maximizing the GPU occupancy and the locality of data. While an optimized algorithm provides significant performance gain, the optimization normally involves code modification and, thus, the executed operations differ, 1) optimizing the code, (2) taking advantage of dedicated architectural solutions, (3) increasing the Degree of Parallelism, or (4) using aggressive optimization flags will also impact the error rate (not necessarily for worse, as we will show). We aim to identify the configurations that provide a performance gain higher than the FIT rate increase, thus improving the MEBF.

impacting both the raw FIT rate (different functional units have different FIT rates) and the PVF (the probability for a fault to reach the program output). The optimizations normally increase the computation density, i.e., the number of operations executed in parallel and, then, increase also the FIT rate. This is confirmed by the experimental data in Section V-A.

Degree of Parallelism (DOP): when a workload needs to be executed on a GPU, the higher the number of parallel threads in which the workload is divided, the higher the DOP. A higher DOP improves performance but a higher number of parallel threads means that a larger area of the GPU is used for computation and, thus, the FIT rate increases.

Architectural Solutions: The GPU architecture includes dedicated functional units for strategic operations. This is the case of tensor cores [START_REF] Jia | [END_REF] and mixed-precision [START_REF][END_REF]. A tensor core is a hardware unit that multiplies 4x4 matrices in 1 clock cycle [START_REF] Jia | [END_REF]. Obviously, a tensor core is much larger than a multiplier or an adder [START_REF] Jia | [END_REF] and, thus, it also has a higher error rate than a single ADD or MUL. The use of mixedprecision (32bit, 16bit, or lower precision) cores is a peculiar optimization as, besides being faster, lower precision cores are also smaller and, thus, are expected to reduce not just the execution time but also the FIT rate of the code.

Compiler Optimizations: the process of translating the C++ CUDA code into GPU executable is very complex. Even if the source code does not change, the compiler optimization modifies the machine code to be executed. The compiler modifies the instructions distribution and memory utilization, which changes the way the software is mapped in the underlying hardware. As we show in Section V-D, generally, the higher the optimization, the more reliable the execution.

IV. EVALUATION METHODOLOGY

To evaluate the FIT/performances trade-off we consider several representative applications from different computing domains executed on three different GPU architectures.

A. Codes and Devices

The codes we test are: Matrix Multiplication (MxM and GEMM), LavaMD, Fast Fourier Transform (FFT), Sorting (Quicksort and Mergesort), Needleman-Wunsch (NW), and an object detection convolutional neural network (YOLO).

We consider three GPUs architectures: Fermi (GTX480), Kepler (Tesla K20), and Volta (Tesla V100) NVIDIA GPUs. Fermi, belonging to an older GPU generation, has computing cores that support only float precision (FP32), while Kepler supports float and double (FP64) precisions and Volta GPUs support three IEEE754 float point precisions (FP64, FP32, and FP16) plus eight tensor cores.

In some configurations, we also compare the results when the available Single Error Correction Double Error Detection (SECDED) Error Correcting Code (ECC) is turned ON or OFF. We only compare different configurations on one device.

B. Neutrons Beam Experiments

To measure the FIT rate, we take advantage of controlled neutron beam experiments performed in the past few years. Beam experiments were performed at ChipIR facility of the Rutherford Appleton Laboratory (RAL) in Didcot, UK, and at the LANSCE facility of the Los Alamos National Laboratory (LANL) in Los Alamos, US. The available neutron flux at ChipIR and LANSCE are about 10 6 n/(cm 2 /s) and 10 5 n/(cm 2 /s), respectively, i.e. about 7-8 orders of magnitude higher than the terrestrial flux (13n/(cm 2 × h) at sea level [7]). Since the terrestrial neutron flux is low, in a realistic application, it is highly unlikely to observe more than a single corruption during the program execution. We have carefully designed the experiments to maintain this property (observed error rates were lower than 1 error per 1,000 executions). Experimental data, then, can be scaled to the natural radioactive environment without introducing artifacts.

Figure 2 shows a portion of the setup installed at ChipIR during one of the multiple test campaigns we performed. While ChipIR and LANSCE have been shown to provide similar FIT rates, we only compare configurations tested at the same facility to reduce data uncertainty. We always compare configurations that differ of only one characteristic. Thus, we can derive that the observed trend is strictly related to the characteristic that we have modified. Moreover, thanks to the application profile and analysis, we can provide a justification for the obtained results.

V. PERFORMANCE AND ERROR RATE COMPARISON

In this section, we characterize four different solutions to improve GPU performance: (a) algorithm optimizations, (b) degree of parallelism, to fully use the GPU parallel architecture, (c) architectural solutions, such as tensor core and mixedprecision, and (d) compiler optimization. For each evaluated solution, we consider the benefit in terms of performance and the impact on the error rate.

To compare the performance and the error rate of the different codes implementations, we plot the relative FIT rate and execution time of the tested configurations. That is, we divide the FIT rate (execution time) of each tested configuration by the lower FIT rate (execution time) we measured. All experimentally measured relative FIT rates and MEBF are shown with 95% confidence intervals.

A. Algorithmic Solutions

The same algorithm can be implemented in various ways and, on GPUs, it is also necessary to find the parallel implementation that delivers better performances. We have selected four common algorithms for GPUs: 512x512 FFT, 2Kx2K MxM, 16K NW, and 32M Sorting. For each, we have considered a naive implementation, i.e., not crafted to fully exploit the GPU, and the optimized implementation.

In the naive implementation of MxM, each thread is responsible for computing one output element, performing 2048 sums and additions, thus saturating the number of registers and caches available in the SM. The efficient implementation, also known as GEMM, is taken from the NVIDIA cuBLAS library and, by dividing the matrices in custom tiles sizes and applying specific algorithm policies, ensures that each thread works only on local data and increases the parallel occupancy of the SM [START_REF] Volkov | Proceedings of the 2008 ACM/IEEE Conference on Supercomputing[END_REF]. For FFT, in the naive implementation, each thread computes an FFT on 512 points independently from the others, forcing most data to be stored in the main memory, increasing the memory latency. The efficient FFT implementation, developed by Volkov and Kazian [START_REF] Volkov | [END_REF], divides the 512 FFT computation among 64 threads (2 warps), using only shared data and avoiding synchronizations. For NW (sequencing of DNA), in the naive implementation, each algorithm step, performing a diagonal search, is done in a dedicated kernel while the optimized implementation organizes the input in matrices of 32x32 cells, and each group is assigned to a warp, reducing the number of kernels and improving the efficiency. For sorting, the naive implementation is quicksort, and the optimized one is mergesort. Quicksort suffers from long memory latencies as the whole input array does not fit into the caches, while mergesort, being recursive, better exploits GPU characteristics.

Figure 3 shows the relative FIT rate, execution time, and MEBF of the naive and optimized implementation of the four algorithms executed on Kepler K20 GPUs, with ECC OFF. As shown, the naive implementation is from 0.5 × (for NW) to 6 × (for MxM) slower than the optimized one. The FIT rate has the opposite trend for all the codes, with the optimized implementation having an up to 2 × higher error rate than the naive one (for FFT). The higher FIT rate is caused by the denser computation imposed by the optimizations.

A promising result is that, for all the considered optimizations, the performances benefit is higher than the FIT rate increase. The MEBF (number of executions completed between failures) is always higher for the optimized version, being up to 5× higher for MxM. A optimized algorithm is then to be preferred, if the goal is to increase the amount of data correctly produced by the GPU. On average, the optimized version of the algorithm allows to correctly complete about 2.2 × more executions than the naive one.

B. Degree of Parallelism

To maximize the performances of a code on a GPU, it is fundamental to increase as much as possible the Degree of Parallelism (DOP). We consider two applications, 2048x2048 matrix multiplication (MxM) and a 512x512 FFT, implemented with increasing DOP. We keep the workload constant, gradually increasing the number of threads in which the workload is divided (reducing the number of operations each thread has to perform). We have used a Fermi GPU with the ECC turned OFF. Figure 4 shows how increasing the DOP influences the FIT, execution time, and MEBF of MxM and FFT. On the left side of Figure 4 we keep constant the number of threads per block (256 for MxM and 512 for FFT), increasing just the number of blocks. On the right side of Figure 4 we keep the number of blocks constant (16K for MxM and 512 for FFT), increasing just the number of threads per block. The number of blocks and threads per block has been chosen to fit the algorithm better and ease its coding.

1) Constant number of threads per block: On the left part of Figure 4 the DOP is increased keeping the number of threads per block constant to 256 for MxM and 512 for FFT. In both algorithms, the execution time and the FIT rate decrease when we increase the number of blocks. This should not surprise, given that even the configuration with the lowest number of blocks fully exploits the GPU hardware, and each block has 256 or 512 threads, ensuring that all 32 CUDA cores in each SM of the Fermi GPU are used. Thus, the DOP is increased without increasing the GPU exposed computing area but actually reducing the SM memory footprint (which influences the error rate as ECC is OFF). In fact, at most 15 blocks can be executed in parallel on the Fermi and, as we increase the DOP, each thread has less operations to perform and smaller memory requirements. As a result, a higher DOP with constant block size implies fewer operations executed and less memory stored in each SM, slightly reducing the FIT rate. Concurrently, the execution time is also reduced as a higher DOP better fits the GPU architecture. This leads to higher MEBF, with the configuration delivering higher performances having a MEBF that is 2 × and 1.6 × higher than the slower one for MxM and FFT, respectively.

2) Constant number of blocks: a different trend is observed when the DOP is increased keeping the number of blocks constant (16K for MxM and 512 for FFT) and increasing the number of threads per block (right figures). For both MxM and FFT, when the number of threads per block is increased, the FIT rate increases, and the execution time decreases. This is because increasing the number of threads per block imposes higher memory and scheduler requirements in each SM (which increases the FIT) and reduces the number of sequential operations of each thread (which reduces the execution time). As shown in Figure 4, using fewer computing cores than the 32 available per SM reduces the FIT rate but jeopardizes performances. Passing from 16 to 64 threads per block for MxM and from 1 to 32 for FFT, in fact, increases the FIT rate of 1.5 × and 3.2 × , respectively, but reduces the execution time to 1/2 and 1/5, respectively. Intuitively, as 32 CUDA cores are available per SM, instantiating less than 32 threads per block underutilizes the GPU hardware reducing the FIT rate but compromising the GPU parallel efficiency. As shown in Figure 4, underutilizing the GPU is not a good reliability solution as the FIT reduction is not sufficient to compensate for the performance degradation.

When the number of threads per block saturates the available cores, as shown in the right side of Figure 4, the FIT rate is basically constant for MxM (number of threads per block from 64 to 256) and, for FFT, the FIT increase is of about 0.5 × (threads per block are from 32 to 512). The (slight) FIT increase is caused by the higher strain in the scheduler (more threads to manage) and the different memory distribution. For FFT, the FIT increases faster as threads interact with each other (butterfly modules), while in MxM, as there is no interaction between threads, it is easier for the scheduler to dispatch parallel threads.

What we can derive from the proposed DOP analysis is that underutilizing the GPU is never a good idea as, even if the FIT rate is reduced, the performances are jeopardized. When the GPU resources are saturated, the FIT rate normally grows slower than the execution time reduction, increasing the MEBF. Particular attention needs to be given to the block size. Having too many threads in a block can lead to memory latencies or scheduler overcharges that reduce the benefit brought by the faster execution. As a general rule, we can state that the number of threads per block should be engineered to fully exploit the GPU parallelism without saturating the SM register/cache.

C. Architectural Solutions

The GPUs architects are making available dedicated hardware resources to improve the efficiency and performances of key operations. This is the case of tensor core and mixedprecision functional units. Lower precision units are smaller, faster, and more efficient. The user can select the precision of the operations to tune the hardware utilization (and the execution time) with the application's needs. A tensor core is a hardware unit that performs a 4×4 FP16 matrix multiplication in one clock cycle.

Figure 5 shows, for the Volta V100 GPU, the relative FIT rate, execution time, MEBF, and MWBF for LavaMD, and Matrix Multiplication (MxM) implemented in FP64, FP32, and FP16. MxM is also executed with ECC ON and, in FP16, using tensor core. We also consider YOLOV3 executed in FP32 and FP16. To have a fair comparison, YOLO was not retrained. We simply cast the FP32 data to FP16 and execute the convolutions. The cast operations reduces of about 23% the benefit, in terms of performances, of executing the CNN in lower precision. In Figure 5, we also show the Mean Work Between Failure (MWBF) to have a fair comparison between the different precisions (an FP64 execution produces 4× more work than an FP16 execution).

From Figure 5, it is evident that both reduced precision and tensor core improve performances. On average, the FP64 implementation has a 3.94 × higher execution time than the FP16 implementation. As mentioned, for YOLO, the execution time of FP32 and PF16 is similar due to data casting. A purely FP16 execution would have a 23% lower execution time. Interestingly, the FIT rate follows the same trend of the execution time, with the FP64 implementation having, on average, a 4 × higher FIT rate compared to the FP16 implementation. Reducing the precision of the operations significantly improves the MEBF, and even if we consider the MWBF (4 × more bits are output for FP64 than for FP16), the reduced precision execution outperforms the other configurations. When ECC is enabled (only shown for MxM), the benefit in the FIT rate is maintained, indicating that the lower error rate is not solely attributed to the lower amount of memory (with ECC ON memory faults are masked) but also to the smaller sizes of the lower precision functional units. Whenever possible, the use of lower precision units is definitely recommended on GPUs.

Tensor core, as shown for GEMM with ECC ON in Figure 5, has a 63% lower FIT rate compared to the software implementation of FP16 matrix multiplication. The tensor core circuit, despite being bigger than an adder and multiplier, is slightly more reliable than the combination of adders, multiplications, and the loop control variables required to implement MxM in software. Tensor core is particularly efficient in improving the GPU performances, providing (in the FP16 implementation) a 2 × higher MEBF (and MWBF) than the FP16 software MxM. Whenever possible, then, the use of a tensor core is highly recommended.

D. Compiler Optimizations

Figure 6 shows the relative execution time, FIT rate, and the MEBF for the Matrix Multiplication (MXM) compiled with four different optimizations, O0, O1, O3+MinRf (O3 restricting the thread register usage to the minimum), and O3. We show the value for each metric relative to the lower value between the tested configurations. That is, the lower value for each metric in each graph is always 1. Results are obtained on a Kepler K20 GPU with ECC ON. We choose to test the compiler optimization with ECC enabled to focus on the computation errors.

The static compiler optimizations, by changing the number and kind of machine instructions of each thread, can impact the reliability of a code running on a GPU. The O0 flag produces the least optimized machine code, including many instructions that could be simplified or reordered. The O0 version has the highest number of instructions, being 290% higher than the O3 version, the O0 version also includes 21× more memory movements instructions than the O3 and 9× more integer instructions than the O3 version. As shown in Figure 6, these inefficient set of instructions increases both the execution time and the error rate, consequently reducing the MEBF.

The execution time decreases significantly when the first optimization level is applied (O1), but the FIT rate remains similar to O0. As discussed in Section II, the FIT rate does not depend on the execution time but on how the code uses the resources. O1 significantly reduces the number and organization of machine instructions compared to O0. In particular, the stalls caused by instruction dependence on the code generated with O1 are still higher than the code generated with O3 (19% higher). O1 only organizes the instructions differently, reaching a similar execution time than O3. As the execution time for O1 is almost half the O0 one, and the FIT remains similar, the MEBF for O1 is ∼ 2x the O0 one.

O3 is the most aggressive optimization that can be applied. We consider the O3 both when the register per thread is limited to 16 (45% fewer registers than only O3 code) and when all registers can be used. The execution time for the two versions are similar (O3+MinRF has a 5% higher execution time), and O3+MinRF has 5.5% more instructions than the O3 version, necessary to perform the register spill to the memory. These additional instructions increase the FIT rate. We recall that, as the ECC is enabled, the extra register footprint of O3 compared to O3+MinRF does not increase the FIT rate.

O3 compiled code provides the optimal resource utilization, with the smallest number of machine instructions to reach the code solution (in O3, only 13% of the stalls are caused by instruction dependency). Consequently, O3 has the lowest FIT rate compared to the other configurations, and the MEBF of O3 can be 6.1× higher than O0. Compiler optimizations, then, can significantly impact the reliability of a code, and better resources utilization can reduce the FIT rate while improving performances. Higher compiler optimizations, then, are definitely to be preferred.

VI. CONCLUSIONS

In this paper, we have evaluated the impact in the final code reliability of the most common optimizations available to GPU developers and architects. We have considered software optimizations, Degree Of Parallelism, compiler optimizations, and architectural solutions to speed up a set of applications executed on GPUs belonging to three different generations. The results we have presented, based on extensive neutron beam experiment campaigns, are highly encouraging, as the faster execution does not necessarily increase the GPU error rate and, even if it does, the execution time is reduced much more than the FIT rate increase. To improve reliability, then, it is not mandatory to add extra hardware or specific hardening solutions. It is actually not even necessary to sacrifice performance as, in most configurations, we can have it all: a faster configuration that can also produce more correct data. In particular, the use of dedicated hardware functional units and compiler optimizations are particularly efficient in improving the reliability-performance trade-ff. Finally, we have only provided relative FIT rates. However, we can mention that the reliability of the tested GPUs have significantly improved through the generations. Overall, the same code on a Kepler GPU has an error rate about one order of magnitude lower than the Fermi, and the same code on the Volta would be almost one order of magnitude more reliable than on the Kepler. The improved reliability is much higher than the one reported for similar technologies and even higher than the reported one for the memory structures of NVIDIA GPUs. If we also consider the improvement in performances and efficiency between the different GPU generations, the number of correctly completed executions are likely to increase of almost 2 orders of magnitudes. The intense research to improve GPUs' reliability has borne fruit and paves a promising path for the future.

Fig. 1 :

 1 Fig.1: Overview of the proposed study. Improving the performance of codes executed on GPUs by (1) optimizing the code, (2) taking advantage of dedicated architectural solutions, (3) increasing the Degree of Parallelism, or (4) using aggressive optimization flags will also impact the error rate (not necessarily for worse, as we will show). We aim to identify the configurations that provide a performance gain higher than the FIT rate increase, thus improving the MEBF.

Fig. 2 :

 2 Fig. 2: Portion of the neutron beam setup at ChipIR.

FITFig. 3 :

 3 Fig. 3: Relative effects of changing the algorithm implementation on the Failure In Time rate, Execution time, and Mean Execution Between Failures for FFT, MxM, NW, and Sort executed on the Kepler K20 with ECC OFF.

Fig. 4 :

 4 Fig.4: Relative effects of varying the Degree of Parallelism on the FIT rate, execution time, and MEBF of MxM and FFT executed on the Fermi GPU with ECC OFF. We first maintain the number of threads per block constant (256 for MxM and 512 for FFT), increasing the number of blocks, and then keep the number of blocks constant (16K for MxM and 512 for FFT), increasing the number of threads per block. The workload is constant for all configurations.

Fig. 5 :

 5 Fig. 5: Relative impact of using mixed-precision or tensor core in the Failure In Time, Execution Time, and Mean Executions or Work Between Failures of the tested Volta V100 GPU. For the FP16 YOLOv3 we have highlighted the time wasted in the casting of data.

Fig. 6 :

 6 Fig. 6: Relative differences of compiler optimizations of Failure In Time, Execution Time, and Mean Executions Between Failures for MxM executed on Kepler GPU with ECC ON.