# Adaptive Robust Blind Watermarking Scheme Improved by Entropy-Based SVM and Optimized Quantum Genetic Algorithm

Authors: Jun Zhang; Xiaoyi Zhou; Jilin Yang; Chunjie Cao; Jixin Ma
Journal: Mathematical Problems in Engineering
 (2019)

---

## Abstract

With the intensive study of machine learning in digital watermarking, its ability to balance the robustness and transparency of watermarking technology has attracted researchers’ attention. Therefore, quantum genetic algorithm, which serves as an intelligent optimized scheme combined with biological genetic mechanism and quantum computing, is widely used in various fields. In this study, an adaptive robust blind watermarking algorithm by means of optimized quantum genetics (OQGA) and entropy classification-based SVM (support vector machine) is proposed. The host image was divided into two parts according to the odd and even rows of the host image. One part was transformed by DCT (discrete cosine transform), and then the embedding intensity and position were separately trained by entropy-based SVM and OQGA; the other part was by DWT (discrete wavelet transform), in which the key fusion was achieved by an ergodic matrix to embed the watermark. Simulation results indicate the proposed algorithm ensures the watermark scheme transparency as well as having better resistance to common attacks such as lossy JPEG compression, image darken, Gaussian low-pass filtering, contrast decreasing, salt-pepper noise, and geometric attacks such as rotation and cropping.


---

## Body

## 1. Introduction

A well-designed watermarking scheme shall have the following two characteristics: ① the watermark is imperceptible in the carrier data (transparency); ② the watermark is difficult to be destroyed by unauthorized parties (robustness). Normally, we can increase the watermark embedding intensity to improve robustness, but it will lose imperceptibility and vice versa. Moreover, the advancement of analysis techniques for digital watermarking, especially the extensive application of deep learning in watermark removal technology [1, 2], has made watermarking techniques unsecure. Hence, how to maintain the balance between robustness and imperceptibility as well as the security of watermarking is a key problem for a complex watermarking system [3]. To maintain transparency by changing the transform domain of watermark embedding [4] and to balance robustness and embedding strength by adaptive selection of embedding strength [5] are our main research so far. With intensive study, the application of statistical learning and machine learning provides new methods for optimization and solution to this problem.Genetic algorithm, which is based on biological natural selection and genetic mechanism, is widely used in digital watermarking due to its powerful robustness and ability to solve optimization problems of a complex system. Huang et al. [6–8] proposed a genetic-based watermarking algorithm for progressive transmission. The GA with an appropriate fitness function is applied to the scheme. The watermark embedding and extraction was implemented in the transform domain, and the JPEG spectrum selection mode was applied to complete the scalable transmission of the watermark image. The algorithm has a large watermark capacity and strong robustness, but the imperceptibility is unsatisfying, and the PSNR (peak signal-to-noise ratio) value is limited to the range of 34 to 35. Zhou et al. [9] used an SVM model to obtain the optimal embedding intensity, as well as the genetic algorithm to find the optimal embedding position in the DCT domain, thus to realize the variable embedding of the watermark. The watermarked image and the original host image have higher PSNR and better anti-attack ability to JPEG compression, Gaussian noise, rotation, low-pass filtering, histogram equalization, etc., but the algorithm complexity is high and the resistance to geometric attacks is weak. In particular, the original watermark cannot be well recovered if the rotation attack is more than 1 degree.Although genetic algorithm has better global search ability and strong robustness, it has the same problems of low search efficiency, poor convergence, and easy to fall into local extremum. As a consequence, quantum genetic algorithm (QGA), a new intelligent optimization algorithm combining quantum computing and genetic algorithm, is introduced to solve combinatorial optimization and complex function maximum value problem [10–13]. In the case of the same number of iterations, the population in QGA tends to be much smaller than that of GA and thus can achieve rapid convergence. In addition, quantum rotation gate adjustment, quantum mutation, and other operations greatly reduce the possibility of falling into the local extremum. However, although quantum genetic algorithms are widely used in various fields, there are few cases of applying them into digital watermarking.Xu et al. [14] proposed a real coded quantum genetic algorithm (RQGA) to replace the binary code with real code and the rotation gate with the approximation operator. The stability, convergence, and global optimization of RQGA are proven by function optimization and 0-1 knapsack problem. Li et al. [15] proposed a double-chain quantum genetic algorithm (DCQGA). The chromosome consists of qubits of which the probability amplitude comprises gene chains, eliminating the need for cumbersome decoding operation and improving the convergence efficiency by parallel search. Quantum chromosome is evolved by quantum rotation gates, and the quantum nongate mutation ensures that the optimal solutions are searched in the global scope. However, since the chromosomes are updated by the quantum rotation gate, their offspring individuals may have been very close to the optimal solution. If the quantum nongate mutations are executed in this case, the quantum bits will tend to change the state “0” to the state “π/2,” which results in the qubits being updated in the opposite direction, and it highly tends to cause population oscillations, thus losing the excellent information. In response to this problem, Xu et al. [16] proposed a quantum genetic algorithm based on Hadamard gate variation. The algorithm used the Hadamard gate mutation operation to apply slight fluctuation to the chromosome, which effectively makes the population jumping out of the current optimal solution, increases the diversity of the population, searches in multiple directions, and maintains the stability of the population, therefore ensuring that the excellent information in the offspring population will not be lost. However, the traditional quantum rotation gate cannot achieve the unification of convergence speed and precision. Huan et al. [17] proposed an improved quantum rotation gate adjustment scheme and introduced a limiting correction operator to make the quantum genetic algorithm have better global optimization ability. Based on Darwin’s natural selection mechanism, the rotation angle was adaptively selected for global search according to the ratio of individual fitness to the optimal fitness. Simultaneously, the stagnation algebra was combined with the local search, and the probability amplitude was corrected in the iterative process. The experimental results of typical continuous multimodal functions show that the proposed algorithm has a great improvement in convergence speed and approximation accuracy compared with the existing quantum rotation gate adjustment methods.In view of the above, quantum genetic algorithm has better characteristics than genetic algorithm in many aspects. In this study, a watermarking algorithm was proposed. Firstly, the SVM was used to train the images to obtain the adaptive embedding intensity, and then the quantum genetic algorithm was used to obtain the optimal embedding position in the DCT domain. However, experiments have found that even though the processing speed, efficiency, and image visual effects are greatly improved, the watermarked image is weak in resistance to geometric attacks. By analyzing the robust watermarking algorithms proposed by Zhou et al. [18] and Kricha et al. [19], we found that embedding watermark information in the DWT domain can effectively improve the anti-geometric attack ability for the watermark scheme by experiments. Consequently, an adaptive robust blind watermarking algorithm was proposed based on optimized quantum genetic (OQGA) and improved by entropy-based SVM. In this study, the host image was divided into two parts according to the odd and even pixel rows, and the watermark was divided into half in order. One part of the host image was performed on the DCT domain. To obtain adaptive embedding intensity, it firstly used entropy classification in SVM training for dividing the texture of the original host image into three categories, namely, smooth area, edge area, and texture area, and then further dividing the image into nine categories according to brightness and variance of the host image. After that, the combination of NC (normalized correlation coefficient) value and PSNR value were taken as the fitness function and made use of niche coevolution to complete the population catastrophe, thus to complete the population mutation by Hadamard gate and optimize the QGA by the adaptive quantum rotation gate adjustment strategy which was combined with global search and local search. Accordingly, one half of the watermark embedding was completed after the optimal embedding position. The other part of the host image was performed on the DWT domain. The ergodic matrix [20, 21] was used to generate the key, and the embedding rules were used to embed the other half of the watermark. The watermarking algorithm further improves the ability of the watermark to resist geometric attacks while keeping the balance of the imperceptibility and robustness of the watermark.The rest of this paper is organized as follows. The second section introduces the basic knowledge of watermarking. The third section describes the proposed image digital watermarking model. The fourth section gives the experimental results and discussion. Finally, the fifth section outlines the conclusions and future work.


## 2. Basic Knowledge

In the proposed watermarking scheme, the optimized quantum genetic algorithm is used as the main algorithm to train the embedding position in the DCT domain embedding watermarking process, and the ergodic matrix is used as an important element to generate the key to complete the fusion operation in the DWT domain embedding watermarking process. Consequently, they will be highlighted in the basics section.


### 2.1. Ergodic Matrix

The ergodic matrix was proposed by Zhao et al. [20–22], and its definition and theorem are as follows [18, 23].Definition 1.
GivenQ∈Fn×nq, if ∀v∈Fn×1q\0, Qv,Q2v,…,Qqn−1v just exhausts ∀v∈Fn×1q\0, then Q is called an ergodic matrix over finite field Fq (here 0 = [0 0 … 0]T).Theorem 1.
Q
∈
F
n
×
n
q is an ergodic matrix if and only if the multiplication period of Q on Fq is (qn−1).Theorem 2.
GivenQ∈Fn×nq is an ergodic matrix, there are φ (qn−1) ergodic matrices in 〈Q〉 (φ(x) is an Euler function). These ergodic matrices are generated by the same spanning set.
It can be seen from the above theorem that alln × n ergodic matrices over the finite field Fq have the same number of elements and are larger than any other n × n ergodic matrices. Taking an arbitrarily selected ergodic matrix Q∈F50×50q256 as an example, the grayscale image and histogram of the ergodic matrix are shown in Figure 1.
Figure1(b) is the plot of the histogram for Figure 1(a). Figure 1 shows that the values of the ergodic matrix spread almost evenly. Therefore, it can be used to construct the key for the embedding rules in the DWT domain.Figure 1
Ergodic matrix and its corresponding histogram: (a) randomly selected ergodic matrixQ∈F50×50256; (b) the corresponding histogram.
(a)
(b)


### 2.2. Optimized Quantum Genetic Algorithm (OQGA)

As a combination of quantum computation and genetic algorithm, QGA uses the bit probability amplitude in quantum properties to represent chromosome coding for increasing the state information of chromosomes and uses quantum rotation gate to update populations to achieve the goal of function optimization. In order to reduce the local optimal solution rate and increase the efficiency of finding the optimal solution to better balance the imperceptibility and robustness of the watermarking scheme, this research optimizes the quantum genetic algorithm as follows.


#### 2.2.1. Quantum Bit Coding

