Performance Evaluation of Independent Component Analysis-Based Fault Detection Using Measurements Corrupted with Noise

The detection of sensor faults has proven to be easier through data-driven methods which rely on historical data collected from sensors that are placed at various locations in a process plant. Since the distribution of industrial process variables is random and non-Gaussian, the independent component analysis (ICA) method has been better suited for fault detection (FD) problems. Whenever data comes with any level of noise, there is difficulty in separating useful information, which hence degrades the monitoring quality of an FD strategy. In this paper, the robustness of FD strategies is assessed for different noise realizations of sensor data using stochastic simulations. The main objective of this work is to demonstrate that ICA-based FD strategies are more robust for different noise levels in comparison with principal component analysis (PCA). The ICA modeling algorithm is improved to avoid random initialization of a de-mixing orthogonal matrix during computation of independent components. Two case studies are considered for evaluating the robustness of FD strategies: a simulated quadruple tank process and a simulated distillation column process. Comparisons have been carried out between ICA, dynamic ICA, modified ICA and PCA strategies for different sensor noise levels. The simulation results reveal that ICA-based FD strategies over-perform PCA FD strategy in monitoring sensor faults for different levels of noise.


Introduction
An active field of research in the process control domain, process monitoring-fault detection (PM-FD) aims to quickly detect sensor faults in chemical process plants (Venkatasubramanian et al. 2003;Ramakrishna Kini and Madakyaru 2020). Until recently, the task of detecting process deviations was dependent on human supervision, which eventually lead to more cost due to the increased complexity in modern process plants. To deal with process deviations or faults, model-based and process history methods have been applied in practice. The model-based PM-FD methods rely on comparing information of measured variables with an explicit known mathematical model of the process. A few commonly used methods include parity relations, expert systems and Kalman filters (Sanjula and Zukui 2019). In contrast, process history PM-FD schemes rely purely on available historical data. Owing to the tremendous advances in the field of computer technology, acquiring data from sensors has become an easy task. Once large data is available, the current status of the process plant can be studied, and this would aid in diagnosing faults.
Process history methods are further classified into two categories: univariate and multivariate techniques. While univariate techniques monitor a single variable at a time, multivariate methods excel in monitoring multiple variables. In recent years, a few of the commonly applied multivariate methods include principal component analysis (PCA), principal component regression (PCR) and partial least squares (PLS) (Yin et al. 2014). PCA, a classic dimensionality reduction method, is one of the oldest multivariate methods used in fault detection problems (Nor et al. 2020;Alauddin et al. 2018). In PCA, large-dimension data is transformed to new dimensional subspace with an aim of retaining the original data information in the transformed variables (Harkat et al. 2020). Despite being the core method in fault detection problems, PCA FD strategy carries with it a few disadvantages. First, modeling in PCA is carried out by considering second-order statistics of mean and variance. Secondly, PCA assumes that the latent variables adhere to Gaussian distribution, which is not often true with industrial processes.
More recently, a PM-FD method based on independent component analysis(ICA) has been proposed for capturing non-Gaussian characteristics of data through higher-order statistical representation. The ICA method, initially developed for blind source separation problems, has progressed very well to find applications in PM-FD problems (Tong et al. 2017;Li and Hongbo 2010). In ICA, large data is decomposed into linear combination of independent components (ICs) that are independent with no orthogonality constraint (Hyvarinen 2013). ICA modeling involves computation of ICs by using an iterative algorithm in which separating the matrix is randomly initialized, and this usually results in different solutions for ICA algorithms (Lee et al. 2006). In this work, the de-mixing orthogonal matrix V is initialized to identity the matrix instead of a random matrix, as this helps in providing a constant solution. Once the ICA model is developed, it is used to monitor new data through three fault indicators, namely, I 2 d , I 2 e and squared prediction error (SPE). Many extensions of conventional ICA technique have been developed in recent years. This includes dynamic ICA (DICA), which considers the dynamics of the data via lagged variables, modified ICA (MICA) (Tong et al. 2017;Yingwei and Yang 2010), where dominant ICs are extracted from data to reduce computational complexity, kernel ICA (KICA) (Tian et al. 2009), which is used to handle both nonlinear and non-Gaussian industrial data, and multi-block ICA where ICA process monitoring is carried out in different blocks (Yingwei and Chi 2012). A multi-way kernel "FastICA" based on feature samples was developed, which was used for diagnosing batch processes, while nonlinear detection using integrated kernel PCA and locality preserving projections (LPPs) was also developed (Lianfang and Xuemin 2014).
Though ICA and its variant methods have been applied successfully in the process monitoring domain, they still face a stiff challenge of extracting information from highdimensional noisy data available from industrial plants. In this work, the robustness of ICA and PCA fault detection strategies is validated by their ability in monitoring process data in the presence of different levels of noise. To test the robustness of an FD strategy towards different noise realizations of sensor data, stochastic simulations have been carried out. The main aim is to show that the ICA-based models are more robust for different levels of noise as compared to conventional PCA method. Stochastic simulations involving 1000 runs will be performed for noise with a signal-to-noise ratio (SNR) of 20, 10 and 5. The SNR is defined as the ratio of the variance of the signal to variance of the noise. SNR = 20 corresponds to quality data with less noise, SNR = 10 corresponds to a medium level of noise, while SNR = 5 indicates a high amount of noise in the data. Comparisons have been carried out between ICA-, DICA-, MICA-and PCA-based fault detection strategies for different noise levels.
The organization of this paper is as follows: Section 2 provides a detailed overview of ICA, DICA, MICA and PCA modeling methods along with the fault indicators. Next, robust fault detection based on ICA involving stochastic simulation with different noise realization is presented through a block diagram representation in Sect. 3. To test the effectiveness of different FD strategies against noise realizations, two case studies have been considered: a simulated quadruple tank process (QTP) and a simulated distillation column (DC) process in Sect. 4. The paper terminates with a conclusion section in Sect. 5.

