#  Deep Scattering Spectrum germaneness to Fault Detection and Diagnosis for
Component-level Prognostics and Health Management (PHM)

Ali Rohan<EMAIL_ADDRESS><EMAIL_ADDRESS>Department of
Mechanical, Robotics and Energy Engineering, Dongguk University, 30 Pildong 1
Gil, Jung gu, 04620, Seoul, Republic of Korea  Faculty of Medicine and Health
Sciences, University of Nottingham, Sutton Bonington, Loughborough, LE12 5RD,
UK

###### Abstract

In fault detection and diagnosis of prognostics and health management (PHM)
systems, most of the methodologies utilizes machine learning (ML) or deep
learning (DL) through which either some features are extracted beforehand (in
case of ML) or filters are used to extract features autonomously (in case of
DL) to perform the critical classification task. Particularly in the fault
detection and diagnosis of industrial robots where electric current, vibration
or acoustic emissions signals are the primary source of information, a feature
domain that can map the signals into their constituent components with
compressed information at different levels can reduce the complexities and
size of typical ML and DL-based frameworks. The Deep Scattering Spectrum (DSS)
is one of the strategies that use the Wavelet Transform (WT) analogy to
separate and extract the information encoded in a signal’s various temporal
and frequency domains. As a result, the focus of this work is on the study of
the DSS’s relevance to fault detection and diagnosis for mechanical components
of industrial robots. We used multiple industrial robots and distinct
mechanical faults to build an approach for classifying the faults using low-
variance features extracted from the input signals. The presented approach was
implemented on the practical test benches and demonstrated satisfactory
performance in fault detection and diagnosis for simple and complex
classification problems with a classification accuracy of 99.7% and 88.1%,
respectively.

###### keywords:

Deep Scattering Spectrum, Wavelet Scattering Network, Prognostics and Health
Management (PHM), Fault Detection and Diagnosis

##  1 Introduction

A significant proportion of tasks in modern industry are performed by devices
such as robots. These robots are composed of a variety of components that work
together to perform a specially designed operation. Due to their continuous
operation, these components are prone to deterioration which most of the time
can lead to fatal damage. Therefore, it is necessary that over a certain
period, proper assessment and maintenance strategies are in place. To help
tackle this, Prognostics and Health Management (PHM) has emerged as an
appealing method for developing methodologies for system health monitoring,
diagnostics, RUL prediction, and prognostics. PHM is regarded as an efficient
method capable of delivering complete, accurate, and customized solutions for
system health management [1]. PHM has three primary functions. The first and
foremost is the early detection of the fault triggers in a system, followed by
the separation and recognition of faults and their respective source, and
finally predicting the Remaining Useful Life (RUL) of a specific component or
a system. Figure 1 depicts the fundamental tasks performed by a PHM system.
PHM can be used at the component, system, or both levels. PHM at the component
level aims to create health monitoring methodologies for electromechanical
components such as electric motors, electronic devices, bearings, gear
reducers, etc. It evaluates if the monitored component’s health has
deteriorated over time because of numerous environmental, operational, and
performance-related parameters. [2, 3]. However, PHM at the system level
evaluates detailed system health, considering system operation, design, and
process-related parameters [4].

Figure 1: Basic tasks involved in a typical PHM system.

The most essential of the three PHM tasks is fault detection and diagnosis.
The sustainability of a PHM system is heavily reliant on it. If the fault is
not adequately diagnosed, the entire PHM system may fail. There has been a
great development in fault detection and diagnosis of a PHM system in recent
years, with both mathematical model-based [5, 6] and data-driven [7, 8]
techniques yielding encouraging results. Data-driven techniques, in
particular, are gaining prominence due to their capacity to adapt and
transform in real-time under diverse conditions. Furthermore, there have been
significant advancements in the computing capability of devices with enhanced
sensor technologies that allow for effective data acquisition. These
techniques are mostly based on machine learning (ML) and deep learning (DL),
in which handcrafted features or deep neural networks with numerous layers of
back-to-back filters are used to learn important information in a given
dataset. To develop a PHM system using current data-driven methodologies,
considerable information and additional real-time data, such as vibration,
acoustic emission, laser displacement, temperature, speed, and electric
current [9, 10, 11, 12, 13] are required.

Previously, several methods have been able to retrieve statistical
characteristics from audio data and categorize them using standard machine
learning algorithms. The Support Vector Machine (SVM) and Artificial Neural
Network (ANN) classifiers were trained using six different statistical
features such as range, mean, standard deviation, kurtosis, skewness, and
crest factor [14]. A method for autonomous mining operations based on
vibration signals was developed in [15]. The authors employed a low pass
Butterworth filter to process and analyze vibration signals before extracting
08 time-domain features, 08 frequency-domain features, and 05 Morlet wavelet
features. A technique using Mel-frequency Cepstrum Coefficients (MFCCs) from
audio signals with SVM for classifications was presented in [16]. The use of
MFCCs features and their comparison with a library of features to determine if
a bike engine is healthy or malfunctioning was proposed in [17]. In contrast,
[18] developed an ML-based fault identification and diagnosis technique for
electric current signals based on a distinct feature selection, extraction,
and infusion procedure.

For DL, [19] created a diagnosis technique for the health status
characterization of the centrifugal pump faults based on sensory data fusion
and a Convolutional Neural Network (CNN). [20] described a similar methodology
for sensor fault diagnosis, in which a computer software was used to generate
defective sensor data, and the signal recognition task was turned into an
image processing task using CWT once the incorrect data was generated. CNN was
used to identify the sensor faults. Recent researches have achieved impressive
results when employing visual representations of audio signals to train
cutting-edge deep learning architectures such as CNNs for the problem of
acoustic emission machine faults [21, 22, 23, 24].

