随着计算机网络的发展，网络安全越来越受到关注。在网络信息中，图像信息因其生动、形象、直观的特点而被广泛应用，成为网络安全的重点关注领域。针对图像信息在网络传输、存储过程中存在的安全隐患，各类图像加密算法不断被提出，特别是混沌系统因其具有对初值的敏感性、伪随机性、无周期性和混沌序列的遍历性成为图像加密的主要方法^[1-5]。

纵观在图像加密领域的混沌密码算法，仍是基于香农对密码设计提出的混淆和扩散原则，图像加密过程中像素位置置乱与像素值替代结合是最普遍的加密方式。图像像素位置的置乱方法诸如：徐亚等人^[6]采用的Arnold变换、郭建胜等人^[7]采用的Baker变换、王静等人^[8]采用的混沌序列排序构造替换表等；像素值的替代方法诸如：Song等人^[9]采用的秘钥序列与明文进行异或运算，考虑到安全问题，随后Zhang等人^[10]引入密文反馈机制。也有部分算法如Zahmoul等人^[11]采用mod运算，廖琪男等人^[12]将mod运算与异或运算相结合，并引入反馈机制。总结当前各类混沌加密算法存在如下安全性问题：1)加密算法中秘钥序列完全依赖于秘钥算法，对于不同的明文进行加密，秘钥序列是一样的，使得攻击者可以通过选择明文攻击破解^[13-16]。2)只进行了一轮扩散替代操作，使得攻击者可以通过明文和密文推导出秘钥从而破解密文^[17]。3)由于单一混沌的频谱、误差函数等特性明确，在EFA(error function attack)攻击下，即使高维混沌也不够安全。

针对以上基于混沌系统的一系列图像加密算法所存在的安全性问题，提出了一种新的整合神经网络置乱图像的动态自反馈混沌系统图像加密算法。该算法通过1维Logistic混沌、chebyshev混沌和自定义$m\left( x \right)$运算构造了一种动态自反馈混沌系统，通过频数检测、序列分布图、平衡度分析、Lyapunov指数等详细验证了混沌系统的随机性；算法对动态自反馈混沌系统的秘钥序列进行了均匀化处理，通过序列分布图、期望和方差给予了均匀性的证明；混沌系统每一组秘钥序列的产生与明文图像自身的特性紧密相连，不同的明文图像产生完全不同的混沌序列，使得选择明文(密文)攻击方法无效。该算法有别于其他加密算法的是：采用新的置乱方法，引入神经网络，产生与原图像大小相等的置乱矩阵，实现图像的全局置乱，在初次置乱的基础之上再进行二轮的比特位置乱。在文献[18]的思想之上，本文算法在像素值替代的过程中使用两组不同的秘钥序列，使得明文、中间密文、密文呈现复杂的非线性关系，且明文、密文没有直接关系，能够有效抵抗选择明文攻击。实验仿真表明，本文加密算法同时具有较好的抗差分攻击能力和抗统计特性分析能力，还具有较大的密钥空间和对初始条件敏感性等密码学特征，体现了良好的安全性。

定义动态自反馈混沌系统为$h\left( x \right)$，$h\left( x \right)$由$f\left( x \right)$和$g\left( x \right)$两个子系统构成，其中$f\left( x \right)$和$g\left( x \right)$都是定义在[0, 1]区间上的1维离散混沌动力系统。