Independent Component Analysis
ICA is a multivariate technique with the sole objective of extracting hidden independent components (ICs) from non-Gaussian data. Available training data, X a = [x 1 , x 2 . . . , x n ] T ∈ m×n is represented as a linear mixture of p(≤ m) unknown ICs s 1 , s 2 . . . s p and a mixing or a separation matrix. The original data and ICs are related which is described as: where n and m represent number of observations and variables, respectively, and A = [a 1 . . . a p ] T ∈ m× p is a deterministic mixing matrix, S = [s 1 , s 2 . . . s n ] ∈ p×n is the IC matrix and E ∈ m×n is a residual matrix. Since two unknown matrices S and A cannot be computed, the aim is focused in finding a separating matrix W such that the reconstructed matrix is given by: The ICA model development stage consists of four steps, namely, centering, whitening, preprocessing and iterative computing. Initially, the data X a is preprocessed to have a mean of zero. Next, whitening is performed by singular value decomposition (SVD) on the covariance of X a which generates the transformation Z = QX a with Q = −1 U T , where is a matrix involving diagonal eigen values and U is an eigen vector matrix obtained using the covariance of X a . The transformation post whitening stage is: where V is an orthogonal matrix. Next, the preprocessing step involves developing important relationships between original data and the separating matrix W, which is useful for determining the ICs. From Eq.
(3), S can be estimated as follows: Comparing Eqs.
(2) and (4), a relation between W and V is obtained, given by: To avoid random solutions in each iteration, instead of random initialization of V, it is initialized to an identity matrix I ∈ R m×m . Next, each column v i of V is updated one by one through an iterative algorithm. The computation of V requires all ICs to be independent, and criteria for independence is based on non-Gaussianity which is measured through kurtosis and negentropy. Kurtosis, the fourth-order moment, is sensitive to outliers which makes negentropy a better choice. For a variable having density function f(Y), H is represented as: Any difference in entropy between the given distribution and normal distribution having the same mean-variance is computed by negentropy which is mathematically given by: Y guass is a random Gaussian variable with a similar covariance matrix as that of Y. For estimating negentropy, a simpler approximation was proposed (Aapo and Erkki 2000): where Y and y g are assumed to have a mean of zero and unity variance, and G(.) is a non-quadratic function which can be evaluated as (Aapo and Erkki 2000): where c ≈ 1. By maximizing the objective function, extraction of independent components takes place sequentially: Each element of V corresponds to one IC and can be estimated by the use of a deflation scheme. By maximizing the sum of the N-unit cost function, one unit cost function is generalized for determining V. The resulting optimization problem can be described as: This undergoes maximization at certain optima of E{g(v T x a )} under the constraint E{(v T x a ) 2 } = ||v|| 2 = 1 which is described by: and v 0 is the optimum value of v. Solving the above by using the Newton method yields the fixed-point algorithm (Ali and Mahdi 2012): The above algorithm extracts only one independent component at a time. This fixed-point algorithm has to be run m times to estimate the remaining ICs. After each iteration, the vectors have to undergo orthogonalization to make sure that they do not converge to the same maxima. Hence, performing the deflationary orthogonalization through the Gram-Schmidt method: where j is the number of previously estimated vectors. Once the ICs are computed using the fixed-point algorithm, the next task is selection of dominant ICs which is needed for two reasons: robust performance and reduction of complexity. A few methods have been commonly preferred for selection of dominant ICs (Lee et al. 2006). However, the cumulative percentage variance (CPV) method has been applied in the current work to determine optimum ICs. First, the L 2 norm is used to sort the rows of separating W from descending to ascending order, λ = ||W || 2 . Next, CPV method is applied to determine p-dominant ICs which is expressed as: Once a model is developed, determining the prediction capability is crucial for its acceptance. Generally, the prediction capability of a model is reduced when large noise is present in the data. In this work, the R 2 test has been used to evaluate the goodness of the developed PCA as well as the ICA model. The coefficient of determination, denoted as R 2 , is the ratio of the sum of squares regression to the total sum of squares denoted mathematically as: SSQ is the sum of squares which determines how far the predicted value is from the mean of the original data, and TSSQ determines how the original data points vary from their mean. Once an ICA model is developed from data exhibiting normal process operation, the future behavior of variables is compared with the present behavior through the fault detection indicators. The ICA model may be represented as: The ICA model consists of three parts: The systematic space represents a model with respect to p optimum ICs, the excluded space represents a model with respect to ignored m-p ICs, and the residual space. To monitor the three parts, three fault indicators are employed namely, I 2 d to monitor systematic space, I 2 e to monitor excluded space and SPE to monitor residual space. The reconstructed vectorsf newp = W p X newp andf newm−p = W m−p X newm−p are computed for a new data X new where W p is a matrix corresponding to p rows, and W m−p corresponds to excluded m-p rows of the separating matrix. The fault detection indicators of ICA at sample i are represented as (Yuan et al. 2019): Kernel density estimation(KDE) technique is employed for computing the threshold limits for fault indicators (Lee et al. 2004).