The techniques described above are effective for performing fault
classification tasks for various fault kinds in a variety of settings.
However, the majority of these techniques have some distinct disadvantages.
Traditional handcrafted feature extraction approaches need a specialist
understanding of the type of input data. In general, the input data for PHM
systems consists of vibration or acoustic emission signals [25]. These signals
are non-stationary and might differ from one fault diagnostics application to
another. As a result, generalizing a feature space becomes a highly
complicated issue. In most cases, the features are only applicable to a single
component’s failure. The DL-based algorithms, on the other hand, provide
leverage over the manual and labor-intensive processes of feature extraction.
At various layers, feature extraction is accomplished autonomously using
filters with varying kernel sizes. However, developing the architecture of an
appropriate deep neural network model remains a time-consuming effort. The
manual hit and trial approach is typically used to determine the number of
filters, kernel size, number of layers, and hyper-parameters tuning. In
certain instances, a certain number of layers may be sufficient, albeit it may
be as computationally expensive for another task where fewer layers might
yield as promising results. Especially when it comes to the classification of
faults where input data is composed of non-stationary signals, an intermediate
feature extraction tool such as Fast Fourier Transform (FFT) [26, 27] focusing
on time continuum signals and frequency analysis, Short-time Fourier Transform
(STFT) [28] focusing on time-frequency domain analysis or Wavelet Transform
(WT) [29, 30] is intensively used. These tools are used to generate a map of
time-frequency domain features in a scalogram or spectrogram image, which is
then inputted into a deep neural network like CNN. As a result, the
computational challenges are exacerbated by one step. To overcome the
aforementioned disadvantages, we investigate the use of the Deep Scattering
Spectrum (DSS), especially the scattering transform, for fault detection and
diagnosis of mechanical components of industrial robots in this work. DSS is
specially developed for signal classification. It builds a deep network using
the WT idea by performing a scattering transform with fixed settings of the
dilated filter. It combines the strength of traditional signal processing
methods with the depth of a deep neural network. According to our
understanding, scattering transforms in the literature are largely focused on
audio applications, but a generic scattering representation for classification
that applies to numerous signal modalities other than audio remains
understudied [31, 32]. For component-level PHM, we propose the application of
the scattering transform to electric current signals, i.e. Motor Current
Signature Analysis (MCSA), rather than audio, acoustic emission, or vibration
data. MCSA outperforms vibration and acoustic emission analysis in various
ways. MCSA uses the inherent current signal of the motor control unit, needing
no extra sensors, that reduces cost and system complexities. Furthermore, the
current signals are distinct and not easily influenced by the adjacent working
conditions. Lastly, earlier solutions have used balanced datasets that are not
easily available in the context of industrial robots. Expanding this work, we
implement DSS with simple classifiers for an Imbalanced Scarce and Multi-
domain (ISMD) dataset.

The details of this study are presented in the following sections. Section 2
defines the materials and methods, including the experimental test bench, and
descriptions of the suggested technique. Section 3 consists of the results and
discussions, and Section 4 presents the conclusion.

##  2 Materials and Methods

Figure 2: Basic overview of the workflow.

Figure 2 depicts the fundamental process of the approach used in this work.
Multiple industrial robot-based experimental setups were utilized to obtain
real-time data by inducing mechanical component faults. The acquired data were
pre-processed to remove ambient noise, and signal segmentation was used to
divide the raw recorded signal into successive segments based on the robots’
motion patterns. The segmented signals were then processed with the scattering
transform to extract the scattering and scalogram coefficients, resulting in a
feature vector that was utilized to identify the faults. The specifics of the
procedures are provided in the following subsections.

###  2.1 Experimental Setups

We concentrated on the practical implementation of the fault detection and
diagnosis system in this study, taking into account the experimental setups
depicted in Figures 3 and 4. We validated the applicability of the given
methodology for diverse industrial equipment settings using experimental
setups for two distinct types of robots. The first experimental setup is
comprised of three major parts: an industrial robot, a Programmable Logic
Controller (PLC), and a command generating module. Robostar Co. manufactured
the robot utilized in this setup, and the model number is R004. The second
experimental setup, like the first, is made up of three major components: an
industrial robot, a controller, and a Personal Computer (PC). Hyundai Robotics
Co. manufactured the robot utilized in this setup, and the model number is
YS080. It can carry a maximum payload of 80kgf. This robot is much larger than
Robostar. The robots in both experimental setups have six axes or joints, each
of which is designed with a different type of electric motor, allowing robots
to move freely around each axis. To increase or decrease rotation speed, the
motors are connected to reducers at each axis. The Hyundai robot is controlled
by sending commands to the controller through a PC, but the Robostar employs a
manual command generating module that delivers operation commands to the PLC,
which then controls electric motors to generate a certain motion. On each axis
of each robot, three-phase servo motors are employed. The motors’ power is set
dependent on the amount of mechanical load on each axis. Figure 5 shows the
details of the robots, which are as follow: (a) Robostar R004, and (b) the
Hyundai Robot YS080.

Figure 3: Experimental Setup for Robostar R004. Figure 4: Experimental Setup
for Hyundai Robot YS080. Figure 5: (a) Robostar R004, and (b) Hyundai Robot
YS080.

###  2.2 Data Acquisition

We recorded single-phase electric current signals at each axis of the robot
for the Robostar experimental setup, with each of the six-axes motors has
current sensors installed. Hall Effect Base Linear Current Sensors WCS6800
were employed as the current sensor. Data were collected for two classes:
normal and faulty. The fault replicated in this experimental setup was related
to the strain wave gear reducer, which is a sort of mechanical gear system
that employs a pliable spline with external teeth that is distorted by a
spinning elliptical plug to connect with the internal gear teeth of an outer
spline. Because of its compactness, lightweight feature, high gear ratio, and
high torque capabilities, it is widely employed in robotic systems. The fault
was mimicked by deforming the teeth of the internal gears in the gear reducer
of the third axis motor. The strain wave gear reducer is shown in Figure 6 for
both the normal and faulty  situations.

Figure 6: Strain wave gear reducer for Robostar Axis 3: normal and faulty.

In the Hyundai robot experimental setup, motor current data were collected
using current sensors for each of the three phases of the electric motor. The
current sensors are mounted on each phase of the electric motors, resulting in
a total of 18 current sensors to acquire 18 current signals for six electric
motors. NI DAQ 9230 devices are used to record the current signals for each
axis motor. This data acquisition device transmits the collected data to a PC
executing LabView. The acquired signals are processed and a comprehensive
database is created, containing the sensed data for each axis. Data are
collected concurrently for each motor in various faulty conditions. In the
first case, a Rotate Vector (RV) reducer eccentricity bearing fault was
introduced into the reducer that was linked to the 4th axis motor. The fault
was introduced in the second case by substituting the RV reducer with an
aged/damaged one. The data were collected for three different classes: normal,
faulty (RV reducer eccentric bearing fault), and faulty aged (RV reducer aging
fault). The faulty  and faulty aged RV reducers can be seen in Figure 7. Both
robots were designed to work in all orientations along every axis of rotation
for several cycles. One cycle denotes the conclusion of a single range of
motion across a single axis. For each axis, data were collected for ten
cycles. Figure 8 illustrates the fundamental block diagram of the data
collection procedure for a single-axis motor. It should be noted that the same
procedure is performed for Robostar, except only electric current signals for
a single-phase are collected to study the influence of different phases on
fault classification.

