Stability and Synchronization of Delayed Quaternion-Valued Neural Networks under Multi-Disturbances

Abstract

This paper discusses a type of mixed-delay quaternion-valued neural networks (QVNNs) under impulsive and stochastic disturbances. The considered QVNNs model are treated as a whole, rather than as complex-valued neural networks (NNs) or four real-valued NNs. Using the vector Lyapunov function method, some criteria are provided for securing the mean-square exponential stability of the mixed-delay QVNNs under impulsive and stochastic disturbances. Furthermore, a type of chaotic QVNNs under stochastic and impulsive disturbances is considered using a previously established stability analysis method. After the completion of designing the linear feedback control law, some sufficient conditions are obtained using the vector Lyapunov function method for determining the mean-square exponential synchronization of drive–response systems. Finally, two examples are provided to demonstrate the correctness and feasibility of the main findings and one example is provided to validate the use of QVNNs for image associative memory.

1. Introduction

Neural networks (NNs) have been extensively used in medicine, biology, communication, and transportation because of their characteristics, such as nonlinear mapping, associative memory, classification and recognition, optimal computing, and others. A quaternion is a four-element vector in the class of hypercomplex numbers first described by Hamilton in 1843 [1]. Because a quaternion with one real and three imaginary parts can carry large amounts of information, it can handle multidimensional issues, such as associative memory for colored images [2], trajectory tracking control of a rotor missile [3], attitude and altitude tracking of a quadrotor unmanned aerial vehicle [4], and cross-modal matching [5] more effectively. Considering image associative memory as an example, a color image can be depicted by a quaternion matrix with pure imaginary numbers of suitable dimensionality, where the three imaginary parts correspond to the three basic monochromatic colors of red (R), green (G), and blue (B) in the RGB color space. This demonstrates that a quaternion matrix can hold more data than a real or complex matrix. Accordingly, a quaternion-valued neural network (QVNN) model is proposed. QVNNs combine the benefits of NNs and quaternions and outperform real-valued NNs and complex-valued NNs in solving three-dimensional and four-dimensional data problems; therefore, they have attracted considerable research interest [2,4,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29].

During the implementation of NNs, the delay phenomenon is inevitable owing to the limitation of the switching speed of the amplifier or the signal transmission speed. The delay phenomenon may reduce the convergence rate to equilibria of NNs and even cause considerable damage to the stability of systems. The delay phenomena of NN models have been examined in the literature, including fixed delays [30], variable delays [9,10,11,12,13,21,22,23,24,25,26,31,32,33,34,35,36], infinitely distributed delays [21,33,34,36], and neutral delays [12]. Uncertain interference factors are also unpreventable in NNs. Impulsive and stochastic disturbances are two types of interference factors that are commonly found in practical systems, and they considerably degrade the system’s stability. Therefore, incorporating them into NN models when examining the dynamic behavior of NNs is imperative. In [13], an impulsive disturbed QVNN with time-varying delays was built. By decomposing the QVNNs into the equivalent real-valued NNs, the exponential stability conditions of the delayed system were derived using generalized norms. After stochastic disturbances were introduced into these models, and without decomposing QVNNs, the researchers of [12] obtained relevant criteria to ensure the mean-square stability of equilibria of a type of QVNNs with variable delays and neutral delays using the fixed-point theorem.

Because synchronization of nonlinear complex systems can be applied to image encryption, secure communication, associative memory, and many other fields [8], research on the synchronization of these systems has always been a priority. Synchronization has attracted considerable research interest in recent years because it is an important dynamic behavior of coupled chaotic NNs. Similar to the equilibria stability of NNs, the delay phenomenon [6,7,8,15,16,17,19,20,28,29,37,38,39,40,41,42,43,44,45,46,47,48], impulsive disturbances [6,14,27], and stochastic disturbances [38,39,40,43,44,45,46,48,49,50,51,52,53] should be considered when examining the synchronization control of chaotic NNs. Currently, most research results [37,38,39,40,41,42,43,44,45,46,47,48,49,50,51,53,54,55] on synchronization control are primarily applicable to NNs in the real number domain. However, research [6,7,8,14,15,16,17,18,19,20,27,28,29] on the synchronization control of NNs in the quaternion domain is scarce. Studies [8,14,27,56] have only considered impulsive disturbances when analyzing the synchronization of QVNNs. Other related studies [15,16,17,18,19,20,28,29] have not considered the influence of uncertain interference factors on synchronization control.

Based on the abovementioned analysis of the existing literature, it is concluded that:

(1)

Considering the existing uncertain QVNNs, studies [8,14,27,56] only considered impulsive disturbances, and the researchers of [12,57] only investigated stochastic disturbances. Furthermore, the stochastic disturbance introduced in the QVNNs in [12] does not exert any influence on connecting neurons.

(2)

Decomposition is a common method when dealing with multivalued NNs. Studies [7,9,10,11,13,14,19,20,21,23,24,25,26,27] decomposed QVNNs and proposed relevant criteria to ensure system stability or synchronization. However, decomposable quaternion activation functions are rare, and the assumptions associated with decomposable activation functions are strict.

(3)

The scalar Lyapunov function method is the most commonly used method for investigating the stability and synchronization of NNs [6,7,8,10,11,13,14,15,16,17,18,19,31,32,33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,48,49,50,51,52,53,54,55,57]. Compared with vector Lyapunov functions, scalar Lyapunov functions are difficult to construct using this method.

(4)

The synchronization analysis of chaotic NNs based on the concept of the drive-response system is equivalent to the stability problem of a synchronization error system. Only a few studies have simultaneously investigated the stability and synchronization of QVNNs.

(5)

Based on certain important conclusions on stability and synchronization reported in the literature [6,8,9,10,11,12,13,14,15,16,17,18,19,21,22,23,24,25,27,28,29,31], examples are provided to validate the obtained results. However, no examples are provided to demonstrate the use of QVNNs.

To sum up, additional research on the dynamic behavior of delayed QVNNs under stochastic and impulsive disturbances is required. Therefore, we further investigate the aforementioned issues. The main contributions and advantages of this study are as follows:

(1)

We will consider both impulsive and stochastic disturbances in mixed-delay QVNNs and examine the interaction of connected neurons after introducing stochastic disturbances.

(2)

The non-decomposition method, which retains the coupling characteristics of each part of the quaternion, is used to conduct the research.

(3)

To avoid the difficulty of constructing the scalar Lyapunov functions, this study combines the vector Lyapunov function method with mathematical induction and differential integration theory and proposes sufficient criteria to ensure the mean-squared exponential stability of the equilibria of the system.

(4)

For a type of chaotic QVNNs with mixed delays as well as impulsive and stochastic disturbances, judgment conditions for ensuring the mean-square exponential synchronization of the drive-response system with a linear feedback controller are obtained using the established mean-square exponential stability conditions.

(5)

In this study, the correctness and feasibility of the stability and synchronization criteria are demonstrated in examples 1 and 2, respectively. Furthermore, the applicability of the criterion for realizing image associative memory is validated via example 3.

This paper comprises six main sections. Section 2 presents a description of the model and the associated assumptions. Section 3 and Section 4 present theorems and corollaries related to stability and synchronization judgments of QVNNs. Section 5 presents three simulation examples. Finally, conclusions and an outlook for future work are presented in Section 6.

Table 1. Necessary notations and means throughout the paper.

Notations	Means	Notations	Means
$ℝ$	real number domain	$|z|$	$|z|=\sqrt{{(z^{(0)})}^2+(z^{(1)}{)}^2+(z^{(2)}{)}^2+(z^{(3)}{)}^2}$
$ℂ$	complex number domain	$\mathrm{Re}(z)$	$\mathrm{Re}(z)=z^{(0)}$
$\mathbb{Q}$	skew field of quaternions	$\mathrm{I}\mathrm{m}(z)$	$\mathrm{I}\mathrm{m}(z)=z^{(1)} i+z^{(2)} \iota +z^{(3)} \kappa$
$\mathbb{N}$	natural number set	$\overline{z}$	conjugate of $z\in \mathbb{Q}$ defined as $\overline{z}=z^{(0)}-z^{(1)} i-z^{(2)} \iota -z^{(3)} \kappa$
$\mathsf{\Lambda }$	set defined as $\{1, 2, \dots , n\}$	$|z|$	modulus of vector $z$ defined as ${(|z_1|, |z_2|,\dots , |z_n|)}^{\mathrm{T}}$
${(\cdot )}^{\mathrm{T}}$	transpose of a vector or matrix	$\Vert z\Vert$	norm of $z$ defined as $\sqrt{{\displaystyle {\sum }_{m=1}^n|z_m{|}^2}}$
$i, \iota , k$	imaginary unit	$|B|$	modulus of the matrix $B\in {\mathbb{Q}}^{n\times n}$ defined as ${(|b_{mr}|)}_{n\times n}\in {ℝ}^{n\times n}$
${\rm E}(.)$	expectation function	$(\mathsf{\Omega }, \mathrm{F}, \mathsf{{\rm P}})$	complete probability space
$z\in \mathbb{Q}$	$z=z^{(0)}+z^{(1)}i+z^{(2)}\iota +z^{(3)}\kappa$

lists some symbols used in this paper.

2. Model Description and Preliminaries

A type of QVNNs with mixed delays as well as impulsive and stochastic disturbances is presented in this section, which can be expressed as follows:

$$\{\begin{array}{l}\mathrm{d}{z}_m(t)=\{-{\omega }_m{z}_m(t)+{\displaystyle \sum _{r=1}^n[a_{mr}{f}_r(z_r(t))+b_{mr}{g}_r(z_r(t-{\tau }_{mr}(t)))+}\\ \hspace{1em}\hspace{1em}\hspace{1em} c_{mr}{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s)h_r(z_r(s))}\mathrm{d}s]\}\mathrm{d}t+\\ \hspace{1em}\hspace{1em}\hspace{1em} {\displaystyle \sum _{r=1}^n{p}_{mr}(z_r(t),z_r(t-{\tau }_{mr}(t)))}\mathrm{d}{w}_r(t)+J_m(t),\hspace{1em}t\ne t_k,\\ \mathsf{\Delta }{z}_m(t_k)=z_m(t_k^{+})-z_m(t_k^{-}),\hspace{1em}t=t_k,\end{array}$$

(1)