$ \left\{ \begin{array}{l} f\left( {{x_i}} \right) = {f_q}\left( {{x_{i - 1}}} \right)\;\;\;\;q = 0,1,2 \cdots ,k\\ g\left( {{x_i}} \right) = {g_q}\left( {{x_{i - 1}}} \right)\;\;\;\;q = 0,1,2 \cdots ,k \end{array} \right. $

(1)

本文算法所构造的动态自反馈混沌系统其子系统$f\left( x \right)$和$g\left( x \right)$分别选取1维Logistic混沌系统和进行绝对值处理的Chebyshew混沌系统，故$h\left( x \right)$的可以表示为

$ h\left( x \right) = f\left( {m\left( {f\left( x \right),g\left( x \right)} \right)} \right) $

(2)

模仿有限域上的乘法运算操作交换群概念，将($m\left( x \right)$)定义为[0, 1]上的模1加法运算$\left( {{\mathit{\boldsymbol{G}}_q}, { \oplus _q}} \right)$，${ \oplus _q}$定义为

$ \forall a,b \in {\mathit{\boldsymbol{G}}_q},a{ \oplus _q}b = \bmod \left( {a + b,1} \right) $

(3)

那么动态自反馈混沌系统就可以具体表示为

$ \left\{ \begin{array}{l} h\left( x \right) = f\left( {\bmod \left( {f\left( x \right) + g\left( x \right),1} \right)} \right)\\ f\left( {{x_i}} \right) = u{x_{i - 1}}\left( {1 - {x_{i - 1}}} \right)\\ g\left( {{x_i}} \right) = {\rm{abs}}\left( {\cos \left( {k\arccos {x_{i - 1}}} \right)} \right) \end{array} \right. $

(4)

式中，abs为取绝对值函数，为了检验动态自反馈混沌系统的性能，将其混沌序列$\left\{ {{a_1}, {a_2} \ldots {a_n}} \right\}$转化为二值序列$\left\{ {{b_1}, {b_2} \ldots {b_n}} \right\}$，再对二值序列进行分析，转换公式为

$ {b_i} = \left\{ \begin{array}{l} 0\;\;\;0 \le {a_i} < 0.5\\ 1\;\;\;\;0.5 \le {a_i} \le 1 \end{array} \right. $

(5)

频数检验通过检测二值序列中0和1的个数，用以判断序列的随机情况。若二值序列具有较强随机性，其序列中0和1的个数必然大致相等，设${n_1}$，${n_0}$分别为二值序列中1和0的个数，$ n$为序列的长度，那么频数检验计算公式就可表示为

$ \varepsilon _1^2 = \frac{{{{\left( {{n_1} - {n_0}} \right)}^2}}}{{{n^2}}} $

(6)

对应5 %的显著性水平，与1自由度的$\varepsilon _1^2$分布比较，可得$\varepsilon _1^2$的值为3.84，因此，只要计算本文构造的动态自反馈混沌系统的$\varepsilon _1^2$值大小，与3.84比较就可判断序列的随机性是否良好。本文分别以0.1为间隔从0~1取值生成长度为1 000的随机数列，得到的频数检验结果如表 1所示。

实验结果表明，动态自反馈混沌系统生成的随机序列可以通过频数检测，并且$\varepsilon _1^2$的值远远小于3.84，可以判断序列的随机性良好。

根据本文定义的动态自反馈混沌系统，设定初值$x$ =0.4，参数u =4, k =3, 生成长度为1 000的随机序列，其分布情况如图 1所示。

从图 1可以看出，动态自反馈混沌系统产生的混沌序列在[0, 1]区间上满足近似均匀的随机分布。

同频数检验一样，平衡度也是通过二值序列中1和0的个数是否相近来表征二值序列随机性是否良好，平衡度的值越小，随机性越好，反之，随机性越差。设长为$L$的二值序列，${n_1}$，${n_0}$分别为序列中1和0的个数，则二值序列的平衡度表达式为

$ E\left( L \right) = \frac{{\left| {{n_1} - {n_0}} \right|}}{L} $

(7)

采取抽取样本的方式，分别从0~1以0.1为间隔取值生成长为度1 000的随机数列，计算11组数据得到平衡度情况如图 2所示。

实验结果表明，动态自反馈混沌系统选取不同的初值得到的平衡度曲线接近于0，平衡度取值接近于0，序列的随机性良好。

序列的自相关性曲线直接体现了序列的随机程度，自相关性系数${a^c}$与序列的步长$m$相关，若${a^c}$值随$m$变化越小，则表示序列的随机性越好，反之，随机性越差。假设序列${a_n}$的长度为$N$，那么自相关系数${a^c}$可以通过下式来计算，即

$ {a^c}\left( m \right) = \frac{1}{N}\sum\limits_{i = 1}^{N - m} {{a_i}{a_{i + m}}} $

(8)

在同等条件下，本文得到1维Logistic混沌系统、Chebyshew混沌系统及动态自反馈混沌系统的自相关性如图 3所示。

从图 3可以看出，动态自反馈混沌系统相较于1维Logistic混沌和chebyshev混沌，自相关性系数曲线变化范围更小，绝大部分在0.02的区间内波动，而Logistic混沌系统自相关性系数曲线整体波动范围较大，chebyshev混沌系统波动范围则在0.2区间内，可以判断动态自反馈混沌系统自相关性更低，更符合加密要求。

Lyapunov指数反映了混沌系统轨线局部发散或收缩的状况以及对初值微小变化的敏感性，它体现了轨迹的不可预测性和随机性。1维混沌系统的Lyapunov指数的计算公式为

$ LE = \mathop {\lim }\limits_{n \to \infty } \frac{1}{n}\sum\limits_{i = 0}^{n - 1} {\ln \left| {f'\left( {{x_i}} \right)} \right|} $

(9)

由第1节可知本文定义的动态自反馈混沌系统其实质为：将运算$m\left( x \right)$的输出值作为$f\left( x \right)$的输入值，$h\left( x \right)$实际由两个子系统组合而成，系统1为$f\left( x \right)$，系统2为$m\left( x \right)$。根据Lyapunov指数的定义公式显然可以得到$h\left( x \right)$的Lyapunov指数为

$ {L_E} = \mathop {\lim }\limits_{n \to \infty } \frac{1}{n}\sum\limits_{i = 0}^{n - 1} {\ln \left| {{{\left[ {f\left( {m\left( {{x_i}} \right)} \right)} \right]}^\prime }} \right|} $

(10)

对上式进一步推导可得

$ \begin{array}{*{20}{c}} {{L_E} = \mathop {\lim }\limits_{n \to \infty } \frac{1}{n}\sum\limits_{i = 0}^{n - 1} {\ln \left| {m'\left( {{x_i}} \right) \cdot f'\left[ {m\left( {{x_i}} \right)} \right]} \right|} = }\\ {\mathop {\lim }\limits_{n \to \infty } \frac{1}{n}\sum\limits_{i = 0}^{n - 1} {\ln \left| {m'\left( {{x_i}} \right)} \right|} + }\\ {\mathop {\lim }\limits_{n \to \infty } \frac{1}{n}\sum\limits_{i = 0}^{n - 1} {\ln \left| {f'\left( {m\left( {{x_i}} \right)} \right)} \right|} } \end{array} $

(11)

由此可以看出$h\left( x \right)$的${L_E}$值可以拆分为两部分，即

$ \mathop {\lim }\limits_{n \to \infty } \frac{1}{n}\sum\limits_{i = 0}^{n - 1} {\ln \left| {m'\left( {{x_i}} \right)} \right|} = {L_{E1}} $

(12)

而${L_{E1}}$实际为$m\left( x \right)$的Lyapunov指数, 又由于

$ {x_{2i}} = m\left( {{x_i}} \right),f'\left( {{x_{2i}}} \right) = \frac{{{\rm{d}}f}}{{{\rm{d}}{x_{2i}}}} $

(13)

$ \mathop {\lim }\limits_{n \to \infty } \frac{1}{n}\sum\limits_{i = 0}^{n - 1} {\ln \left| {f'\left( {m\left( {{x_i}} \right)} \right)} \right|} = {L_{E2}} $

(14)

${L_{E2}}$实际为$f\left( x \right)$的Lyapunov指数。由混沌理论可知：若一个系统为混沌状态其Lyapunov指数必然大于0。因此，对于混沌系统$f\left( x \right)$和$m\left( x \right)$，其Lyapunov指数存在下列关系

$ \left\{ \begin{array}{l} L{E_1} > 0\\ L{E_2} > 0 \end{array} \right. $

(15)

$ {\rm{故可得}}\;\;\;\;LE = L{E_1} + L{E_2} > 0 $

(16)

由此可以判断，当$h\left( x \right)$处于混沌状态时，其Lyapunov指数即为$h\left( x \right)$和$m\left( x \right)$ 2个子系统Lyapunov指数之和。对Lyapunov指数推导结果验证如下：分别计算Logistic混沌系统、自定义$m\left( x \right)$运算、动态自反馈混沌系统${L_E}$值，结果如表 2所示。

从表 2可以看出，相较于$f\left( x \right)$，自定义$m\left( x \right)$运算的${L_E}$值改变不大，即$m\left( x \right)$运算没有提升系统的混沌特性，但是将$m\left( x \right)$和$f\left( x \right)$组合后的动态自反馈混沌系统${L_E}$值远远大于$f\left( x \right)$和$m\left( x \right)$混沌系统，近似于$f\left( x \right)$和$m\left( x \right)$的Lyapunov指数之和，验证了对Lyapunov指数的推导结果。也由此可见，动态自反馈混沌系统Lyapunov指数更大，具有更好的随机性，更适合于加密。

在同等条件下，分别对Logistic混沌系统、chebyshev混沌系统、$m\left( x \right)$混沌系统、动态自反馈混沌系统的序列取1 000个样本进行统计分析，直方图如图 4所示。

通过1.1—1.4节验证可知动态自反馈混沌系统随机性得到了增强，但是从取样直方图可以看出$m\left( x \right)$混沌系统的均匀性强于原1维混沌系统，而动态自反馈混沌序列的均匀性没有得到改变，保留了Logistic混沌的分布特征。通过引入概率均匀性指数，对混沌序列的均匀性验证如下：

设0~1二值随机变量为$\beta $，$P\left( {\beta = 0} \right)$和$P\left( {\beta = 1} \right)$分别表示${\beta = 0}$和${\beta = 1}$的概率，当$P\left( {\beta = 0} \right)$= $P\left( {\beta = 1} \right)$时，$\beta $均匀分布，$\beta $的概率均匀性指数为$\phi $，$\phi $的值越大则$\beta $的均匀性越好，其计算公式为

$ \phi = \frac{{\min \left( {P\left( {\beta = 0} \right),P\left( {\beta = 1} \right)} \right)}}{{\max \left( {P\left( {\beta = 0} \right),P\left( {\beta = 1} \right)} \right)}},\;\;\;\;\;\;\phi \in \left( {0,1} \right] $

(17)

设$\alpha $, $\beta $为1维Logistic混沌和chebyshev混沌的0~1二值随机变量，$\gamma $为$m\left( x \right)$混沌的二值随机变量，$\gamma $ =mod($\alpha $ + $\beta $, 2)。设

$ \left\{ \begin{array}{l} P\left( {\alpha = 0} \right) = a\\ P\left( {\beta = 0} \right) = b \end{array} \right. $

(18)

可得

$ \left\{ \begin{array}{l} P\left( {\alpha = 1} \right) = 1 - a\\ P\left( {\beta = 1} \right) = 1 - b \end{array} \right. $

(19)

故

$ \begin{array}{*{20}{c}} {P\left( {\gamma = 0} \right) = P\left( {\alpha = 0,\beta = 0} \right) + P\left( {\alpha = 1,} \right.}\\ {\left. {\beta = 1} \right) = a + \left( {1 - b} \right)\left( {1 - 2a} \right)} \end{array} $

(20)

$ \begin{array}{*{20}{c}} {P\left( {\gamma = 1} \right) = P\left( {\alpha = 0,\beta = 1} \right) + P\left( {\alpha = 1,} \right.}\\ {\left. {\beta = 0} \right) = a + b\left( {1 - 2a} \right)} \end{array} $

(21)

$ \left\{ \begin{array}{l} a < a + \left( {1 - b} \right)\left( {1 - 2a} \right) < 1 - a\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;a < 1/2\;或\;a > 1/2\\ a = 1 - a\;\;\;\;\;a = 1/2\\ a < a + b\left( {1 - 2a} \right) < 1 - a\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;a < 1/2\;或\;a > 1/2 \end{array} \right. $

(22)

根据式(20)可知

$ \left\{ \begin{array}{l} \min \left( {a,1 - a} \right) \le \max \left( {a,1 - a} \right)\\ \min \left( {b,1 - b} \right) \le \max \left( {b,1 - b} \right) \end{array} \right. $

(23)

结合式(15)可得

$ 1 \ge {\phi _\gamma } \ge {\phi _\alpha }\;且\;1 \ge {\phi _\gamma } \ge {\phi _\beta } $

(24)

由此可以推断$m\left( x \right)$混沌序列均匀性优于1维Logistic混沌系统和chebyshev混沌系统，但当将$m\left( x \right)$输入$f\left( x \right)$后，序列继承了$f\left( x \right)$混沌系统的特性，也就是说和$f\left( x \right)$具有同样的分布特征。因此本文采用文献[19]方式对动态自反馈混沌系统进行均匀化处理：

$ y = \frac{2}{{\rm{ \mathsf{ π} }}}\arcsin \sqrt x $

(25)

均匀化处理后，混沌序列在[0, 1]区间均匀分布的证明如下：由Logistic混沌定义可知其实际为定义在[0, 1]区间，值域为[0, 1]的关于随便变量$x$和$y$的单调函数，$x$在[0, 1]上的概率密度函数为$\rho \left( x \right)$，计算公式为

$ \begin{array}{*{20}{c}} {\rho \left( x \right) = y' = \frac{2}{{\rm{ \mathsf{ π} }}} \times \frac{1}{{\sqrt {1 - x} }} \times \frac{1}{{2\sqrt x }} = }\\ {\frac{1}{{{\rm{ \mathsf{ π} }}\sqrt {x\left( {1 - x} \right)} }}} \end{array} $

(26)

那么$x$的概率密度函数就可表示为

$ \rho \left( x \right) = \left\{ \begin{array}{l} \frac{1}{{{\rm{ \mathsf{ π} }}\sqrt {x\left( {1 - x} \right)} }}\;\;\;0 < x < 1\\ 0\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;其他 \end{array} \right. $

(27)

设动态自反馈混沌系统迭代序列均匀化处理后为[c₁, c₂, c₃... c_n]，将混沌吸引域[0, 1]划分为由$N$个点组成的子域${\mathit{\boldsymbol{k}}_i}$。

$ \left\{ \begin{array}{l} {\mathit{\boldsymbol{k}}_i} = \left[ {{S_i},{S_{i + 1}}} \right);i = 0,1, \cdots ,N - 2\\ {\mathit{\boldsymbol{k}}_{N - 1}} = \left[ {{S_{N - 1}},{S_N}} \right] \end{array} \right. $

(28)

${S_i}$可以定义为

$ {S_i} = \sin \left( {\frac{i}{N} \cdot \frac{{\rm{ \mathsf{ π} }}}{2}} \right);i = 0,1, \cdots ,N - 1 $

(29)

根据式(25)(27)可得序列中元素$P\left( i \right)$出现的概率为

$ \begin{array}{*{20}{c}} {P\left( i \right) = \int_{{S_i}}^{{S_{i + 1}}} {\rho \left( x \right){\rm{d}}x} = \int_{{S_i}}^{{S_{i + 1}}} {\frac{1}{{{\rm{ \mathsf{ π} }}\sqrt {x\left( {1 - x} \right)} }}{\rm{d}}x} = }\\ {\frac{2}{{\rm{ \mathsf{ π} }}}\arcsin \sqrt x \left| {_{{S_i}}^{{S_{i + 1}}}} \right. = \frac{1}{N}} \end{array} $

(30)

由此可知均匀化处理后的序列在区间[0, 1]上服从均匀分布。处理之后的序列直方图如图 4(d)所示，可以看出，处理之后的序列分布特点发生了改变，整体变得均匀，不再集中在两端。

由概率论知识可知，若一个序列随机且分布均匀，则满足

$ f\left( x \right) = \left\{ \begin{array}{l} \frac{1}{{b - a}}\;\;\;a < x < b\\ 0\;\;\;\;\;\;\;\;\;其他 \end{array} \right. $

(31)

其期望值$E\left( X \right)$和方差$D\left( X \right)$计算公式分别为

$ \begin{array}{*{20}{c}} {E\left( X \right) = \int\limits_{ - \infty }^{ + \infty } {xf\left( x \right){\rm{d}}x} = \int\limits_a^b {\frac{x}{{b - a}}{\rm{d}}x} = }\\ {\frac{{{b^2} - {a^2}}}{{2\left( {b - a} \right)}} = \frac{{a + b}}{2}} \end{array} $

(32)

$ \begin{array}{*{20}{c}} {D\left( X \right) = E\left( {{X^2}} \right) - E\left( {{X^2}} \right) = }\\ {\int\limits_{ - \infty }^{ + \infty } {{x^2}f\left( x \right){\rm{d}}x} - {{\left( {\frac{{a + b}}{2}} \right)}^2} = \frac{{{{\left( {b - a} \right)}^2}}}{{12}}} \end{array} $

(33)

对于混沌系统，显然其序列值属于[0, 1]。取$b$ =1, $a$ =0，代入计算序列$\mathit{\boldsymbol{x}}\left( n \right)$的理论期望值和方差值，再对动态自反馈混沌系统进行均匀化处理，生成序列$\mathit{\boldsymbol{x}}\left( n \right)$，分别计算实际期望值和方差值，得到结果如表 3所示。

由表 3可以看出，动态自反馈混沌系统均匀化处理后生成的混沌序列，其实际期望值和方差值与理论值相差很小，由此可以判断动态自反馈混沌系统满足均匀分布。

当前图像像素矩阵的置乱普遍采用简单的像素矩阵行列互换、Arnold变换、Baker变换、3D猫映射、以及混沌序列排序构造替换表、bit位的分块置乱等，本文则是采用了新的方法来实现像素全局置乱。构造一个3层的神经网络，神经网络模型如图 5所示。

神经网络的输入层为$\mathit{\boldsymbol{I}}$，从动态自反馈混沌序列${\mathit{\boldsymbol{S}}_x}\left( \mathit{\boldsymbol{i}} \right)$中随机选取；输出层为$\mathit{\boldsymbol{F}}$，由${F_1}$和${F_2}$构成；传输函数由函数$p\left( x \right)$、$r\left( x \right)$和$t\left( x \right)$构成，实现对输入数据的处理。

$ \mathit{\boldsymbol{I}} = {\left[ {{X_1},{X_2},{X_3},{X_4}} \right]^{\rm{T}}} $

(34)

此时$\mathit{\boldsymbol{I}}$的取值均在(0, 1)之间，通过函数$p\left( x \right)$将$\mathit{\boldsymbol{I}}$的部分值扩散到(0, 1)外，变换后的矩阵为${\mathit{\boldsymbol{I'}}}$，即

$ \mathit{\boldsymbol{I'}} = p\left( x \right) = \mathit{\boldsymbol{G}} \times \mathit{\boldsymbol{I}} = {\left[ {{{I'}_1},{{I'}_2},{{I'}_3},{{I'}_4}} \right]^{\rm{T}}} $

(35)

式中，${\mathit{\boldsymbol{G}}}$为4×4方阵。${\mathit{\boldsymbol{I'}}}$通过传输函数$r\left( x \right)$和$t\left( x \right)$得到输出值${\mathit{\boldsymbol{F}}}$，其公式为

$ \mathit{\boldsymbol{F}} = r\left( {t\left( {\mathit{\boldsymbol{WI'}} + \mathit{\boldsymbol{B}}} \right),\mathit{\boldsymbol{Q}}} \right) $

(36)

式中，${\mathit{\boldsymbol{W}}}$为2行4列矩阵，${\mathit{\boldsymbol{B}}}$为2行1列矩阵，均从序列${\mathit{\boldsymbol{S}}_x}\left( \mathit{\boldsymbol{i}} \right)$中随机选取，${\mathit{\boldsymbol{W}}}$和${\mathit{\boldsymbol{B}}}$分别是神经元的权值和阙值，${\mathit{\boldsymbol{Q}}}$为2行1列矩阵，是函数$r\left( x \right)$的控制参数，函数$r\left( x \right)$和$t\left( x \right)$的表达式分别为

$ \begin{array}{l} r\left( {x\left( k \right),q} \right) = \\ \left\{ \begin{array}{l} x\left( k \right)/q\;\;\;\;\;\;\;\;\;\;\;0 < x\left( k \right) \le q\\ 1 - x\left( k \right)/1 - q\;\;\;q < x\left( k \right) < 1 \end{array} \right. \end{array} $

(37)

$ t\left( x \right) = 1/1 + {{\rm{e}}^{ - ax}} $

(38)

式中，$a $为函数$t\left( x \right)$的参数，函数$t\left( x \right)$将数值映射到区间(0, 1)上，使得传输函数$r\left( x \right)$刚好与$t\left( x \right)$区间吻合，$x\left( k \right)$是由$t\left( x \right)$传输过来的状态，$q $为控制参数，并且0 < $x\left( k \right)$ < 1, 0 < $q $ < 1。

在每一次输出后更新$\mathit{\boldsymbol{Q}}$的值为

$ \mathit{\boldsymbol{Q}} = \left( {0.2 \times \mathit{\boldsymbol{F}}} \right) + 0.4 $

(39)

每运行一次神经网络，都会产生两个随机值${F_1}$，${F_2}$，并更新$\mathit{\boldsymbol{Q}}$值。采用更新后的$\mathit{\boldsymbol{Q}}$重复运行神经网络得到${F_1}$，${F_2}$，并将${F_1}$，${F_2}$处理到[0, $M$]和[0, $N$]区间内，$M$，$N$分别为原始图像矩阵的行和列，若已产生过相同的一对值则舍去，直至产生与明文像素大小相等的矩阵。由产生的置乱矩阵下标索引，将明文图像矩阵的像素值依次放到新矩阵中，全局置乱的过程完成。

由BP神经网络知识可知，神经网络由一个输入层、一个输出层和一个或多个隐含层组成。神经网络通过信号的正向传播与误差的反向传播完成学习过程。当进行正向传播时，输入样本从输入层传入，经隐含层处理后传向输出层，如果输出层输出与期望输出不一致，就跳转为误差反向传播。误差反向传播则是采用梯度下降法将输出误差以某种形式通过隐层向输入层逐层反传，逐层迭代修改各神经元之间的网络权值和各神经元阈值。本文采用神经网络置乱的实质是用来产生随机矩阵，因此在每一轮的迭代中，默认随机选取了扩散矩阵$\mathit{\boldsymbol{G}}$、权值$\mathit{\boldsymbol{W}}$、阙值$\mathit{\boldsymbol{B}}$、初始控制参数矩阵$\mathit{\boldsymbol{Q}}$的神经网络是已经训练好的神经网络，也就是说忽略学习效率对神经网络计算的影响。神经网络计算效率就可简化为每轮运算所需要的时间，以256×256灰度图像为例，选取初值，生成随机序列${\mathit{\boldsymbol{S}}_x}\left( \mathit{\boldsymbol{i}} \right)$，得到每轮计算所需时间、完成$M$×$N$轮计算所需时间及最终获得完整置乱矩阵所需时间如表 4所示。

由置乱算法设计思想可知，算法需要生成一个无重复的$M$×$N$大小的矩阵，也就是说最终需要神经网络产生$M$×$N$对完全不同的数组，以256×256像素灰度图像为例，选取初值，得到神经网络迭代轮次和结果遍历性如表 5所示。

从表 5可以看出，在迭代轮次较少时，结果的重复性较低，当迭代的轮次增多时，结果的重复性增高，直至将迭代轮次设置至350万次，才得到最终遍历所有数组的结果，产生完整的置乱矩阵。由此可以判断，只要迭代轮次足够大，计算结果可以遍历所需大小数组，但是也因此降低了算法的计算效率，当图像像素矩阵越大时，所需的计算时间就越久。

本文算法原理：在加密过程中与明文紧密相关，不同的密文将产生不同的密钥，有效防止了攻击者对等效秘钥的破解；采用神经网络产生置乱矩阵进行全局置乱，再从bit位进行置乱；在进行像素替代时采用两组秘钥进行了两轮替代扩散，明文、中间密文、密文三者非线性关系更为复杂，有效的防止了选择明文攻击。本文加密算法设计图如图 6所示。

读取明文图像像素矩阵，L = $M$ × $N$，$M$，$N$分别为像素矩阵的行和列。计算Chebyshev混沌系统的参数$k$和动态自反馈混沌系统的初值${x_0}$，计算公式为

$ \left\{ \begin{array}{l} k = \bmod \left( {{\rm{sum}},4} \right)\\ {x_0} = {x_0} \times q/256 \end{array} \right. $

(40)

式中，$sum$为图像所有像素值和，$q$为所有像素按位异或的值。

选取不同的初值$x$，$x = x \times {x_0}$，运行动态自反馈混沌系统，迭代$N$₀+L次消除暂态影响，选取后L个值为秘钥序列，产生4个混沌序列：${\mathit{\boldsymbol{S}}_x}\left( \mathit{\boldsymbol{i}} \right)$, $\mathit{\boldsymbol{x}}\left( i \right)$，$\mathit{\boldsymbol{y}}\left( i \right)$，$\mathit{\boldsymbol{z}}\left( i \right)$，i =1, 2, 3，…，L。

1) 像素级全局置乱。像素位置置乱操作是为了实现像素位置的全局置乱，破坏原来相邻像素之间的相关性。神经网络置乱像素的具体步骤为：

(1) 扫描明文图像生成明文图像矩阵$\mathit{\boldsymbol{B}}$，设置乱后的图像矩阵为$\mathit{\boldsymbol{A}}$。

(2) 从${\mathit{\boldsymbol{S}}_x}\left( \mathit{\boldsymbol{i}} \right)$中随机选取数值构成输入层$\mathit{\boldsymbol{I}}$、扩散矩阵$\mathit{\boldsymbol{G}}$、权值$\mathit{\boldsymbol{W}}$、阙值$\mathit{\boldsymbol{B}}$、初始控制参数矩阵$\mathit{\boldsymbol{Q}}$。

(3) 运行神经网络得到${F_1}$, ${F_2}$，将${F_1}$，${F_2}$分别改造为[1, $M$]和[1, $N$]的整数${{F'}_1}$，${{F'}_2}$，并且更新$\mathit{\boldsymbol{Q}}$值，$M$, $N$为明文像素矩阵的行和列。

(4) 设置数组$\mathit{\boldsymbol{R}}$和$\mathit{\boldsymbol{C}}$分别存储${{F'}_1}$，${{F'}_2}$的值。

(5) 重复步骤(1)—(4)，如果${{F'}_1}$，${{F'}_2}$值同时已经生成过则舍去这组数值。将$\mathit{\boldsymbol{R}}$和$\mathit{\boldsymbol{C}}$对应为元素的行和列生成置换矩阵，直至产生与明文图像大小相等的矩阵。

(6) 完成图像置乱，计算公式为

$ A\left( i \right) = B\left( {R\left( i \right),C\left( i \right)} \right) $

(41)

2) Bit位置乱。Bit置乱可以增强算法的安全性，使选择明文攻击无效，其具体步骤描述为：

(1) 读取秘钥序列$\mathit{\boldsymbol{x}}\left( i \right)$，$x = \left\{ {{x_1}, {x_2}, {x_3}, \ldots {x_L}} \right\}$。

(2) 选取$\mathit{\boldsymbol{x}}\left( i \right)$，$\mathit{\boldsymbol{x}}\left( i \right)$为混沌序列$x$中第$i$个实数值，提取实数值$\mathit{\boldsymbol{x}}\left( i \right)$的小数点后8位数字，并由这8个数字组成数组$\mathit{\boldsymbol{M}} = \left\{ {{m_1}, {m_2}, {m_3}, \ldots {m_8}} \right\}$。

(3) 对步骤(2)中生成的数组$\mathit{\boldsymbol{M}} $按照由小到大进行排序，得到有序数组$\mathit{\boldsymbol{M’}} $；然后利用数组$\mathit{\boldsymbol{W}} $保存$\mathit{\boldsymbol{M’}} $各元素在$\mathit{\boldsymbol{W}} $中的位置。

(4) 将置乱图像矩阵$\mathit{\boldsymbol{A}} $中第$i$点的像素值${P_i}$转换成二进制形式，生成数组$\mathit{\boldsymbol{P}}_i^1 = \left\{ {{B_1}, {B_2}, {B_3}, \ldots {B_8}} \right\}$。

(5) 利用数组$\mathit{\boldsymbol{W}} $置乱数组$\mathit{\boldsymbol{P}}_i^1 = \left\{ {{B_1}, {B_2}, {B_3}, \ldots {B_8}} \right\}$，得到像素值比特置乱后重新排列的新数组$\mathit{\boldsymbol{P}}{{_{i}^{1}}^{\prime }}=\left\{ {{B}_{1}}^{\prime }, {{B}_{2}}^{\prime }, {{B}_{3}}^{\prime }, \ldots {{B}_{8}}^{\prime } \right\} $。

(6) 将二进制数组$\mathit{\boldsymbol{P}}{{_{i}^{1}}^{\prime }}$转换成十进制数，得到bit位置乱后的密文$\mathit{\boldsymbol{C}}\left( i \right)$。

(7) 重复执行，直到所有点转为密文，密文矩阵为$\mathit{\boldsymbol{C}}$。

3) 像素的替代扩散。像素的替代扩散分为两轮进行，所用秘钥参数为${C_0}$，秘钥组为${\mathit{\boldsymbol{y}}\left( i \right)}$，${\mathit{\boldsymbol{z}}\left( i \right)}$，此时，${\mathit{\boldsymbol{y}}\left( i \right)}$, ${\mathit{\boldsymbol{z}}\left( i \right)}$均为小数，通过${\mathit{\boldsymbol{x}}\left( i \right)}$=mod(round ${\mathit{\boldsymbol{x}}\left( i \right)}$×10¹⁴, 256);将其改造为[0, 255]区间的整数，得到新的整数形势的序列${\mathit{\boldsymbol{y}}\left( i \right)}$，${\mathit{\boldsymbol{z}}\left( i \right)}$，round为四舍五入函数从置乱后的图像像素中随机选取两个分别添加到秘钥序列${\mathit{\boldsymbol{y}}\left( i \right)}$和${\mathit{\boldsymbol{z}}\left( i \right)}$中，使其长度变为L+1。

(1) 第1轮扩散替代加密

当 $i = 1$ 时

$ \left\{ \begin{array}{l} M1 = \bmod \left( {C\left( 1 \right) + y\left( 2 \right),256} \right)\\ D\left( 1 \right) = M1 \oplus \bmod \left( {y\left( 1 \right) + {C_0},256} \right) \end{array} \right. $

(42)

当 $1 < i \le L$ 时

$ \left\{ \begin{array}{l} M1 = \bmod \left( {C\left( i \right) + y\left( {i + 1} \right),256} \right)\\ D\left( i \right) = M1 \oplus \bmod \left( {y\left( i \right) + D\left( {i - 1} \right),256} \right) \end{array} \right. $

(43)

(2) 第2轮扩散替代加密

当$i = 1$时

$ \left\{ \begin{array}{l} M2 = \bmod \left( {D\left( 1 \right) + z\left( 2 \right),256} \right)\\ D1\left( 1 \right) = M2 \oplus \bmod \left( {z\left( 1 \right) + D\left( L \right),256} \right) \end{array} \right. $

(44)

当$1 < i \le L$时

$ \left\{ \begin{array}{l} M2 = \bmod \left( {D\left( i \right) + z\left( {i + 1} \right),256} \right)\\ D1\left( i \right) = M2 \oplus \bmod \left( {z\left( i \right) + D1\left( {i - 1} \right),256} \right) \end{array} \right. $

(45)

(3) 重复步骤(1)(2)，直到$i$ = L，得到密文图像${\mathit{\boldsymbol{D}}_{\rm{1}}}$，加密过程结束。

解密过程是加密过程的逆过程，分别对应替代扩散、bit位置乱和全局置乱进行反操作。

1) 首先由下式得到置乱后的密文${\mathit{\boldsymbol{C}}}$：

(1) 第1轮解密

当$1 < i \le L$时

$ \left\{ \begin{array}{l} M2 = D1\left( i \right) \oplus \bmod \left( {z\left( i \right) + D1\left( {i - 1} \right),256} \right)\\ D\left( i \right) = \bmod \left( {M2 + 256 - z\left( {i + 1} \right),256} \right) \end{array} \right. $

(46)

当$i = 1$时

$ \left\{ \begin{array}{l} M2 = D1\left( 1 \right) \oplus \bmod \left( {z\left( 1 \right) + D\left( L \right),256} \right)\\ D\left( 1 \right) = \bmod \left( {M2 + 256 - z\left( 2 \right),256} \right) \end{array} \right. $

(47)

(2) 第2轮解密

当$1 < i \le L$时

$ \left\{ \begin{array}{l} M1 = D\left( i \right) \oplus \bmod \left( {y\left( i \right) + D\left( {i - 1} \right),256} \right)\\ J\left( i \right) = \bmod \left( {M1 + 256 - y\left( {i + 1} \right),256} \right) \end{array} \right. $

(48)

当$i = 1$时

$ \left\{ \begin{array}{l} M1 = J\left( 1 \right) \oplus \bmod \left( {y\left( 1 \right) + {C_0},256} \right)\\ J\left( 1 \right) = \bmod \left( {M1 + 256 - y\left( 2 \right),256} \right) \end{array} \right. $

(49)

${\mathit{\boldsymbol{J}}}$即为反扩散替代得到的解密矩阵，$\mathit{\boldsymbol{J = C}}$。

2) 把密文${\mathit{\boldsymbol{J}}}$转化为2进制，对应置乱数组${\mathit{\boldsymbol{W}}}$，反置乱得到原始的比特序列，再转换为十进制数，得到位置乱前的置乱矩阵${\mathit{\boldsymbol{A}}}$。

3) 根据全局置乱后的矩阵${\mathit{\boldsymbol{A}}}$，得到原始矩阵${\mathit{\boldsymbol{B}}}$，至此全部解密过程完成。

$ B\left( {R\left( i \right),C\left( i \right)} \right) = A\left( i \right) $

(50)

本文算法使用256级灰度图像lena进行测试，仿真环境为Windows 7操作系统，Corei5 CPU，4 GB RAM，Matlab 2010b。经过仿真得到的像素置乱图像、像素替代图像、解密图像分别如图 7所示。

从仿真结果可以看出，经过加密后的图像已经完全改变了明文图像的特征，为了评估该算法的总体性能，下面将分别从秘钥空间、抗统计攻击、抗差分攻击等方面进行分析。

由本文算法可知，动态自反馈混沌系统的初始秘钥有4个，为$x$₁, $x$₂, $x$₃, $x$₄。另有预迭代次数$N$₀，加密引入参数C₀，不妨假设$N$₀的取值个数为500，C_0，的取值个数为255，$x$₁, $x$₂, $x$₃, $x$₄均为15位小数，则可计算秘钥空间为：500×255×10¹⁵×10¹⁵×10¹⁵×10¹⁵≈2²¹⁶，从计算结果可以看出该加密算法秘钥空间足够抵抗穷举攻击。

实验得到明文与密文图像的直方图如图 8所示，从图 8中可以看出：明文图像直方图的像素值的分布范围与数量都不均匀，经过本文算法加密的图像像素值在0~ 255范围内分布均匀并且各像素值的出现概率基本相同，明文像素的统计特性发生了根本改变。由此可以判断，该算法能够有效的抵抗基于统计分析的攻击。

相邻像素相关系数用来反映相邻像素之间的相关程度，检验加密效果，相关系数越接近0则加密效果越好。相关系数通过选取图像中$N$对相邻像素，计算公式为

$ \left\{ \begin{array}{l} \bar x = \frac{1}{N}\sum\limits_{i = 1}^N {{X_i}} \\ D\left( x \right) = \frac{1}{N}\sum\limits_{i = 1}^N {{{\left( {{x_i} - \bar x} \right)}^2}} \\ Conv\left( {x,y} \right) = \frac{1}{N}\sum\limits_{i = 1}^N {\left( {{x_i} - \bar x} \right)\left( {{y_i} - \bar y} \right)} \\ {\gamma _{xy}} = \frac{{Conv\left( {x,y} \right)}}{{\sqrt {D\left( x \right)} \sqrt {D\left( y \right)} }} \end{array} \right. $

(51)

式中，$x$, $y$分别表示图像中两个相邻像素的像素值，${\bar x}$和${\bar y}$是所有$x$和$y$的平均值，$N$为选取的像素组数，${\gamma _{xy}}$为相关系数。分别计算加密后的图像水平、垂直、对角方向的相关性系数结果及其他文献的数据如表 6所示。从表 6可以看出, 明文图像相邻像素水平、垂直、对角3个方向的相关性都大于0.9，相邻像素相关程度高，经过本文算法加密后的图像相邻像素相关性趋近于0，达到了10^-4数量级。分别从水平、垂直、对角方向选取原始图像和加密图像的2 000组像素点，其相邻像素灰度值关系如图 9所示。从图 9中可以看出，明文图像水平、垂直、对角的相邻像素灰度值围绕$y = x$分布，而加密后的图像相邻像素灰度值在0~255之间随机分布。

信息熵用来反应图像中信息分布的随机性强弱，灰度值分布的越均匀，信息熵越大，反之则越小，信息熵$H$可表示为

$ H\left( S \right) = - \sum\limits_{i = 0}^{N - 1} {P\left( {{S_i}} \right)\log \left[ {P\left( {{S_i}} \right)} \right]} $

(52)

式中，$N$表示图像中信息有$N$种不同的取值，取值的集合为$\left\{ {{s_0}, {s_1}, \ldots {s_{N - 1}}} \right\}$, $P\left( {{S_i}} \right)$表示${S_i}$在图像S中出现的概率。本文选取8位Lena图像，其理想值$ H\left( S \right)$ =8。$ H\left( S \right)$=8表示这是一种随机信息源图像，也就是说对于加密算法所得的密文，信息熵越接近8，越接近随机图像。采用本文算法得到的信息熵与相关文献信息熵如表 7所示。

从数据可以看出，明文图像信息熵较小，经过本文算法加密后的信息熵达到了7.998 3, 非常接近于8，可以判断该算法可以极大的改变明文图像像素值特征。

由4.2节可知本文算法可以组成一个6位的秘钥组，即$\left( {{N_0}, {C_0}, {x_1}, {x_2}, {x_3}, {x_4}} \right)$，不妨设置对应的初始秘钥${\mathit{\boldsymbol{s}}_1}$为(500, 235, 0.442 1, 0.386 5, 0.693 2, 0.843 7)，对Lena图像进行加密，然后采用只改变了秘钥${x_1}$的秘钥组${\mathit{\boldsymbol{s}}_2}$ (500, 235, 0.442 2, 0.386 5, 0.693 2, 0.843 7)进行解密。得到解密后的图像如图 10所示，从实验结果可以看出，即使只改变一个秘钥，图像也无法正确解密，解密的图像与原始图像也区别很大。实验改变秘钥组其他秘钥进行测试，都得到类似的结果，说明秘钥具有充分的敏感性。

为观察秘钥改变后密文的变化，本文进一步进行了实验，分别将秘钥组中的${x_1}$, ${x_2}$改变10^-4，比较密文图像中前100个像素点的像素值差异，结果如图 11所示。

实验结果表明：秘钥即使只改变了10^-4，得到的密文图像每个像素值都发生了变化，实验也验证了密文对秘钥组其他秘钥${N_0}, {C_0}, {x_3}, {x_4}$同样敏感。

密文对明文的敏感性强弱体现了算法抵抗差分攻击能力的强弱，实验证实本文加密算法对明文具有充分的敏感性，能够抵抗差分攻击。实验采取计算$NPCR$ (像素变化率)和$UACI$ (归一化像素值平均改变强度)验证，计算公式为

$ \left\{ \begin{array}{l} NPCR = \frac{1}{{M \times N}}\sum\limits_{i = 1}^M {\sum\limits_{j = 1}^N {D\left( {i,j} \right) \times 100\% } } \\ UACI = \frac{1}{{M \times N}}\sum\limits_{i = 1}^M {\sum\limits_{j = 1}^N {\frac{{\left| {{C_1}\left( {i,j} \right) - {C_2}\left( {i,j} \right)} \right|}}{{255}} \times 100\% } } \end{array} \right. $

(53)

式中，$M$、$N$表示图像像素的行数和列数，${C_1}\left( {i, j} \right)$和${C_2}\left( {i, j} \right)$分别表示两幅只有一个像素点不同的明文图像加密后对应的图像中第$\left( {i, j} \right)$点的像素值，如果${C_1}\left( {i, j} \right)$ = ${C_2}\left( {i, j} \right)$，$D\left( {i, j} \right)$ =0, 否则，$D\left( {i, j} \right)$ =1。

$NPCR$和$UACI$的理想期望值为

$ \left\{ \begin{array}{l} NPC{R_E} = \left( {1 - {2^{ - n}}} \right) \times 100\% \\ UAC{I_E} = \frac{1}{{{2^{2n}}}}\frac{{\sum\nolimits_{i = 1}^{{2^n} - 1} {i\left( {i + 1} \right)} }}{{{2^n} - 1}} \times 100\% \end{array} \right. $

(54)

式中，$n$表示图像颜色位深，对于8位灰度图像，$NPCR$和$UACI$理想期望值为0.996 1和0.334 6。

实验选取8位Lena图进行加密，每2个像素值为一组，每组中一个为原始图像，另一个是对原始图像像素加2，随机选取100组图像进行加密，得到的$NPCR$和$UACI$值结果如图 12所示，另实验分别改变明文图像中第10、16 000、32 000、48 000、65 530个像素点的值，计算得到$NPCR$和$UACI$的具体值如表 8所示。

从实验结果可以看出，$NPCR$值在(0.995 8~0.996 6)区间，接近理想值0.996 1，$UACI$值在(0.333~0.338)区间，接近理想值0.334 6，可以表明本文算法具有抵抗差分攻击的能力。

鉴于当前众多图像加密算法存在的普遍安全问题，本文提出了整合神经网络置乱图像的动态自反馈混沌系统图像加密算法，该算法体现了良好的加密效果和性能，且加密过程中一次一密，多重组合抵抗选择明文攻击，通过了抗差分攻击测试、抗统计攻击测试，有效提升了图像加密的安全性。优点方面：该算法创新的构造了一种动态自反馈混沌系统作为秘钥生成源；区别于当前常用的置乱方式，引入了神经网络产生置乱矩阵；并且进行了两轮的像素替代来提升密文安全性。但是在算法实现过程中也存在2个明显的缺点：对于神经网络的输入层$\mathit{\boldsymbol{I}}$、扩散矩阵$\mathit{\boldsymbol{G}}$、权值$\mathit{\boldsymbol{W}}$、阙值$\mathit{\boldsymbol{B}}$、初始控制参数矩阵$\mathit{\boldsymbol{Q}}$从序列${\mathit{\boldsymbol{x}}\left( i \right)}$中随机选取时不能很好的实现遍历性；由于神经网络的重复运算，以及去除重复值的过程，神经网络产生完整的置乱矩阵时间较长。总体而言，该算法可以被广泛应用在灰度图像加密中乃至扩展到彩色图像加密中，对网络传输、存储中的图像信息起到隐私保护的作用。下一步将致力于解决神经网络的输入值及参数值选取的遍历性问题及神经网络计算的时间复杂度问题。