- About this Journal ·
- Abstracting and Indexing ·
- Aims and Scope ·
- Annual Issues ·
- Article Processing Charges ·
- Author Guidelines ·
- Bibliographic Information ·
- Citations to this Journal ·
- Contact Information ·
- Editorial Board ·
- Editorial Workflow ·
- Free eTOC Alerts ·
- Publication Ethics ·
- Recently Accepted Articles ·
- Reviewers Acknowledgment ·
- Submit a Manuscript ·
- Subscription Information ·
- Table of Contents

Computational and Mathematical Methods in Medicine

Volume 2012 (2012), Article ID 674243, 13 pages

http://dx.doi.org/10.1155/2012/674243

## Theorems and Application of Local Activity of CNN with Five State Variables and One Port

^{1}State Key Laboratory of Management and Control for Complex Systems, Institute of Automation, Chinese Academy of Sciences, Beijing 100190, China^{2}Department of Information Science and Electronics Engineering, Zhejiang University, Hangzhou 310027, China^{3}The Department of Electrical, Computer, Software, and Systems Engineering, Embry-Riddle Aeronautical University, Daytona Beach, FL 32114, USA

Received 22 September 2011; Accepted 15 November 2011

Academic Editor: Sheng-yong Chen

Copyright © 2012 Gang Xiong et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

Coupled nonlinear dynamical systems have been widely studied recently. However, the dynamical properties of these systems are difficult to deal with. The local activity of cellular neural network (CNN) has provided a powerful tool for studying the emergence of complex patterns in a homogeneous lattice, which is composed of coupled cells. In this paper, the analytical criteria for the local activity in reaction-diffusion CNN with five state variables and one port are presented, which consists of four theorems, including a serial of inequalities involving CNN parameters. These theorems can be used for calculating the bifurcation diagram to determine or analyze the emergence of complex dynamic patterns, such as chaos. As a case study, a reaction-diffusion CNN of hepatitis B Virus (HBV) mutation-selection model is analyzed and simulated, the bifurcation diagram is calculated. Using the diagram, numerical simulations of this CNN model provide reasonable explanations of complex mutant phenomena during therapy. Therefore, it is demonstrated that the local activity of CNN provides a practical tool for the complex dynamics study of some coupled nonlinear systems.

#### 1. Introduction

Coupled nonlinear dynamical systems have been widely studied in recent years. However, the dynamical properties of these systems are difficult to deal with. Although the research on emergence and complexity has gained much attention during the past decades, the determination, prediction, and control of the complex patterns generated from high-dimensional coupled nonlinear systems are still far from perfect. Nature abounds with complex patterns and structures emerging from homogeneous media, and the local activity is the origin of these complexities [1, 2]. The cellular neural network (CNN), firstly introduced by Chua and Yang [3] as an implementable alternative to fully connected Hopfield neural network, has been widely studied for image processing, robotic, biological versions, and higher brain functions, and so on [3]. Many of the coupled nonlinear systems can be modeled and studied via the CNN paradigm [4]. The local activity proposed by Chua asserts that a wide spectrum of complex behaviors may exist if the cell parameters of the corresponding CNN are chosen in or nearby the edge of chaos [2, 4]. There have been quite a few new methods developed for complex systems [5–8], and local activity has attracted the attention of many researchers. Now, local activity has been successfully applied to the research of complex patterns generated from several CNNs in physical, biological, and chemical domains, such as Fitzhugh-Nagumo equation [9], Brusselator equation [10], Gierer-Meinhart equation [11], Oregonator equation [12], Hodgkin-Huxley equation [13], Van Der Pol equation [14], the biochemical model [15], coupled excitable cell model [16], tumor growth and immune model [17], Lorenz model [18], advanced image processing [19], Rossler equation [20], images analysis [21, 22], data prediction [23], neutron transport equation [24], vision safety [25], retinomorphic model [26], and theory research [27–30], and so forth.

Although Chua presents the main theorem of local activity at a cell equilibrium point [1, 2], it is actually difficult to “test” directly the complex patterns of the high-dimensional coupled nonlinear systems, since the theorem contains no recipe for finding whether a variable actually exists or not. It is necessary to develop some mathematical criteria according to the numbers of the variables and ports; that is the topic addressed in this paper.