where: $z_m\in \mathbb{Q}$ denotes the state of the m-th neuron, $m\in \mathsf{\Lambda }$, and n denotes the number of neurons. $\mathsf{\Delta }{z}_m(t_k)$ represents an abrupt change in the system state at a discrete time $t_k$, $k\in \mathbb{N}$. The discrete set $\left\{t_k\right\}$ satisfies $0\le t_0<t_1<\dots <t_k<\dots$, and $t_k\to \infty$ when $k\to \infty$. $A={(a_{mr})}_{n\times n}$.$\in {\mathbb{Q}}^{n\times n}$ and $B={(b_{mr})}_{n\times n}\in {\mathbb{Q}}^{n\times n}$ $C={(c_{mr})}_{n\times n}\in {\mathbb{Q}}^{n\times n}$ are the weight matrices of the system (1) for no delay, variable delays, and infinite distributed delays between neurons, respectively. The function $P={(p_{mr})}_{n\times n}\in {\mathbb{Q}}^{n\times n}$ signifies the intensity of stochastic disturbances among neurons. $W(t)={(w_1(t),w_2(t), \dots ,w_n(t))}^{\mathrm{T}}$ represents the Brown motion defined $(\mathsf{\Omega }, \mathrm{F}, \mathsf{{\rm P}})$ with natural filtration ${\left\{{\mathrm{F}}_t\right\}}_{t\ge 0}$. $J={(J_1(t),J_2(t),\dots ,J_n(t))}^{\mathrm{T}}\in {\mathbb{Q}}^n$ indicates the external input of the system. $f(z)={(f_1(z_1),f_2(z_2),\dots ,f_n(z_n))}^{\mathrm{T}}$, $g(z)=(g_1(z_1),g_2(z_2),$$\dots ,g_n(z_n){)}^{\mathrm{T}}$ and $h(z)=(h_1(z_1),$ $h_2(z_2),\dots ,h_n(z_n){)}^{\mathrm{T}}$ represent activation functions of the system (1) for no delays, variable delays, and infinite distributed delays among neurons, respectively. The positive matrix $\omega =\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}({\omega }_1,{\omega }_2,\dots ,{\omega }_n)\in {ℝ}^{n\times n}$ denotes the self-feedback interconnected matrix.

Multiple assumptions are expressed for the system (1).

Assumption 1.

The following assumptions are made for the delay part of the system (1).

(i) The variable time delay ${\tau }_{mr}(t)$ in model (1) is a bounded function $\tau ={\max }_{m,r\in \Lambda }{\sup }_{t\ge 0}{\tau }_{mr}(t)$.

(ii) The kernel function ${\theta }_{mr}$:$[0,+\infty )$$\to$$[0,+\infty )$ in the infinite distributed delay is a piecewise continuous function, and satisfies

$${\displaystyle {\int }_0^{+\infty } e^{\beta s}{\theta }_{mr}(s)ds={\mu }_{mr}(\beta )}, m,r\in \mathsf{\Lambda },$$

(2)

Here ${\mu }_{mr}(\beta )$ is a continuous function on $[0, \delta )$ with ${\mu }_{mr}(0)=1$ and $\delta >0$.

Assumption 2.

It is presumed that activation functions $f_m(.)$ with $f_m(0)=0$ in system (1) satisfy the Lipschitz condition, indicating that $l_m^f>0$ such that the inequalities $|f_m(z_m)-f_m(v_m)|\le l_m^f|z_m-v_m|$ are valid, here $z_m,v_m\in \mathbb{Q}$, $m\in \Lambda$.

Remark 1.

Note that the activation functions $g_m(.)$ and  $h_m(.)$ with $g_m(0)=h_m(0)=0$ in the system (1) satisfy the Lipschitz condition. Only the notation f in the Assumption 2 is substituted by $g$ or h.

Here, let $L^f=\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(l_1^f,l_2^f,\dots ,l_n^f)$, $L^g=\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(l_1^g,l_2^g,\dots ,l_n^g)$, and $L^h=\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(l_1^h,l_2^h,\dots ,l_n^h)$.

Assumption 3.

Suppose $z_m(t_k)=z_m(t_k^{+})$ and $z_m(t_k^{-})={\lim }_{t\to t_k^{-}}z(t)$, $m\in \Lambda$, $k\in \mathbb{N}$. Let $\mathsf{\Delta }{z}_m(t_k)=$$z_m(t_k^{+})-z_m(t_k^{-})$$=I_m^k(z_m(t_k^{-}))$, where $I_m^k(.)$ with $I_m^k(0)=0$ implying an impulsive function. We assume that the inequalities $|z_m(t_k^{-})+I_m^k(z_m(t_k^{-})){|}^2\le {\eta }_m^{(k)}|z_m(t_k^{-}){|}^2$ are valid for all $m\in \mathsf{\Lambda }$, $k\in \mathbb{N}$. Let ${\eta }^{(k)}={\max }_{m\in \mathsf{\Lambda }}\{1, {\eta }_m^{(k)}\}$, $k\in \mathbb{N}$.

Assumption 4.

There exist positive numbers ${\pi }_{mr}$ and ${\pi }_{mr}^{(\tau )}$ such that all $z_r\ne q_r$ the weighted functions of the stochastic term $p_{mr}$ with $p_{mr}^{ }(0,0)=0$ ($m,r\in \mathsf{\Lambda }$) satisfy the follow-up inequalities:

$$\begin{gathered} \hspace{1em}{\left(p_{mr}(z_r(t),z_r(t-{\tau }_{mr}(t)))\right)}^2-{\left(p_{mr}(q_r(t),q_r(t-{\tau }_{mr}(t)))\right)}^2 \\ \le {\pi }_{mr}{|z_r(t)-q_r(t)|}^2+{\pi }_{mr}^{(\tau )}{|z_r(t-{\tau }_{mr}(t))-q_r(t-{\tau }_{mr}(t))|}^2 \end{gathered}$$

(3)

Here, we denote $\Pi ={({\pi }_{mr}^{ })}_{n\times n}$ and ${\Pi }^{(\tau )}={({\pi }_{mr}^{(\tau )})}_{n\times n}$.

Remark 2.

In studies [12,31,33,36,38,39,40,44,45,48,49,50,51,53,57], stochastic disturbances were introduced into several NNs. As per the interaction of stochastic disturbance on each neuron, the function $p_{mr}$ in (1) can be split into strong-coupling forms [31,33,38,48,49,50,53] and weak-coupling types [36,40,44,45,51]. The research objects of the studies [31,33,36,38,39,40,44,45,48,49,50,51,53] are real-valued NNs. Reference [12] investigated the stability of a class of QVNNs with firmly coupled stochastic perturbations, but without considering the delay in the stochastic term. In this analysis, the assumption conditions for stochastic disturbances involve not only strong coupling but also delays. Obviously, Assumption 4 in this paper includes Hypothesis H2 in [12].

Let ${\phi }_m(s)$ be the continuous function mapping from $(-\infty ,0]$ to $\mathbb{Q}$, and $z_m(s)={\phi }_m(s)$ be the original condition of the system (1), where $s\in (-\infty ,0]$ and $m\in \mathsf{\Lambda }$.

