# Detection of Weak Signal Based on Parameter Identification of Delay Differential System with Noise Disturbance

Authors: QiuBao Wang; Yuejuan Yang; Xing Zhang
Journal: Mathematical Problems in Engineering
 (2020)

---

## Abstract

Extraction of the weak signals is crucial for fault prognostics in which case features are often very weak and masked by noise. In the time domain, the detection of weak signals depends on the identification of nonlinear parameters. A new signal detection and estimation method based on the incremental harmonic balance (IHB) is developed from the stochastic Van der Pol–Duffing equation with delayed feedback under parametric excitation. This is the first time that the IHB has been applied for the identification of parameters in stochastic delay differential equations (SDDEs). Compared to the method of intermittency transition between order and chaos to detect weak signals, this new method is more direct and the calculation result is what we want to obtain. This new method is suitable for the generalization and application of SDDEs.


---

## Body

## 1. Introduction

It is well known that a self-excited system can exhibit a remarkable complex dynamical behaviour. The earliest mathematical model has been proposed by Van der Pol to describe the vacuum tube circuit [1], and subsequently, many researchers have studied the complex phenomena observed from experiments and computer simulations in the fields of mechanical, electric, ecological, and chemical systems and so on. For rotating machinery, undesirable bifurcations, high-amplitude vibrations caused by resonance, quasiperiodic motion, and chaotic behaviour may occur and cause degradation or catastrophic failure of the system.Bearing failure is one of the foremost causes of breakdowns in rotating machinery, and such failure can be catastrophic, resulting in huge economic losses. To prevent these kinds of failures from happening, various bearing condition monitoring techniques have been developed. Among them, vibration analysis has been used extensively due to its intrinsic advantage of revealing bearing failure [2–4]. Detecting the dangerous vibrations is very crucial for engineering science, and various papers have been dedicated for the control of resonantly forced systems in various applicative fields; Atay [5] studied the effect of delayed position feedback on the response of a Van der Pol oscillator, showing that the delay can change the amplitude of limit cycle oscillations or suppress them altogether. Lai and Leng [6] considered a two-dimensional Duffing oscillator and proposed the weak-signal detection approach based on the GPASR model. A filtering technique through chaotic regime-switching is proposed by Denis et al. [7]. The inherent deficiency of the measuring mechanism leads to a mass of various noises to the effective signal inevitably. The signal of a defective bearing is spread across a wide frequency band and hence can become masked by noise and low-frequency signals with ease [8]. Prognostics is achieved by detecting the failure and defect at its initial stage and alerting maintenance personnel before it develops into a sharp loss. Nevertheless, in the early failure stage, it is difficult to detect the weak rolling bearing failure signal.To solve the problem of detecting effective weak failure signals under strong noise, numerous scholars have proposed various methods: (1) apply the transform to achieve a goal of signal demodulation. Pang et al. [9] present an improved Hilbert time transform (IHTT), and the wavelet transform has been used for signal demodulation [10]. (2) Construct filters based on noise type and application to conduct the denoising [11, 12]. (3) Utilize the characteristics of chaos system which is immune to noise and sensitive to weak signals at the same time to detect the weak effective failure signals. Zhao [13] detects the weak signal submerged in strong noise successfully by the transition from the chaotic state to a large-scale periodic state. (4) Use noise to enhance the amplitude of weak failure signal. Yang and Cao [14] uses the principle of primary resonance and stochastic resonance to study time delay-controlled SD oscillator with stiffness nonlinearities, and it shows that the weak harmonic signal can be enhanced with appropriate noise and time delay. (5) In recent years, detecting weak signals based on the deep network has become a new way in the background of machine learning artificial intelligence that can automatically extract bearing signals and identify faults [15, 16]. Xia et al. [17] propose an ensemble framework based on convolutional bidirectional long short-term memory with multiple time windows (MTW CNN-BLSTM ensemble) for accurately predicting remaining useful life.It is worth noting that stochastic resonance is a nonlinear system involving noise and weak periodic signals; under the right conditions, the noise enhances the periodic signal of the weak signal where Kramers rate [18] is the key. And, the chaotic threshold is of great importance for detecting harmonic weak signal in the way of transition from the chaos to large-scale period, where the method of Melnikov function [19] can solve this problem very well. Moreover, a mass of data is needed to train the network layer with deep network detecting weak signals. In addition, Dong et al. [20] integrate a Bayesian updating prognostic model using sensor-based degradation information for computing each machine’s time-to-failure (TTFs), with an opportunistic maintenance policy handling flexible system structures for optimizing the maintenance scheduling. In this paper, a new method for detecting weak signals is proposed based on SDDEs.In the time domain, the basis of detecting weak signals is parameter identification, starting from the 1970s, showing a variety of system parameter recognition methods. Nayfeh [21] proposed a nonlinear system parameter recognition method, which is based on the response of nonlinear systems and the comparison between the system hypothesis model.A method of recognition parameters based on the wavelet transform principle is proposed by Staszewski [22], which can identify parameters in multi-DOF nonlinear systems, but this method can only deal with free vibration system. Chen et al. [23] carried out parameter identification of nonlinear systems based on neural networks. In recent years, parameter identification is applied in various fields, such as automatic identification and medical research [24–26].The incremental harmonic balance (IHB) method was originally developed by Lau et al. [27] for treating periodic structural vibrations and which is a semianalytical method; moreover, it is remarkably effective in computer implementation for obtaining both the stable and unstable solutions. Then, through the efforts of many scholars, the IHB method has been successfully applied to many fields of nonlinear vibration, and Shen [28] successfully applied the IHB method to study the bifurcation and the route-to-chaos analyses. Dou and Ye [29, 30] raised the method of IHB to rectify the harmonic coefficients in frequency response calculating, which may help the identification algorithm to correct the noise-induced distortion in the experiment, and compared with the harmonic balance (HB) method, the accuracy of parameter identification has obvious advantages. Moreover, these methods are applied in ordinary differential equations (ODEs), and it is necessary to study identifying parameters with the delay differential equations (DDEs). The study of time-delay system is more consistent with the real vibration because a state in the past will affect the vibration state now [31]. Recently, Zhang [32] applies the IHB method to solve the superharmonic resonance solution of nonlinear time-delay differential equations. The research methods of parameter estimation in time-delay chaotic systems are endless, for example, neural network method [33] and self-synchronization method [34].In this paper, for the first time, the IHB method is applied to SDDEs for the identification of weak signals, and from the result of parameter identification, it is more direct than the chaos detection method. We take stochastic delay Van der Pol–Duffing equation, for example, to identify parameters. The delayτ is reduced by scaling, and the results show that it does not affect the identification. Calculating the explicit expressions K,P,andR is the important point of this paper, and then the explicit expressions of K,P,andR are substituted into equation (9) and iterated until the R is small enough; then, the result is the value of the parameter to be identified, which is clearly expressed. Finally, the ability of stochastic delay Van der Pol–Duffing system to detect weak signals is evaluated through the signal to noise ratio (SNR).