Figure 7: RV reducer gear for Hyundai Robot: faulty and faulty aged. Figure 8:
Basic overview of the data acquisition process.

Subsequently, the motors were executed at various speeds ranging from 10% to
100% of their rated speed to see how the speed variation affected the PHM
system. It is also done to produce an imbalanced and multi-domain dataset.
Figures 9 and 10 display the equipment utilized in the data collection
procedure for the Robostar and Hyundai robots. The electric current are
collected for each axis motor, even if the fault is only put into the gear
reducer of a single axis because a fault in just one axis may impact the
performance and effectiveness of some other axes motors owing to mechanical
linkage.

Figure 9: Data Acquisition system for Hyundai Robot. Figure 10: Data
Acquisition system for Robostar.

###  2.3 Data Pre-Processing

Unlike mathematical-model-based approaches, data-driven approaches are heavily
impacted by how data is fed to a feature extractor, whether manual or
automatic feature extraction. The effectivenss of a machine learning (ML) or
deep learning (DL) method may be greatly enhanced if the data are pre-
processed on the fly to eliminate superfluous information, redundancy, and
undesirable features. Data preparation may not be essential in the case of
image data since the image has a compact two-dimensional structure rather than
a one-dimensional signal space. Because of the participation of multiple
practical components existing in a single location, the signals, whether audio
or electric current, tend to be noisy. The recorded signal data is extremely
susceptible to ambient noise, redundancy, and superfluous information produced
by the recording process, especially in the case of industrial robots with
many electromechanical components conjoined to perform a particular task.
Furthermore, the information must be organized in such a way that the
computational complexity is maintained to a bare minimum at all times. We
conducted various operations to eliminate ambient noise and superfluous
information from the raw recorded electric current signals for this purpose.
Unlike acoustic emission methods, which are commonly used to detect mechanical
component faults, the Motor Current Signature Analysis (MCSA) provides
leverage in safeguarding the original form of the signal. When compared to
audio signals, electric current signals are less susceptible to noise.

Figure 11: An example of the raw recorded current signal for Robostar at 10%
speed of rotation. Figure 12: General block architecture of the segmentation
method.

Figure 11 shows a raw recorded electric current signal for Robostar when it is
continuously operating along the faulty strain wave gear reducer axis. The
fault is mimicked on axis 3 in the case of Robostar. The electric current
signals were recorded for 110 seconds at a sampling rate of 2048 Hz at various
speed ranges, from 10% to 100%. The signal shown in Figure 11 is for a 10%
speed scenario and class: normal. The raw input signal was pre-processed by
first eliminating the ambient noise with Savitzky-Golay filtering. Digital
smoothing polynomial filters or least-squares smoothing filters are used in
this case. Savitzky-Golay filters outperform typical averaging FIR filters,
which tend to remove high-frequency components with noise and are optimum in
the notion that they reduce the least-squares error in matching a polynomial
to frames of the noisy dataset [33, 34]. Following the elimination of silent
sections from the signal, the resulting signal is divided into cycles of
rotations. The rotation cycle is the movement of the robot in a direction
along a single axis. A cycle is defined as the two clockwise movements: from
origin to a maximum range of motion, and two anticlockwise movements back to
the origin of the robot along one axis. The purpose of signal segmentation is
to reduce the computational complexity of an ML or DL model. We used the idea
of image dilation [35, 36] rather than the traditional envelope detection
approach to detect the envelope of the recorded signal since the traditional
envelope detection technique failed to detect the changing envelope of the
recorded signals in real time. Using image-dilation for envelope detection is
beneficial as it successfully produces a spectrum of the boundary around the
peaks of the signal regardless of the change in shape or sampling frequency.
Figure 11 identifies the detected envelope in red color. For signal
segmentation, we defined a threshold value (shown in Figure 11 in green
color), in this case, the median value, and compared it to the signal envelope
to find the points of actual recurrence of the cycle. Based on these points,
the cycles were segmented and placed in a concatenated database. Figure 12
shows the general block architecture of the segmentation method, and Figure 13
features an example of segmented cycles for the signal shown in Figure 11.
Note that the original recorded current signal contained 10 cycles whereas
Figure 13 just shows 5 cycles as an example. The same procedure was adopted
for the case of the Hyundai robot for data pre-processing for the classes of
normal, faulty, and faulty aged.

Figure 13: An example of segmented cycles for the signal shown in Figure 11.

###  2.4 Deep Scattering Spectrum

Figure 14: Wavelet ScatteringTransform with two levels architecture for 1D
signal.

Deep Scattering Spectrum (DSS) has a similar architecture to a Deep Neural
Network (DNN) where different levels of numerous filters are utilized to
extract features from input data. Unlike the DNN, which uses a customizable
dilation method to extract information from filters, DSS uses a fixed set of
dilated filters. These filters are pre-defined, adhering to the Wavelet
Transform’s (WT) characteristics. The Scattering Transform (ST) or Wavelet
Scattering Transform (WST) is the DSS’s core component. The WST is based on
latent wavelets, providing consistent and stable informative signals from the
input data set. Distortion tolerant, WST preserves class discriminability,
making it suitable for classification. To acquire the signal characteristics,
a sequential multi-wavelet decomposition based on modular arithmetic and local
averaging is used. Multi-scale complex wavelets are used to acquire the low-
frequency characteristics of the signals. High-frequency coefficients, on the
other hand, are calculated, while relatively stable frequency characteristics
are obtained by local averaging. The second high-frequency convolution wavelet
is used to recover high-frequency information that has been lost as a result
of the local averaging required to ensure the robustness of higher frequency
terms. WST works with data in phases, which means the output of one phase is
used as the input for the next. Each phase consists of three operations: (1)
convolution, (2) non-linearity (modulus), and (3) averaging (scaling
function), which is comparable to the Convolutional Neural Network approach
(CNN). Similar procedures are done at each layer of the network in CNN,
however, they are referred to as (1) Convolution, (2) ReLU (Rectification
Linear Unit), and (3) Max Pooling. WST has a clear advantage over CNN because
the former employs a fixed set of dilated filters or wavelets, and all
features are typically compacted in the form of scattering coefficients within
the network’s two levels or layers, owing to the wavelet filters’ high
effectiveness in decomposing frequency domain information into substituent
levels. Through the different scattering pathways, the scattering coefficients
balance the data invariance and discriminations. WST resembles physiological
models of the cochlea and the auditory system [37] [38] which are also
employed in audio processing [39]. Furthermore, the DSS architecture is
typically designed based on the properties of the input signal, such as the
length of the signal, the scaling or averaging function, which can also be
determined by determining the frequency spectrum of the signal, and the number
of wavelet filters, which typically range from 1 filter per octave to a
maximum of 32 or 64 filters per octave. The wavelet decomposition and WST
provide an obvious advantage over conventional DNN for any application of
signals where the essential information is in the frequency domain. As
mentioned previously in the introduction section, WST in the literature is
largely focused on audio applications, but a generic scattering representation
for classification that applies to numerous signal modalities other than audio
is still understudied, especially the applications where the signal type is
electric current. In contrast to audio, electric current signal contains a
repetitive pattern in both the time and frequency domains. The information in
the signal is highly correlated in the frequency spectrum, hence, making it an
ideal signal type for the application of WST.