Similar to information bits (0 or 1), qubits are also the smallest unit of information. A qubitφ is a linear combination of two ground states 0 and 1. It represents not only 0 or 1 but also any superposition of the two:(1)φ=α0+β1,α2+β2=1,where α and β are complex numbers, corresponding to the probability amplitude of the two ground states.QGA encodes chromosomes in the form of probability amplitudes, and a qubit is given as a vector (2):(2)αiβi.The encoding of a quantum chromosomeS of length n (n qubits) can be denoted as a matrix (3):(3)S=α1α2α3⋯αnβ1β2β3⋯βn,where αi and βi are the probability amplitudes of the two ground states of the i-th qubit. In quantum measurement, by extracting a random number between [0, 1] and comparing it with the size of αi2, the extraction of information encoded in the quantum state is achieved.In this research, the concept of cluster center in the cluster analysis is applied to record the number of possible embedding positions, which was denoted ascenterNum, and the decoding step is calculated according to equation (4).(4)step=log2centerNum.


#### 2.2.2. Fitness Function

Similar to GAs, to determine whether to keep these chromosomes or discard them, QGAs also need to define a fitness function to calculate the fitness of each chromosome. Hence, in order to better balance the relationship between watermark robustness and imperceptibility, the scheme defines the fitness function as follows:(5)F=C×PSNR+1−C×NC,where PSNR is the peak signal-to-noise ratio of the unattacked original host image and the watermarked host image; NC is the normalized correlation coefficient of the extracted watermark and the original watermark after attack; and C is a constant coefficient which is to balance the Transparency and Robustness and can be modified according the experimental results, and it is 0.0217 in the experiments.PSNR is calculated by using the following equation:(6)PSNR=20×log10MAXIMSEdB,MSE=1Ih×Iw∑x=1Ih∑y=1IwIx,y−I′x,y2,where MSE represents the mean square error between the original image and the processed image; Ih and Iw represent the length and width of the original image, respectively; Ix,y and I′x,y represent the pixel values of the original image and the watermarked image, respectively; and MAXI represents the maximum value of the original image color. If there are 8 sampling points, MAXI is 255.NC is calculated by using the following equation:(7)NC=1Wh×Ww∑x=1Wh∑y=1WwWx,y×W′x,y,where Wh and Ww represent the length and width of the original watermark image, respectively, and Wx,y and W′x,y represent the pixel values of the original watermark image and the extracted watermark image, respectively.


#### 2.2.3. Quantum Rotation Gate Adjustment Strategy

In QGA, the population is updated by the quantum rotation gate, and the probability amplitude of each qubit on the chromosome is converged to 0 or 1, thus obtaining the optimal solution.The traditional quantum rotation gate is expressed as follows:(8)UΔθ=cosΔθ−sinΔθsinΔθcosΔθ.Thei-th qubit adjusted after the traditional quantum rotation gate is given in scalar product (9):(9)αiβi=cosϕsinϕ,αi′βi′=cosΔϕ−sinϕsinϕcosΔϕαiβi=cosϕ+Δθsinϕ+Δθ,where αiβi′ represents the i-th qubit in the chromosome, αi2+βi2=1, Δθ is the rotation angle, and φ is the amplitude of the original qubit.At present, the adjustment strategies of quantum rotation gate can be divided into two types: fixed rotation angle adjustment [24] and adaptive dynamic adjustment [25, 26]. The former adjustment strategy sacrifices the accuracy of the algorithm and increases the convergence speed while the latter has the disadvantage that the rotation angle is too small, and hence it is easy to fall into the local optimal solution. Moreover, because the evolution rules of the population are not considered, the convergence is difficult and slow. To overcome the above imperfections, the proposed scheme optimizes the quantum rotation gate adjustment strategy based on the laws of biological evolution. In general, the optimized strategy is completed by two parts: global search and local search.Due to the rich individual diversity in the early stage of population evolution, the global search can be quickly completed through adaptive adjustment, so that the optimal solution of the current generation can constantly approach the global optimal solution. Hence, at the beginning, the algorithm is to implement the global search according to the following equations:(10)Δθ=π1−FbFre−g/gmax,(11)αi′=αicosΔθ−βisinΔθ,(12)θi=sgnbi−αiαi′−αiΔθ,where b is the chromosome coding sequence with the highest fitness of the current algebra, r is the chromosome coding sequence in the current population, F is the fitness function, g is the current algebra, gmax is the total number of iterations, and sgn is the Heaviside function.The evolution rate of the population increases while the individual fitness decreases. When the fitness is close to the optimal value of the current generation, given the rotation angleΔθ of each qubit is large, the individual will be near the optimal solution and will not converge, and then it will easily fall into the local minimum. At this point, a small fixed rotation angle shall be adopted to achieve a local search. The fitness stagnation algebra is used to determine whether the fitness is close to the optimal fitness, and if so, the local search is implemented according to the following equation:(13)θi=sgnFr−Fbαi′−bi×αiβiΔθe−g/gmax.


#### 2.2.4. Quantum Mutation

In order to avoid the premature phenomenon and the local optimal solution, we replace the traditional quantum nongate with the Hadamard gate to perform mutation on the basis of the rotation angleθi after the quantum rotation gate is adjusted. The mutation process is realized by the following equation:(14)121212−12cosθisinθi=cosπ4−θisinπ4−θi.


#### 2.2.5. Quantum Catastrophe

Although QGA has strong adaptability in solving practical problems, its evolutionary process is irreversible. As a result, when the evolution is stagnant, the population is greatly disturbed by the quantum catastrophe, and thus it effectively makes the population escape from the local extremum. The quantum catastrophic operation of the algorithm in this study is to reinitialize all individuals in the population, and each qubit is reassigned according to the niche coevolution strategy of the following equation:(15)αkβk=i/N1−i/N,where N represents the population size, i represents the i-th population, and k represents the k-th qubit. It shall be noted that in the proposed scheme, the qubits were initialized to 1/2 when the population was initialized for the first time. And if the evolution is stagnant, the niche coevolution strategy will be carried out for population catastrophe to find the best chromosome.


## 2.1. Ergodic Matrix

The ergodic matrix was proposed by Zhao et al. [20–22], and its definition and theorem are as follows [18, 23].Definition 1.
GivenQ∈Fn×nq, if ∀v∈Fn×1q\0, Qv,Q2v,…,Qqn−1v just exhausts ∀v∈Fn×1q\0, then Q is called an ergodic matrix over finite field Fq (here 0 = [0 0 … 0]T).Theorem 1.
Q
∈
F
n
×
n
q is an ergodic matrix if and only if the multiplication period of Q on Fq is (qn−1).Theorem 2.
GivenQ∈Fn×nq is an ergodic matrix, there are φ (qn−1) ergodic matrices in 〈Q〉 (φ(x) is an Euler function). These ergodic matrices are generated by the same spanning set.
It can be seen from the above theorem that alln × n ergodic matrices over the finite field Fq have the same number of elements and are larger than any other n × n ergodic matrices. Taking an arbitrarily selected ergodic matrix Q∈F50×50q256 as an example, the grayscale image and histogram of the ergodic matrix are shown in Figure 1.
Figure1(b) is the plot of the histogram for Figure 1(a). Figure 1 shows that the values of the ergodic matrix spread almost evenly. Therefore, it can be used to construct the key for the embedding rules in the DWT domain.Figure 1
Ergodic matrix and its corresponding histogram: (a) randomly selected ergodic matrixQ∈F50×50256; (b) the corresponding histogram.
(a)
(b)


## 2.2. Optimized Quantum Genetic Algorithm (OQGA)

As a combination of quantum computation and genetic algorithm, QGA uses the bit probability amplitude in quantum properties to represent chromosome coding for increasing the state information of chromosomes and uses quantum rotation gate to update populations to achieve the goal of function optimization. In order to reduce the local optimal solution rate and increase the efficiency of finding the optimal solution to better balance the imperceptibility and robustness of the watermarking scheme, this research optimizes the quantum genetic algorithm as follows.


### 2.2.1. Quantum Bit Coding

Similar to information bits (0 or 1), qubits are also the smallest unit of information. A qubitφ is a linear combination of two ground states 0 and 1. It represents not only 0 or 1 but also any superposition of the two:(1)φ=α0+β1,α2+β2=1,where α and β are complex numbers, corresponding to the probability amplitude of the two ground states.QGA encodes chromosomes in the form of probability amplitudes, and a qubit is given as a vector (2):(2)αiβi.The encoding of a quantum chromosomeS of length n (n qubits) can be denoted as a matrix (3):(3)S=α1α2α3⋯αnβ1β2β3⋯βn,where αi and βi are the probability amplitudes of the two ground states of the i-th qubit. In quantum measurement, by extracting a random number between [0, 1] and comparing it with the size of αi2, the extraction of information encoded in the quantum state is achieved.In this research, the concept of cluster center in the cluster analysis is applied to record the number of possible embedding positions, which was denoted ascenterNum, and the decoding step is calculated according to equation (4).(4)step=log2centerNum.


### 2.2.2. Fitness Function

Similar to GAs, to determine whether to keep these chromosomes or discard them, QGAs also need to define a fitness function to calculate the fitness of each chromosome. Hence, in order to better balance the relationship between watermark robustness and imperceptibility, the scheme defines the fitness function as follows:(5)F=C×PSNR+1−C×NC,where PSNR is the peak signal-to-noise ratio of the unattacked original host image and the watermarked host image; NC is the normalized correlation coefficient of the extracted watermark and the original watermark after attack; and C is a constant coefficient which is to balance the Transparency and Robustness and can be modified according the experimental results, and it is 0.0217 in the experiments.PSNR is calculated by using the following equation:(6)PSNR=20×log10MAXIMSEdB,MSE=1Ih×Iw∑x=1Ih∑y=1IwIx,y−I′x,y2,where MSE represents the mean square error between the original image and the processed image; Ih and Iw represent the length and width of the original image, respectively; Ix,y and I′x,y represent the pixel values of the original image and the watermarked image, respectively; and MAXI represents the maximum value of the original image color. If there are 8 sampling points, MAXI is 255.NC is calculated by using the following equation:(7)NC=1Wh×Ww∑x=1Wh∑y=1WwWx,y×W′x,y,where Wh and Ww represent the length and width of the original watermark image, respectively, and Wx,y and W′x,y represent the pixel values of the original watermark image and the extracted watermark image, respectively.


### 2.2.3. Quantum Rotation Gate Adjustment Strategy