## 2. Stochastic Van der Pol–Duffing Oscillator under Time-Delayed Feedback

The Van der Pol–Duffing oscillator is a representative kind of system in the nonlinear system, and the system can show rich nonlinear dynamic behaviour with the change in strength in the period force [19, 35, 36]. As a classic chaotic system, Van der Pol–Duffing oscillator is often used in the modeling of dynamical systems; nowadays, there are many nonlinear problems in the fields of material resources, biology, engineering, psychiatry, and economics, which can be discussed and analyzed based on the modeling of nonlinear systems [37]. With the further study of vibration equation and the process of modeling, the stochastic Van der Pol–Duffing equation under the time-delayed feedback is more consistent with the process of actual measurement and signal acquisition, which has rich dynamic behaviour. The characteristic of being sensitive to a weak signal with strong noise perturbation makes it possible to be used in the field of weak signal detection. For stochastic Van der Pol oscillator with time-delayed feedback, some scholars have studied its response and stability analysis [38].The stochastic Van der Pol–Duffing oscillator under time-delayed feedback is in the form of the following:(1)x¨t+αx2t−τ−1x˙t+βxt=fcosΩt+σξt,where α and β are the damping coefficient and stiffness coefficient, respectively; fcosΩt is the external excitation; f and Ω are the amplitude and angular frequency of external excitation, respectively; τ is the time delay; x2t−τx˙t indicates that the square of displacement at t−τ will affect the velocity at t, namely, displacement feedback; and ξt is the Gaussian white noise with a constant power spectral density K˜ and autocorrelation function Rυ=Eξtξt+υ=δυ, where υ is autocorrelation time and δ is Diarc function; moreover, E⋅ denotes the expectation operator and σ=2πK˜ denotes noise intensity.For various parametersf,τ, and Ω, the system has rich nonlinear dynamic behaviours, seen from Figures 1–3.Figure 1
The bifurcation diagram off with certain fixed parameters τ=2, α=5, β=1, and Ω=1 shows that chaotic state ⟶ quasiperiodic motion ⟶ period-doubling bifurcation ⟶ periodic motion with the increasing f.Figure 2
The bifurcation diagram ofτ with certain fixed parameters α=0.5, β=1, Ω=1, and f=1.Figure 3
The bifurcation diagram ofΩ with certain fixed parameters α=0.45, β=1, τ=2, and f=1.