The equilibria of model (1) are denoted as $z^{\#}={(z_1^{\#},z_2^{\#},\dots ,z_n^{\#})}^{\mathrm{T}}$.

Definition 1. 

If there exist constants $\lambda >0$ and $\Gamma >0$ such that $E(||z(t)-z^{\#}|{|}^2)$$\le \Gamma {\sup }_{s\in (-\infty ,0]}E(||\phi (s)-z^{\#}|{|}^2)e^{-\lambda (t-t_0)}$ hold for external input $J\in {\mathbb{Q}}^n$, the equilibria $z^{\#}={(z_1^{\#},z_2^{\#},\dots ,z_n^{\#})}^{\mathrm{T}}$ of system (1) are mean-square exponential stability,

• where:

  $$\mathrm{E}(||z(t)-z^{\#}|{|}^2)={[\mathrm{E}(|z_1(t)-z_1{ }^{\#}{|}^2),\mathrm{E}(|z_2(t)-z_2{ }^{\#}{|}^2),\dots ,\mathrm{E}(|z_n(t)-z_n{ }^{\#}{|}^2)]}^{\mathrm{T}}, t>t_0\ge 0,$$

  $$\mathrm{E}(||\phi (s)-z^{\#}|{|}^2)={[\mathrm{E}(|{\phi }_1(s)-z_1{ }^{\#}{|}^2),\mathrm{E}(|{\phi }_2(s)-z_2{ }^{\#}{|}^2),\dots ,\mathrm{E}(|{\phi }_n(s)-z_n{ }^{\#}{|}^2)]}^{\mathrm{T}}, s\in (-\infty ,0].$$

Lemma 1

([58]). For a differentiable quaternionic function $z(t)\in \mathbb{Q}$, the formula $\frac{\mathrm{d}{|z(t)|}^2}{\mathrm{d}t}=2\mathrm{Re}\left(\overline{z(t)}\frac{\mathrm{d}z(t)}{\mathrm{d}t}\right)$ is valid.

Lemma 2

([59]). It is assumed that $z(t)$ is an $It\widehat{o}$ process described by

$$\mathrm{d}z(t)=f(t)\mathrm{d}t+p(t)\mathrm{d}w(t).$$

(4)

$V(t,z(t)):\left[0,\infty \right)\times \mathbb{Q}\to ℝ$ is supposed an $It\widehat{o}$ process with second-order derivative respect to z(t), then

$$\mathrm{d}V(t,z(t))=\frac{\partial V(t,z(t))}{\partial t}\mathrm{d}t+\frac{\partial V(t,z(t))}{\partial z(t)}\mathrm{d}z(t)+\frac{1}{2}\frac{{\partial }^{2}V(t,z(t))}{\partial z^2(t)}{\left(\mathrm{d}z(t)\right)}^2.$$

(5)

Here ${\left(\mathrm{d}z(t)\right)}^2=\left(\mathrm{d}z(t)\right)\left(\mathrm{d}z(t)\right)$ is estimated using $\mathrm{d}t$. $\mathrm{d}t$ = $\mathrm{d}t$. $\mathrm{d}w(t)$ = $\mathrm{d}w(t).\mathrm{d}t$ = 0 and $\mathrm{d}w(t)$. $\mathrm{d}w(t)=\mathrm{d}t$.

Remark 3. 

After substituting (4) into (5) and using Lemma 2, Equation (6) can be derived, which will be employed in the subsequent proof of the theorems,

$$\begin{gathered} \mathrm{d}V(t,z(t))=\left(\frac{\partial V(t,z(t))}{\partial t}+\frac{\partial V(t,z(t))}{\partial z(t)}f(t)+\frac{1}{2}\frac{{\partial }^{2}V(t,z(t))}{\partial z^2(t)}{(p(t))}^2\right)\mathrm{d}t+ \\ \frac{\partial V(t,z(t))}{\partial z(t)}p(t)\mathrm{d}w(t). \end{gathered}$$

(6)

Lemma 3

([21]). Let $B={(b_{mr})}_{n\times n}\in {ℝ}^{n\times n}$ be a matrix $b_{mr}\le 0$ $(m\ne r, m,r\in \mathsf{\Lambda })$. If the original parts of all eigenvalues in B are positive, then the conclusion that B is an M-matrix is consistent with the conclusion that there is a positive vector $\varsigma \in {ℝ}^n$ such that $B\varsigma >0$.

3. Main Results of Stability

Let $q(t)=z(t)-z^{\#}$. Accordingly, the system (1) is rewritten as follows:

$$\{\begin{array}{l}\mathrm{d}{q}_m(t)=\left\{-{\omega }_m{q}_m(t)+{\displaystyle \sum _{r=1}^n[a_{mr}{\tilde{f}}_r(q_r(t))+b_{mr}{\tilde{g}}_r(q_r(t-{\tau }_{mr}(t)))+c_{mr}{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s){\tilde{h}}_r(q_r(s))}\mathrm{d}s]}\right\}\mathrm{d}t\\ \hspace{1em}\hspace{1em}\hspace{1em} +{\displaystyle \sum _{r=1}^n{\tilde{p}}_{mr}(q_r(t),q_r(t-{\tau }_{mr}(t)))}\mathrm{d}{w}_r(t),\hspace{1em}t\ne t_k,\\ \mathsf{\Delta }{q}_m(t_k)=q_m(t_k^{+})-q_m(t_k^{-}),\hspace{1em}t=t_k,\end{array}$$

(7)

where:

$${\tilde{f}}_r(q_r(t))=f_r(q_r(t)+z_r^{\#})-f_r(z_r^{\#}), {\tilde{h}}_r(q_r(s))=h_r(q_r(s)+z_r^{\#})-h_r(z_r^{\#}),$$

$${\tilde{g}}_r(q_r(t-{\tau }_{mr}(t)))={\tilde{g}}_r(q_r(t-{\tau }_{mr}(t))+z_r^{\#})-g_r(z_r^{\#}),$$

$${\tilde{p}}_{mr}(q_r(t),q_r(t-{\tau }_{mr}(t)))=p_{mr}(q_r(t)+z_r^{\#},q_r(t-{\tau }_{mr}(t))+z_r^{\#})-p_{mr}(z_r^{\#},z_r^{\#}), m,r\in \mathsf{\Lambda }.$$

The initial condition of the system (7) is ${\psi }_m(s)={\phi }_m(s)-z_m^{\#}$, $m\in \mathsf{\Lambda }$, $-\infty <s\le 0$.

Apparently, the mean-square exponential stability of the zero solution of system (7) indicates that the equilibria of system (1) show mean-square exponential stability. In what follows, a theorem is suggested to determine the mean-square exponential stability of the equilibria $z^{\#}$ of (1).

Theorem 1. 

Suppose that Assumptions 1–4 hold. Let $\eta ={\lim }_{k\to \infty }\left(\sup \frac{\ln {\eta }^{(k)}}{t_k-t_{k-1}}\right)$. If there is a positive constant $\lambda >0$ with $\lambda >\eta >0$ and a positive number ${\zeta }_m$ such that for all $J\in {\mathbb{Q}}^n$ and $m\in \mathsf{\Lambda }$, the following inequalities (8) hold; accordingly, the equilibria $z^{\#}$ of system (1) are mean-square exponential stability, and the converge rate is $\lambda -\eta$,

$$\begin{gathered} {\vartheta }_m(\lambda )=\left[(\lambda -w_m)+0.5{\displaystyle \sum _{r=1}^n(l_r^f|a_{mr}|+l_r^g|b_{mr}|+l_r^h|c_{mr}|)}\right]{\zeta }_m+ \\ {\displaystyle \sum _{r=1}^n[(0.5l_r^f|a_{mr}|+{\pi }_{mr})+(0.5l_r^g|b_{mr}|+{\pi }_{mr}^{(\tau )})e^{\lambda \tau }}+0.5l_r^h|c_{mr}|{\mu }_{mr}(\lambda )]{\zeta }_r<0. \end{gathered}$$

(8)

Proof. 

Selecting the candidate vector the Lyapunov function is as follows:

$$V_m(t)=e^{\lambda t}{|q_m(t)|}^2, m\in \mathsf{\Lambda }.$$

When $t\in (t_0,t_1)$, according to Lemmas 1 and 2, we obtain:

$$\begin{array}{l}d{V}_m\left(t\right)\\ =\lambda e^{\lambda t}{\left|q_m\left(t\right)\right|}^2+e^{\lambda t}d{\left|q_m\left(t\right)\right|}^2\\ =\lambda e^{\lambda t}{\left|q_m\left(t\right)\right|}^2+e^{\lambda t}Re\left(\overline{q_m\left(t\right)}\mathrm{d}{q}_m(t)\right)\\ =\lambda e^{\lambda t}{\left|q_m\left(t\right)\right|}^2+e^{\lambda t}\left\{-{\omega }_m{\left|q_m\left(t\right)\right|}^2+\mathrm{R}\mathrm{e}\left(\overline{q_m\left(t\right)}{\displaystyle \sum _{r=1}^n}\left[a_{mr}{\tilde{f}}_r\left(q_r\left(t\right)\right)+b_{mr}{\tilde{g}}_r\left(q_r\left(t-{\tau }_{mr}\left(t\right)\right)\right)+\right.\right.\right.\\ \left.\left.\left.c_{mr}{\int }_{-\mathrm{\infty }}^t{\theta }_{mr}(t-s){\tilde{h}}_r\left(q_r\left(s\right)\right)\mathrm{d}s\right]\right)+{\displaystyle \sum _{r=1}^n}{\left({\tilde{p}}_{mr}\left(q_r\left(t\right),q_r\left(t-{\tau }_{mr}\left(t\right)\right)\right)\right)}^2\mathrm{d}{w}_r(t)\right\}dt+\\ e^{\lambda t}Re\left(\overline{q_m\left(t\right)}{\displaystyle \sum _{j=1}^n}{\tilde{p}}_{mr}\left(q_r\left(t\right),q_r\left(t-{\tau }_{mr}\left(t\right)\right)\right)\mathrm{d}{w}_r(t)\right),m\in \Lambda \\ \end{array}$$

(9)

In what follows, the $It\widehat{o}$ differential form of (9) will be modified to the $It\widehat{o}$ integral form. Based on the theory of stochastic differential equations [59], after integrating both sides of (9) are from t to $t+\Delta t$ for any $\mathsf{\Delta }t>0$, and calculating the mathematical expectation of (9), we derive:

$$\begin{array}{l}\mathrm{E}(V_m(t+\mathsf{\Delta }t))-\mathrm{E}(V_m(t))\\ =e^{\lambda t}\mathrm{E}({\displaystyle {\int }_t^{t+\mathsf{\Delta }t}\{\lambda |q_m(t){|}^2-{\omega }_m|q_m(t){|}^2+\mathrm{Re}(\overline{q_m(t)}{\displaystyle \sum _{r=1}^n[a_{mr}{\tilde{f}}_r(q_r(t))+b_{mr}{\tilde{g}}_r(q_r(t-{\tau }_{mr}(t)))+}}\\ c_{mr}{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s){\tilde{h}}_r(q_r(s))}\mathrm{d}s])+{\left({\displaystyle \sum _{r=1}^n{\tilde{p}}_{mr}(q_r(t),q_r(t-{\tau }_{mr}(t)))}\mathrm{d}{w}_r(t)\right)}^2\}\mathrm{d}t)+\\ e^{\lambda t}\mathrm{E}\left({\displaystyle {\int }_t^{t+\mathsf{\Delta }t}\mathrm{Re}\left(\overline{q_m(t)}{\displaystyle \sum _{r=1}^n{\tilde{p}}_{mr}(q_r(t),q_r(t-{\tau }_{mr}(t)))}\right)\mathrm{d}{w}_r(t)}\right), m\in \mathsf{\Lambda }.\end{array}$$

(10)

Given that the mathematical expectation of $It\widehat{o}$ process has continuity, Equation (10) can be determined using the properties of $It\widehat{o}$ integral [59]. Furthermore, by using Holder inequality and Young inequality, we have:

$$\begin{array}{l}{\mathrm{D}}^{+}\mathrm{E}\left(V_m(t)\right)\\ =e^{\lambda t}\mathrm{E}(\left(\lambda -{\omega }_m\right)|q_m(t){|}^2+\mathrm{Re}(\overline{q_m(t)}{\displaystyle \sum _{r=1}^n[a_{mr}{\tilde{f}}_r(q_r(t))+b_{mr}{\tilde{g}}_r(q_r(t-{\tau }_{mr}(t)))+}\\ c_{mr}{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s){\tilde{h}}_r(q_r(s))}\mathrm{d}s])+{\displaystyle \sum _{r=1}^n{\left({\tilde{p}}_{mr}(q_r(t),q_r(t-{\tau }_{mr}(t)))\right)}^2})\\ \le e^{\lambda t}\mathrm{E}((\lambda -{\omega }_m)|q_m(t){|}^2+|q_m(t)|{\displaystyle \sum _{r=1}^n[l_r^f|a_{mr}||q_r(t)|+l_r^g|b_{mr}||q_r(t-{\tau }_{mr}(t))|+}\\ l_r^h|c_{mr}|{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s)|q_r(s)|}\mathrm{d}s])+{\displaystyle \sum _{r=1}^n\left({\pi }_{mr}|q_r(t){|}^2+{\pi }_{mr}^{(\tau )}|q_r(t-{\tau }_{mr}(t)){|}^2\right)})\\ \le e^{\lambda t}\mathrm{E}(\left(\lambda -{\omega }_m\right)|q_m(t){|}^2+0.5{\displaystyle \sum _{r=1}^n[l_r^f|a_{mr}|\left(|q_m(t){|}^2+|q_r(t){|}^2\right)+}\\ l_r^g|b_{mr}|\left(|q_m(t){|}^2+|q_r(t-{\tau }_{mr}(t)){|}^2\right)+\\ l_r^h|c_{mr}|{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s)\left(|q_m(t){|}^2+|q_r(s){|}^2\right)}\mathrm{d}s])+\\ {\displaystyle \sum _{r=1}^n\left({\pi }_{mr}|q_r(t){|}^2+{\pi }_{mr}^{(\tau )}|q_r(t-{\tau }_{mr}(t)){|}^2\right)})\\ \le \left[\left(\lambda -{\omega }_m\right)+0.5{\displaystyle \sum _{r=1}^n\left(l_r^f|a_{mr}|+l_r^g|b_{mr}|+l_r^h|c_{mr}|\right)}\right]\mathrm{E}\left(e^{\lambda t}{|q_m(t)|}^2\right)+\\ {\displaystyle \sum _{r=1}^n[\left(0.5l_r^f|a_{mr}|+{\pi }_{mr}\right)\mathrm{E}\left(e^{\lambda t}{|q_r(t)|}^2\right)+\left(0.5l_r^g|b_{mr}|+{\pi }_{mr}^{(\tau )}\right)e^{\lambda \tau }\mathrm{E}\left(e^{\lambda \left(t-{\tau }_{mr}(t)\right)}{|q_r(t-{\tau }_{mr}(t))|}^2\right)}+\\ 0.5l_r^h|c_{mr}|{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s)e^{\lambda \left(t-s\right)}\mathrm{E}\left(e^{\lambda s}{|q_r(s)|}^2\right)}\mathrm{d}s]\\ \le \left[\left(\lambda -{\omega }_m\right)+0.5{\displaystyle \sum _{r=1}^n\left(l_r^f|a_{mr}|+l_r^g|b_{mr}|+l_r^h|c_{mr}|\right)}\right]\mathrm{E}\left(V_m(t)\right)+\\ {\displaystyle \sum _{r=1}^n[\left(0.5l_r^f|a_{mr}|+{\pi }_{mr}\right)\mathrm{E}\left(V_r(t)\right)+\left(0.5l_r^g|b_{mr}|+{\pi }_{mr}^{(\tau )}\right)e^{\lambda \tau }\mathrm{E}\left(V_r(t-{\tau }_{mr}(t))\right)}+\\ 0.5l_r^h|c_{mr}|{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s)e^{\lambda \left(t-s\right)}\mathrm{E}\left(V_r\left(s\right)\right)}\mathrm{d}s], t_0\le t<t_1, m\in \mathsf{\Lambda }.\end{array}$$

(11)

Define curve

$$\zeta =\{\omega (\chi ):{\omega }_m={\zeta }_m\chi , \chi >0, m\in \mathsf{\Lambda }\}$$

and set

$$\Omega (\omega )=\{v:0\le v\le \omega , \omega \in \xi \}.$$

When $\chi >{\chi }^{\prime }$ we have $\Omega (\omega (\chi ))\supset \Omega (\omega ({\chi }^{\prime }))$. Let ${\zeta }_{\max }={\max }_{m\in \mathsf{\Lambda }}\left\{{\zeta }_m\right\}$, ${\zeta }_{\min }={\min }_{m\in \mathsf{\Lambda }}\left\{{\zeta }_m\right\}$, ${\chi }_0=\frac{(1+\delta )\mathsf{{\rm E}}(|{\psi }_m(s){|}^2)}{{\zeta }_{\min }}$, where $\delta >0$ is a constant. Then

$$\left\{V(s):V(s)=\mathrm{E}\left(e^{\lambda s}||\psi (s)|{|}^2\right), -\infty <s\le 0\right\}\subset \Omega ({\omega }_0({\chi }_0)),$$

which means $\mathsf{{\rm E}}\left(V_m(s)\right)=\mathsf{{\rm E}}\left(e^{\lambda s}{|{\psi }_m(s)|}^2\right)<{\zeta }_m{\chi }_0$, $-\infty <s\le 0$, $m\in \mathsf{\Lambda }$.

Furthermore, inequalities $\mathrm{E}\left(V_m(t)\right)<{\zeta }_m{\chi }_0$ are necessarily valid, here $m\in \mathsf{\Lambda }$ and $0<t<t_1$. If $\mathrm{E}\left(V_m(t)\right)<{\zeta }_m{\chi }_0$ does not hold, then there will be some $m$ at the time instant $t*$$(0<t*<t_1)$, such that $\mathrm{E}\left(V_m(t*)\right)={\zeta }_m{\chi }_0$, $D^{+}\mathrm{E}\left(V_m(t*)\right)\ge 0$ and $\mathrm{E}\left(V_r(t*)\right)\le {\zeta }_r{\chi }_0$, here $m, r\in \mathsf{\Lambda }$. By substituting these assumptions in (9) and considering (7), the following is obtained:

$$\begin{array}{l}{\mathrm{D}}^{+}\mathrm{E}\left(V_m(t*)\right)\\ \le \left[\left(\lambda -{\omega }_m\right)+0.5{\displaystyle \sum _{r=1}^n\left(l_r^f|a_{mr}|+l_r^g|b_{mr}|+l_r^h|c_{mr}|\right)}\right]\mathrm{E}\left(V_m(t*)\right)+\\ {\displaystyle \sum _{r=1}^n[\left(0.5l_r^f|a_{mr}|+{\pi }_{mr}\right)\mathrm{E}\left(V_r(t*)\right)+\left(0.5l_r^g|b_{mr}|+{\pi }_{mr}^{(\tau )}\right)e^{\lambda \tau }\mathrm{E}\left(V_r(t*-{\tau }_{mr}(t*))\right)}+\\ 0.5l_r^h|c_{mr}|{\displaystyle {\int }_{-\infty }^{t*}{\theta }_{mr}(t*-s)e^{\lambda \left(t*-s\right)}\mathrm{E}\left(V_r(s)\right)}\mathrm{d}s]\\ \le \left[\left(\lambda -{\omega }_m\right)+0.5{\displaystyle \sum _{r=1}^n\left(l_r^f|a_{mr}|+l_r^g|b_{mr}|+l_r^h|c_{mr}|\right)}\right]{\zeta }_m{\chi }_0+\\ {\displaystyle \sum _{r=1}^n[\left(0.5l_r^f|a_{mr}|+{\pi }_{mr}\right)+\left(0.5l_r^g|b_{mr}|+{\pi }_{mr}^{(\tau )}\right)e^{\lambda \tau }}+0.5l_r^h|c_{mr}|{\mu }_{mr}(\lambda )]{\zeta }_r{\chi }_0\\ <0, m\in \mathsf{\Lambda }.\end{array}$$

The preceding analysis ${\mathrm{D}}^{+}\mathrm{E}\left(V_m(t*)\right)<0$ disproves the assumption ${\mathrm{D}}^{+}\mathrm{E}\left(V_m(t*)\right)\ge 0$. Furthermore, we have $\mathrm{E}\left(V_m(t)\right)<{\zeta }_m{\chi }_0$, which equals to

$$\mathrm{E}\left(|q_m(t){|}^2\right)<e^{-\lambda t}{\zeta }_m{\chi }_0, m\in \mathsf{\Lambda }, 0<t<t_1.$$

(12)

Next, according to mathematical induction, we prove the follow-up inequalities (13) hold:

$$\mathrm{E}(|q_m(t){|}^2)<{\eta }_{ }^{(0)}{\eta }_{ }^{(1)}{\eta }_{ }^{(2)}\dots {\eta }_{ }^{(k-1)}{\zeta }_m{\chi }_0{e}^{-\lambda t}, m\in \mathsf{\Lambda }, t_{k-1}\le t<t_k, k\in \mathbb{N}.$$

(13)

When $k=1$, it follows from (12) that $\mathrm{E}(|q_m(t){|}^2)<{\eta }_m^{(0)}{\zeta }_m{\chi }_0{e}^{-\lambda t}$ holds for $m\in \mathsf{\Lambda }$, $t_0\le t<t_1$, where ${\eta }^{(0)}=1$.

Suppose the follow-up inequalities (14) are true, then

$$\mathrm{E}(|q_m(t){|}^2)<{\eta }_{ }^{(0)}{\eta }_{ }^{(1)}{\eta }_{ }^{(2)}\dots {\eta }_{ }^{(v-1)}{\zeta }_m{\chi }_0{e}^{-\lambda t}, m\in \mathsf{\Lambda }, t_{v-1}\le t<t_v \mathrm{and} v=1,2,\dots ,k.$$

(14)

When $t=t_k$ using Assumption 2, we derive:

$$|q_m(t_k^{+}){|}^2=|q_m(t_k^{-})+I_m^k(q_m(t_k^{-})){|}^2\le {\eta }_m^{(k)}|q_m(t_k^{-}){|}^2\le {\eta }_{ }^{(k)}|q_m(t_k^{-}){|}^2.$$

(15)

Because ${\eta }_{ }^{(k)}\ge 1$ inequalities (14) are modified as follows:

$$\mathrm{E}(|q_m(t){|}^2)<{\eta }_{ }^{(0)}{\eta }_{ }^{(1)}{\eta }_{ }^{(2)}\dots {\eta }_{ }^{(k-1)}{\zeta }_m{\chi }_0{e}^{-\lambda t}, m\in \mathsf{\Lambda }, t_{k-1}-\tau \le t\le t_k, k\in \mathbb{N}.$$

(16)

Therefore, it is also possible that the following inequalities exist:

$$\mathrm{E}(|q_m(t){|}^2)<{\eta }_{ }^{(0)}{\eta }_{ }^{(1)}{\eta }_{ }^{(2)}\dots {\eta }_{ }^{(k-1)}{\eta }_{ }^{(k)}{\zeta }_m{\chi }_0{e}^{-\lambda t}, m\in \mathsf{\Lambda }, t_k\le t<t_{k+1}, k\in \mathbb{N}.$$

(17)

The proof by contradiction method will be used to demonstrate that inequalities (17) are true. If inequalities (17) do not hold, there will be some $m^{\prime }$ and time instant $t^{\prime }$, in which ${\mathrm{D}}^{+}\mathrm{E}\left(V_{m^{\prime }}(t^{\prime })\right)\ge 0$, equality (18) and inequalities (19) hold,

$$\mathrm{E}(|q_{m^{\prime }}(t^{\prime }){|}^2)={\eta }_{ }^{(0)}{\eta }_{ }^{(1)}{\eta }_{ }^{(2)}\dots {\eta }_{ }^{(k-1)}{\eta }_{ }^{(k)}{\zeta }_{m^{\prime }}{\chi }_0{e}^{-\lambda t^{\prime }}, t_k\le t^{\prime }<t_{k+1},$$

(18)

$$\mathrm{E}(|q_r(t){|}^2)<{\eta }_{ }^{(0)}{\eta }_{ }^{(1)}{\eta }_{ }^{(2)}\dots {\eta }_{ }^{(k-1)}{\eta }_{ }^{(k)}{\zeta }_r{\chi }_0{e}^{-\lambda t}, t_k-\tau <t\le t^{\prime }, r\in \mathsf{\Lambda }.$$

(19)

Substituting equality (18) and inequalities (19) into inequalities (11), and considering inequalities (9), we obtain

$$\begin{array}{l}{\mathrm{D}}^{+}\mathrm{E}(V_{m^{\prime }}(t^{\prime }))\\ \le \left[(\lambda -{\omega }_{m^{\prime }})+0.5{\displaystyle \sum _{r=1}^n(l_r^f|a_{m^{\prime }r}|+l_r^g|b_{m^{\prime }r}|+l_r^h|c_{m^{\prime }r}|)}\right]\mathrm{E}(V_{m^{\prime }}(t^{\prime }))+\\ {\displaystyle \sum _{r=1}^n[(0.5l_r^f|a_{m^{\prime }r}|+{\pi }_{m^{\prime }r})\mathrm{E}(V_r(t^{\prime }))+(0.5l_r^g|b_{m^{\prime }r}|+{\pi }_{m^{\prime }r}^{(\tau )})e^{\lambda \tau }\mathrm{E}(V_r(t^{\prime }-{\tau }_{m^{\prime }r}(t^{\prime })))}+\\ 0.5l_r^h|c_{m^{\prime }r}|{\displaystyle {\int }_{-\infty }^{t^{\prime }}{\theta }_{m^{\prime }r}(t^{\prime }-s)e^{\lambda \left(t^{\prime }-s\right)}\mathrm{E}(V_r(s))}\mathrm{d}s]\\ \le \left[(\lambda -{\omega }_m)+0.5{\displaystyle \sum _{r=1}^n(l_r^f|a_{m^{\prime }r}|+l_r^g|b_{m^{\prime }r}|+l_r^h|c_{m^{\prime }r}|)}\right]{\eta }_{ }^{(0)}{\eta }_{ }^{(1)}{\eta }_{ }^{(2)}\dots {\eta }_{ }^{(k-1)}{\eta }_{ }^{(k)}{\zeta }_{m^{\prime }}{\chi }_0{e}^{-\lambda t^{\prime }}+\\ {\displaystyle \sum _{r=1}^n[(0.5l_r^f|a_{m^{\prime }r}|+{\pi }_{m^{\prime }r})+(0.5l_r^g|b_{m^{\prime }r}|+{\pi }_{m^{\prime }r}^{(\tau )})e^{\lambda \tau }}+0.5l_r^h|c_{m^{\prime }r}|{\mu }_{m^{\prime }r}(\lambda )]\times \\ {\eta }_{ }^{(0)}{\eta }_{ }^{(1)}{\eta }_{ }^{(2)}\dots {\eta }_{ }^{(k-1)}{\eta }_{ }^{(k)}{\zeta }_r{\chi }_0{e}^{-\lambda t^{\prime }}\\ <0, m^{\prime }\in \mathsf{\Lambda }, t_k\le t^{\prime }<t_{k+1}, k\in \mathbb{N}.\end{array}$$

(20)

This conclusion ${\mathrm{D}}^{+}\mathrm{E}(V_{m^{\prime }}(t^{\prime }))<0$ is inconsistent with the assumption ${\mathrm{D}}^{+}{V}_{m^{\prime }}(t^{\prime })\ge 0$. Therefore, we can deduce that inequalities (17) are always true for all $m\in \mathsf{\Lambda }$, $t_k\le t<t_{k+1}$, $k\in \mathbb{N}$.

According to the mathematical induction method, we derive that the inequalities (21) hold,

$$\mathrm{E}(|q_m(t){|}^2)<{\eta }_{ }^{(0)}{\eta }_{ }^{(1)}{\eta }_{ }^{(2)}\dots {\eta }_{ }^{(k-1)}{\zeta }_m{\chi }_0{e}^{-\lambda t}, m\in \mathsf{\Lambda }, t_{k-1}\le t<t_k, \mathrm{and} k\in \mathbb{N}.$$

(21)

It follows from the condition $\eta ={\lim }_{k\to \infty }\left(\sup \frac{\ln {\eta }^{(k)}}{t_k-t_{k-1}}\right)$ in Theorem 1 that ${\eta }^{(k)}\le e^{\eta (t_k-t_{k-1})}$, $k\in \mathbb{N}$. Substituting it into inequalities (21), we obtain

$$\mathrm{E}(|q_m(t){|}^2)$$

$$<{\eta }_{ }^{(0)}{\eta }_{ }^{(1)}{\eta }_{ }^{(2)}\dots {\eta }_{ }^{(k-1)}{\zeta }_m{\chi }_0{e}^{-\lambda t}\le e^{\eta (t_1-t_0)}{e}^{\eta (t_2-t_1)}\cdots e^{\eta (t_{k-1}-t_{k-2})}{\zeta }_m{\chi }_0{e}^{-\lambda t}$$

$$\le {\zeta }_m{\chi }_0{e}^{-(\lambda -\eta )(t-t_0)}, t_{k-1}\le t<t_k, k\in \mathbb{N}.$$

Furthermore, for all $t_{k-1}\le t<t_k$, $m\in \mathsf{\Lambda }$, $k\in \mathbb{N}$, we obtain

$$\mathrm{E}(|q_m(t){|}^2)<\frac{(1+\delta )\mathrm{E}(|{\psi }_m(s){|}^2)}{{\zeta }_{\min }}{\zeta }_m{e}^{-(\lambda -\eta )(t-t_0)}=\Gamma \mathrm{E}(|{\psi }_m(s){|}^2)e^{-(\lambda -\eta )(t-t_0)},$$

which denotes $\mathrm{E}(||q(t)|{|}^2)<\frac{(1+\delta )\mathrm{E}(||\psi (s)|{|}^2)}{{\zeta }_{\min }}{\zeta }_m{e}^{-(\lambda -\eta )(t-t_0)}=\Gamma \mathrm{E}(||\psi (s)|{|}^2)e^{-(\lambda -\eta )(t-t_0)}$ hold, where $\Gamma =\frac{(1+\delta ){\zeta }_{\max }}{{\zeta }_{\min }}$.

Based on Definition 1, we can infer that the solution $q=0$ of system (7) demonstrates mean-square exponential stability. That implies that the equilibria $z^{\#}$ of system (1) are exponentially stable in the mean square. □

When only stochastic disturbances are evaluated in the model (1), and no impulsive disturbances are considered, the corresponding corollary can be derived as follows.

Corollary 1. 

Suppose that Assumptions 1, 2, and 4 hold. If for all $J\in {\mathbb{Q}}^n$ and $m\in \mathsf{\Lambda }$, the matrix $\Psi$ is an M-matrix, then the equilibria $z^{\#}$ of system (1) are mean-square exponential stability,

$${\Psi }_{mm}={\omega }_m-0.5{\displaystyle \sum _{r=1}^n(l_r^f|a_{mr}|+l_r^g|b_{mr}|+l_r^h|c_{mr}|)},$$

$${\Psi }_{mr}={\displaystyle \sum _{r=1}^n[(0.5l_r^f|a_{mr}|+{\pi }_{mr})+(0.5l_r^g|b_{mr}|+{\pi }_{mr}^{(\tau )})}+0.5l_r^h|c_{mr}|].$$

Proof. 

Because $\Psi$ is an M-matrix, it follows from Lemma 3 that there is a positive vector $\zeta$ such that all inequalities (22) hold,

$$\begin{array}{l}\left[-{\omega }_m+0.5{\displaystyle \sum _{r=1}^n(l_r^f|a_{mr}|+l_r^g|b_{mr}|+l_r^h|c_{mr}|)}\right]{\zeta }_m+\\ {\displaystyle \sum _{r=1}^n[(0.5l_r^f|a_{mr}|+{\pi }_{mr})+(0.5l_r^g|b_{mr}|+{\pi }_{mr}^{(\tau )})}+0.5l_r^h|c_{mr}|]{\zeta }_r<0, m\in \mathsf{\Lambda }.\end{array}$$

(22)

The functions related to the system (1) are constructed as follows:

$$\begin{array}{l}{F}_m(\epsilon )\\ =\left[(\epsilon -{\omega }_m)+0.5{\displaystyle \sum _{r=1}^n(l_r^f|a_{mr}|+l_r^g|b_{mr}|+l_r^h|c_{mr}|)}\right]{\zeta }_m+\\ {\displaystyle \sum _{r=1}^n[(0.5l_r^f|a_{mr}|+{\pi }_{mr})+(0.5l_r^g|b_{mr}|+{\pi }_{mr}^{(\tau )})e^{\lambda \epsilon }}+0.5l_r^h|c_{mr}|{\mu }_{mr}(\epsilon )]{\zeta }_r, m\in \mathsf{\Lambda }.\end{array}$$

(23)

It is evident that $F_m^{ }(\epsilon )$ is a continuous function on $\epsilon$, $F_m^{ }(\epsilon )<0$, $m\in \mathsf{\Lambda }$. Therefore, there is a constant $\lambda >0$ such that $F_m^{ }(\lambda )<0$ ($m\in \mathsf{\Lambda }$) holds. The left part of the proof for Corollary 1 can be obtained directly using the same procedure as that used for proving Theorem 1. □

When only impulsive disturbances are incorporated in model (1) and no stochastic disturbances are considered, the corresponding corollary is as follows.

Corollary 2. 

Suppose that Assumptions 1–3 hold. Let $\eta ={\lim }_{k\to \infty }\left(\sup \frac{\ln {\eta }^{(k)}}{t_k-t_{k-1}}\right)$, $k\in \mathbb{N}$. If there is a positive constant $\lambda >0$ with $\lambda >\eta >0$ and a positive number ${\zeta }_m$ such that for all $J\in {\mathbb{Q}}^n$ and $m\in \mathsf{\Lambda }$, the following inequalities hold, then the equilibria $z^{\#}$ of system (1) are exponentially stable, and the converge rate is $\lambda -\eta$,

$${\vartheta }_m(\lambda )=(\lambda -{\omega }_m){\zeta }_m+{\displaystyle \sum _{r=1}^n[l_r^f|a_{mr}|+l_r^g|b_{mr}|e^{\lambda \tau }}+l_r^h|c_{mr}|{\mu }_{mr}(\lambda )]{\zeta }_r<0.$$

(24)

Because Corollary 2 can be obtained directly by proving Theorem 1 in a similar manner, we skip the proof of Corollary 2.

When there is no disturbance in model (1), the corollary is as follows.

Corollary 3. 

It is supposed that Assumptions 1 and 2 hold. If the matrix $\Psi$ is an M-matrix, where ${\Psi }_{mm}={\omega }_m$, ${\Psi }_{mr}={\displaystyle \sum _{r=1}^n[(l_r^f|a_{mr}|+l_r^g|b_{mr}|}+l_r^h|c_{mr}|]$, $m,r\in \mathsf{\Lambda }$, then the equilibria $z^{\#}$ of system (1) show exponential stability for all $J\in {\mathbb{Q}}^n$.

Remark 4. 

Selecting an activation function is one of the most difficult problems in studying the dynamic behavior of NNs. Currently, boundedness and the Lipschitz conditions are the two most common assumptions. In studies [31,32,35,37,42,45,46,47,51] and [10], the activation functions are assumed to be in the real number field and quaternion field, respectively. In studies [30,33,36,38,39,41,43,44,48,49,53] and [14,15,16,17,18,19,20,21,22,26,29,34], the activation functions are assumed to meet the Lipschitz condition in the real field and quaternion field, respectively.

Table 2. Literature summary of assumption for quaternion activation functions and study method.

Types of Assumption	Decomposition Method	Non-Decomposition Method
Lipschitz condition	9, 11, 13, 14, 19–21, 26	8, 12, 15–18, 22
Boundedness	10, 24	none

summarizes the review of existing work. According to this table, the existing results that the activation function satisfies the boundedness are very rare among the research results on the dynamic behavior of QVNNs.

Table 2 illustrates that almost no study has been conducted on the dynamic behavior of QVNNs using the non-decomposition method under the assumption that the quaternion activation function satisfies boundedness. Subsequently, Assumption 5 is used when the activation function in model (1) satisfies boundedness.

Assumption 5. 

The neuron activation functions $f_m(.)$ with $f_m(0)=0$ in the system (1) satisfy the condition ${\underset{\_}{l}}_m^f\le \frac{f_m(z_m)-f_m(v_m)}{z_m-v_m}\le {\overline{l}}_m^f$ for any $z_m,v_m\in \mathbb{Q}$, $m\in \Lambda$, where the ${\underset{\_}{l}}_m^f$ and ${\overline{l}}_m^f$ are constants.

Remark 1 also applies to Assumption 5; therefore, we exclude it herein.

Remark 5. 

Note that the constants ${\underset{\_}{l}}_m^f$ and ${\overline{l}}_m^f$ in Assumption 5 are allowed to be positive, negative, or zero. If we take $l_m^f=\max \left\{|{\underset{\_}{l}}_m^f|, |{\overline{l}}_m^f|\right\}$, Assumption 5 implies that $|f_m(z_m)-f_m(v_m)|\le l_m^f|z_m-v_m|$. The previously employed Lipschitz condition can be considered as a unique case of Assumption 5. According to this, we can directly derive the stability criteria when the activation function satisfies Assumption 5. The corresponding stability criteria will not be replicated here.

4. Main Results of Synchronization

As stated earlier, the synchronization control problem of chaotic systems can be modified into a stability problem of synchronization error systems to be analyzed. In this section, we will implement the stability study method presented in Section 3 to assess the mean square exponential synchronization problem by developing a linear feedback controller for a class of chaotic QVNNs with mixed time delays, stochastic disturbances, and impulsive disturbances.

First, the models are provided for chaotic QVNNs based on the drive-response conception.

The drive systems of chaotic QVNNs are defined as follows:

$$\{\begin{array}{l}\mathrm{d}{x}_m(t)=\left\{-{\omega }_m{x}_m(t)+{\displaystyle \sum _{r=1}^n\left[a_{mr}{f}_r(x_r(t))+b_{mr}{g}_r(x_r(t-{\tau }_{mr}(t)))+c_{mr}{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s)h_r(x_r(s))}\mathrm{d}s\right]}\right\}\mathrm{d}t+\\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}{\displaystyle \sum _{r=1}^n{p}_{mr}(x_r(t),x_r(t-{\tau }_{mr}(t)))}\mathrm{d}{w}_r(t)+J_m(t), t\ne t_k, m\in \mathsf{\Lambda } , k\in \mathbb{N},\\ \mathsf{\Delta }{x}_m(t_k)=x_m(t_k^{+})-x_m(t_k^{-}), t=t_k, m\in \mathsf{\Lambda } , k\in \mathbb{N}.\end{array}$$

(25)

The response systems of chaotic QVNNs are described as follows:

$$\{\begin{array}{l}\mathrm{d}{y}_m(t)=\left\{-{\omega }_m{y}_m(t)+{\displaystyle \sum _{r=1}^n\left[a_{mr}{f}_r(y_r(t))+b_{mr}{g}_r(y_r(t-{\tau }_{mr}(t)))+c_{mr}{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s)h_r(y_r(s))}\mathrm{d}s\right]}\right\}\mathrm{d}t+\\ \hspace{1em}\hspace{1em}\hspace{1em} {\displaystyle \sum _{r=1}^n{p}_{m r}(y_r(t),y_r(t-{\tau }_{mr}(t)))}\mathrm{d}{w}_r(t) +J_m(t)-u_m(t), t\ne t_k, m\in \mathsf{\Lambda } , k\in \mathbb{N},\\ \mathsf{\Delta }{y}_m(t_k)=y_m(t_k^{+})-y_m(t_k^{-}), t=t_k, m\in \mathsf{\Lambda }.\end{array}$$

(26)

In systems (25) and (26), $x_m(t)$ and $y_m(t)$ denote the neuronal states of the drive system and the response system, respectively, and $u_m(t)$ denotes the driving signal or synchronous control signal of the response system. The definitions of other symbols and functions are entirely coherent with the relevant definitions in system (1), which are excluded here.

The synchronization error signal is characterized as ${\epsilon }_m(t)=x_m(t)-y_m(t)$, $m\in \mathsf{\Lambda }$. The synchronization errors of drive-response systems are expressed as (27),

$$\{\begin{array}{l}\mathrm{d}{\epsilon }_m(t)=\left\{-{\omega }_m{\epsilon }_m(t)+{\displaystyle \sum _{r=1}^n[a_{mr}{\widehat{f}}_r({\epsilon }_r(t))+b_{mr}{\widehat{g}}_r({\epsilon }_r(t-{\tau }_{mr}(t)))+c_{mr}{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s){\widehat{h}}_r({\epsilon }_r(s))}\mathrm{d}s]}\right\}\mathrm{d}t\\ \hspace{1em}\hspace{1em}\hspace{1em} + {\displaystyle \sum _{r=1}^n{\widehat{p}}_{mr}({\epsilon }_r(t),{\epsilon }_r(t-{\tau }_{mr}(t)))}\mathrm{d}{w}_r(t)+u_m(t), t\ne t_k, k\in \mathbb{N},\\ \mathsf{\Delta }{\epsilon }_m(t_k)={\epsilon }_m(t_k^{+})-{\epsilon }_m(t_k^{-}), t=t_k, k\in \mathbb{N},\end{array}$$

(27)

where:

$${\widehat{f}}_r({\epsilon }_r(t))=f_r(x_r(t))-f_r(y_r(t)), {\widehat{h}}_r({\epsilon }_r(s))=h_r(x_r(s))-h_r(y_r(s)),$$

$${\widehat{g}}_r({\epsilon }_r(t-{\tau }_{mr}(t)))=g_r(x_r(t-{\tau }_{mr}(t)))-g_r(y_r(t-{\tau }_{mr}(t))),$$

$${\widehat{p}}_{mr}({\epsilon }_r(t),{\epsilon }_r(t-{\tau }_{mr}(t)))=p_{mr}(x_r(t),x_r(t-{\tau }_{mr}(t)))-p_{mr}(y_r(t),y_r(t-{\tau }_{mr}(t))), m,r\in \mathsf{\Lambda }.$$

Here, we use the coupling of each neuron output signal as the synchronization control signal in the following form:

$$u_m(t)={\displaystyle \sum _{r=1}^n{\varpi }_{mr}\left[f_r(x_r(t))-f_r(y_r(t))\right]}, m\in \mathsf{\Lambda }.$$

(28)

Let:

$$x(t)={[x_1(t),x_2(t),\dots ,x_n(t)]}^{\mathrm{T}}, y(t)={[y_1(t),y_2(t),\dots ,y_n(t)]}^{\mathrm{T}}, u(t)={[u_1(t),u_2(t),\dots ,u_n(t)]}^{\mathrm{T}}.$$

The definition of mean-square exponential synchronization is stated below.

Definition 2. 

The response state $y(t)$ of system (26) is exponentially synchronized in the mean square sense with the drive state $x(t)$ of system (25) under the control signal $u(t)$ with the form of (28), if there are constants $\lambda >0$ and $\Gamma >0$ such that $E(||x(t)-y(t)|{|}^2)\le \Gamma {\sup }_{s\in (-\infty ,0]}E(||x(s)-y(s)|{|}^2)e^{-\lambda (t-t_0)}$ hold,

where

$$\mathrm{E}(||x(t)-y(t)|{|}^2)={[\mathrm{E}(|x_1(t)-y_1(t){|}^2),\mathrm{E}(|x_2(t)-y_2(t){|}^2),\dots ,\mathrm{E}(|x_n(t)-y_n(t){|}^2)]}^{\mathrm{T}}, t>t_0\ge 0,$$

$$\mathrm{E}(||x(s)-y(s)|{|}^2)={[\mathrm{E}(|x_1(s)-y_1(s){|}^2),\mathrm{E}(|x_2(s)-y_2(s){|}^2),\dots ,\mathrm{E}(|x_n(s)-y_n(s){|}^2)]}^{\mathrm{T}},$$

here $s\in (-\infty ,0]$.

By substituting (28) into (27), the error of the synchronization system is rewritten as follows:

$$\{\begin{array}{l}\mathrm{d}{\epsilon }_m(t)=\{-{\omega }_m{\epsilon }_m(t)+{\displaystyle \sum _{r=1}^n[(a_{mr}+{\varpi }_{mr}){\widehat{f}}_r({\epsilon }_r(t))+b_{mr}{\widehat{g}}_r({\epsilon }_r(t-{\tau }_{mr}(t)))+}\\ \hspace{1em}\hspace{1em}\hspace{1em} c_{mr}{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s){\widehat{h}}_r({\epsilon }_r(s))}\mathrm{d}s]\}\mathrm{d}t+ {\displaystyle \sum _{r=1}^n{\widehat{p}}_{mr}({\epsilon }_r(t),{\epsilon }_r(t-{\tau }_{mr}(t)))}\mathrm{d}{w}_r(t), t\ne t_k, k\in \mathbb{N},\\ \mathsf{\Delta }{\epsilon }_m(t_k)={\epsilon }_m(t_k^{+})-{\epsilon }_m(t_k^{-}), t=t_k, k\in \mathbb{N}.\end{array}$$

(29)

Remark 6. 

As an extended application of Theorem 1, assume that the delay parts of the error system (29) satisfy Assumption 1, the activation functions ${\widehat{f}}_r$, ${\widehat{g}}_r$ and ${\widehat{h}}_r$ satisfy Assumption 2, the impulsive parts $\mathsf{\Delta }{\epsilon }_m(t_k)$ satisfy Assumption 3, and the stochastic coupling functions ${\widehat{p}}_{mr}$ satisfy Assumption 4.

Theorem 2. 

Suppose that Assumptions 1–4 hold. Let $\eta ={\lim }_{k\to \infty }\left(\sup \frac{\ln {\eta }^{(k)}}{t_k-t_{k-1}}\right)$. If there is a positive constant $\lambda$ with $\lambda >\eta >0$ and a positive number ${\zeta }_m$ such that the follow-up inequalities (30) hold, then the response state $y(t)$ of system (26) can determine mean-square exponential synchronization with the drive state $x(t)$ of system (25) under the control signal $u(t)$ designed as (28), and the converge rate of synchronization error system is $\lambda -\eta$,

$$\begin{array}{l}{F}_m(\lambda )\\ =\left\{\left(\lambda -{\omega }_m\right)+0.5{\displaystyle \sum _{r=1}^n\left[l_r^f(|a_{mr}|+|{\varpi }_{mr}|)+l_r^g|b_{mr}|+l_r^h|c_{mr}|\right]}\right\}{\zeta }_m+\\ {\displaystyle \sum _{r=1}^n[\left(0.5l_r^f(|a_{mr}|+|{\varpi }_{mr}|)+{\pi }_{mr}\right)+\left(0.5l_r^g|b_{mr}|+{\pi }_{mr}^{(\tau )}\right)e^{\lambda \tau }}+0.5l_r^h|c_{mr}|{\mu }_{mr}(\lambda )]{\zeta }_r\\ <0, m\in \mathsf{\Lambda }.\end{array}$$

(30)

Proof. 

Selecting vector Lyapunov functions as follows:

$$V_m(t)=e^{\lambda t}{|{\epsilon }_m(t)|}^2, m\in \mathsf{\Lambda }.$$

When $t\in (t_0,t_1)$, according to Lemmas 1 and 2, we obtain:

$$\begin{array}{l}\mathrm{d}{V}_m(t)\\ =\mathrm{d}\left(e^{\lambda t}{|{\epsilon }_m(t)|}^2\right)\\ =\lambda e^{\lambda t}{|{\epsilon }_m(t)|}^2+e^{\lambda t}\{-{\omega }_m{|{\epsilon }_m(t)|}^2+\mathrm{Re}(\overline{{\epsilon }_m(t)}{\displaystyle \sum _{r=1}^n[(a_{mr}+{\varpi }_{mr}){\widehat{f}}_r({\epsilon }_r(t))+b_{mr}{\widehat{g}}_r({\epsilon }_r(t-{\tau }_{mr}(t)))+}\\ c_{mr}{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s){\widehat{h}}_r({\epsilon }_r(s))}\mathrm{d}s])+{\displaystyle \sum _{r=1}^n{\left({\widehat{p}}_{mr}({\epsilon }_r(t),{\epsilon }_r(t-{\tau }_{mr}(t)))\right)}^2}\mathrm{d}{w}_r(t)\}\mathrm{d}t+\\ e^{\lambda t}\mathrm{Re}\left(\overline{{\epsilon }_m(t)}{\displaystyle \sum _{r=1}^n{\widehat{p}}_{mr}({\epsilon }_r(t),{\epsilon }_r(t-{\tau }_{mr}(t)))}\mathrm{d}{w}_r(t)\right), m\in \mathsf{\Lambda }.\end{array}$$

(31)

Subsequently, the $It\widehat{o}$ differential form of (31) is transposed to the $It\widehat{o}$ integral form. Based on the theory of stochastic differential equations [59], after integrating both sides of (31) are from t to $t+\Delta t$ for any $\mathsf{\Delta }t>0$ and calculating the expectation of (31), we obtain:

$$\begin{array}{l}\mathrm{E}(V_m(t+\Delta t))-\mathrm{E}(V_m(t))\\ =e^{\lambda t}\mathrm{E}({\displaystyle {\int }_t^{t+\Delta t}\{(\lambda -{\omega }_m)|{\epsilon }_m(t){|}^2+\mathrm{Re}(\overline{{\epsilon }_m(t)}{\displaystyle \sum _{r=1}^n[(a_{mr}+{\varpi }_{mr}){\widehat{f}}_r({\epsilon }_r(t))+b_{mr}{\widehat{g}}_r({\epsilon }_r(t-{\tau }_{mr}(t)))+}}\\ c_{mr}{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s){\widehat{h}}_r({\epsilon }_r(s))}\mathrm{d}s])+{({\displaystyle \sum _{r=1}^n{\widehat{p}}_{mr}({\epsilon }_r(t),{\epsilon }_r(t-{\tau }_{mr}(t)))}\mathrm{d}{w}_r(t))}^2\}\mathrm{d}t+\\ e^{\lambda t}\mathrm{E}\left({\displaystyle {\int }_t^{t+\Delta t}\mathrm{Re}\left(\overline{{\epsilon }_m(t)}{\displaystyle \sum _{r=1}^n{\widehat{p}}_{mr}({\epsilon }_r(t),{\epsilon }_r(t-{\tau }_{mr}(t)))}\right)\mathrm{d}{w}_r(t)}\right),m\in \mathsf{\Lambda }.\end{array}$$

(32)

Given that the math expectation of $It\widehat{o}$ process has continuity, Equation (32) can be determined using the properties of $It\widehat{o}$ integral [59]. Furthermore, by using Holder inequality and Young inequality, we have:

$$\begin{array}{l}{\mathrm{D}}^{+}\mathrm{E}\left(V_m(t)\right)\\ =e^{\lambda t}\mathrm{E}((\lambda -{\omega }_m)|{\epsilon }_m(t){|}^2+\mathrm{Re}(\overline{{\epsilon }_m(t)}{\displaystyle \sum _{r=1}^n[(a_{mr}+{\varpi }_{mr}){\widehat{f}}_r({\epsilon }_r(t))+b_{mr}{\widehat{g}}_r({\epsilon }_r(t-{\tau }_{mr}(t)))+}\\ c_{mr}{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s){\widehat{h}}_r({\epsilon }_r(s))}\mathrm{d}s])+{\left({\displaystyle \sum _{r=1}^n{\widehat{p}}_{mr}({\epsilon }_r(t),{\epsilon }_r(t-{\tau }_{mr}(t)))}\mathrm{d}{w}_r(t)\right)}^2)\\ \le \left\{\left(\lambda -{\omega }_m\right)+0.5{\displaystyle \sum _{r=1}^n\left[l_r^f\left(|a_{mr}|+|{\varpi }_{mr}|\right)+l_r^g|b_{mr}|+l_r^h|c_{mr}|\right]}\right\}\mathrm{E}\left(V_m(t)\right)+\\ {\displaystyle \sum _{r=1}^n[\left(0.5l_r^f\left(|a_{mr}|+|{\varpi }_{mr}|\right)+{\pi }_{mr}\right)\mathrm{E}\left(V_r(t)\right)+\left(0.5l_r^g|b_{mr}|+{\pi }_{mr}^{(\tau )}\right)e^{\lambda \tau }\mathrm{E}\left(V_r(t-{\tau }_{mr}(t))\right)}+\\ 0.5l_r^h|c_{mr}|{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s)e^{\lambda \left(t-s\right)}\mathrm{E}\left(V_r(s)\right)}\mathrm{d}s], t_0\le t<t_1, m\in \mathsf{\Lambda }.\end{array}$$

(33)

The other part of the proof of Theorem 2 conforms with Theorem 1; therefore, it is omitted herein.

Furthermore, we can determine that the follow-up conclusion covers all $t_{k-1}\le t<t_k$, $m\in \mathsf{\Lambda }$, $k\in \mathbb{N}$,

$$\mathrm{E}(||\epsilon (t)|{|}^2)<\frac{(1+\delta )\mathrm{E}\left(||\epsilon (s)|{|}^2\right)}{{\zeta }_{\min }}{\zeta }_m{e}^{-(\lambda -\eta )(t-t_0)}=\Gamma \mathrm{E}\left(||\epsilon (s)|{|}^2\right)e^{-(\lambda -\eta )(t-t_0)}, m\in \mathsf{\Lambda },$$

which implies

$$\mathrm{E}(||x(t)-y(t)|{|}^2)<\frac{(1+\delta )\mathrm{E}\left(||x(s)-y(s)|{|}^2\right)}{{\zeta }_{\min }}{\zeta }_m{e}^{-(\lambda -\eta )(t-t_0)}=\Gamma \mathrm{E}\left(||x(s)-y(s)|{|}^2\right)e^{-(\lambda -\eta )(t-t_0)},$$

where $\Gamma =\frac{(1+\delta ){\zeta }_{\max }}{{\zeta }_{\min }}$ and $\delta >0$.

Definition 2 infers that the zero solution $\epsilon (t)=0$ of the error system (29) demonstrates mean-square exponential stability, which indicates that the neuron $\mathrm{E}(x(t))$ of drive system (25) can be exponentially synchronized in the mean-square sense with the neuron $\mathrm{E}(y(t))$ of response system (26). □

Because Corollary 1–3 regarding the stability discussion in section three applies to the synchronization discussion as well, it will not be restated here.

Notably, most researchers [6,7,8,10,11,13,14,15,16,17,18,19,31,32,33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,48,49,50,51,52,53,54,55,57] have used the scalar Lyapunov function method for the stability analysis of the synchronization error of chaotic NNs. By observing the scalar Lyapunov functions constructed in literature [6,7,8,10,11,13,14,15,16,17,18,19,31,32,33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,48,49,50,51,52,53,54,55,57], it is clear that the form of the scalar Lyapunov function is rather complicated. In addition, when the system being analyzed is an infinite dimensional system, the convergence of the scalar Lyapunov function must be analyzed first, which leads to a complicated research process. The vector Lyapunov function method with simple forms is employed herein to investigate infinite dimensional error systems (27) in a manner that avoids the convergence analysis of this function.

5. Examples

5.1. Example 1

Consider the following system:

$$\{\begin{array}{l}\mathrm{d}{z}_m(t)=\{-9z_m(t)+{\displaystyle \sum _{r=1}^2[a_{mr}\tanh (z_r(t))+b_{mr}\tanh (z_r(t-{\tau }_{mr}(t)))+}\\ \hspace{1em}\hspace{1em}\hspace{1em} c_{mr}{\displaystyle {\int }_{-\infty }^t\exp (-(t-s))\tanh (z_r(s))}\mathrm{d}s]\}\mathrm{d}t+\\ \hspace{1em}\hspace{1em}\hspace{1em} {\displaystyle \sum _{r=1}^2{p}_{mr}(z_r(t),z_r(t-{\tau }_{mr}(t)))}\mathrm{d}{w}_r(t), t\ne t_k, m,r=1,2,\\ z_m(t_k^{ })=z_m(t_k^{+})-z_m(t_k^{-})=z_m(t_k^{-}), t_k=2k, k\in \mathbb{N}.\end{array}$$

(34)

It is assumed that the quaternion matrices are as follows:

$$A=\left[\begin{array}{cc}0.1-0.8i-0.6\iota +0.4\kappa & 0.8-0.2i-0.7\iota +0.1\kappa \\ 0.9-0.3i+0.6\iota -0.7\kappa & 0.1-0.6i+0.5\iota -0.2\kappa \end{array}\right],$$

$$B=\left[\begin{array}{cc}0.6-0.4i-0.1\iota +0.3\kappa & 0.2-0.5i-0.7\iota +0.1\kappa \\ 0.3-0.6i+0.5\iota +0.9\kappa & 0.7-0.4i+0.2\iota +0.1\kappa \end{array}\right],$$

$$C=\left[\begin{array}{cc}0.7-0.2i+0.7\iota -0.5\kappa & 0.4-0.6i+0.4\iota +0.2\kappa \\ 0.1-0.3i-0.6\iota +0.2\kappa & 0.2-0.3i-0.4\iota +0.7\kappa \end{array}\right].$$

Here, let the delays in the model (34) be ${\tau }_{1r}(t)=0.3+0.1\sin t$, ${\tau }_{2r}(t)=0.2+0.1\cos t$, $r=1,2$, and $t\ge 0$. Let $\lambda =0.5$ and $\xi ={[1, 1]}^{\mathrm{T}}$.

Let:

$$p_{1r}(z_{1r}(t),z_{1r}(t-{\tau }_{1r}(t)))=0.3z_{1r}(t)+0.15z_{1r}(t-{\tau }_{1r}(t)),$$

$$p_{2r}(z_{2r}(t),z_{2r}(t-{\tau }_{2r}(t)))=0.2z_{2r}(t)-0.03z_{2r}(t-{\tau }_{2r}(t)).$$

After computation, we acquire $L^f=L^g=L^h=\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(1, 1)$, $\pi =\left[0.18, 0.18; 0.08, 0.08\right]$. ${\pi }^{(\tau )}=\left[0.045, 0.045; 0.0018, 0.0018\right]$, $\eta =0.35$, $\tau =0.4$.

After substituting the above-assumed conditions into (8), we acquire ${\vartheta }_1=-1.08<0$, ${\vartheta }_1=-1.54<0$, which implies that Theorem 1 holds. Based on Theorem 1, the zero solution ${(0,0)}^{\mathrm{T}}$ of the model (34) shows mean-square exponential stability, and the convergence rate is $\lambda -\eta =0.15$.

Furthermore, we estimate that the matrix outlined in Corollary 3 is $\Psi =\left[-9, 5.82; 5.79, -9\right]$. Evidently, the matrix $\Psi$ is an M-matrix, and thus Corollary 3 holds.

As the external inputs in the model (34) are $J_1=J_2=0$, it follows from Theorem 1 that the equilibria of the system (34) are zero with mean-square exponential stability. The state trajectories of $z_1^{ }(t)$ and $z_2^{ }(t)$ in model (34) without and with disturbances are depicted in Figure 1 and Figure 2 with the same initial conditions ${\psi }_1(s)=-1.7+1.6i-2.6\iota +2.2\kappa$ and ${\psi }_2(s)=-1.1+1.9i+$$2.7\iota -2.8\kappa$, respectively. Figure 1 shows that the equilibria of the system (34) without disturbances are zero solution and unique, and that these equilibria are convergent. When there are impulsive disturbances and stochastic disturbances in the system (34), it can be seen from Figure 2 that the equilibria of system (34) converge to zero, but the convergence speed becomes slower relative to that shown in Figure 1. The simulation results prove the correctness of Theorem 1.

5.2. Example 2

Considering the following chaotic QVNNs as a drive system:

$$\{\begin{array}{l}\mathrm{d}{x}_m(t)=\{-{\omega }_m{x}_m(t)+{\displaystyle \sum _{r=1}^n[a_{mr}\tanh (x_r(t))+b_{mr}\tanh (x_r(t-{\tau }_{mr}(t)))}\\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}+c_{mr}{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s)\tanh (x_r(s))}\mathrm{d}s]\}\mathrm{d}t\\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}+{\displaystyle \sum _{r=1}^n\left[0.01x_{mr}(t)+0.02x_{mr}(t-{\tau }_{mr}(t))\right]}\mathrm{d}{w}_r(t), t\ne t_k,\\ \Delta x_m(t_k^{+})=x_m(t_k^{+})-x_m(t_k^{-})=0.4x_m(t_k^{-}), t_k=3, 6, 9, 12\dots \end{array}$$

(35)

Considering the following chaotic QVNNs as a response system:

$$\{\begin{array}{l}\mathrm{d}{y}_m(t)=\{-{\omega }_m{y}_m(t)+{\displaystyle \sum _{r=1}^n[a_{mr}\tanh (y_r(t))+b_{mr}\tanh (y_r(t-{\tau }_{mr}(t)))+}\\ \hspace{1em}\hspace{1em}\hspace{1em} c_{mr}{\displaystyle {\int }_{-\infty }^t{\theta }_{mr}(t-s)\tanh (y_r(s))}\mathrm{d}s]\}\mathrm{d}t\\ \hspace{1em}\hspace{1em}\hspace{1em} +{\displaystyle \sum _{r=1}^n\left[0.01x_{mr}(t)+0.02x_{mr}(t-{\tau }_{mr}(t))\right]}\mathrm{d}{w}_r(t) +{\displaystyle \sum _{r=1}^n{\varpi }_{mr}\left[x_r(t)-y_r(t)\right]}, t\ne t_k,\\ \Delta y_m(t_k^{ })=y_m(t_k^{+})-y_m(t_k^{-})=0.4y_m(t_k^{-}), t_k=3, 6, 9, 12\dots \end{array}$$

(36)

Let $W=\mathrm{diag}(1.1, 1.5)$. The following are believed to be the interconnected matrices and control parameters:

$$A=\left[\begin{array}{l}- 1.32+0.06i-0.67\iota -0.06\kappa -0.03- 0.88i-0.42\iota -0.59\kappa \\ - 1.34- 0.51i-1.65\iota -1.57\kappa -0.52+0.51i-1.25\iota +0.54\kappa \end{array}\right],$$

$$B=\left[\begin{array}{l}0.18-0.04i+0.10\iota + 0.16\kappa 0.02-1.25i -1.17\iota +0.04\kappa \\ -0.26+0.40i-0.55\iota -0.51\kappa -0.02+0.65i-0.01\iota -1.11\kappa \end{array}\right],$$

$$C=\left[\begin{array}{l} 0.31-0.47i-0.63\iota -1.19\kappa -1.43-1.26i+0.70\iota -1.53\kappa \\ -0.25-1.08i-0.29\iota -0.74\kappa -1.26-1.40i+0.53\iota -0.02\kappa \end{array}\right],$$

$$\varpi =\left[\begin{array}{l}-0.6400+0.0045i+0.0010\iota +0.0060\kappa -0.0015-0.0320i-0.0335\iota +0.0125\kappa \\ -0.0479-0.0026i+0.0006\iota +0.0003\kappa -0.0105-0.0004i-0.0006\iota -0.0008\kappa \end{array}\right].$$

(37)

Furthermore, the other parameters of the systems (35) and (36) are consistent with those of the system (34).

After substituting the above-assumed conditions into inequality (8), we acquire ${\vartheta }_1=-0.979<0$ and ${\vartheta }_1=-1.458<0$. Based on verification, the conditions of Theorem 2 hold. It follows from Theorem 2 that the mean square exponential synchronization of the drive and response systems can be realized, and the convergence rate is $\lambda -\eta =0.388$.

The initial conditions of the drive system (35) and response system (36) are $x_1(s)=-2.5-1i+ 2.2\iota +1.8\kappa$, $x_2(s)=-1.9+1.1i+1.7\iota -2.4\kappa$, $y_1(s)=1.6+1.3i-1.8\iota -1.9\kappa$ and $y_2(s)=-2.7+ 1.9i+1.7\iota -2.2\kappa$, here $\mathrm{s}\in (-\infty ,0]$.

When there are no disturbances, the phase plane of drive system (35) shows that the states of system (35) are in chaos (Figure 3). The state trajectories of the response and drive systems with impulsive and stochastic disturbances as well as their error curves are depicted in Figure 4 and Figure 5, where the blue and red curves represent the state curves of the drive and response systems, respectively, and the black curves represent their synchronization error curves. The error curves shown in Figure 4 and Figure 5 show that the response system will eventually synchronize with the drive system, which further verifies the correctability of the synchronization determination conditions established in Theorem 2.

5.3. Example 3

Here, a representative Lena image is selected for testing, which is widely used in image processing–related research owing to its complete inclusion of flat areas, shadows, textures and other details. The original Lena image used is shown in Figure 6, which is a color image with dimensions of 128 pixels × 128 pixels. Because of the limitation of the research group’s computer processing power, this study only evaluates time-varying delays and not the infinite distributed delays when designing QVNNs to implement color image associative memory. The designed model is defined as follows:

$$\{\begin{array}{l}\frac{\mathrm{d}{z}_m(t)}{\mathrm{d}t}=\left\{-z_m(t)+{\displaystyle \sum _{r=1}^{16384}\left[a_{mr}{f}_r(z_r(t))+b_{mr}{g}_r(z_r(t-{\tau }_{mr}(t)))\right]}+J_m\right\}\mathrm{d}t+\\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}{\displaystyle \sum _{r=1}^{16384}{p}_{mr}\left(z_r(t),z_r(t-{\tau }_{mr}(t))\right)}\mathrm{d}{w}_r(t), t\ne t_k, k\in \mathbb{N},\\ \Delta z_m(t_k)=z_m(t_k^{+})-z_m(t_k^{-}), t=t_k, k\in \mathbb{N}.\end{array}$$

(38)

wherein: the activation functions are expressed as $f_r=g_r= |z_r+2|-|z_r+1|$. It is assumed that $W=\mathrm{diag}(-1,-1,\dots ,-1)$, and the interconnected matrices are:

$${\displaystyle a_{mr}=\{\begin{array}{l}0.63-0.17i-0.75\iota +0.87\kappa , m=r\\ 0+0i+0\iota -0\kappa , m\ne r\end{array}, b_{mr}=\{\begin{array}{l}-1.5+0.85i-0.3\iota -0.2\kappa , m=r\\ 0+0i+0\iota -0\kappa , m\ne r\end{array},}$$

where $m, r=1,2,\dots ,16384$.

As the associative memory process based on NNs is simulated on electronic devices, which is characterized by a complex structure as well as large amounts of data and calculation, the associative memory process will inevitably be affected by the impulsive disturbances of electronic devices, which will reduce the efficiency of the associative memory of the color image and even lead to erroneous associative memory results directly. Therefore, impulsive disturbances are considered in system (37).

Figure 6 comprises 128 pixels × 128 pixels = 16,384 pixels in total. When QVNNs are used for associative memory, each pixel corresponds to a neuron. The value of RGB three channels of each pixel is assigned to the three imaginary parts of the quaternion in turn, as the initial state of each neuron of QVNNs. For easy observation, we arbitrarily select the state curves of five neurons in model (37) to be shown, i.e., $z_{1802}(t)$$z_{3722}(t)$, $z_{5642}(t)$, $z_{7562}(t)$ and $z_{9482}(t)$. For subsequent comparison, Figure 7 shows the state curves of each part of the five neurons selected without considering any disturbances. As illustrated in Figure 7, the states of each neuron eventually shift to their equilibria, respectively. Notably, because the RGB requires only three parts of the quaternion, the three imaginary parts of the quaternion neuron are used in the simulation.

To compare the effects comprehensively, we conducted simulations for both cases of weak stochastic disturbances and strong stochastic disturbances. Here, it is assumed that $z_m(t_k^{+})=[1, 1.5]z_m(t_k^{-})$ where the intensity of the impulsive disturbance to each neuron is taken randomly in the interval of [1, 1.5], where $m=1,2,\dots ,16384$ and $t_k=\{3s,6s,\dots ,18s\}$.

Case 1. Weak stochastic disturbances:

Taking $p_{mr}(z_r(t),z_r(t-{\tau }_{mr}(t)))=0.12z_{mr}(t)+0.09z_{mr}(t-{\tau }_{mr}(t))$, $m,r=1,2,\dots ,16384$.

Case 2. Strong stochastic disturbances:

Taking $p_{mr}(z_r(t),z_r(t-{\tau }_{mr}(t)))=z_{mr}(t)+z_{mr}(t-{\tau }_{mr}(t))$, $m, r=1,2,\dots ,16384$.

The simulations are depicted in Figure 8, Figure 9, Figure 10 and Figure 11, where Figure 8 (under the weak stochastic disturbances) and Figure 10 (under the strong stochastic disturbances) show the state curves of five randomly selected neurons, and Figure 9 (under the weak stochastic disturbances) and Figure 11 (under the strong stochastic disturbances) depict the image association memory process based on the designed QVNNs implemented in (37).

For case 1 (weak interference) and case 2 (strong interference), several simulation experiments demonstrate that the average PSNR value of case 2 is 11.1% lower than that of case 1 when the simulation time is taken as 20 s and the simulation step size is 0.1 s. This indicates that the greater the intensity of the stochastic disturbance, the greater the impact of QVNNs on the image associative memory. These simulation results further validate the correctness of the conclusions of this paper.

Remark 7. 

The magnitude of stochastic disturbance strength exerts a direct effect on the undulation intensity of gray value data at each pixel in the image. The stronger the stochastic disturbance, the stronger the fluctuation of gray value data of each pixel in the image. However, it is difficult to directly observe the influence of stochastic disturbance strength on the quality of image restoration from the visual perspective of Figure 9 and Figure 11. To measure the effect of stochastic disturbance strength on the efficiency of image recovery during neural network associative memory, the peak signal noise ratio (PSNR), a frequently employed image evaluation metric, was selected to examine the image quality. The higher the value of PSNR, the closer the current image is to the original image.

Table 3. PSNR values of images in three cases.

Three Cases	PSNR Value (dB) at Different Times (s)
	t = 0	t = 1	t = 3	t = 5	t = 10	t = 20
case 0 (only impulsive disturbances)	8.662	17.348	34.720	52.091	95.521	182.380
case 1 (weak stochastic disturbances)	8.662	17.266	16.880	31.904	24.492	31.141
case 2 (strong stochastic disturbances)	8.662	16.302	15.969	21.953	20.914	21.573

demonstrates the recovery quality of the images at the corresponding moments in Figure 9 and Figure 11 for different stochastic disturbance strength through PSNR data. By comparing the PSNR data in Table 3, we conclude that: (1) the stability of the neural network in any of the three cases is not destroyed, and (2) the PSNR values in the strong stochastic disturbance case are lower than those in the weak stochastic disturbance case at the same moments.

6. Conclusions and Future Work

For a type of mixed-delay QVNNs with stochastic and impulsive disturbances, several conditions are proposed to determine the mean-square exponential stability of the system in light of the nondecomposition method. This research is being expanded to include the synchronization control of chaotic QVNNs using stochastic and impulsive disturbances. By designing a linear feedback controller and using the previously established stability analysis method, some sufficient conditions were obtained for realizing/determining the mean-square exponential synchronization of the drive–response system. The correctness and feasibility of the main results were validated using two numerical examples. Furthermore, the associative memory of color images by designing an appropriate QVNNs was achieved. Motivated by the studies [60,61,62], our approach in this paper will be applied to further discuss the stability of discrete QVNNs with variable coefficients and Markovian jumping parameters.