In QGA, the population is updated by the quantum rotation gate, and the probability amplitude of each qubit on the chromosome is converged to 0 or 1, thus obtaining the optimal solution.The traditional quantum rotation gate is expressed as follows:(8)UΔθ=cosΔθ−sinΔθsinΔθcosΔθ.Thei-th qubit adjusted after the traditional quantum rotation gate is given in scalar product (9):(9)αiβi=cosϕsinϕ,αi′βi′=cosΔϕ−sinϕsinϕcosΔϕαiβi=cosϕ+Δθsinϕ+Δθ,where αiβi′ represents the i-th qubit in the chromosome, αi2+βi2=1, Δθ is the rotation angle, and φ is the amplitude of the original qubit.At present, the adjustment strategies of quantum rotation gate can be divided into two types: fixed rotation angle adjustment [24] and adaptive dynamic adjustment [25, 26]. The former adjustment strategy sacrifices the accuracy of the algorithm and increases the convergence speed while the latter has the disadvantage that the rotation angle is too small, and hence it is easy to fall into the local optimal solution. Moreover, because the evolution rules of the population are not considered, the convergence is difficult and slow. To overcome the above imperfections, the proposed scheme optimizes the quantum rotation gate adjustment strategy based on the laws of biological evolution. In general, the optimized strategy is completed by two parts: global search and local search.Due to the rich individual diversity in the early stage of population evolution, the global search can be quickly completed through adaptive adjustment, so that the optimal solution of the current generation can constantly approach the global optimal solution. Hence, at the beginning, the algorithm is to implement the global search according to the following equations:(10)Δθ=π1−FbFre−g/gmax,(11)αi′=αicosΔθ−βisinΔθ,(12)θi=sgnbi−αiαi′−αiΔθ,where b is the chromosome coding sequence with the highest fitness of the current algebra, r is the chromosome coding sequence in the current population, F is the fitness function, g is the current algebra, gmax is the total number of iterations, and sgn is the Heaviside function.The evolution rate of the population increases while the individual fitness decreases. When the fitness is close to the optimal value of the current generation, given the rotation angleΔθ of each qubit is large, the individual will be near the optimal solution and will not converge, and then it will easily fall into the local minimum. At this point, a small fixed rotation angle shall be adopted to achieve a local search. The fitness stagnation algebra is used to determine whether the fitness is close to the optimal fitness, and if so, the local search is implemented according to the following equation:(13)θi=sgnFr−Fbαi′−bi×αiβiΔθe−g/gmax.


### 2.2.4. Quantum Mutation

In order to avoid the premature phenomenon and the local optimal solution, we replace the traditional quantum nongate with the Hadamard gate to perform mutation on the basis of the rotation angleθi after the quantum rotation gate is adjusted. The mutation process is realized by the following equation:(14)121212−12cosθisinθi=cosπ4−θisinπ4−θi.


### 2.2.5. Quantum Catastrophe

Although QGA has strong adaptability in solving practical problems, its evolutionary process is irreversible. As a result, when the evolution is stagnant, the population is greatly disturbed by the quantum catastrophe, and thus it effectively makes the population escape from the local extremum. The quantum catastrophic operation of the algorithm in this study is to reinitialize all individuals in the population, and each qubit is reassigned according to the niche coevolution strategy of the following equation:(15)αkβk=i/N1−i/N,where N represents the population size, i represents the i-th population, and k represents the k-th qubit. It shall be noted that in the proposed scheme, the qubits were initialized to 1/2 when the population was initialized for the first time. And if the evolution is stagnant, the niche coevolution strategy will be carried out for population catastrophe to find the best chromosome.


## 2.2.1. Quantum Bit Coding

Similar to information bits (0 or 1), qubits are also the smallest unit of information. A qubitφ is a linear combination of two ground states 0 and 1. It represents not only 0 or 1 but also any superposition of the two:(1)φ=α0+β1,α2+β2=1,where α and β are complex numbers, corresponding to the probability amplitude of the two ground states.QGA encodes chromosomes in the form of probability amplitudes, and a qubit is given as a vector (2):(2)αiβi.The encoding of a quantum chromosomeS of length n (n qubits) can be denoted as a matrix (3):(3)S=α1α2α3⋯αnβ1β2β3⋯βn,where αi and βi are the probability amplitudes of the two ground states of the i-th qubit. In quantum measurement, by extracting a random number between [0, 1] and comparing it with the size of αi2, the extraction of information encoded in the quantum state is achieved.In this research, the concept of cluster center in the cluster analysis is applied to record the number of possible embedding positions, which was denoted ascenterNum, and the decoding step is calculated according to equation (4).(4)step=log2centerNum.


## 2.2.2. Fitness Function

Similar to GAs, to determine whether to keep these chromosomes or discard them, QGAs also need to define a fitness function to calculate the fitness of each chromosome. Hence, in order to better balance the relationship between watermark robustness and imperceptibility, the scheme defines the fitness function as follows:(5)F=C×PSNR+1−C×NC,where PSNR is the peak signal-to-noise ratio of the unattacked original host image and the watermarked host image; NC is the normalized correlation coefficient of the extracted watermark and the original watermark after attack; and C is a constant coefficient which is to balance the Transparency and Robustness and can be modified according the experimental results, and it is 0.0217 in the experiments.PSNR is calculated by using the following equation:(6)PSNR=20×log10MAXIMSEdB,MSE=1Ih×Iw∑x=1Ih∑y=1IwIx,y−I′x,y2,where MSE represents the mean square error between the original image and the processed image; Ih and Iw represent the length and width of the original image, respectively; Ix,y and I′x,y represent the pixel values of the original image and the watermarked image, respectively; and MAXI represents the maximum value of the original image color. If there are 8 sampling points, MAXI is 255.NC is calculated by using the following equation:(7)NC=1Wh×Ww∑x=1Wh∑y=1WwWx,y×W′x,y,where Wh and Ww represent the length and width of the original watermark image, respectively, and Wx,y and W′x,y represent the pixel values of the original watermark image and the extracted watermark image, respectively.


## 2.2.3. Quantum Rotation Gate Adjustment Strategy

In QGA, the population is updated by the quantum rotation gate, and the probability amplitude of each qubit on the chromosome is converged to 0 or 1, thus obtaining the optimal solution.The traditional quantum rotation gate is expressed as follows:(8)UΔθ=cosΔθ−sinΔθsinΔθcosΔθ.Thei-th qubit adjusted after the traditional quantum rotation gate is given in scalar product (9):(9)αiβi=cosϕsinϕ,αi′βi′=cosΔϕ−sinϕsinϕcosΔϕαiβi=cosϕ+Δθsinϕ+Δθ,where αiβi′ represents the i-th qubit in the chromosome, αi2+βi2=1, Δθ is the rotation angle, and φ is the amplitude of the original qubit.At present, the adjustment strategies of quantum rotation gate can be divided into two types: fixed rotation angle adjustment [24] and adaptive dynamic adjustment [25, 26]. The former adjustment strategy sacrifices the accuracy of the algorithm and increases the convergence speed while the latter has the disadvantage that the rotation angle is too small, and hence it is easy to fall into the local optimal solution. Moreover, because the evolution rules of the population are not considered, the convergence is difficult and slow. To overcome the above imperfections, the proposed scheme optimizes the quantum rotation gate adjustment strategy based on the laws of biological evolution. In general, the optimized strategy is completed by two parts: global search and local search.Due to the rich individual diversity in the early stage of population evolution, the global search can be quickly completed through adaptive adjustment, so that the optimal solution of the current generation can constantly approach the global optimal solution. Hence, at the beginning, the algorithm is to implement the global search according to the following equations:(10)Δθ=π1−FbFre−g/gmax,(11)αi′=αicosΔθ−βisinΔθ,(12)θi=sgnbi−αiαi′−αiΔθ,where b is the chromosome coding sequence with the highest fitness of the current algebra, r is the chromosome coding sequence in the current population, F is the fitness function, g is the current algebra, gmax is the total number of iterations, and sgn is the Heaviside function.The evolution rate of the population increases while the individual fitness decreases. When the fitness is close to the optimal value of the current generation, given the rotation angleΔθ of each qubit is large, the individual will be near the optimal solution and will not converge, and then it will easily fall into the local minimum. At this point, a small fixed rotation angle shall be adopted to achieve a local search. The fitness stagnation algebra is used to determine whether the fitness is close to the optimal fitness, and if so, the local search is implemented according to the following equation:(13)θi=sgnFr−Fbαi′−bi×αiβiΔθe−g/gmax.


## 2.2.4. Quantum Mutation

In order to avoid the premature phenomenon and the local optimal solution, we replace the traditional quantum nongate with the Hadamard gate to perform mutation on the basis of the rotation angleθi after the quantum rotation gate is adjusted. The mutation process is realized by the following equation:(14)121212−12cosθisinθi=cosπ4−θisinπ4−θi.


## 2.2.5. Quantum Catastrophe

Although QGA has strong adaptability in solving practical problems, its evolutionary process is irreversible. As a result, when the evolution is stagnant, the population is greatly disturbed by the quantum catastrophe, and thus it effectively makes the population escape from the local extremum. The quantum catastrophic operation of the algorithm in this study is to reinitialize all individuals in the population, and each qubit is reassigned according to the niche coevolution strategy of the following equation:(15)αkβk=i/N1−i/N,where N represents the population size, i represents the i-th population, and k represents the k-th qubit. It shall be noted that in the proposed scheme, the qubits were initialized to 1/2 when the population was initialized for the first time. And if the evolution is stagnant, the niche coevolution strategy will be carried out for population catastrophe to find the best chromosome.


## 3. Methodology

This section mainly discusses the methodology of how a watermark is embedded and extracted. The meaning of several important notations is shown in Table1.Table 1
The meaning of several important notations.
Notation
Meaning
I
Original grayscale image
W
Watermark image
R
i
DCT coefficient;i is i-th pixel after the block was read in zigzag order
LL
2
Low-frequency subbands of 2-level DWT
I
w
Complete watermarked image
W
e
″
Extracted watermark


### 3.1. Embedding Intensity Trained by SVM Based on Entropy Classification

Entropy-based classification plays a vital role in embedding intensity trained by SVM. Therefore, the concept of entropy and its application in SVM is firstly introduced in this section.


#### 3.1.1. Entropy

The application of information entropy to an image is called image entropy, which can reflect not only the richness of image information but also the comprehensive characteristics of gray level information of image pixel position and gray level distribution in pixel neighborhood. Consequently, this paper applies the image entropy as a measure of the complexity of image texture to the image block classification standard in SVM, which was defined as the entropy accumulation of different gray levels, denoted ase:(i)
Entropye:(16)e=−∑i=0255Pzilog2Pzi,wherezi denotes gray level and Pzi represents the probability of gray level in a subgraph block. When the Pzi changes little, the value of e is large, indicating that the image has a uniform gray distribution and rich texture; when the Pzi difference is large, the e value is small and the image gray distribution is discrete and has few textures.There are three main steps for the embedding intensity trained by SVM based on entropy classification in the following.


#### 3.1.2. Select Training Pictures