## 3. The Process of Parameter Identification

The influence of noise on the result of parameter identification can be fully considered in the process of numerical simulation. We can write equation (1) without noise as the equation of state, as follows:(2)y˙1t=y2t,y˙2t=fcosΩt−αy12t−τy2t+αy2t−βy1t.Then, take scale transformationt⟶t/τ and write it as follows for convenience:(3)y¨t−ατ2y˙t+βτ2yt+∑i=1Mbigiyt,y˙t,yt−1,τ,t=0,where M=2,b1=τ2α,b2=−τ2f,g1=y2t−1y˙t,g2=cosΩτt, and b1 b2 are nonlinear parameters to be identified.


## 4. Construct Increment Process

y
,
y
t
−
1, and bii=1,2 are as follows:(4)yt=y0t+Δyt,yt−1=y0t−1+Δyt−1,bi=bi0+Δbi,i=1,2.Working with initial valuesbi0=0i=1,2. Determine the initial values of the displacement y0t and y0t−1 based on the experimental data, that is, depending on the displacement response on the t0,t0+T when the system enters the steady-state period. Let ω=Ω/μ (μ is an integer), period T=2π/ω. The initial values of displacement are expanded into the Fourier series in plural form, and they are written in the following matrix form, where N is the degree of expansion of harmonic term:(5)y0t=∑k=−NNejkωt⋅ak=Y1⋅A0,(6)y0t−1=∑k=−NNejkωt−1⋅ak=Y2⋅A0,where(7)Y1=ej−Nωt…1…ejNωt,Y2=ej−Nωt−1…1…ejNωt−1,A0=a−N…a0…aNT.Each element inA0 is the Fourier coefficient of each order, where the k-order Fourier coefficient is obtained by numerical integral transformation of the data series of the displacement response:(8)ak1=1T∫0Ty˜0t⋅e−jkωtdt,ak2=1T∫0Ty˜0t−1⋅e−jkωt−1dt,where y˜0t and y˜0t−1 are a sequence of displacement response measured and collected in experiments. To the extent that the error allows ak1≈ak2, ak=1/T∫0Ty˜0t⋅e−jkωtdt. Working with the incremental processes. Write Δyt and Δyt−1 in the form of Fourier series in the plural form:(9)Δyt=Y1⋅ΔA,Δyt−1=Y2⋅ΔA,ΔA=Δa−N…Δa0…aNT.The increments of the parameters are Δb=Δb1Δb2T.Formulas in (4) are substituted into equation (3), by the means of Taylor unfold, and omitting the incremental two times and more than two high-order items, we can obtain(10)Δy¨t−ατ2Δy˙t+βτ2Δyt+2b10y0t−1y˙0tΔyt−1+b10y02t−1Δy˙t+y02t−1y˙tΔb1+cosΩτtΔb2=−y¨0t+ατy˙0t−βτ2y0t−b10y02t−1y˙0t−b20cosΩτt.


## 5. Construct Increment Equation