The remaining of this paper is organized as follows. The local activity of CNN is introduced in Section 2. A set of theorems for testing the local activity of reaction-diffusion CNN with five state variables and one port are set up in Section 3. As an application of the theorems, a coupled reaction-diffusion CNN of hepatitis B Virus (HBV) mutation-selection model is introduced, aiming at describing HBV mutation in the therapeutic process. The bifurcation diagrams of this CNN are developed and some numerical simulations are presented in Section 4. Concluding remarks are given in Section 5.

#### 2. Local Activity Theory of CNN

The CNN architecture is composed of a two-dimensional array of cells. Each cell is denoted by , where . The dynamics of each cell is given by the equation: where are the state, output, and input variables of the cell, respectively. are the elements of the A-template, the B-template, and threshold, respectively. is the radius of influence sphere. The output is the piece-wise linear function given by Clearly, CNN with different template elements may have different functions.

A vast majority of active homogeneous media that are known to exhibit complexity are modeled by a reaction-diffusion partial differential equation (PDE): where is state variables, is spatial coordinates, is a coupled nonlinear vector function called the kinetic term, and are constants called diffusion coefficients. Replacing the Laplace in above formulation by its discrete version yields where Chua et al. have introduced reaction-diffusion CNN equations: where , denotes the state variable located at a point in three-dimensional space with spatial coordinates. Chua refers to the process of transforming a PDE into a reaction-diffusion CNN [2].

From Chua and his collaborators’ point, PDEs are merely mathematical abstractions of nature, and the concept of a continuum is in fact an idealization of reality. Even the collection of all electrons in a solid does not form a continuum, because much volume separating the electrons from the nucleus represents a vast empty space [2]. Reaction-diffusion CNNs have been used to model some phenomena with important practical backgrounds, which were described by PDEs.

Generally speaking, in a reaction-diffusion CNN, every cell has state variables, but only () state variables couple directly to their nearest neighbors via “reaction-diffusion”. Consequently, each cell has the following state equations: where

The cell equilibrium point of equation (7) can be determined via

The Jacobian matrix at the equilibrium point has the following form: where are called cell parameters and The local state equations at the cell equilibrium point are defined via

*Definition 1. *
is called the admittance matrix at the cell equilibrium point .

Lemma 2. *A reaction-diffusion CNN cell is called locally active at the equilibrium point if and only if, its admittance matrix at satisfies at least one of the following four conditions [4].*(1)* has a pole in .*(2)* for some , where is any real number.*(3)* has a simple pole on the imaginary axis, where its associated residue matrix:
is either a complex number or a negative real number.*(4)* has a multiple pole on the imaginary axis.*

*Definition 3. *The cell equilibrium point is called stable if and only if, all the real parts of eigenvalue of *Jacobian *matrix at the equilibrium point are negative [2].

*Definition 4. *A “reaction-diffusion” CNN with state variables and ports is said to be operating on the “edge of chaos” with respect to an equilibrium point if and only if, is both locally active and stable when [4].

Using the above lemma and definitions, the bifurcation of CNN with respect to an equilibrium point can be divided into three parts: the edge of chaos domains (the locally active and stable domains), the locally active and unstable domains, and the locally passive domains. Numerical simulations indicated that many complex dynamical behaviors, such as oscillatory patterns, chaotic patterns, or divergent patterns, may emerge if the selected cell parameters are located in or nearby the edge of chaos domains.

#### 3. Analytical Criteria for Local Activity of CNN with Five State Variables and One Port

For the reaction-diffusion CNN with five state variables and one port, its local state equations have the form where

The corresponding CNN cell admittance matrix is given by [1]. where are the parameters of ’s.

Theorem 5. *A necessary and sufficient condition for to satisfy condition (1) in Lemma 2 is that , such that , and any one of the following conditions holds.*(1)*.*(2)*, and , where is and orders zero point of and , respectively, where .*

*Proof. *Obviously proved.

Denote

Theorem 6. *Let the following parameters be defined as in Theorem 5, then for some if any one of the following conditions holds.*(1)*. *(2)*. *(3)*. *(4)*. *(5)*,
*(6)*. *(7)*. *(8)*. *(9)*and , for or 2.*(10)*. *(11)*, and , for or 3.*(12)*. *(13)*and , , for or 2.*