Select 50 pairs of pictures as training pictures and divide them into 8 × 8 image blocks.


#### 3.1.3. Determine the Types for the Blocks

Equations (17) and (18) respectively calculate the entropy(e) and variance(v) of the image block. In this way, the texture complexity of the image block was determined, and the image blocks were classified into 3 categories: the smooth area, edge area, and texture area.(i)
Variancev:(17)v=∑i=0255zi−u2Pzi,(18)u=∑i=0255ziPzi,wherezi represents the gray level and Pzi represents the probability of the gray level in the image blocks.On the basis of entropy classification, the brightness (b) of the image block was calculated according to equation (19) and then subdivided again based on brightness.(i)
Brightnessb:(19)b=∑m=18∑n=18Am,n−1282β,whereAm,n is the pixel value of the 8-by-8 image block and β is the scaling factor used to balance the value of b; in the experiment, β = 500.Finally, in the light of the value of entropye, variance v, and brightness b, image blocks were divided into 9 categories according to Table 2.Table 2
Types of the image blocks.
Texture
Brightness
Insensitive
Sensitive
More sensitive
Smooth
T1
T2
T3
Edge
T4
T5
T6
Texture
T7
T8
T9


#### 3.1.4. Train the Embedding Intensity by SVM

The training was performed by thesvmtrain( ) and svmpredict( ) functions in Matlab’s SVM toolbox, and the image block was divided into 9 types according to the criteria and classification types (refer to Section 3.1.2). Through several experiments, the embedding intensity of the 9 image block types is 21.1905, 15.4096, 17.6974, 18.1180, 25.1315, 18.4404, 18.4904, 19.8725, and 18.6187.


### 3.2. Embedding Position Trained by OQGA

Based on the idea in [9], the embedding position in the proposed scheme was improved by quantum genetic algorithm optimized with the adjustment strategy of quantum rotation gate.(1)
Determine algorithm parameters, such as population size, maximum number of iterations, mutation probability, and number of cluster centers(2)
Initialize the population:①
Initialize each qubit of the chromosome to1/2②
Compare the qubitαi2 with the generated random number between [0, 1] to extract the binary coded information in the quantum state③
Change the binary code into decimal code for the position embedding, which was calculated by using equation (4)(3)
Calculate the fitness of each chromosome in the first generation population according to the fitness function equation (5), and store the optimal chromosome, that is, the optimal solution(4)
Update the rotation angle with the optimized quantum rotation door according to equations (10)–(13)(5)
Perform mutation operations to the quantum Hadamard gate by using equation (14), thus to update the population(6)
Calculate the fitness value of the updated population and find the current optimal solution; if it is better than the currently stored optimal solution, replace the optimal solution with the original stored optimal solution and update the global optimal solution(7)
Calculate the stagnation algebra according to the fitness value (refer to equation (15)), and determine whether to carry out the population catastrophe(8)
Determine whether the maximum number of iterations has reached; if yes, jump out of the loop and output the optimal solution; otherwise, go to step 5. to continue


### 3.3. Watermark Embedding

LetI be the original grayscale image with the size of 256 × 256, and W be a binary watermark image with the size of 32 × 32. The watermark embedding model is shown in Figure 2. The detailed procedure is as follows.Figure 2
Watermark embedding process.


#### 3.3.1. Preprocess the Watermark Image

Watermark imageW was scrambled by Arnold algorithm and converted into a sequence W′. Denote the first half of W′ as W1′ and embed it in the DCT transform domain. The last half was then divided into two parts in order and converted into two matrices W2′ and W3′. Both matrices were 1/4 scrambled watermark images to be embedded in the DWT domain. W1′, W2′, and W3′ are denoted by equations (20) and (21).(20)W′=w1,w2,w3,…,w512,W2′,W3′=w11…w1n⋮⋱⋮wm1⋯wmn,where m=n=16.


#### 3.3.2. Preprocess the Host Image

The original grayscale imageI was divided into 8 × 8 nonoverlap blocks represented by a 32 × 32 matrix, totaling 1,024 blocks. According to whether the row number is odd or even, the blocks was divided into two parts. There were 512 blocks with even row number (for example, row 2, 4, 6, …, 32) were taken into the DCT transformation, and each block was sequentially read in zigzag order to form a 512 × 64 matrix I1′ (which is a half of the scrambled host image). I1′ was denoted by the following equation:(21)I1′=I1′i,j1≤i≤512,1≤j≤64,where i is the i-th block and j is the j-th pixel after the block was read in zigzag order.Divide the other 512 blocks with odd row number (for example, row 1, 3, 5, …, 31) into two groups and resynthesize each of them into 128 × 128 matrix (namely, 1/4 scrambled host image). Denote the matrices asI2′ and I3′, respectively.


#### 3.3.3. Embed the Watermark in the DCT Domain

According toI1′, calculate the correlation value of the DCT coefficient according to equation (22), denoted as Ri:(22)Ri=8×8M×N∑k=1512I1′k,iI1′k,1,2≤i≤64,where M×N is the size of the host image. In the experiment, the size was 256 × 256.From equation (23), according to the chromosome generated by OQGA, 512 different blocks were selected from I1′. Each block was embedded with only one watermark bit. The selected block was denoted as Bi presented in the following equation:(23)Bi=I1′i,j1≤i≤512,1≤j≤64,where i is the i-th selected block and j is the j-th position where the selected block embeds the watermark, and it is the phenotype of each gene of the chromosome after decoding. The value of Bi is the embedding position in the i-th block.C
i is used in the watermark extraction, as shown in the following equation:(24)Ci=I1′i,1×Rj1≤i≤512,2≤j≤64,where i is the i-th selected block. It is consistent with that in Bi. Rj is the correlation value of the DCT coefficient.Thek-th watermark information was embedded according to the relationship between Bi and Ci. The embedding rule is as follows:(25)B′k=Ck−AkifCk−Bk<Ak,W1′k=0,Ck+AkifCk−Bk>Ak,W1′k=1,Bk,else,1≤k≤512,where Ak is the embedding intensity corresponding to the k-th selected block trained by the entropy classification-based SVM, Bk is the coefficient sequence after the watermark is embedded, and W1′k is the watermark sequence.When all the watermark information was embedded, restored the newBk to I1′ and denoted as I1w′.


#### 3.3.4. Embed Watermark in the DWT Domain

Perform two levels of 2D DWT on the 1/4 host imageI2′ and I3′, respectively. As a consequence, there are two different low-frequency subbands LL2.According to equation (26), each of the two subbands was divided into blocks of size 2 × 2, and they are denoted as Sk:(26)Sk=LL22i−1,2j−1,LL22i−1,2j,LL22i,2j−1,LL22i,2j,1≤i,j≤16;1≤k≤256.The maximum value ofSk are denoted as LL2max. Therefore, we have two sets constructed by LL2max. For convenience, we denote the sets as Smax(1) and Smax(2), respectively.According to the embedding equation (27), respectively embed 1/4 watermark W2′ and W3′ into Smax(1) and Smax(2) to obtain LL2max′:(27)LL2max′=LL2max−modLL2max,T+0.75T,W∗′=0,LL2max−modLL2max,T+0.25T,W∗′=1,where W∗′ is W2′ or W3′ and T is a 16 × 16 ergodic matrix over the finite field F13.After the watermark was embedded, all 2 × 2 subblocks were recombined to obtain a new LL2, which was performed two-level inverse DWT to obtain the final watermarked imagesI2w′ and I3w′.


#### 3.3.5. Combine the Image

I
1
w
′ was performed inverse zigzag and inverse DCT and then combined with I2w′ and I3w′ to form the complete watermarked image Iw.


#### 3.3.6. Generate the Final Watermarked Image

Last but not least, through OQGA to find the embedding method with the highest fitness value and then obtain the final watermarked imageIw.


### 3.4. Watermark Extraction

Watermark extraction is the inverse operation of watermark embedding. The watermark extraction model is shown in Figure3. The process consisted of two parts: performing extraction in the DCT domain and in the DWT domain, then finally combining the two parts into a complete watermark image.Figure 3
Watermark extraction process.


#### 3.4.1. Preprocess Watermarked Image

The watermarked imageIw was divided into nonoverlap 8 × 8 blocks according to equation (2). A total of 512 blocks, of which the row index was even, were selected for the DCT transformation. Then, each block was read in the order of zigzag to obtain a 512 × 64 matrix Iw1′ (i.e., 1/2 scrambled watermarked image). Iw1′ was denoted as follows:(28)Iw1′=Iw1′i,j∣1≤i≤512,1≤i≤64,where i is the i-th block and j is the j-th pixel after the block was read in zigzag order.A total of 512 blocks of odd rows were divided into 8 × 8 subblocks and respectively combined into two 128 × 128 matrices (namely, 1/4 scrambled watermarked image)I2′ and I3′.


#### 3.4.2. Extract the Watermark after DCT

According toIw1′, the correlation value of the DCT coefficient was calculated by equation (29), which was denoted as Ri′:(29)Ri′=8×8M×N∑k=1512Iw1′k,iIw1′k,1,2≤i≤64,where M × N is the size of the watermarked image.From equation (30), 512 different block coefficient sequences were selected from Iw1′ for watermark extraction from the DCT domain. Only one watermark bit was extracted from each block. The selected block is denoted as Bi′:(30)Bi′=Iw1′i,j1≤i≤512,1≤j≤64,where i denotes the selected i-th block, j denotes the j-th adaptable chromosome generated by OQGA and the value of Bi′ indicates the i-th embedding position in the block.C
i
′ is defined to extract watermark, as shown in the following equation:(31)Ci′=Iw1′i,1×Rj′∣1≤i≤512,2≤j≤64,where i denotes sequentially selected blocks, the order is consistent with the sequence Bi′, and Ri′ denotes the correlation value of the DCT coefficient.The watermark sequenceW1″ was extracted according to the following equation:(32)W1″k=1,ifBk′≥Ck′,0,else.


#### 3.4.3. Watermark Extraction on the DWT Transform Domain

The 1/4 watermarked imagesIw2′ and Iw3′ of size 128 × 128 were, respectively, performed by two-level 2D DWT to obtain two different low-frequency subbands LL2.The two low frequency subbands were further divided into 2 × 2 blocks, which are denoted asSk′, as shown in the following equation:(33)Sk′=LL22i−1,2j−1,LL22i−1,2j,LL22i,2j−1,LL22i,2j,1≤i,j≤16;1≤k≤256.Take the maximum value ofSk′ and denote it as LL2max.Extract watermarksW2″ and W3″ according to the following equation:(34)W∗″=1,ifmodLL2max,T<0.5T,0,else,where W∗″ is W2″ or W3″.