Dynamic ICA
In chemical process plants, variables do not remain at a fixed state, as they have a tendency to move around a nominal operating range, and this will result in variables having strong autocorrelation. The conventional process monitoring methods perform modeling assuming that the present sampling instants are independent from the past sampling instants. However, due to the autocorrelation exhibited by variables, present sampling instants have dependency on past sampling instants, and this information can be included in the modeling stage through dynamic ICA strategy. In dynamic ICA strategy, the variables are supplemented with past observations to generate an augmented matrix which is represented by (Lee et al. 2004): where x T t is an observation vector at sampling instant t, and l denotes the number of lagged measurements. The selection of a proper value for time lag l is crucial for capturing dynamic information. A few methods available in the literature for choosing the right number of lags include a subspace identification criterion and the Akaike information criterion. However, it may be noted that selecting l = 1 or l = 2 yields good results for dynamic multivariate techniques (Lee et al. 2004). Once a reference dynamic ICA model is constructed using fault-free data, the fault indicators described in Eqs. (19)-(21) are used for detecting faults in a new data set.

Modified ICA
The modified ICA (MICA) algorithm was proposed as an improvement of the ICA algorithm in order to have structured ordering of the independent components (ICs). In the MICA algorithm, the initial estimate of ICs is done using normalized principal components (PCs) from the PCA model, and then, the FastICA algorithm is utilized to update the retained ICs with same variance being maintained as that of PCs. The starting step of MICA includes applying PCA to extract PCs from data X a = [x 1 , x 2 . . . , x n ] T ∈ m×n .
where T and P are score and loading vector matrices, respectively, generated from covariance of X a . The last few eigenvalues of = diag(λ 1 , λ 2 , . . . λ m ) are close to zero, and they can be left out. The retained eigenvalues or PCs undergo whitening which is expressed as Z = QX a with Q = −1 P T . The whitened components of Z serve as initial estimates for the ICs. The main objective of MICA is to ensure that elements of S become as independent as possible, and the MICA model is defined as: where C∈ m×p such that C T C = D and D = diag(λ 1 , λ 2 , . . . λ p ) considering all the dominant eigenvalues. The variance of each element of S is the same as the variance of the score elements of PCA since the following condition is satisfied: This is followed by normalization of S such that: where E(SS T ) = I. Hence, the objective of MICA is reduced in determining C no . The elements of the C no matrix are computed by a fixed-point algorithm using the negentropy approximation (Lee et al. 2007). Before the fixed-point algorithm is executed, the following initialization is ensured: where I m is a m-dimensional identity matrix and 0 is an m × (p − m) zero matrix. After the computation of C no , W and A matrices are computed using: The next step is to determine dominant ICs. As the variance of each IC in S is the same as that of PCs, any standard criterion for determining optimum PCs in the PCA strategy can be employed. The CPV technique is considered for determining optimum ICs in this study. Two fault indices are defined for the MICA strategy: the T 2 and SPE for modeling systematic and residual parts, respectively. For a new data X new , they are defined as follows: where s = WX new andX = As. If the fault indices exceed a predetermined reference threshold, then a fault is declared (Lee et al. 2006).