*Proof. *, so to satisfy condition (2) in Lemma 2 equals to ,
(1)If , then when is large enough (See (1) of Theorem 6).(2)If , then
Let ,(I) If , then when is large enough (See (2) of Theorem 6).(II)If , then(i)If , then when is large enough (See (3) of Theorem 6).(ii)If ,(a)If , such that (See (4) of Theorem 6).(b)If , solve , we can get .Then when , such that (See (5) of Theorem 6).(iii)If , then .(a)If , then when is large enough (See (6) of Theorem 6).(b)If , then , such that (See (7) of Theorem 6).(III)If , let ,(i)If , then , such that (See (8) of Theorem 6).(ii)If , solve , we can get
Then, for , if , , then , such that (See (9) of Theorem 6).(3)If , let , then . Let , then the above becomes , then are the roots of , are the roots of . If any one of the (10)–(13) of Theorem 6 holds, we can get .So, if any one of conditions (1)–(13) holds, . Satisfies condition (2) in Lemma 2. This completes the proof.

Theorem 7. *For , or 2, let
**
Then satisfies condition (3) of Lemma 2, if any one of the following conditions holds.*(I)*, and any one of the following conditions holds.(1) . (2). (3). (4). *(II)

*, and any one of the following conditions holds.(1)*(III)

*.*(2)*.*(3)*.**.*(IV)

*, orand or and any one of the following conditions holds for , or 2.(1)*

*.*(2)*.**Proof. *Let , obviously, is not a single pole of on the imaginary axis.

If has a simple pole on the imaginary axis, where its associated residue
is either a complex number or a negative real number, then , so , which implies that is not a zero point of is not a removed pole of .(I) If has four poles on the imaginary axis. In this case,. Hence we obtain , , . Then, we can get , , which implies that , . So,
Then, when condition () in Theorem 7 holds, is a complex number or a negative real number. satisfies condition (3) in Lemma 2.(II)If has a simple pole and two conjugate poles on the imaginary axis, and another pole is .In this case, it follows that , and has the form:
which implies that , Therefore,
(1)The residue of at is
Then, we conclude that if , is a negative real number. satisfies condition (3) in Lemma 2 (See (1) of () in Theorem 7).(2)The residue of at isConsequently, we conclude that if (2) or (3) in () in Theorem 7 holds, is either an imaginary number or a negative real number. satisfies condition (3) in Lemma 2.(III)If has a simple pole on the imaginary axis, and the other poles are , it follows that , and has the form
Therefore we obtain that , hence the reside of at is
Then, when , is a negative real number. satisfies condition (3) in Lemma 2 (See () of Theorem 7).(IV)If has two conjugate poles on the imaginary axis, and the other poles are . In this case, has the form
Therefore, we obtain that . Then, . Solving it, we have
It implies that or , and . Then, the residue of at is
Hence, if condition () in Theorem 7 holds, or it is an imaginary number. satisfies condition (3) in Lemma 2.

Therefore, when any one of conditions ()–() holds, satisfies condition (3) in Lemma 2. This completes the proof.

Theorem 8. * has a multiple pole on the imaginary axis if any one of the following conditions holds.*(I)*. *(II)*. *(III)*, and or or .*(IV)*, and any one of the following conditions holds.(1) . (2). *

*Proof. *Let . Obviously, when the conditions ()–() hold, 0 is the multiply poles of .

If has two multiply nonzero poles , then has the form:
which implies that .

If are the multiply poles of , then where

Obviously, when any one of (1)-(2) of () in Theorem 8 holds, , any one of is a multiply pole of .

So, when any one of condition ()–() holds, satisfies condition (4) in Lemma 2. This completes the proof.

When any one of Theorems 5–8 holds, satisfies Lemma 2, which implies that the reaction-diffusion CNN with five state variables and one port at the equilibrium point is active.

These theorems can be implemented by a computer program for calculating the bifurcation diagram of the general corresponding CNN to determine emergence of complex dynamic patterns of the corresponding CNN.

#### 4. Analysis and Simulations of Reaction-Diffusion CNN of HBV Mutation-Selection Model

Life systems consist of locally coupled homogeneous media. Mostly, dynamics of life systems are suitable to be described via locally connected reaction-diffusion CNNs. It may be expected that reaction-diffusion CNN will become a promising candidate for modeling life phenomena.