With the method of Galerkin and the regularΔyt=Y1ΔA=ΔATY1T=ΔAY1T=ΔytT, we can obtain(11)δΔAT∫t0t0+TY1TY¨1−ατY˙1+βτ2Y1+2b10y0t−1y˙0tY2+b10y02t−1Y˙1dt⋅ΔA+δΔAT∫t0t0+TY1Ty02t−1y˙tdt⋅Δb1+δΔAT∫t0t0+TY1TcosΩτtdt⋅Δb2=δΔAT∫t0t0+TY1T−y¨0t+ατy˙0t−βτ2y0t−b10y02t−1y˙0t−b20cosΩτtdt.With arbitrary of δΔA, we can obtain the incremental equation as follows:(12)KΔA+PΔb=R.The elements ofK, P, and R are(13)K=1T∫t0t0+TY1TY¨1−ατY˙1+βτ2Y1+2b10y0t−1y˙0tY2+b10y02t−1Y˙1dt,P1=1T∫t0t0+TY1Ty02t−1y˙tdt,P2=1TY1TcosΩτtdt,R=1T∫t0t0+TY1T−y¨0t+ατy˙0t−βτ2y0t−b10y02t−1y˙0t−b20cosΩτtdt.To know the explicit expression of elements inK, P, and R, according to the plural form of Fourier series of the initial displacement expansion, the coefficients of the Fourier series in the plural form corresponding to y˙0t,y02t−1,y02t−1y˙0t,y¨0t,cosΩt are obtained:(14)y˙0t=∑k=−NNak3⋅ejkωt,ak3=am⋅jmωδmk,y02t−1=∑k=−NNak4⋅ejkωt,ak4=amal⋅e−jωm+lδm+lk,y0t−1y˙0t=∑k=−NNak5⋅ejkωt,ak5=aman3e−jmωδm+nk,y02t−1y˙0t=∑k=−NNak6⋅ejkωt,ak6=am4al3δm+lk,y¨0t=∑k=−NNak7⋅ejkωt,ak7=−m2ω2amδmk.Calculate the Fourier coefficientqk of the harmonic excitation item cosΩτt based on the experimental data:(15)cosΩτt=∑k=−NNqk⋅ejkωt,qk=1T∫0TF˜tρe−jkωtdt,where F˜t is the excitation force signal collected experimentally and ρ is shaking force amplitude. Let m¯=m+N+1,n¯=n+N+1, then we can write the explicit expression of elements in K, P, and R:(16)Km¯,n¯=−n2ω2−ατjnω+βτ2δm+n0+∑k=−NN2b10ak2+b10ak3jnωe−jnωδm+n+k0,P1m¯,1=ak4δm+k0P2m¯,2=qkδm+k0,Rm¯,1=−ak5+ατak1−βτ2ak−b10ak4−b20qkδm+k0.Equation (16) is substituted into equation (12) iteratively. It is clear that equation (12) has 2N+1+2 unknowns Δak,k=−N,…,N;Δb1,Δb2, but there are 2N+1 equations. So the increment of the harmonic item coefficients ΔA and Δb as the tuning amount is obtained by equation (12) on condition that no less than 2 principal harmonic coefficients increment Δak zero. Δy, y, and b1,b2can be obtained from equations (4) and (9), respectively, and then replace the initial values of y0 and b10, b20 to (16) and enter the next iteration, and so on. We do not stop iterating until R is small enough. When the iteration stops, the results of b1 and b2 are the value of the parameter to be identified.


## 6. The Results of Identification