#### 3.4.4. Watermark Image Synthesis

Combine the subwatermark imagesW1″, W2″, and W3″ into the complete watermark W″. Then, perform inverse Arnold transformation to obtain the original watermark We″.


## 3.1. Embedding Intensity Trained by SVM Based on Entropy Classification

Entropy-based classification plays a vital role in embedding intensity trained by SVM. Therefore, the concept of entropy and its application in SVM is firstly introduced in this section.


### 3.1.1. Entropy

The application of information entropy to an image is called image entropy, which can reflect not only the richness of image information but also the comprehensive characteristics of gray level information of image pixel position and gray level distribution in pixel neighborhood. Consequently, this paper applies the image entropy as a measure of the complexity of image texture to the image block classification standard in SVM, which was defined as the entropy accumulation of different gray levels, denoted ase:(i)
Entropye:(16)e=−∑i=0255Pzilog2Pzi,wherezi denotes gray level and Pzi represents the probability of gray level in a subgraph block. When the Pzi changes little, the value of e is large, indicating that the image has a uniform gray distribution and rich texture; when the Pzi difference is large, the e value is small and the image gray distribution is discrete and has few textures.There are three main steps for the embedding intensity trained by SVM based on entropy classification in the following.


### 3.1.2. Select Training Pictures

Select 50 pairs of pictures as training pictures and divide them into 8 × 8 image blocks.


### 3.1.3. Determine the Types for the Blocks

Equations (17) and (18) respectively calculate the entropy(e) and variance(v) of the image block. In this way, the texture complexity of the image block was determined, and the image blocks were classified into 3 categories: the smooth area, edge area, and texture area.(i)
Variancev:(17)v=∑i=0255zi−u2Pzi,(18)u=∑i=0255ziPzi,wherezi represents the gray level and Pzi represents the probability of the gray level in the image blocks.On the basis of entropy classification, the brightness (b) of the image block was calculated according to equation (19) and then subdivided again based on brightness.(i)
Brightnessb:(19)b=∑m=18∑n=18Am,n−1282β,whereAm,n is the pixel value of the 8-by-8 image block and β is the scaling factor used to balance the value of b; in the experiment, β = 500.Finally, in the light of the value of entropye, variance v, and brightness b, image blocks were divided into 9 categories according to Table 2.Table 2
Types of the image blocks.
Texture
Brightness
Insensitive
Sensitive
More sensitive
Smooth
T1
T2
T3
Edge
T4
T5
T6
Texture
T7
T8
T9


### 3.1.4. Train the Embedding Intensity by SVM

The training was performed by thesvmtrain( ) and svmpredict( ) functions in Matlab’s SVM toolbox, and the image block was divided into 9 types according to the criteria and classification types (refer to Section 3.1.2). Through several experiments, the embedding intensity of the 9 image block types is 21.1905, 15.4096, 17.6974, 18.1180, 25.1315, 18.4404, 18.4904, 19.8725, and 18.6187.


## 3.1.1. Entropy

The application of information entropy to an image is called image entropy, which can reflect not only the richness of image information but also the comprehensive characteristics of gray level information of image pixel position and gray level distribution in pixel neighborhood. Consequently, this paper applies the image entropy as a measure of the complexity of image texture to the image block classification standard in SVM, which was defined as the entropy accumulation of different gray levels, denoted ase:(i)
Entropye:(16)e=−∑i=0255Pzilog2Pzi,wherezi denotes gray level and Pzi represents the probability of gray level in a subgraph block. When the Pzi changes little, the value of e is large, indicating that the image has a uniform gray distribution and rich texture; when the Pzi difference is large, the e value is small and the image gray distribution is discrete and has few textures.There are three main steps for the embedding intensity trained by SVM based on entropy classification in the following.


## 3.1.2. Select Training Pictures

Select 50 pairs of pictures as training pictures and divide them into 8 × 8 image blocks.


## 3.1.3. Determine the Types for the Blocks

Equations (17) and (18) respectively calculate the entropy(e) and variance(v) of the image block. In this way, the texture complexity of the image block was determined, and the image blocks were classified into 3 categories: the smooth area, edge area, and texture area.(i)
Variancev:(17)v=∑i=0255zi−u2Pzi,(18)u=∑i=0255ziPzi,wherezi represents the gray level and Pzi represents the probability of the gray level in the image blocks.On the basis of entropy classification, the brightness (b) of the image block was calculated according to equation (19) and then subdivided again based on brightness.(i)
Brightnessb:(19)b=∑m=18∑n=18Am,n−1282β,whereAm,n is the pixel value of the 8-by-8 image block and β is the scaling factor used to balance the value of b; in the experiment, β = 500.Finally, in the light of the value of entropye, variance v, and brightness b, image blocks were divided into 9 categories according to Table 2.Table 2
Types of the image blocks.
Texture
Brightness
Insensitive
Sensitive
More sensitive
Smooth
T1
T2
T3
Edge
T4
T5
T6
Texture
T7
T8
T9


## 3.1.4. Train the Embedding Intensity by SVM

The training was performed by thesvmtrain( ) and svmpredict( ) functions in Matlab’s SVM toolbox, and the image block was divided into 9 types according to the criteria and classification types (refer to Section 3.1.2). Through several experiments, the embedding intensity of the 9 image block types is 21.1905, 15.4096, 17.6974, 18.1180, 25.1315, 18.4404, 18.4904, 19.8725, and 18.6187.


## 3.2. Embedding Position Trained by OQGA

Based on the idea in [9], the embedding position in the proposed scheme was improved by quantum genetic algorithm optimized with the adjustment strategy of quantum rotation gate.(1)
Determine algorithm parameters, such as population size, maximum number of iterations, mutation probability, and number of cluster centers(2)
Initialize the population:①
Initialize each qubit of the chromosome to1/2②
Compare the qubitαi2 with the generated random number between [0, 1] to extract the binary coded information in the quantum state③
Change the binary code into decimal code for the position embedding, which was calculated by using equation (4)(3)
Calculate the fitness of each chromosome in the first generation population according to the fitness function equation (5), and store the optimal chromosome, that is, the optimal solution(4)
Update the rotation angle with the optimized quantum rotation door according to equations (10)–(13)(5)
Perform mutation operations to the quantum Hadamard gate by using equation (14), thus to update the population(6)
Calculate the fitness value of the updated population and find the current optimal solution; if it is better than the currently stored optimal solution, replace the optimal solution with the original stored optimal solution and update the global optimal solution(7)
Calculate the stagnation algebra according to the fitness value (refer to equation (15)), and determine whether to carry out the population catastrophe(8)
Determine whether the maximum number of iterations has reached; if yes, jump out of the loop and output the optimal solution; otherwise, go to step 5. to continue


## 3.3. Watermark Embedding

LetI be the original grayscale image with the size of 256 × 256, and W be a binary watermark image with the size of 32 × 32. The watermark embedding model is shown in Figure 2. The detailed procedure is as follows.Figure 2
Watermark embedding process.


### 3.3.1. Preprocess the Watermark Image

Watermark imageW was scrambled by Arnold algorithm and converted into a sequence W′. Denote the first half of W′ as W1′ and embed it in the DCT transform domain. The last half was then divided into two parts in order and converted into two matrices W2′ and W3′. Both matrices were 1/4 scrambled watermark images to be embedded in the DWT domain. W1′, W2′, and W3′ are denoted by equations (20) and (21).(20)W′=w1,w2,w3,…,w512,W2′,W3′=w11…w1n⋮⋱⋮wm1⋯wmn,where m=n=16.


### 3.3.2. Preprocess the Host Image

The original grayscale imageI was divided into 8 × 8 nonoverlap blocks represented by a 32 × 32 matrix, totaling 1,024 blocks. According to whether the row number is odd or even, the blocks was divided into two parts. There were 512 blocks with even row number (for example, row 2, 4, 6, …, 32) were taken into the DCT transformation, and each block was sequentially read in zigzag order to form a 512 × 64 matrix I1′ (which is a half of the scrambled host image). I1′ was denoted by the following equation:(21)I1′=I1′i,j1≤i≤512,1≤j≤64,where i is the i-th block and j is the j-th pixel after the block was read in zigzag order.Divide the other 512 blocks with odd row number (for example, row 1, 3, 5, …, 31) into two groups and resynthesize each of them into 128 × 128 matrix (namely, 1/4 scrambled host image). Denote the matrices asI2′ and I3′, respectively.


### 3.3.3. Embed the Watermark in the DCT Domain

According toI1′, calculate the correlation value of the DCT coefficient according to equation (22), denoted as Ri:(22)Ri=8×8M×N∑k=1512I1′k,iI1′k,1,2≤i≤64,where M×N is the size of the host image. In the experiment, the size was 256 × 256.From equation (23), according to the chromosome generated by OQGA, 512 different blocks were selected from I1′. Each block was embedded with only one watermark bit. The selected block was denoted as Bi presented in the following equation:(23)Bi=I1′i,j1≤i≤512,1≤j≤64,where i is the i-th selected block and j is the j-th position where the selected block embeds the watermark, and it is the phenotype of each gene of the chromosome after decoding. The value of Bi is the embedding position in the i-th block.C
i is used in the watermark extraction, as shown in the following equation:(24)Ci=I1′i,1×Rj1≤i≤512,2≤j≤64,where i is the i-th selected block. It is consistent with that in Bi. Rj is the correlation value of the DCT coefficient.Thek-th watermark information was embedded according to the relationship between Bi and Ci. The embedding rule is as follows:(25)B′k=Ck−AkifCk−Bk<Ak,W1′k=0,Ck+AkifCk−Bk>Ak,W1′k=1,Bk,else,1≤k≤512,where Ak is the embedding intensity corresponding to the k-th selected block trained by the entropy classification-based SVM, Bk is the coefficient sequence after the watermark is embedded, and W1′k is the watermark sequence.When all the watermark information was embedded, restored the newBk to I1′ and denoted as I1w′.


### 3.3.4. Embed Watermark in the DWT Domain

Perform two levels of 2D DWT on the 1/4 host imageI2′ and I3′, respectively. As a consequence, there are two different low-frequency subbands LL2.According to equation (26), each of the two subbands was divided into blocks of size 2 × 2, and they are denoted as Sk:(26)Sk=LL22i−1,2j−1,LL22i−1,2j,LL22i,2j−1,LL22i,2j,1≤i,j≤16;1≤k≤256.The maximum value ofSk are denoted as LL2max. Therefore, we have two sets constructed by LL2max. For convenience, we denote the sets as Smax(1) and Smax(2), respectively.According to the embedding equation (27), respectively embed 1/4 watermark W2′ and W3′ into Smax(1) and Smax(2) to obtain LL2max′:(27)LL2max′=LL2max−modLL2max,T+0.75T,W∗′=0,LL2max−modLL2max,T+0.25T,W∗′=1,where W∗′ is W2′ or W3′ and T is a 16 × 16 ergodic matrix over the finite field F13.After the watermark was embedded, all 2 × 2 subblocks were recombined to obtain a new LL2, which was performed two-level inverse DWT to obtain the final watermarked imagesI2w′ and I3w′.