In Chapter 11 “Timing the emergence of resistance” (Page 110) of the book “Virus dynamic: mathematical principles of immunology and virology” (Oxford university press), Nowak et al. proposed a mathematical model which describes the mutation selection of HBV infection during the therapy [31]: where the five variables— represent the numbers of uninfected cells, infected cells infected by normal virus, normal virus, infected cells infected by mutated virus, and mutant viruses, respectively. is the rate of reproduction of uninfected cells. Uninfected cells die at rate and become infected at rate by normal virus and infected at rate by mutated virus. Infected cells infected by normal and mutated virus are removed at rate and , respectively. Normal virus is produced at rate and removed at rate , mutated virus is produced at rate and removed at rate . is the rate constant describing the probability of mutation of virus (usual ), are positive constants. The model was briefly analyzed in Nowak’s book.

The reaction-diffusion CNN of HBV mutation selection of model has the form: where .

Let equation (38) be zeros () and solve it, we can get the two equilibrium points: where and , stand for the patient’s complete recovery and HBV persistent infection, respectively.

Consequently, the Jacobian matrix at the equilibrium point is

Taking as variables, and , , , , , , and , using Theorems 5–8, we can calculate the bifurcation of the reaction-diffusion CNN model equation (38) at the equilibrium points and at , see Figures 1 and 2.

In Figures 1 and 2, the domains are coded as follows: edge of chaos (locally active and stable) domain (shown red), locally active and unstable domain (shown green) and locally passive domain (shown blue). From Figure 1(a), we can see that the bifurcation at equilibrium point does not exist at the edge of chaos domain.

Take , , , , , ,, and, , we model the dynamic trajectories of equation (37) using MATLAB, see Table 1.

In the following discussions, we select some parameters in No. 8-14 and , . The simulation results are shown in Figures 3, 4, 5, 6, and 7. During the simulation, we reached a new conclusion.

From Table 1 and Figures 3–7, we can conclude that(I) when is smaller (less than 5),(1)these parameters are located in the green domain (the local and unstable domains);(2)regardless of the value of , the dynamic pattern of equation (37) is convergent or divergent depending on initial values;(3)the No. 1, No. 4, and No. 5 variables in equation (37) increase and the No. 2 and No. 3 variables in equation (37) decrease to 0. This means the numbers of the mutant virus and of infected cells infected by mutant virus both increase, and the numbers of normal virus and of the infected cells infected by normal virus both decrease, even to near zero. Also, the No. 1 variable in equation (37), that is, the number of uninfected cells increases as compared to the initial number. All these indicate that the potency is perfect except for some virus mutation. The potency is ideal.(II)When is larger (greater than 5) but less than a threshold value (according to initial values and parameters, for example 12.5 in Figure 6), we can conclude the following.(1)These parameters are located in the red domain (edge of chaos).(2)Regardless of the value of , the dynamic pattern of equation (37) is convergent.(3)The No. 1, No. 2, and No. 3 variables in equation (37) increase and No. 4 and No. 5 variables decrease to 0. This means that the number of uninfected cells, the numbers of the normal virus, and of the cells infected by normal virus all increase. Meanwhile, the numbers of the mutant virus and of the cells infected by mutant virus both decrease, even to near zero. All these imply that the drug cannot clean the normal virus, but can destroy the mutant virus and increase the infection cells. The potency is also ideal.(III)When and greater than a threshold value (according to initial values and parameters),(1)these parameters are located in the red domain (edge of chaos);(2)regardless of the value, the dynamic pattern of equation (37) is convergent;(3)The No. 2 and No. 3 variables in equation (37) increase and No. 1, No. 4, and No. 5 variables decrease, which means the numbers of the mutant virus and of the uninfected cells decrease, and the number of normal virus increases. These imply that although the drug can prevent the mutation of the HBV effectively, it may destroy uninfected cells and the liver. The potency is not ideal.

#### 5. Conclusions

The local activity of CNN has provided a powerful tool for studying the emergence of complex patterns in a homogeneous lattice formed by coupled cells. Based on the local activity principle, the analytic criteria for the local activity in reaction-diffusion CNN with five state variables and one port are set up. The analytical criteria include four theorems, which provide the inequalities involving the parameters of the CNN. The inequalities can be used for calculating the bifurcation diagram to determine emergence of complex dynamic patterns of the reaction-diffusion CNN. As an application example, a reaction-diffusion CNN of HBV mutation-selection model is analyzed and simulated, and the bifurcation diagrams are calculated. Numerical simulations show this CNN model may explain certain complex mutant conditions during the therapy. We conclude that the local activity theory provides a practical tool for the study of the complex dynamics of certain coupled nonlinear systems.