Consider the parametersb1=ατ2 and b2=−fτ2 to be identified, and given parameter initial value b10=b20=0; constructing incremental equations, select a0,a±1,a±2 as the main harmonic coefficients, that is, Δa0=Δa±1=Δa±2=0; solve the remaining harmonic item increments and parameter increments. and iterate over and over until the algorithm converges steadily. The identification results are shown below (Table 1).Table 1
Parameter recognition results of incremental harmonic equilibrium method in different states of time-delay system, whereb1∗=ατ and b2∗=−τ2f are exact values.
Parameters
σ
0%
1%
5%
N
f
=
6
f
=
7
f
=
8
f
=
6
f
=
7
f
=
8
f
=
6
f
=
7
f
=
8
b
1
10
10.00968
10.00962
10.00515
10.00962
10.00962
10.00512
10.00975
10.00965
10.00515
15
10.00956
10.00957
10.00497
10.00961
10.00951
10.00511
10.00957
10.00950
10.00503
20
10.00959
10.00954
10.00503
10.00961
10.00951
10.00506
10.00972
10.00952
10.00503
25
10.00961
10.00951
10.00501
10.00957
10.00950
10.00510
10.00963
10.00954
10.00502
30
10.00959
10.00949
10.00503
10.00973
10.00950
10.00500
10.00951
10.00963
10.00501
b
2
10
24.00924
28.00882
32.00723
24.00920
28.00860
32.00723
24.00903
28.00873
32.00723
15
24.00941
28.00886
32.00769
24.00950
28.00871
32.00710
24.00933
28.00871
32.00810
20
24.00937
28.00889
32.00769
24.00934
28.00882
32.00820
24.00932
28.00872
32.00841
25
24.00940
28.00882
32.00774
24.00931
28.00881
32.00720
24.00891
28.00881
32.00790
30
24.00932
28.00886
32.00769
24.00943
28.00881
32.00770
24.00937
28.00877
32.00757Figures4–6 show that the recognition results of parameters b1 and b2 when the system is in different states, such as period, period-2, and chaos. We can get from the figures that the results of parameter recognition are stable when the system is in different states. It is worth noting that the difficulty of parameter recognition using the IHB method is finding the appropriate period initial displacement value from the steady-state response to recognize the parameters.Figure 4
The influence of noise on parameter recognition in a periodic state withf=8andτ=2.Figure 5
The influence of noise on parameter recognition in period-2 state withf=7andτ=2.Figure 6
The influence of noise on parameter recognition in a chaotic state withf=6andτ=2.


## 7. Identify the Amplitude of Weak Signals

Add the weak signal under strong noiserτ2cosΩτt+στ2ξt to equation (3) with r=0.001 and τ=2, and we can obtain(17)y¨−ατ2y˙t+βτ2yt+ατ2y2t−1y˙t=f+rτ2cosΩτt+στ2ξt,where the parameter to be identified is b = −(f + r)τ2. From the results of identifying Table 2 and Figure 7 with certain parameters, we can obtain that b≈−32.004 and then the amplitude of the weak harmonic signal is r=0.001, which is consistent with the added signal amplitude so that the weak signal can be identified well.Table 2
Results of parameters identification with the method of IHB.
Parameter
σ
0%
1%
5%
N
Results of parameter identification
b
b
∗
=
−
τ
2
f
+
r
10
−32.020632
−32.020633
−32.020623
15
−32.020632
−32.020637
−32.020638
20
−32.020633
−32.020643
−32.020641
25
−32.020641
−32.020631
−32.020630
30
−32.020619
−32.020620
−32.020625Figure 7
The influence of noise on parameter recognition with fixed parametersf=8,r=0.001,α=5,β=1,Ω=1,andτ=2.Consider the SNR to judge how well the system detects weak signal under the noise intensityσ=0.05:(18)SNR=10lg0.5r2σ2=10lg0.50.00120.052=−36.99dB.The amplitude of the weak signal can be detected under the condition of low SNR, and it is proved that this method is worthy of application and popularization in the field of DDEs.According to the system response power spectral density function under different noises, we can obtain the power spectral density (PSD) images, as shown in Figure8; the frequency at the highest energy is fre=1/T=0.3180 obviously. The angular frequency of a weak signal is Ωτ=2π/T, so T=π, and then, the frequency of a weak signal is 1/T=0.3182. We can detect weak signals from strong noise.Figure 8
PSD under different noise intensityσ=1%andσ=5%.


## 8. Conclusions

In this paper, there is a rich dynamic behaviour of the stochastic Van der Pol–Duffing oscillator with noise interference, and the parameters can be identified by the method of IHB in different states. For a very low SNR, the detection results are very good and this method can be extended to detect weak signals in nonlinear DDEs containing the certain noise.



---

*Source: 2047952-2020-04-14.xml*