Principal Component Analysis
PCA is the most widely applied technique for capturing correlation between the process data variables. After normalization to zero mean and unity variance, the scaled data X ∈ n×m is used to estimate principal components (PCs) through singular value decomposition: where T and P are the PC and loading vector matrix, respectively (Harrou et al. 2013). The loading vector matrix is related with covariance of X, which is described mathematically as: where = diag(λ 1 , λ 2 , . . . λ m ) is a matrix comprising eigenvalues of 6 arranged in diagonal format in descending order. Once optimum PCs p are computed via CPV technique, the PCA model can be represented as: TheX andT represent model parameters for p retained PCs. The residual matrix F represent model parameters for m-p ignored PCs.
The T 2 statistics monitors the modeled part of the developed PCA model. It is represented in the following manner for a new data X new : The matrixˆ contains diagonal eigenvalues of retained PCs. The Q statistics monitors the residual part of the developed PCA model, and it is expressed as: For X new , a fault may be declared if the value of T 2 and Q fault indicators is greater than the threshold (Sanjula and Zukui 2019).

ICA-Based Robust Fault Detection Strategy
Once a multivariate model is developed using ICA, evaluating the model for robustness and prediction capability is crucial for its acceptance. The developed multivariate model should be robust enough to perform in harsh and noisy industrial environments. It has been observed that a multivariate FD strategy which performs well in less noisy environments may not be able to replicate the same performance in a noisier environment. Hence, the presence of noise does have its effect and can degrade the fault handling capacity of an FD strategy. In this paper, we present a multivariate FD strategy based on ICA and present the behavior of the developed ICA model when noise is introduced in the underlying data. The aim of this work is to study the effect of different levels of measurement noise on performance of ICA, DICA and MICA models. After this, similar noise realizations are performed on a PCA model followed by drawing comparisons between PCA-and ICA-based FD methods.