Figure 14 is a visual representation of the WST procedure of this study,
including the architecture of the network for two levels for the 1D signal.
Where $m$ represents the first-order, second-order, and third-order scattering
transforms, $f$ represents the input 1D signal, $\ast$ represents the
convolutional operator, $\phi(t)$ represents the low pass filter,
$J\;$represents the scale, $\psi(t)$ represents the wavelet function, and
$\Lambda\;$the family of wavelet indices. Mathematically, if$f\left(t\right)$
is the signal to be considered for WST, the wavelet function $\psi$ and the
low-pass filter $\phi$ are designed to create filters that cover all of the
frequencies in a signal. The local invariant feature $f$ is generated by the
convolution functions$\;S_{0}f\left(t\right)=f\phi_{J}(t)$, and the high
frequencies can be calculated as:

$\displaystyle\begin{array}[]{@{}l}|W_{1}|f=\\{S_{0}f(t),|f\ast\psi_{j1}(t)|\\}_{j_{1}\in\wedge_{1}}\end{array}$
(1)

By averaging the wavelet coefficient, the first-order scattering coefficient
can be calculated as:

$\displaystyle\begin{array}[]{@{}l}S_{1}f(t)=\\{|f\ast\psi_{j1}|\ast\phi_{j}(t)\\}_{j_{1}\in\wedge_{1}}\end{array}$
(2)

where $S_{1}f\left(t\right)$ can be considered as a low-frequency component of
$|f\ast\psi_{j1}|$, however, the high-frequency component can be extracted as:

$\displaystyle\begin{array}[]{@{}l}|W_{2}||f\ast\psi_{j1}|=\\{S_{1}f(t),\;\left||f\ast\psi_{j1}|\ast\psi_{j_{2}}(t)\right|\\}_{j_{2}\in\wedge_{2}}\end{array}$
(3)

The second-order scattering coefficient can be obtained as:

$\displaystyle\begin{array}[]{@{}l}S_{2}f(t)=\\{\left||f\ast\psi_{j1}|\ast\psi_{j_{2}}\right|\ast\phi_{j}(t)\\}_{j_{i}\in\wedge_{i}},\;i=1,2\end{array}$
(4)

The wavelet modulus convolutions can be obtained by iterating the above
process:

$\displaystyle\begin{array}[]{@{}l}\begin{array}[]{l}U_{n}f(t)=\\{\left|\left||f\ast\psi_{j1}|\ast.....\right|\ast\psi_{j_{n}}\right|\\}_{j_{i}\in\wedge_{i}},\\\
\;\;\;\;\;\;\;\;i=1,2,...,n.\end{array}\end{array}$ (5)

The nthorder scattering coefficient $(S_{n}f\left(t\right))$can be obtained by
averaging all the wavelet modulus convolutional
coefficients${\;(U}_{n}f\left(t\right))$ as:

$\displaystyle\begin{array}[]{@{}l}\begin{array}[]{l}S_{n}f(t)=\\{\left|\left||f\ast\psi_{j1}|\ast.....\right|\ast\psi_{j_{n}}\right|\ast\phi_{j}(t)\\}_{j_{i}\in\wedge_{i}},\\\
\;\;\;\;\;\;\;\;i=1,2,...,n.\end{array}\end{array}$ (6)

The final scattering matrix can be represented as:

$\displaystyle\begin{array}[]{@{}l}Sf(t)=\\{S_{m}f(t)\\}_{0\leq n\leq
l}\end{array}$ (7)

In Equation (7), $l\;$denotes the maximal decomposition order, and the final
scattering matrix is calculated using features collected from each level of
decomposition along several paths. The extracted features, in the form of a
feature vector, are then utilized by a classifier to categorize the classes.
In this paper, we used various classifiers to validate the efficacy of the
WST’s feature domain. The section that follows contains information on the
classifiers that were employed as well as the classification results.

##  3 Results and Discussion

Figures 15 and 16 illustrate the electric current signal cycles for the
Robostar and Hyundai robot experimental setups, respectively. Because each
robot has its payload and operational parameters, the pattern of the current
signals varies. Since the Robostar is a smaller robot, the electric motors
used are of smaller specifications than those used in the Hyundai robot. As a
result, the current amplitude during complete operation along a single axis is
significantly lower than in the Hyundai robot. We present results and
experiments from two types of problems. The first problem is for two classes:
normal and faulty for the Robostar, whereas the second problem is for a
complex classification problem created by the inclusion of the faulty aged
class among the normal and faulty classes of the Hyundai robot due to the
relevance of the faulty aged class with the faulty one. This is because, in
the faulty aged class, the reducer gear is merely replaced with an old worn-
out reducer. Secondly, the Hyundai robot’s size, capacity, and component
structure render classification a considerably more challenging task. This can
also be seen in Figure 16, where the electric current signals for a single
rotational cycle tend to follow a fairly similar pattern across classes. The
difference may still be seen, although not as clearly as in the case of the
Robostar, where there is a noticeable difference in the signal amplitude and
frequency domain.

Figure 15: Electric current signal for one cycle of rotation for Robostar.
Figure 16: Electric current signal for one cycle of rotation for Hyundai
robot.

These current signals were recorded at a sampling frequency of 2048Hz and
down-sampled to produce a comparable signal length to compare the outcomes for
the same wavelet scattering network. The scattering network employed is
comprised of two layers. For the scattering filter banks, the first layer had
a quality factor of 8 while the second layer had a quality factor of 1. The
number of wavelet filters per octave is the quality factor for each filter
bank. The wavelet transform factorizes the scales by deploying the provided
number of wavelet filters. The invariance scale, another important parameter
in the architectural design of a wavelet scattering network, was computed
using the sampling frequency $f_{s}$ and the lenght of the input signal $N$.
It is given in the Equation (8).