### 3.3.5. Combine the Image

I
1
w
′ was performed inverse zigzag and inverse DCT and then combined with I2w′ and I3w′ to form the complete watermarked image Iw.


### 3.3.6. Generate the Final Watermarked Image

Last but not least, through OQGA to find the embedding method with the highest fitness value and then obtain the final watermarked imageIw.


## 3.3.1. Preprocess the Watermark Image

Watermark imageW was scrambled by Arnold algorithm and converted into a sequence W′. Denote the first half of W′ as W1′ and embed it in the DCT transform domain. The last half was then divided into two parts in order and converted into two matrices W2′ and W3′. Both matrices were 1/4 scrambled watermark images to be embedded in the DWT domain. W1′, W2′, and W3′ are denoted by equations (20) and (21).(20)W′=w1,w2,w3,…,w512,W2′,W3′=w11…w1n⋮⋱⋮wm1⋯wmn,where m=n=16.


## 3.3.2. Preprocess the Host Image

The original grayscale imageI was divided into 8 × 8 nonoverlap blocks represented by a 32 × 32 matrix, totaling 1,024 blocks. According to whether the row number is odd or even, the blocks was divided into two parts. There were 512 blocks with even row number (for example, row 2, 4, 6, …, 32) were taken into the DCT transformation, and each block was sequentially read in zigzag order to form a 512 × 64 matrix I1′ (which is a half of the scrambled host image). I1′ was denoted by the following equation:(21)I1′=I1′i,j1≤i≤512,1≤j≤64,where i is the i-th block and j is the j-th pixel after the block was read in zigzag order.Divide the other 512 blocks with odd row number (for example, row 1, 3, 5, …, 31) into two groups and resynthesize each of them into 128 × 128 matrix (namely, 1/4 scrambled host image). Denote the matrices asI2′ and I3′, respectively.


## 3.3.3. Embed the Watermark in the DCT Domain

According toI1′, calculate the correlation value of the DCT coefficient according to equation (22), denoted as Ri:(22)Ri=8×8M×N∑k=1512I1′k,iI1′k,1,2≤i≤64,where M×N is the size of the host image. In the experiment, the size was 256 × 256.From equation (23), according to the chromosome generated by OQGA, 512 different blocks were selected from I1′. Each block was embedded with only one watermark bit. The selected block was denoted as Bi presented in the following equation:(23)Bi=I1′i,j1≤i≤512,1≤j≤64,where i is the i-th selected block and j is the j-th position where the selected block embeds the watermark, and it is the phenotype of each gene of the chromosome after decoding. The value of Bi is the embedding position in the i-th block.C
i is used in the watermark extraction, as shown in the following equation:(24)Ci=I1′i,1×Rj1≤i≤512,2≤j≤64,where i is the i-th selected block. It is consistent with that in Bi. Rj is the correlation value of the DCT coefficient.Thek-th watermark information was embedded according to the relationship between Bi and Ci. The embedding rule is as follows:(25)B′k=Ck−AkifCk−Bk<Ak,W1′k=0,Ck+AkifCk−Bk>Ak,W1′k=1,Bk,else,1≤k≤512,where Ak is the embedding intensity corresponding to the k-th selected block trained by the entropy classification-based SVM, Bk is the coefficient sequence after the watermark is embedded, and W1′k is the watermark sequence.When all the watermark information was embedded, restored the newBk to I1′ and denoted as I1w′.


## 3.3.4. Embed Watermark in the DWT Domain

Perform two levels of 2D DWT on the 1/4 host imageI2′ and I3′, respectively. As a consequence, there are two different low-frequency subbands LL2.According to equation (26), each of the two subbands was divided into blocks of size 2 × 2, and they are denoted as Sk:(26)Sk=LL22i−1,2j−1,LL22i−1,2j,LL22i,2j−1,LL22i,2j,1≤i,j≤16;1≤k≤256.The maximum value ofSk are denoted as LL2max. Therefore, we have two sets constructed by LL2max. For convenience, we denote the sets as Smax(1) and Smax(2), respectively.According to the embedding equation (27), respectively embed 1/4 watermark W2′ and W3′ into Smax(1) and Smax(2) to obtain LL2max′:(27)LL2max′=LL2max−modLL2max,T+0.75T,W∗′=0,LL2max−modLL2max,T+0.25T,W∗′=1,where W∗′ is W2′ or W3′ and T is a 16 × 16 ergodic matrix over the finite field F13.After the watermark was embedded, all 2 × 2 subblocks were recombined to obtain a new LL2, which was performed two-level inverse DWT to obtain the final watermarked imagesI2w′ and I3w′.


## 3.3.5. Combine the Image

I
1
w
′ was performed inverse zigzag and inverse DCT and then combined with I2w′ and I3w′ to form the complete watermarked image Iw.


## 3.3.6. Generate the Final Watermarked Image

Last but not least, through OQGA to find the embedding method with the highest fitness value and then obtain the final watermarked imageIw.


## 3.4. Watermark Extraction

Watermark extraction is the inverse operation of watermark embedding. The watermark extraction model is shown in Figure3. The process consisted of two parts: performing extraction in the DCT domain and in the DWT domain, then finally combining the two parts into a complete watermark image.Figure 3
Watermark extraction process.


### 3.4.1. Preprocess Watermarked Image

The watermarked imageIw was divided into nonoverlap 8 × 8 blocks according to equation (2). A total of 512 blocks, of which the row index was even, were selected for the DCT transformation. Then, each block was read in the order of zigzag to obtain a 512 × 64 matrix Iw1′ (i.e., 1/2 scrambled watermarked image). Iw1′ was denoted as follows:(28)Iw1′=Iw1′i,j∣1≤i≤512,1≤i≤64,where i is the i-th block and j is the j-th pixel after the block was read in zigzag order.A total of 512 blocks of odd rows were divided into 8 × 8 subblocks and respectively combined into two 128 × 128 matrices (namely, 1/4 scrambled watermarked image)I2′ and I3′.


### 3.4.2. Extract the Watermark after DCT

According toIw1′, the correlation value of the DCT coefficient was calculated by equation (29), which was denoted as Ri′:(29)Ri′=8×8M×N∑k=1512Iw1′k,iIw1′k,1,2≤i≤64,where M × N is the size of the watermarked image.From equation (30), 512 different block coefficient sequences were selected from Iw1′ for watermark extraction from the DCT domain. Only one watermark bit was extracted from each block. The selected block is denoted as Bi′:(30)Bi′=Iw1′i,j1≤i≤512,1≤j≤64,where i denotes the selected i-th block, j denotes the j-th adaptable chromosome generated by OQGA and the value of Bi′ indicates the i-th embedding position in the block.C
i
′ is defined to extract watermark, as shown in the following equation:(31)Ci′=Iw1′i,1×Rj′∣1≤i≤512,2≤j≤64,where i denotes sequentially selected blocks, the order is consistent with the sequence Bi′, and Ri′ denotes the correlation value of the DCT coefficient.The watermark sequenceW1″ was extracted according to the following equation:(32)W1″k=1,ifBk′≥Ck′,0,else.


### 3.4.3. Watermark Extraction on the DWT Transform Domain

The 1/4 watermarked imagesIw2′ and Iw3′ of size 128 × 128 were, respectively, performed by two-level 2D DWT to obtain two different low-frequency subbands LL2.The two low frequency subbands were further divided into 2 × 2 blocks, which are denoted asSk′, as shown in the following equation:(33)Sk′=LL22i−1,2j−1,LL22i−1,2j,LL22i,2j−1,LL22i,2j,1≤i,j≤16;1≤k≤256.Take the maximum value ofSk′ and denote it as LL2max.Extract watermarksW2″ and W3″ according to the following equation:(34)W∗″=1,ifmodLL2max,T<0.5T,0,else,where W∗″ is W2″ or W3″.


### 3.4.4. Watermark Image Synthesis

Combine the subwatermark imagesW1″, W2″, and W3″ into the complete watermark W″. Then, perform inverse Arnold transformation to obtain the original watermark We″.


## 3.4.1. Preprocess Watermarked Image

The watermarked imageIw was divided into nonoverlap 8 × 8 blocks according to equation (2). A total of 512 blocks, of which the row index was even, were selected for the DCT transformation. Then, each block was read in the order of zigzag to obtain a 512 × 64 matrix Iw1′ (i.e., 1/2 scrambled watermarked image). Iw1′ was denoted as follows:(28)Iw1′=Iw1′i,j∣1≤i≤512,1≤i≤64,where i is the i-th block and j is the j-th pixel after the block was read in zigzag order.A total of 512 blocks of odd rows were divided into 8 × 8 subblocks and respectively combined into two 128 × 128 matrices (namely, 1/4 scrambled watermarked image)I2′ and I3′.


## 3.4.2. Extract the Watermark after DCT

According toIw1′, the correlation value of the DCT coefficient was calculated by equation (29), which was denoted as Ri′:(29)Ri′=8×8M×N∑k=1512Iw1′k,iIw1′k,1,2≤i≤64,where M × N is the size of the watermarked image.From equation (30), 512 different block coefficient sequences were selected from Iw1′ for watermark extraction from the DCT domain. Only one watermark bit was extracted from each block. The selected block is denoted as Bi′:(30)Bi′=Iw1′i,j1≤i≤512,1≤j≤64,where i denotes the selected i-th block, j denotes the j-th adaptable chromosome generated by OQGA and the value of Bi′ indicates the i-th embedding position in the block.C
i
′ is defined to extract watermark, as shown in the following equation:(31)Ci′=Iw1′i,1×Rj′∣1≤i≤512,2≤j≤64,where i denotes sequentially selected blocks, the order is consistent with the sequence Bi′, and Ri′ denotes the correlation value of the DCT coefficient.The watermark sequenceW1″ was extracted according to the following equation:(32)W1″k=1,ifBk′≥Ck′,0,else.


## 3.4.3. Watermark Extraction on the DWT Transform Domain