Algorithm 1 : Robust ICA fault detection algorithm
1. For a given training and testing data set, perform scaling to have a mean of zero.
2. Whiten training data to remove cross-correlation and develop relations between Q, V and W to generate W = V T Q.
3. Initialize V to identity matrix I ∈ R m×m and calculate the number of ICs using a fixed-point algorithm. It is unfair to come to a conclusion regarding the robustness of a fault detection scheme based on few rounds of simulations. In the field of system identification and signal processing, when faced with significant uncertainty in a process model involving making a estimation, stochastic simulation has proved to be a good solution. Thus, in this study, stochastic simulations have been carried out with 1000 realizations. The different noise realizations that would be carried out involve introducing noise with a signal-to-noise ratio (SNR) to the data. In the present task, noise with a defined SNR, i.e., SNR = 20, SNR = 10 and SNR = 5, is added to the data by performing stochastic simulation. The data with SNR = 20 corresponds to quality data with less noise, the data with SNR = 10 corresponds to a medium level of noise and SNR = 5 indicates very noisy data. The algorithm of the robust fault detection strategy is described in Algorithm 1. Figure 1 describes a block diagram representation of the proposed strategy. The available process data is initially split into training and testing groups. The training group data is preprocessed to a zero mean which is then used to develop the ICA model. After determining the optimum ICs, KDE technique is employed to compute the threshold for fault indicators. For the testing data group, stochastic simulation of 1000 runs is performed for noise with an SNR that is fixed at 20, 10 and 5. Next, the fault indicators for ICA are computed and compared with the threshold. During the absence of a fault, the value of the fault indicator will be less than the control limit value. A fault is declared when the fault indicator exceeds the control limit value. Similar noise realizations are also performed on DICA, MICA and PCA strategies for having comparison. In process monitoring problems, the performance of an FD strategy is assessed using false alarm rate (FAR) and missed detection rate (MDR). The FAR is defined as the ratio of the total number of false alarms to the total number of data in a normal operating range. The MDR is defined as the ratio of the total number of faults which do not exceed the control limits to the total number of data in the faulty range. A fault detection strategy is deemed to be robust if FAR and MDR values are minimum for different levels of noise.

Case Studies
In this section, two case studies are considered: simulated quadruple tank process and simulated distillation column (DC) process to demonstrate the effectiveness of ICA-based strategies over PCA-FD strategy.

Types of Sensor Faults
A fault is defined as a continuous series of events where a variable deviates from its regular operating range. Consider a process having a large number of variables. A measurement M(t) associated with a process variable may be represented as: where M n (t) is the true value of the measurement M(t), and F(t) is a possible fault at time instant t. The value of F(t) in the equation would be zero when the process variable is operating normally and it takes a non-zero value in the presence of a fault. A bias fault is recognized by sudden jump of a variable from normal operating value to higher or lower abnormal value. In this fault, measurement shifts from its normal value of M n (t) to a new value M n (t) + F. This can be represented in time series format as: where F is a bias term representing the rate of increase (%), and t a is a time instant where the fault appears. The aging or drift fault occurs due to aging of a sensor and can be represented in time series format as:

Simulation Study on a Quadruple Tank Process
In the last few years, the quadruple tank system has been applied in many control as well as process problems since it demonstrates very useful multivariable phenomena between the four tanks. In this process, voltages of a pump v 1 and v 2 In the above equations, a 1 · · · a 4 represent area of the outlet pipes; A 1 · · · A 4 represent the cross-sectional area of the tanks; q 1 and q 2 are valve ratios; k 1 and k 2 denote pump constants; and g is the gravity of the earth. In this case, dynamic simulations are carried out to generate a data set for which standard nominal operating points are considered (Karl 2000). Perturbation of inputs is done at their normal operation range using a pseudo-random binary signal in the range of frequency [0 0.03 ω n ], where ω n = π/T represents Nyquist frequency. The idinput function in the System Identification Toolbox of MATLAB is used to generate 2048 samples of data, which were later split equally into training and testing data. For ICA and modified ICA(MICA) strategies, three optimum ICs are selected, whereas three PCs are selected for the PCA strategy. For the dynamic ICA strategy, a lag of l = 1 is considered, and six optimum ICs are selected for the same. The model performance of ICA is assessed using R 2 values which are 0.975, 0.941 and 0.891 for SNR = 20, SNR = 10 and SNR = 5, respectively. A value of R 2 = 0.975 suggests that 97% of the total sum of squares in testing data is described by the ICA model and only 3% in the residuals.  First, a sensor bias fault with a magnitude of +1 is introduced in the tank 1 measurement at sample 175 of the testing data. The fault indicators of different strategies evaluated through FAR as well as MDR for different noise realizations is presented in Table 1. While PCA-Q statistics has no false alarms, the PCA-T 2 statistics has many more false alarms than the fault indicators of ICA, DICA and MICA methods. The PCA strategy has very high value of MDR with increased noise level in comparison to other strategies. The performance of ICA, MICA and DICA methods is acceptable even though MDR value increases with added levels of noise. For providing clarity to the reader, monitoring plots of PCA, ICA and DICA strategies in monitoring bias fault for the case of SNR = 5 are presented in Figs. 3, 4 and 5, respectively.
In the next case, a fault with a rate of 0.008 resembling a drift sensor fault is introduced in sensor data-3 (i.e tank 3 height measurement) from sampling time instant 150 till the end of testing data. The performance of different fault