Figure 17: Scaling function and the coarsest-scale wavelet plot for the
wavelet scattering network.

$\displaystyle\begin{array}[]{@{}l}Invariance\;Scale\;=\;\displaystyle\frac{\displaystyle
N/f_{s}}{2}\end{array}$ (8)

Figure 18: Wavelet scattering network’s filter banks at each layer.

When the invariance scales are specified, the network becomes invariant to
translations up to the invariance scale. The extent of the invariance in time
or space is determined by the scaling function’s support. The scaling function
and the coarsest-scale wavelet plot for the wavelet scattering network are
shown in Figure 17. The invariance scale is computed using equation 8, and
based on signals temporal domain the scale of 1.6 seconds was chosen. The
distances of the central frequencies of the wavelets in the filter banks are
likewise affected by the invariance scale. However, in a scattering network, a
wavelet’s temporal support cannot exceed beyond the invariance scale. The
coarsest-scale wavelet plot exemplifies this feature. Frequencies less than
the invariant scale are linearly separated, with the scale maintained
constant, such that the invariant’s size is not surpassed. Figure 18 depicts
the wavelet scattering network’s filter banks at each layer. The first layer
features eight filters per octave, whereas the second uses a single filter per
octave. These filters are spread over the sampling frequency range to convolve
and extract the features for the whole frequency spectrum of the signal.

Figure 19 illustrates the network path design based on the invariance scale,
signal input length, and sampling frequency. As previously stated, the network
contains two filter banks. The scattering paths are comprised of wavelets from
both the first and second filter banks. The wavelet number and filter bank
level of each wavelet filter on at least one path are indicated on the
associated node. The network architecture consists of 329 paths, each
containing 53 wavelets in the first layer and 8 wavelets in the second layer.
The features at each layer are extracted using Equation 7, and a final feature
vector consisting of 336 data features for each dataset for the Robostar and
Hyundai robot is generated. These feature vectors were then used to classify
faults using various sorts of classifiers.

Figure 19: The wavelet scattering network’s paths based on the invariance
scale, signal input length, and sampling frequency.

We present classification task results for the six best and most popular
classification algorithms. These algorithms /classifiers are Support Vector
Machine (SVM), K-Nearest Neighbors (KNN), Decision Tree (DT), Ensemble
Learning, Nave Bayes (NB), and Discriminant Analysis (DA). The following
metrics were used to assess the performance of these classifiers. Some of
these metrics are commonly used to evaluate the performance of an ML/DL model
that has been trained. We utilized Accuracy, Sensitivity, Specificity,
Precision, and F-score between these metrics. The most essential metric is
accuracy, which indicates how many samples are properly categorized out of all
the samples. It is commonly stated as the proportion of True Positives (TPs)
to True Negatives (TNs) divided by the number of TPs, TNs, False Positives
(FPs), and False Negatives (FNs) (FNs). A TP or TN is a data sample correctly
classified as true or false by the algorithm. An FP or FN, on the other hand,
is a data sample that the algorithm incorrectly classifies. This metric is
represented by Equation (9):

$\displaystyle\begin{array}[]{@{}l}Accuracy\;=\;\displaystyle\frac{TP+TN}{TP+TN+FP+FN}\end{array}$
(9)

The recall metric is another name for the sensitivity metric. It is defined as
the number of correctly categorized positive samples, implying how many
samples of the positive classifications are correctly identified. It is given
by Equation (10):

$\displaystyle\begin{array}[]{@{}l}Sensitivity\;/\;Recall\;=\;\displaystyle\frac{TP}{TP+FN}\end{array}$
(10)

The specificity of the given class pertains to the forecast of the possibility
of a negative label turning true. It can be represented as in Equation (11):

$\displaystyle\begin{array}[]{@{}l}Specificity\;=\;\displaystyle\frac{TN}{TN+FP}\end{array}$
(11)

Equation 12 gives precision as an amount of TP divided by the number of TPs
plus the number of FPs. This metric is all about regularity. In other words,
it assesses the algorithm’s prediction accuracy. It determines how exactly a
model is based on what is predicted to be positive, and how many of them are
truly positive:

$\displaystyle\begin{array}[]{@{}l}Precision\;=\;\displaystyle\frac{TP}{TP+FP}\end{array}$
(12)

Finally, Equation (13), which is described as the relative average of
precision and recall, provides the F-score. It is dependent on favorable class
evaluation. This parameter’s high value indicates that the model works better
in the positive class:

$\displaystyle\begin{array}[]{@{}l}F-score\;=\;2\times(\displaystyle\frac{Precision\times
Recall}{Precision+Recall})\end{array}$ (13)

Table 1 shows the results obtained using the six classification methods stated
above, depending on the performance metrics of the Robostar experimental
setup. It delivers that the overall performance metrics for the Robostar fault
detection show extremely promising results and the specified method works well
for fault classification.

Table 1: Classifiers performance score for Robostar fault detection.  |  |  |  |  |  |
---|---|---|---|---|---|---
Robostar Fault Detection |  Metrics
Classifier |  Number of Classes |  Accuracy (%) |  Sensitivity (%) |  Specificity (%) |  Precision (%) |  F-score (%)
SVM |  2 |  99.9 |  99.806 |  99.943 |  99.829 |  99.796
Decision Tree |  99 |  98.906 |  99.043 |  98.929 |  98.896
Ensemble |  99.6 |  99.506 |  99.643 |  99.529 |  99.496
KNN |  99.4 |  99.306 |  99.443 |  99.329 |  99.296
Discriminant Analysis |  97.7 |  97.606 |  97.743 |  97.629 |  97.596
Naïve Bayes |  82 |  81.906 |  82.043 |  81.929 |  81.896
Average Performance Score |  96.26 |  96.17 |  96.30 |  96.19 |  96.16

It’s important to note that there are two classes in this situation: normal
and faulty. The SVM classifier yielded the best accuracy. The NB had the
lowest accuracy score. The SVM’s overall performance in distinguishing between
positive and negative observations was likewise quite good. All of the
findings given in the paper were achieved using 5-fold cross-validation to
avoid overfitting in the classification model. The dataset is split into 5
distinct sets of disjoints, each comprising observations from each of the
classes, in 5-fold cross-validation. During the training process, one set is
kept aside for testing while the others are utilized for training. This
procedure is repeated for all folds, and the average accuracy is computed
based on the testing results. To validate the results, we performed the
training of the models for the repetition of 10 with a different range of
random numbers. In the end, the average performance score was calculated by
taking the average of the performance metrics for each classifier. The score
shows that regardless of the type of the classifier, on an average 96.26%
accuracy, 96.17% sensitivity, 96.30% specificity, 96.19% precision, and 96.16%
F-score can be achieved using the wavelet scattering network. Further, the
parallel coordinate plot of the 10 most prominent features are presented in
Figures 20 and 21 highlighting the two-class confusion matrix for the highest
scoring classifier, i.e. SVM.