The 1/4 watermarked imagesIw2′ and Iw3′ of size 128 × 128 were, respectively, performed by two-level 2D DWT to obtain two different low-frequency subbands LL2.The two low frequency subbands were further divided into 2 × 2 blocks, which are denoted asSk′, as shown in the following equation:(33)Sk′=LL22i−1,2j−1,LL22i−1,2j,LL22i,2j−1,LL22i,2j,1≤i,j≤16;1≤k≤256.Take the maximum value ofSk′ and denote it as LL2max.Extract watermarksW2″ and W3″ according to the following equation:(34)W∗″=1,ifmodLL2max,T<0.5T,0,else,where W∗″ is W2″ or W3″.


## 3.4.4. Watermark Image Synthesis

Combine the subwatermark imagesW1″, W2″, and W3″ into the complete watermark W″. Then, perform inverse Arnold transformation to obtain the original watermark We″.


## 4. Experimental Results and Analysis

In this paper, the peak signal-to-noise ratio (PSNR) and the normalized correlation coefficient (NC) were used to evaluate the similarity between the original image and the embedded watermark image. The robustness of the watermark scheme was tested after various attacks, while the imperceptibility of the watermark was evaluated by the PSNR value. The higher PSNR value indicates the better imperceptibility of the watermark. The higher the NC value indicates the better extraction effect of the watermark, which means the higher similarity between the extracted watermark and the original watermark.In order to fully verify the validity and feasibility of the proposed watermarking scheme, the three aspects of the watermark, namely, imperceptibility, robustness, and time complexity, are calculated and analyzed in this section.The population size of OQGA was 5 and the number of iterations was 100 in the following experimental analysis.


### 4.1. Imperceptibility Analysis

From the experiments, we found that the imperceptibility of the watermarking scheme is related to the ergodic matrix and the attack types. As a result, before analyzing the imperceptibility of the watermarking scheme, the ranges of the elements in the ergodic matrix were determined, and the influence of different attacks on the imperceptibility of the watermarking scheme was explained.In the experiment, a 256-by-256 pixel gray image was selected as the host image and a 32-by-32 pixel binary image as the original watermark image; see Figure4 as an example.Figure 4
Original images: (a) cover image; (b) watermark.
(a)
(b)In the DWT-based watermark embedding process, the embedding rules of the watermark were related to the ranges of the ergodic matrix, and the imperceptibility and robustness of the watermarking scheme were related to the embedding rules. Therefore, in order to maintain a good balance between robustness and imperceptibility, Figures5 and 6 are plotted, respectively, for the different values of the ergodic matrix T between [1, 100] regarding the following three changes of the parameters: the PSNR value of the watermarked image, the NC value of the extracted watermark under nonattacks, and the NC value of the extracted watermark under rotation attack.Figure 5
NC and PSNR values under nonattacks if the ergodic matrix changes.Figure 6
NC values under rotation attack if the ergodic matrix changes.From Figure5, it is easy to see that the elements in the ergodic matrix are in the range of [10, 66]. The PSNR value of the watermarked image and the NC value of the extracted watermark can be well balanced under nonattacks. The PSNR values are all above 40.83, and the NC values are all above 0.98. It can be seen from Figure 6 that the elements of the ergodic matrix is in the range of [35, 100] and the NC value of the watermark is extracted between [0.90, 0.93] after the rotation attack. After many experiments, we found that if the elements of the ergodic matrix are between [36, 66], the imperceptibility and robustness of the watermarking scheme can be better balanced. Here, in the experimental analysis section, the elements of the ergodic matrix are between [36, 66].In addition, in the DWT-based watermark embedding process, the fitness function of the OQGA (see equation (5)) is related to the NC value of the extracted watermark under various attacks and the PSNR value of the watermarked image under nonattacks. This results in different PSNR values of watermarked images under different attacks. Therefore, in this study, Figure 7 is drawn for illustration of the changes in fitness values caused by various attacks, which leads to the changes of the PSNR value of the watermarked image and the NC value of the extracted watermark even there are no attacks.Figure 7
NC and PSNR values under nonattacks after the fitness value were pretrained under 18 attacks.According to Figure7, the PSNR of the watermarked image is above 41 and the NC value of the extracted watermark is 0.99181 under nonattacks, which means that the watermarking scheme has a good imperceptibility.After determining the range of the elements in the ergodic matrix and explaining the impact of different attacks on the watermark scheme, the performance of the scheme is shown in Table3 for different host images.Table 3
Experimental results on different host images.
Name of host image
Lena
Baboon
Yacht
Gold hill
Barbara
Fruits
Host image
Watermarked image
Extracted watermark
PSNR
41.1284
41.2590
41.9569
42.3813
41.8349
40.5838
NC
0.99181
0.99243
0.98626
0.98565
0.98687
0.98080Table3 shows that if the ergodic matrix and the attack type are determined, the change of the different host images has less influence on the imperceptibility of the watermark scheme. Based on the ranges of the ergodic matrix, the attack types, and the host image, it can be verified that the watermarking scheme proposed in this study performs better than the other schemes in transparency.


### 4.2. Robustness Analysis

As mentioned above, if the ergodic matrix and attack type are determined, different host images have little influence on the value of PSNR and NC. At the same time, all the referenced papers take Lena as the host image. In order to make the analysis more contrast and achieve accuracy, the watermark image with Figure4(a) as the host image was selected in the robustness test. Firstly, 18 attacks were used to perform on the watermarked image. Then, the NC and PSNR calculated by the extracted watermark and the original watermark were obtained after the attacks, as shown in Figure 8.Figure 8
NC and PSNR values after 18 attacks.As can be seen from Figure8, the proposed algorithm is better to resist attacks such as lossy JPEG compression, Gaussian low-pass filtering, image darken, contrast decreasing, salt-pepper noise, rotation, cropping, and image resizing. Among them, the performance of resisting lossy JPEG compression and some geometric attacks is the best.After some attacks, the watermarked images and extracted watermarks are shown in Table4.Table 4
Attacked image and extracted watermark.
Attack type
JPEG (QF = 20)
JPEG (QF = 50)
JPEG (QF = 70)
JPEG (QF = 90)
Gaussian low-pass filtering (3, 0.5)
Attacked image
Extracted watermark
Attack type
Histogram
equalization
Image darken
Contrast decreasing
Guassian noise (0.003)
Crop (1/4)
Attacked image
Extracted watermark
Attack type
Salt-pepper noise (0.005)
Salt-pepper noise (0.01)
Rotation (−0.25)
Rotation (0.25)
Image resize (10)
Attacked image
Extracted watermarkThe NC value of the extracted watermark changes with various attacks is shown in Figure9Figure 9
NC values of extracted watermark changes from the watermarked Lena image under the attacks of various parameters.
(a)
(b)
(c)
(d)To further verify the feasibility of the proposed watermarking scheme, some of the anti-attack capabilities are compared with other algorithms in [9, 27, 28, 29]. It can be found that the proposed algorithm is superior in many aspects.Figure10 shows that a major breakthrough in the proposed algorithm is that under most attacks, especially under JPEG compression and other geometric attacks, the NC values of the extracted watermark are higher than those in [9, 29].Figure 10
Comparison with the papers [9] and [29].From Figure11, it is clear to see that the NC value of the extracted watermark by the proposed algorithm after JPEG attack, rotation attack, and cropping attack is higher and more stable than those in [27]. Therefore, when facing these attacks, the proposed algorithm performs better. Figure 11 clearly shows that under the rotation attack with a small rotation angle, the performance of the proposed algorithm is better than that in [27]. Although the performance of the algorithm is little weak when the rotation angle slightly changes larger, since the image is obviously deformed, the difference between the original image and the attacked image can be reduced by artificial adjustment before watermark extracting. However, subtle attacks are often undetectable by the human eyes, which is the advantage of the proposed scheme. Therefore, by comparison, it can be seen that our scheme has better robustness against geometric attacks.Figure 11
Robustness comparison with the lossy JPEG compression, rotation attack, and cropping attack in [27].Figure12 shows that when rotation angle of rotation attack is between 0 and 45 degrees, the proposed algorithm is superior to that in [28] in most angles. According to the experiment, the algorithm in [28] has the best robustness at 30 degrees, and the NC value of the extracted watermark reaches 0.8915. It can be seen from Figure 13 that the proposed algorithm has the best robustness between the rotation angles 0 and 0.27, and the NC values of the extracted watermarks are all above 0.90. Meanwhile, the performances at other angles are also stable, and the NC values are maintained at around 0.5.Figure 12
The robustness comparison with the rotation attack in [28].Figure 13
Robustness of the proposed algorithm under rotation attack with rotation angle changed by 0.01.


### 4.3. Time Complexity Analysis

In order to determine the time complexity of the watermark embedding process, Table5 has been made according to the profile summary of Matlab. The experiment was run on macOS, and the PC was configured with Intel Core i5 CPU at 3.3 GHz and 8 GB of memory.Table 5
Time cost of some important functions.
Function name
Calls
Total time
% time
Attack
500
96.490 s
55.282
Qgate
495
48.332 s
27.690
Synthesize
501
25.998
14.895
SVM
1
1.072 s
0.614
dwt2
4
0.043 s
0.024
dct2
512
0.035 s
0.020
idwt2
4
0.028 s
0.016Table5 indicates that the time in the watermark embedding process is mainly spent on the three functions of Attack, Qgate, and Synthesize, which together account for 97.86% of the total running time. These three functions are used to calculate the attack function for fitness, the strategy function for quantum gate adjustment, and the composite function for the images. They are all subprocesses of OQGA, so the time complexity of the embedding process mainly depends on OQGA. The time cost of OQGA mainly depends on the population size Sp, the coding length Lc, and the number of the iteration Ng. Accordingly, the watermark embedding process time complexity can be expressed as OSp×Lc×Ng. It is worth noting that since the time cost of OQGA is on Sp, Lc, and Ng, the time complexity of other subfunctions can be ignored. For the reason that the calls and execution times of the above dwt2, dwt2, and dct2 functions do not belong to OQGA, idct2 is not called outside the Synthesize function, except for OQGA. As a sequence, Table 5 is not listed here to avoid duplication.The watermark extraction process is the opposite of the embedding one, except that it does not require SVM and OQGA. That is, most of the time overhead is concentrated on DCT, IDCT, DWT, and IDWT. DCT and IDCT have the same complexity ofOMlogM [30], and DWT and IDWT have the same complexity of OM2L [31], where M refers to the size of the host image and L refers to the length of the filter which is used to implement the decomposition and synthesis of DWT. In summary, the time complexity of the watermark extraction process is OMlogM+M2L.


## 4.1. Imperceptibility Analysis