Distillation Column Simulation Study
A distillation column (DC) is an energy-consuming unit in any chemical process plant and is used for separating compo-nents from mixtures of components based on the difference in vapor pressure. Proper monitoring of distillation columns is necessary in an industry to avoid any accident and loss of product quality. A schematic of a DC process applied in an industrial setup can be observed in Fig. 6. In this example, the distillation column consists of 32 plates, and 10 resistance temperature detectors (RTDs) are used to monitor temperature at different locations of the column. The flow rates of feed and reflux are used to perturb the distillation column. The Aspen simulator is used for generating distillation column process data. To begin, flow rates of the feed as well as reflux streams are perturbed from their nominal operating ranges (Harrou et al. 2018;Madakyaru et al. 2013). Once the system has reached a steady-state condition, these perturbations are used to generate data. The input variables consist of 10 temperatures corresponding to measurements at various locations of the column along with flow rates of feed and reflux. The variations of the output (i.e., x B and x D ) for changes in the input perturbations are presented in Fig. 7. A total data length of 4096 samples is generated with 14 variables, which are then split equally into training and testing data sets.
The training data is used for developing PCA-and ICAbased models which will then be used for monitoring faults in testing data. For ICA and MICA strategies, 8 and 10 ICs are selected, whereas 8 PCs are selected for the PCA strategy. For the DICA strategy, a lag of l = 1 is considered and 14 ICs are selected for the same. Figure 8 demonstrates a scatter plot of observed values vs. predicted values for the case SNR = 5, and it is observed that predicted values follow observed values perfectly. The goodness of the model in fitting training data is found using R 2 statistics. The value of R 2 for the cases with SNR = 20, SNR = 10 and SNR = 5 is found to be 0.978, 0.955 and 0.9152, respectively. Since R 2 values reduces as more noise is present in the data, it may be concluded that noise has its effect on the prediction capability of the ICA model.
First, the ability of different FD strategies in handling sensor bias fault is being investigated. To simulate this fault, a bias with a magnitude of 20 is introduced in testing data temperature variable 5 at sampling time instant 750. To obtain a valid conclusion from the simulation results, 1000 stochastic simulations of different noise realizations have been carried out, and the results are tabulated in Table 3. From the table, it is observed that despite having no false alarms, the T 2 and Q statistics have a very high MDR which indicates that they are unable to detect bias fault in the presence of noise. In contrast to PCA, the fault indicators of ICA, DICA and MICA strategies have less FAR and MDR values, thus demonstrating enhanced performance in detecting the bias fault even with increase in noise levels.
Next, the capability of different FD strategies in monitoring sensor drift fault is presented in this section. A slow

Conclusion
The aim of this study was to find the effectiveness of ICAbased fault detection strategies against different levels of noise in the measured variables, which was then compared with the PCA-based fault detection strategy. Whenever data comes with any level of noise, there is difficulty in separating useful information, and this leads to poor model prediction. Hence, to check the robustness of the FD strategies, stochastic simulations were performed for noise with SNR = 20, SNR = 10 and SNR = 5. The data with SNR = 20 corresponded to quality data with less noise, the data with SNR = 10 corresponded to a medium level of noise and SNR = 5 indicated very noisy data. The performance of ICA FD, dynamic ICA FD and modified ICA FD strategies was compared with the PCA-FD strategy through two case studies, i.e., the quadruple tank process and the distillation column process. The FAR and MDR metrics were used to assess the monitoring performances of different strategies. It was observed that the value of MDR as well as FAR became increasingly higher with increased level of noise in the data for all the FD strategies. The MDR value for the PCA strategy was very high, and it further increased with increased levels of noise. In contrast, the ICA-based FD strategies had much lower values of MDR for different noise levels. Hence, it can be inferred that ICA, dynamic ICA and modified ICA strategies were more effective to deal with sensor faults in comparison to the PCA strategy for different levels of noise. For the future, the proposed strategy can be extended for considering waveletbased ICA for handling process plants with multi-resolution data.
Funding Open access funding provided by Manipal Academy of Higher Education, Manipal  Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecomm ons.org/licenses/by/4.0/.