Figure 20: Parallel coordinate plot for the highest scoring classifier (SVM)
for Robostar fault classification. Figure 21: Confusion matrix for the highest
scoring classifier (SVM) for Robostar fault classification.

Figure 20 depicts the connection between the retrieved features and the
predictors for class separation. The orange lines represent the normal class,
whereas the blue lines represent the faulty class. The wrongly categorized
classes are shown by the light orange and blue lines. The standard deviation
among the features is also presented along the y-axis of the plot to notice
the difference in the data features. It can be seen that the classifier
performs better when the standard deviation between the feature values of the
classes is considerably different. Only one of the 800 samples/observations
was incorrectly classified as normal  when it was actually faulty. This is
also supported by the model’s confusion matrix (shown in Figure 21).

Table 2: Classifiers performance score for Hundai robot fault detection.  |  |  |  |  |  |
---|---|---|---|---|---|---
Hyundai Robot Fault Detection
Classifier |  Number of classes |  Accuracy (%) |  Sensitivity (%)

|  Specificity (%) |  Precision (%) |  F1 score (%)
3 |  88.1 |  88.006 |  88.143 |  88.029 |  87.996
Discriminant Analysis |  85.6 |  85.506 |  85.643 |  85.529 |  85.496
SVM |  83.2 |  83.106 |  83.243 |  83.129 |  83.096
KNN |  80.3 |  80.206 |  80.343 |  80.229 |  80.196
Decision Tree |  68.3 |  68.206 |  68.343 |  68.229 |  68.196
Naive Bayes |  48.9 |  48.806 |  48.943 |  48.829 |  48.796
Average Performance Score |  75.733 |  75.639 |  75.776 |  75.662 |  75.629
Figure 22: Parallel coordinate plot for the highest scoring classifier
(Ensemble learning) for Hyundai robot fault classification. Figure 23:
Confusion matrix for the highest scoring classifier (Ensemble learning) for
Hyundai robot fault classification.

Table 2 shows the results obtained using the six comparable classification
methods based on the performance metrics for the Hyundai robot. In this case,
three classes are used: normal, faulty, and faulty aged. The Ensemble learning
classification method achieved the highest accuracy, followed by the DA and
SVM. With 88.1% accuracy, 88.006% sensitivity, 88.14% specificity, 88.029%
precision, and 87.99 % F-score, Ensemble learning surpassed the DT and NB in
the classification of faults associated with the Rotate Vector (RV) reducer
for the Hyundai robot. The results demonstrate that selecting a suitable
classification method based on criteria such as input data type, number of
classes, data quality, and data quantity is a critical step in any
classification task. There has been no particular research on the selection of
classification algorithms based on the aforementioned characteristics up to
this point. As a result, to obtain greater classification accuracy,
alternative sets of classifiers must be evaluated. It should also be noted
that the NB and DT had the lowest accuracy in both the Robostar and Hyundai
robot fault classification procedures. One of the reasons for this is because
these two classification algorithms are built on a distinct set of ideas in
contrast to the SVM, KNN, and DA. In terms of class clustering, SVM and KNN
are virtually identical. Furthermore, as compared to the Robostar fault
classification, the classification accuracy for the three-class classification
task was relatively poor. This is because the inclusion of a third class in
the classification problem caused substantial confusion among the normal,
faulty, and faulty aged classes. This can be seen in Figures 22 and 23, which
exhibit a parallel coordinate plot of the ten most prominent features as well
as a confusion matrix. The temporal and frequency domain differences in the
Hyundai robot’s current signal for different classes were likewise not
significant enough. A similar trend can be seen in the wavelet scattering
network-extracted data features. There are several misclassified classes,
particularly between normal and faulty, followed by the normal  and faulty
aged combination. In terms of computational cost, the results are still
considerably more practical than other approaches that employ an extra 2D
Scalogram images-based signal representation with a very deep layered network.
This was verified by using a similar dataset for the Hyundai robot with a
Scalogram image representation plus a CNN model. The highest accuracy achieved
was no more than 80%. The most significant advantage of the wavelet scattering
network is its computational simplicity, which allows for faster and shorter
response time in real-time applications, particularly in applications like
this one. This approach is non-expansive if the wavelets and network
architecture are chosen correctly. Energy dissipates as it passes down through
the network. The energy of the $m^{th}$-order scattering coefficients rapidly
coheres to 0 as the rank of the network’s levels/layers rises [40]. Energy
dissipation has a useful value. It reduces the number of wavelet filter banks
in the network compared to a standard deep neural network while retaining as
much signal energy as possible. The results demonstrate that by adding the
Deep Scattering Spectrum into a PHM application with real-time implementation,
domain knowledge, and 1D signals-based data majored environment, the computing
complexity of the PHM systems can be decreased significantly. Furthermore, in
our future study, we plan to investigate the plausible solution and
enhancement of the classification capacity for complex problems such as the
Hyundai robot to develop a strategy that utilizes the wavelet scattering
network’s fundamental advantages.

##  4 Conclusion

The use of Deep Scattering Spectrum (DSS) for 1D signal domains for the
application of PHM systems, rather than traditional techniques such as
scalogram conversion and image-based deep neural network, has numerous
advantages. Signals, being another type of information with its own set of
properties, necessitate the use of techniques capable of extracting meaningful
information from its non-stationary frequency-domain complicated space.
Wavelets and other transformation techniques have been created and utilized
for this purpose throughout the decades. The DSS blends the power of these
signal processing techniques with the deep learning paradigm. Surprisingly,
the techniques involved are identical to those used in deep neural networks,
but with the added benefit of more accurate information extraction. This
research focuses on the implementation of the core concept of DSS, which is a
wavelet scattering network for fault classification of the mechanical
components of various industrial robots. The utilization of an electric
current signal for such fault classification becomes a difficult task when an
indirect electromechanical relationship of electric motors with their
constituent mechanical components is exploited as the foundation for fault
pattern recognition. Regardless of the intricacy, the results demonstrate that
the DSS has enormous potential in addressing fault classification issues under
various circumstances. Furthermore, the future study would concentrate on the
DSS’s enhancements, limits, and disadvantages for PHM applications, since we
intend to integrate more faults linked to other mechanical components with a
generalized feature space to construct a system-level real-time PHM system.