#### Acknowledgment

This work is supported in part by NSFC 70890084, 60921061, 90920305.

#### References

- L. O. Chua, “CNN: a version of complexity,”
*International Journal of Bifurcation and Chaos*, vol. 7, no. 10, pp. 2219–2425, 1997. - L. O. Chua, “Local activity is the origin of complexity,”
*International Journal of Bifurcation and Chaos*, vol. 15, no. 11, pp. 3435–3456, 2005. View at Publisher · View at Google Scholar · View at Scopus - L. O. Chua and L. Yang, “Cellular neural networks: theory and applications,”
*IEEE Transactions on Circuits and Systems*, vol. 35, no. 10, pp. 1273–1290, 1988. View at Publisher · View at Google Scholar · View at Scopus - L. O. Chua, “Passivity and complexity,”
*IEEE Transactions on Circuits and Systems*, vol. 46, no. 1, pp. 71–82, 1999. View at Scopus - S. Y. Chen, H. Y. Tong, and C. Cattani, “Markov models for image labeling,”
*Mathematical Problems in Engineering*, vol. 2012, Article ID 814356, 18 pages, 2012. View at Publisher · View at Google Scholar - S. Y. Chen, W. Huang, and C. Cattani, “Traffic dynamics on complex networks: a survey,”
*Mathematical Problems in Engineering*, vol. 2012, Article ID 732698, 23 pages, 2012. View at Publisher · View at Google Scholar - S. Y. Chen, M. Zhao, G. Wu, C. Yao, and J. Zhang, “Recent advances in morphological cell image analysis,”
*Computational and Mathematical Methods in Medicine*. In press. - S. Y. Chen, Y. Zheng, C. Cattani, and W. Wang, “Modeling of biological intelligence for SCM system optimization,”
*Computational and Mathematical Methods in Medicine*, vol. 2012, Article ID 769702, 10 pages, 2012. View at Publisher · View at Google Scholar - R. Dogaru and L. O. Chua, “Edge of chaos and local activity domain of Fitzhugh-Nagumo,”
*International Journal of Bifurcation and Chaos*, vol. 8, no. 2, pp. 671–693, 1998. - R. Dogaru and L. O. Chua, “Edge of chaos and local activity domain of the Brusselator CNN,”
*International Journal of Bifurcation and Chaos*, vol. 8, no. 6, pp. 1107–1130, 1998. View at Scopus - R. Dogaru and L. O. Chua, “Edge of chaos and local activity domain of the Gierer-Meinhardt CNN,”
*International Journal of Bifurcation and Chaos*, vol. 8, no. 12, pp. 2321–2340, 1998. View at Scopus - L. Q. Min, K. R. Crounse, and L. O. Chua, “Analytical criteria for local activity and applications to the oregonator CNN,”
*International Journal of Bifurcation and Chaos*, vol. 10, no. 1, pp. 25–71, 2000. View at Scopus - L. Q. Min, K. R. Crounse, and L. O. Chua, “Analytical criteria for local activity of reaction-diffusion CNN with four state variables and applications to the Hodgkin-Huxley equation,”
*International Journal of Bifurcation and Chaos*, vol. 10, no. 6, pp. 1295–1343, 2000. View at Scopus - L. Q. Min and G. Chen, “Local activity of the Van der Pol CNN,”
*International Journal of Bifurcation and Chaos*, vol. 14, no. 7, pp. 2211–2222, 2004. View at Publisher · View at Google Scholar · View at Scopus - L. Q. Min and N. Yu, “Analytical criteria for local activity of CNN with two-port and application to biochemical model,”
*Journal of University of Science and Technology Beijing*, vol. 7, no. 4, pp. 305–314, 2000. View at Scopus - L. Q. Min and N. Yu, “Some analytical criteria for local activity of two-port CNN with three or four state variables: analysis and applications,”
*International Journal of Bifurcation and Chaos*, vol. 12, no. 5, pp. 931–963, 2002. View at Publisher · View at Google Scholar · View at Scopus - L. Q. Min, J. Wang, X. Dong, and G. Chen, “Some analytical criteria for local activity of three-port CNN with four state variables: analysis and applications,”
*International Journal of Bifurcation and Chaos*, vol. 13, no. 8, pp. 2189–2239, 2003. View at Publisher · View at Google Scholar · View at Scopus - L. Q. Min and N. Yu, “Application of local activity theory of cellular neural network with two ports to the coupled Lorenz-cell model,”
*Communications in Theoretical Physics*, vol. 37, no. 6, pp. 759–767, 2002. View at Scopus - M. Itoh and L. O. Chua, “Advanced image processing cellular neural networks,”
*International Journal of Bifurcation and Chaos*, vol. 17, no. 4, pp. 1109–1150, 2007. View at Publisher · View at Google Scholar · View at Scopus - Y. Meng, L. Q. Min, and Y. Ji, “Application of local activity theory to rossler CNN model,” in
*Proceedings of the 9th International Conference for Young Computer Scientists (ICYCS '08)*, pp. 2846–2851, November 2008. View at Publisher · View at Google Scholar · View at Scopus - H. Li, X. Liao, C. Li, H. Huang, and C. Li, “Edge detection of noisy images based on cellular neural networks,”
*Communications in Nonlinear Science and Numerical Simulation*, vol. 16, no. 9, pp. 3746–3759, 2011. View at Publisher · View at Google Scholar - W. Borges Sampaio, E. Moraes Diniz, A. Corrĉa Silva, A. Cardoso de Paiva, and M. Gattass, “Detection of masses in mammogram images using CNN, geostatistic functions and SVM,”
*Computers in Biology and Medicine*, vol. 41, no. 8, pp. 653–664, 2011. View at Publisher · View at Google Scholar - S. U. Alver, B. Cuma, and U. N. Osman, “Application of cellular neural network (CNN) to the prediction of missing air pollutant data,”
*Atmospheric Research*, vol. 101, no. 1-2, pp. 314–326, 2011. View at Publisher · View at Google Scholar - A. Pirouzmand and K. Hadad, “Cellular neural network to the spherical harmonics approximation of neutron transport equation in x-y geometry—part I: modeling and verification for time-independent solution,”
*Annals of Nuclear Energy*, vol. 38, no. 6, pp. 1288–1299, 2011. View at Publisher · View at Google Scholar - A. Grabowski, R. A. Kosiński, and M. Dźwiarek, “Vision safety system based on cellular neural networks,”
*Machine Vision and Applications*, vol. 22, no. 3, pp. 581–590, 2011. View at Publisher · View at Google Scholar - J. J. Martínez, J. Garrigós, J. Toledo, E. Fernández, and J. Manuel Ferrández, “Implementation of a CNN-based retinomorphic model on a high performance reconfigurable computer,”
*Neurocomputing*, vol. 74, no. 8, pp. 1290–1297, 2011. View at Publisher · View at Google Scholar - M. Di Marco, M. Grazzini, and L. Pancioni, “Global robust stability criteria for interval delayed full-range cellular neural networks,”
*IEEE Transactions on Neural Networks*, vol. 22, no. 4, pp. 666–671, 2011. View at Publisher · View at Google Scholar - S. Long and D. Xu, “Stability analysis of stochastic fuzzy cellular neural networks with time-varying delays,”
*Neurocomputing*, vol. 74, no. 14-15, pp. 2385–2391, 2011. View at Publisher · View at Google Scholar - L. Wan and Q. Zhou, “Attractor and ultimate boundedness for stochastic cellular neural networks with delays,”
*Nonlinear Analysis. Real World Applications*, vol. 12, no. 5, pp. 2561–2566, 2011. View at Publisher · View at Google Scholar - M. Syed Ali and P. Balasubramaniam, “Global asymptotic stability of stochastic fuzzy cellular neural networks with multiple discrete and distributed time-varying delays,”
*Communications in Nonlinear Science and Numerical Simulation*, vol. 16, no. 7, pp. 2907–2916, 2011. View at Publisher · View at Google Scholar · View at Scopus - M. A. Nowak, et al.,
*Virus Dynamics—Mathematical Principles of Immunology and Virology*, chapter 11, Oxford university press, 2000.