From the experiments, we found that the imperceptibility of the watermarking scheme is related to the ergodic matrix and the attack types. As a result, before analyzing the imperceptibility of the watermarking scheme, the ranges of the elements in the ergodic matrix were determined, and the influence of different attacks on the imperceptibility of the watermarking scheme was explained.In the experiment, a 256-by-256 pixel gray image was selected as the host image and a 32-by-32 pixel binary image as the original watermark image; see Figure4 as an example.Figure 4
Original images: (a) cover image; (b) watermark.
(a)
(b)In the DWT-based watermark embedding process, the embedding rules of the watermark were related to the ranges of the ergodic matrix, and the imperceptibility and robustness of the watermarking scheme were related to the embedding rules. Therefore, in order to maintain a good balance between robustness and imperceptibility, Figures5 and 6 are plotted, respectively, for the different values of the ergodic matrix T between [1, 100] regarding the following three changes of the parameters: the PSNR value of the watermarked image, the NC value of the extracted watermark under nonattacks, and the NC value of the extracted watermark under rotation attack.Figure 5
NC and PSNR values under nonattacks if the ergodic matrix changes.Figure 6
NC values under rotation attack if the ergodic matrix changes.From Figure5, it is easy to see that the elements in the ergodic matrix are in the range of [10, 66]. The PSNR value of the watermarked image and the NC value of the extracted watermark can be well balanced under nonattacks. The PSNR values are all above 40.83, and the NC values are all above 0.98. It can be seen from Figure 6 that the elements of the ergodic matrix is in the range of [35, 100] and the NC value of the watermark is extracted between [0.90, 0.93] after the rotation attack. After many experiments, we found that if the elements of the ergodic matrix are between [36, 66], the imperceptibility and robustness of the watermarking scheme can be better balanced. Here, in the experimental analysis section, the elements of the ergodic matrix are between [36, 66].In addition, in the DWT-based watermark embedding process, the fitness function of the OQGA (see equation (5)) is related to the NC value of the extracted watermark under various attacks and the PSNR value of the watermarked image under nonattacks. This results in different PSNR values of watermarked images under different attacks. Therefore, in this study, Figure 7 is drawn for illustration of the changes in fitness values caused by various attacks, which leads to the changes of the PSNR value of the watermarked image and the NC value of the extracted watermark even there are no attacks.Figure 7
NC and PSNR values under nonattacks after the fitness value were pretrained under 18 attacks.According to Figure7, the PSNR of the watermarked image is above 41 and the NC value of the extracted watermark is 0.99181 under nonattacks, which means that the watermarking scheme has a good imperceptibility.After determining the range of the elements in the ergodic matrix and explaining the impact of different attacks on the watermark scheme, the performance of the scheme is shown in Table3 for different host images.Table 3
Experimental results on different host images.
Name of host image
Lena
Baboon
Yacht
Gold hill
Barbara
Fruits
Host image
Watermarked image
Extracted watermark
PSNR
41.1284
41.2590
41.9569
42.3813
41.8349
40.5838
NC
0.99181
0.99243
0.98626
0.98565
0.98687
0.98080Table3 shows that if the ergodic matrix and the attack type are determined, the change of the different host images has less influence on the imperceptibility of the watermark scheme. Based on the ranges of the ergodic matrix, the attack types, and the host image, it can be verified that the watermarking scheme proposed in this study performs better than the other schemes in transparency.


## 4.2. Robustness Analysis

As mentioned above, if the ergodic matrix and attack type are determined, different host images have little influence on the value of PSNR and NC. At the same time, all the referenced papers take Lena as the host image. In order to make the analysis more contrast and achieve accuracy, the watermark image with Figure4(a) as the host image was selected in the robustness test. Firstly, 18 attacks were used to perform on the watermarked image. Then, the NC and PSNR calculated by the extracted watermark and the original watermark were obtained after the attacks, as shown in Figure 8.Figure 8
NC and PSNR values after 18 attacks.As can be seen from Figure8, the proposed algorithm is better to resist attacks such as lossy JPEG compression, Gaussian low-pass filtering, image darken, contrast decreasing, salt-pepper noise, rotation, cropping, and image resizing. Among them, the performance of resisting lossy JPEG compression and some geometric attacks is the best.After some attacks, the watermarked images and extracted watermarks are shown in Table4.Table 4
Attacked image and extracted watermark.
Attack type
JPEG (QF = 20)
JPEG (QF = 50)
JPEG (QF = 70)
JPEG (QF = 90)
Gaussian low-pass filtering (3, 0.5)
Attacked image
Extracted watermark
Attack type
Histogram
equalization
Image darken
Contrast decreasing
Guassian noise (0.003)
Crop (1/4)
Attacked image
Extracted watermark
Attack type
Salt-pepper noise (0.005)
Salt-pepper noise (0.01)
Rotation (−0.25)
Rotation (0.25)
Image resize (10)
Attacked image
Extracted watermarkThe NC value of the extracted watermark changes with various attacks is shown in Figure9Figure 9
NC values of extracted watermark changes from the watermarked Lena image under the attacks of various parameters.
(a)
(b)
(c)
(d)To further verify the feasibility of the proposed watermarking scheme, some of the anti-attack capabilities are compared with other algorithms in [9, 27, 28, 29]. It can be found that the proposed algorithm is superior in many aspects.Figure10 shows that a major breakthrough in the proposed algorithm is that under most attacks, especially under JPEG compression and other geometric attacks, the NC values of the extracted watermark are higher than those in [9, 29].Figure 10
Comparison with the papers [9] and [29].From Figure11, it is clear to see that the NC value of the extracted watermark by the proposed algorithm after JPEG attack, rotation attack, and cropping attack is higher and more stable than those in [27]. Therefore, when facing these attacks, the proposed algorithm performs better. Figure 11 clearly shows that under the rotation attack with a small rotation angle, the performance of the proposed algorithm is better than that in [27]. Although the performance of the algorithm is little weak when the rotation angle slightly changes larger, since the image is obviously deformed, the difference between the original image and the attacked image can be reduced by artificial adjustment before watermark extracting. However, subtle attacks are often undetectable by the human eyes, which is the advantage of the proposed scheme. Therefore, by comparison, it can be seen that our scheme has better robustness against geometric attacks.Figure 11
Robustness comparison with the lossy JPEG compression, rotation attack, and cropping attack in [27].Figure12 shows that when rotation angle of rotation attack is between 0 and 45 degrees, the proposed algorithm is superior to that in [28] in most angles. According to the experiment, the algorithm in [28] has the best robustness at 30 degrees, and the NC value of the extracted watermark reaches 0.8915. It can be seen from Figure 13 that the proposed algorithm has the best robustness between the rotation angles 0 and 0.27, and the NC values of the extracted watermarks are all above 0.90. Meanwhile, the performances at other angles are also stable, and the NC values are maintained at around 0.5.Figure 12
The robustness comparison with the rotation attack in [28].Figure 13
Robustness of the proposed algorithm under rotation attack with rotation angle changed by 0.01.


## 4.3. Time Complexity Analysis

In order to determine the time complexity of the watermark embedding process, Table5 has been made according to the profile summary of Matlab. The experiment was run on macOS, and the PC was configured with Intel Core i5 CPU at 3.3 GHz and 8 GB of memory.Table 5
Time cost of some important functions.
Function name
Calls
Total time
% time
Attack
500
96.490 s
55.282
Qgate
495
48.332 s
27.690
Synthesize
501
25.998
14.895
SVM
1
1.072 s
0.614
dwt2
4
0.043 s
0.024
dct2
512
0.035 s
0.020
idwt2
4
0.028 s
0.016Table5 indicates that the time in the watermark embedding process is mainly spent on the three functions of Attack, Qgate, and Synthesize, which together account for 97.86% of the total running time. These three functions are used to calculate the attack function for fitness, the strategy function for quantum gate adjustment, and the composite function for the images. They are all subprocesses of OQGA, so the time complexity of the embedding process mainly depends on OQGA. The time cost of OQGA mainly depends on the population size Sp, the coding length Lc, and the number of the iteration Ng. Accordingly, the watermark embedding process time complexity can be expressed as OSp×Lc×Ng. It is worth noting that since the time cost of OQGA is on Sp, Lc, and Ng, the time complexity of other subfunctions can be ignored. For the reason that the calls and execution times of the above dwt2, dwt2, and dct2 functions do not belong to OQGA, idct2 is not called outside the Synthesize function, except for OQGA. As a sequence, Table 5 is not listed here to avoid duplication.The watermark extraction process is the opposite of the embedding one, except that it does not require SVM and OQGA. That is, most of the time overhead is concentrated on DCT, IDCT, DWT, and IDWT. DCT and IDCT have the same complexity ofOMlogM [30], and DWT and IDWT have the same complexity of OM2L [31], where M refers to the size of the host image and L refers to the length of the filter which is used to implement the decomposition and synthesis of DWT. In summary, the time complexity of the watermark extraction process is OMlogM+M2L.


## 5. Conclusion and Future Work

In order to balance the imperceptibility and robustness of the watermarking scheme, the robustness against various attacks shall be improved as much as possible. In this paper, a new adaptive robust blind watermarking algorithm was proposed not only based on the dual transform domain of DCT and DWT but also combined with optimized QGA and key fusion. On the DCT domain, adaptive embedding intensity was trained by SVM based on the entropy classification, and the optimal embedding positions were trained by optimizing QGA to reduce the possibility of QGA falling into the local optimal solutions, so as to obtain a high quality watermarked image and improve the robustness of the watermark scheme. In the DCT domain, key fusion was achieved by applying the ergodic matrix in the embedded rules to achieve a good balance between watermark robustness and imperceptibility; at the same time, the ability of watermark algorithm to resist geometric attacks is further improved. Through experiments, it is found that in the face of various attacks (e.g., lossy JPEG compression, image darken, histogram equalization, rotation, cropping, and noise), the proposed watermarking algorithm outperforms the most existing algorithms.As a result, the proposed scheme can perform well in the application of covert identification and copyright protection. However, there are some imperfections in the following aspects:(1)
The proposed algorithm performs poorly against attacks such as image brighten, filtering, and Guassian noise.(2)
The proposed OQGA has a large time complexity, which makes it difficult to be applied in practice.In order to remedy and improve these deficiencies, a lot of work needs to be done in the future, such as to optimize the time complexity of the algorithm by using a more concise and efficient coding method in the quantum genetic algorithm optimized at the present stage, to optimize the SVM to achieve true adaptive intensity training to further improve image quality by using several other approaches that could be found from literature to perform better against rotation attacks, especially at large angles, and to improve the performance of the proposed algorithm in those attack types by combining different transform domains with other machine learning algorithms [32].



---

*Source: 7817809-2019-10-28.xml*