## Acknowledgements

This research was financially supported by the Ministry of Trade, Industry,
and Energy (MOTIE) and the Korea Institute for Advancement of Technology
(KIAT), through the International Cooperative R&D program (Project No.
P059500003).The author is grateful for the support provided by the team of SMD
Lab, Dongguk University, South Korea including Heung-soo Kim, Hyewon Lee, and
Izaz Rauf in helping with communication with the industry partner and
providing the resources and funding for the project.

Conflicts of Interest: The authors declare that they have no conflict of
interest.

## References

  * [1] J. Lee, F. Wu, W. Zhao, M. Ghaffari, L. Liao, D. Siegel, Prognostics and health management design for rotary machinery systems—Reviews, methodology and applications, Mechanical Systems and Signal Processing 42 (1-2) (2014) 314–334.  doi:10.1016/j.ymssp.2013.06.004.
URL https://dx.doi.org/10.1016/j.ymssp.2013.06.004

  * [2] P. Lall, R. Lowe, K. Goebel, Prognostics and health monitoring of electronic systems, 12th Intl. Conf. on Thermal, Mechanical & Multi-PhysicsSimulation and Experiments in Microelectronics and Microsystems (2011).
  * [3] C. Bittencourt (2014).  [link].
URL https://doi.org/10.3384/diss.diva-109335

  * [4] B. Abichou, A. Voisin, B. Iung, Bottom-up capacities inference for health indicator fusion within multi-level industrial systems, 2012 IEEE Conference on Prognostics and Health Management (2012) 1–7.
  * [5] J. Fan, K. C. Yung, M. Pecht, Physics-of-Failure-Based Prognostics and Health Management for High-Power White Light-Emitting Diode Lighting, IEEE Transactions on Device and Materials Reliability 11 (3) (2011) 407–416.  doi:10.1109/tdmr.2011.2157695.
URL https://dx.doi.org/10.1109/tdmr.2011.2157695

  * [6] M. Pecht, J. Gu, Physics-of-failure-based prognostics for electronic products, Transactions of the Institute of Measurement and Control 31 (3-4) (2009) 309–322.  doi:10.1177/0142331208092031.
URL https://dx.doi.org/10.1177/0142331208092031

  * [7] S. Yin, S. X. Ding, X. Xie, H. Luo, A Review on Basic Data-Driven Approaches for Industrial Process Monitoring, IEEE Transactions on Industrial Electronics 61 (11) (2014) 6418–6428.  doi:10.1109/tie.2014.2301773.
URL https://dx.doi.org/10.1109/tie.2014.2301773

  * [8] C. Hu, B. D. Youn, P. Wang, J. T. Yoon, Ensemble of data-driven prognostic algorithms for robust prediction of remaining useful life, Reliability Engineering & System Safety 103 (2012) 120–135.  doi:10.1016/j.ress.2012.03.008.
URL https://dx.doi.org/10.1016/j.ress.2012.03.008

  * [9] W. Zhou, T. G. Habetler, R. G. Harley, Bearing Condition Monitoring Methods forElectric Machines: A General Review, IEEE International Symposium onDiagnostics for Electric Machines (2007) 3–6.
  * [10] M. Hamadache, D. Lee, K. C. Veluvolu, Rotor Speed-Based Bearing Fault Diagnosis (RSB-BFD) Under Variable Speed and Constant Load, IEEE Transactions on Industrial Electronics 62 (10) (2015) 6486–6495.  doi:10.1109/tie.2015.2416673.
URL https://dx.doi.org/10.1109/tie.2015.2416673

  * [11] A. Rohan, S. H. Kim, RLC Fault Detection Based on Image Processing and Artificial Neural Network, INTERNATIONAL JOURNAL of FUZZY LOGIC and INTELLIGENT SYSTEMS 19 (2) (2019) 78–87.  doi:10.5391/ijfis.2019.19.2.78.
URL https://dx.doi.org/10.5391/ijfis.2019.19.2.78

  * [12] A. Rohan, S. H. Kim, Fault Detection and Diagnosis System for a Three-Phase Inverter Using a DWT-Based Artificial Neural Network, The International Journal of Fuzzy Logic and Intelligent Systems 16 (4) (2016) 238–245.  doi:10.5391/ijfis.2016.16.4.238.
URL https://dx.doi.org/10.5391/ijfis.2016.16.4.238

  * [13] A. Rohan, M. Rabah, S. H. Kim, An Integrated Fault Detection and Identification System for Permanent Magnet Synchronous Motor in Electric Vehicles, INTERNATIONAL JOURNAL of FUZZY LOGIC and INTELLIGENT SYSTEMS 18 (1) (2018) 20–28.  doi:10.5391/ijfis.2018.18.1.20.
URL https://dx.doi.org/10.5391/ijfis.2018.18.1.20

  * [14] P. K. Kankar, S. C. Sharma, S. P. Harsha, Fault diagnosis of ball bearings using machine learning methods, Expert Systems with Applications 38 (3) (2011) 1876–1886.  doi:10.1016/j.eswa.2010.07.119.
URL https://dx.doi.org/10.1016/j.eswa.2010.07.119

  * [15] J. Kumar, N. K. Ramkumar, S. Verma, Dixit, ‘Detection and classification for faults in drilling process using vibration analysis, Proc. Int. Conf. Prognostics Health Manage (2014) 1–6.
  * [16] J. Lee, H. Choi, D. Park, Y. Chung, H. Y. Kim, S. Yoon, ‘Fault detection and diagnosis of railway point machines by sound analysis, Sensors 16 (4) (2016) 549–549.
  * [17] K. Kemalkar, V. K. Bairagi, ‘Engine fault diagnosis using sound analysis, Proc. Int. Conf. Automat. Control Dyn. Optim. Techn. (ICACDOT) (2016) 943–946.
  * [18] A. Rohan, I. Raouf, H. S. Kim, Rotate Vector (RV) Reducer Fault Detection and Diagnosis System: Towards Component Level Prognostics and Health Management (PHM), Sensors 20 (23) (2020) 6845–6845.  doi:10.3390/s20236845.
URL https://dx.doi.org/10.3390/s20236845

  * [19] S. Li, H. Wang, L. Song, P. Wang, L. Cui, T. Lin, ‘An adaptive data fusion strategy for fault diagnosis based on the convolutional neural network, ” Measurement 165 (2020).
  * [20] L. Gou, H. Li, H. Zheng, H. Li, X. Pei, Aeroengine Control System Sensor Fault Diagnosis Based on CWT and CNN, Mathematical Problems in Engineering 2020 (5357146) (2020) 1–12.  doi:10.1155/2020/5357146.
URL https://dx.doi.org/10.1155/2020/5357146

  * [21] A. Khamparia, D. Gupta, N. G. Nguyen, A. Khanna, B. Pandey, P. Tiwari, Sound Classification Using Convolutional Neural Network and Tensor Deep Stacking Network, IEEE Access 7 (2019) 7717–7727.  doi:10.1109/access.2018.2888882.
URL https://dx.doi.org/10.1109/access.2018.2888882

  * [22] J. Wang, J. Zhuang, L. Duan, W. Cheng, ‘A multi-scale convolution neural network for featureless fault diagnosis, Proc. Int. Symp. Flexible Automat. (ISFA) (2016) 65–70.
  * [23] H. Liu, L. Li, J. M. . T. Tran, J. Lundgren, ‘Rolling bearing fault diagnosis based on STFT-deep learning and sound signals, Drill Fault Diagnosis Based on the Scalogram and Mel Spectrogram of Sound Signals 2016 (2016) 1–12.
  * [24] O. Janssens, V. Slavkovikj, B. Vervisch, K. Stockman, M. Loccufier, S. Verstockt, R. V. de Walle, S. V. Hoecke, Convolutional Neural Network Based Fault Detection for Rotating Machinery, Journal of Sound and Vibration 377 (2016) 331–345.  doi:10.1016/j.jsv.2016.05.027.
URL https://dx.doi.org/10.1016/j.jsv.2016.05.027

  * [25] R. Ziani, A. Hammami, F. Chaari, A. Felkaoui, M. Haddar, Gear fault diagnosis under non-stationary operating mode based on EMD, TKEO, and Shock Detector, Comptes Rendus Mécanique 347 (9) (2019) 663–675.  doi:10.1016/j.crme.2019.08.003.
URL https://dx.doi.org/10.1016/j.crme.2019.08.003

  * [26] A. Sapena-Bano, M. Pineda-Sanchez, R. Puche-Panadero, J. Martinez-Roman, D. Matic, Fault Diagnosis of Rotating Electrical Machines in Transient Regime Using a Single Stator Current’s FFT, IEEE Transactions on Instrumentation and Measurement 64 (11) (2015) 3137–3146.  doi:10.1109/tim.2015.2444240.
URL https://dx.doi.org/10.1109/tim.2015.2444240

  * [27] V. K. Rai, A. R. Mohanty, Bearing fault diagnosis using FFT of intrinsic mode functions in Hilbert–Huang transform, Mechanical Systems and Signal Processing 21 (6) (2007) 2607–2615.  doi:10.1016/j.ymssp.2006.12.004.
URL https://dx.doi.org/10.1016/j.ymssp.2006.12.004

  * [28] M. Portnoff, Time-frequency representation of digital signals and systems based on short-time Fourier analysis, IEEE Transactions on Acoustics, Speech, and Signal Processing 28 (1) (1980) 55–69.  doi:10.1109/tassp.1980.1163359.
URL https://dx.doi.org/10.1109/tassp.1980.1163359

  * [29] T. Tarasiuk, Hybrid Wavelet-Fourier Spectrum Analysis, IEEE Transactions on Power Delivery 19 (3) (2004) 957–964.  doi:10.1109/tpwrd.2004.824398.
URL https://dx.doi.org/10.1109/tpwrd.2004.824398

  * [30] J. A. Anatonio-Daviu, M. Riera-Guasp, J. R. Floch, M. P. Palomares, Validation of a New Method for the Diagnosis of Rotor Bar Failures via Wavelet Transform in Industrial Induction Machines, Vol. 42, IEEE, Piscataway, NJ, USA, 2006.
  * [31] J. Anden, S. Mallat, Deep Scattering Spectrum, IEEE Transactions on Signal Processing 62 (16) (2014) 4114–4128.  doi:10.1109/tsp.2014.2326991.
URL https://dx.doi.org/10.1109/tsp.2014.2326991

  * [32] V. C. cek, J, S. Mallat, P. Abry, M. Doret, Scattering ´transform for intrapartum fetal heart rate characterization and acidosis detection, Proc. IEEE EMBC (2013).
  * [33] S. J. Orfanidis, Introduction to Signal Processing, Prentice Hall, Englewood Cliffs, NJ, 1996.
  * [34] R. Schafer, What Is a Savitzky-Golay Filter? [Lecture Notes], IEEE Signal Processing Magazine 28 (4) (2011) 111–117.  doi:10.1109/msp.2011.941097.
URL https://dx.doi.org/10.1109/msp.2011.941097

  * [35] R. M. Haralick, L. G. Shapiro, Computer and Robot Vision I (1992) 158–205.
  * [36] R. van den Boomgaard, R. van Balen, Methods for fast morphological image transforms using bitmapped binary images, CVGIP: Graphical Models and Image Processing 54 (3) (1992) 252–258.  doi:10.1016/1049-9652(92)90055-3.
URL https://dx.doi.org/10.1016/1049-9652(92)90055-3

  * [37] T. Dau, B. Kollmeier, A. Kohlrausch, Modeling auditory processing of amplitude modulation. I. Detection and masking with narrow-band carriers, The Journal of the Acoustical Society of America 102 (5) (1997) 2892–2905.  doi:10.1121/1.420344.
URL https://dx.doi.org/10.1121/1.420344

  * [38] T. Chi, P. Ru, S. A. Shamma, Multiresolution spectrotemporal analysis of complex sounds, The Journal of the Acoustical Society of America 118 (2) (2005) 887–906.  doi:10.1121/1.1945807.
URL https://dx.doi.org/10.1121/1.1945807

  * [39] N. Mesgarani, M. Slaney, S. A. Shamma, Discrimination of speech from nonspeech based on multiscale spectro-temporal Modulations, IEEE Transactions on Audio, Speech and Language Processing 14 (3) (2006) 920–930.  doi:10.1109/tsa.2005.858055.
URL https://dx.doi.org/10.1109/tsa.2005.858055

  * [40] J. Bruna, S. Mallat, Invariant Scattering Convolution Networks, IEEE Transactions on Pattern Analysis and Machine Intelligence 35 (8) (2013) 1872–1886.  doi:10.1109/tpami.2012.230.
URL https://dx.doi.org/10.1109/tpami.2012.230