1.
Introduction
Reaction-diffusion equations are a fundamental part in modelling the spread of biological populations. These equations were proposed in 1937 in papers by Fisher [1] and Kolmogorov et al. [2]. They are based on the following equation:
This equation represents the change of the amount of cells u=u(x,t) in time t and space x, for a diffusion term D∈R and a proliferation rate ρ∈R. The so-called Fisher's Equation (1.1) has been extensively used in population dynamics studies [3,4] and has been deeply analysed in the literature [5,6] in relation to their solutions and travelling waves. This equation has also been studied in other fields from mathematical biology in general [7,8,9,10]. The Fisher's Equation and its extensions are a family of reaction-diffusion models arising prominently also in cancer modelling [11,12], applications to brain tumour dynamics [13], in the description of propagating crystallisation/polymerisation fronts [14], chemical kinetics [15], geochemistry [16] and many others fields. We have considered then to review the use of a mathematical tool, the so-called Lie symmetries, on mathematical models based on reaction-diffusion equations. Considering the importance and application of several of these equations, we have focused on those applied to cellular dynamics and tumour invasion characteristics.
Given the difficulty in obtaining exact solutions to many of these equations, mathematicians and researchers have often resorted to numerical analysis in order to obtain insights into the dynamical properties of a system. This type of analysis, however, lacks the analytic understanding that exact solutions can provide. This is particularly true for the case of non-linear physical phenomena, which is not as interpretable as linear processes. One of the most extended methods for retrieving exact solutions is the Lie symmetries analysis of differential equations, also called group analysis. The Lie classical method is employed to obtain reductions of a system to ODEs and, if possible, families of exact solutions. It is based on the pioneering work of Lie, and gained popularity in recent decades due to the work of Birkhoff, Sedov and Ovsiannikov [17]. Nowadays, this method is widely employed in several branches of science, mainly belonging to physics and mathematics. Examples of the use of the Lie classical method to find exact solutions can be found in [18,19,20,21].
Processes involving reaction and diffusion may require a generalization of the Fisher's equation in order to be appropriately modeled. These generalizations can also be studied by means of Lie symmetries. For example, the invariance of the generalised Fisher's Equation
was first studied by S. Lie for the case A=1, B=C=0 (classical heat equation), in terms of maximal invariance algebra [22]. The general non-linear heat equation (B=C=0) was classified with Lie symmetries by Ovsiannikov [23]. The case with a source term (B=0) was completely described in [24], and the Lie symmetries of the full equation were later described in [25].
In this work, we focus on generalisations of Fisher's equations with application to biological systems. In particular, we describe the use of Lie symmetry groups to obtain analytic solutions related to tumour dynamics. For instance, Lie symmetries of the density dependent reaction-diffusion equation
were calculated in [26]. The optimal system of one-dimensional subalgebras of the invariant equation was obtained, together with reductions and exact solutions. Here, f(u) is an arbitrary function representing proliferation. The diffusion coefficient g(u) depends on the variable u, with independent variables x and t. Symmetries of the differential equations were also used to obtain non-trivial conservation laws [27]. An extension of this equation to a non-linear multidimensional reaction-diffusion system with variables diffusivities was also considered in [28].
Including an explicit space dependence c(x) in the diffusion coefficient yields a generalised Fisher's equation of the form
Reductions and symmetries of this equation were studied in [29]. This equation arises in a broad range of biological processes [6] and specifically in cancer modeling problems [30]. For example, [11] used this equation to study the motility of cells in the complex geometry of the brain, distinguishing between gray and white matter with different diffusion coefficients. This was also analysed in [13] in order to describe malignancy of gliomas as an invasion of grey matter.
Some variations of these generalised equations may involve the choice of an specific system of coordinates. For example, a particular Fisher's equation with space-dependent diffusion coefficient in cylindrical coordinates is given by
Again, f(u) is an arbitrary function and g(u) represents the diffusion coefficient dependent on variable u, with independent variables x and t. In this case, x is the radial variable with assumption of radial symmetry. Exact solutions for this equation were obtained in [31] by means of symmetry analysis. Particular cases of functions f and g have been considered by Bokhari et al. [32]. This equation also appears in the context of heat conduction problems. An example of the application of Lie symmetries with a power law source term and rectangular, cylindrical or spherical coordinates can be found in [33].
The spatial dependence of the previous equation can be generalised again as follows:
Now the function c(x) accounts for both spatial heterogeneity of the medium and coordinate transformation, with f, g and the independent variables having the same meaning as Eq. (1.5). Lie point symmetries of this equation were studied in [34]. This equation is particularly suited for studies of tumour growth, as shown by the many works that consider particular cases of functions f(u), g(u) and c(x) [7,11,12,13]. Further works related to the general equation include the derivation of non-trivial conservation laws [35] and conservation laws associated to the symmetries for g=kfu and f(u),c(x) arbitrary functions and k∈R [36]. Symmetry reductions and exact solutions obtained with classical and potential symmetries can be found in [37].
One last equation that we will consider here is inspired by a recent proposition that mutations conferring proliferative advantage drive super-exponential growth in tumours [38]. Given that mutations are more likely to occur as tumour size increases, this can be mathematically implemented by including a size-dependent term in the proliferation rate:
where the logistic proliferation function incorporates a new term describing the total mass of the tumour and δ∈R. Since the integral term only depends on t, we can simplify this equation to
with F(t) representing the way in which tumour size influences proliferation. The possibility to derive biologically meaningful exact solutions of this equation was explored in [39].
Considering these particular equations, the structure of this Review is as follows: Firstly, the Lie classical method for the derivation of solutions for differential equations is described in Section 2. Secondly, we apply this method to obtain a group classification for Eqs. (1.3), (1.4), (1.5), (1.6) and (1.8) in Section 3. Finally, in Section 4 we focus on cases with special biological meaning, and then obtain some exact solutions.
2.
Lie symmetries and reductions
Lie classical method is used to determine point symmetries of ordinary and partial differential equations. This group of transformations are able to map solutions of the equation into one another. Lie symmetry of Eqs. (1.3), (1.5), (1.4), (1.8) or (1.6) will be given by generators of the form
These equations would admit a infinitesimal point symmetry whenever
where Δ=Δi, for i=1,...5, are each of the Eqs. in study (1.3), (1.4), (1.5), (1.6) or (1.8):
and pr(2)(v) is the second prolongation of the vector field v:
where
with J=(j1,…,jk),1≤jk≤2 and 1≤k≤2 and u(2) denotes the sets of partial derivatives up to second order [17].
The transformation group associated to the Lie symmetry generator (2.1) with group parameter ϵ is given by
where the identity transformation is
We can solve then the system
with initial conditions
We can then apply the symmetry (2.1) on a solution u(t,x) to any of the Eqs. (2.3). We denote this by u=u(t,x)→u∗=u∗(t,x), this is, solution u is mapped into u∗, with
The so-called characteristic form of the infinitesimal point symmetry (2.1) is defined by
Applying the invariance condition from Eq. (2.2) it yields
for
A system of determining equations for the infinitesimals ξ=ξ(x,t,u), τ=τ(x,t,u) and η=η(x,t,u) is then obtained by means of Eq. (2.12). The corresponding determining system is expanded in the respective papers [26,29,31,34,39]. This method will be used in Section 4 to obtain solutions of each of the Eqs. (1.4), (1.6) and (1.8).
3.
Lie symmetry generators
We focus now in obtaining symmetries from Eqs. (1.3), (1.4), (1.5), (1.6) and (1.8). They admit a Lie point symmetry provided that
where Δ is the equation in study and pr(2)v is the second prolongation of the vector field (2.1). For each Equation, we obtain a set of determining equations for the infinitesimals ξ=ξ(x,t,u), τ=τ(x,t,u) and η=η(x,t,u). Here we present the corresponding symmetries.
3.1. Lie symmetry generators for Eqs. (1.3) and (1.4)
We recall Eq. (1.3) as
whose Lie symmetries were published in [26] as shown in Table 1.
In [29] we presented Eq. (1.4) as a generalisation of Eq. (1.3):
whose corresponding generators for special functions f, g and c are shown in Table 2.
Generators vk from Table 2 stand as follows:
3.2. Lie symmetry generators for Eqs. (1.5) and (1.6)
Considering [31], Eq. (1.5)
yielded the generators present in Table 3.
We also present the corresponding generators for the generalisation of the prior equation, i.e. Eq. (1.6), which we recall as
In [34] the following function α=α(x) is introduced as α(x)=c′(x)c(x), yielding
For arbitrary functions f=f(u), g=g(u), and α, the only symmetry generator admitted by Eq. (3.5) is
Moreover, whenever the function α(x) is constant, Eq. (3.5) also admits the symmetry generator
Considering the case whenever g is not arbitrary, other symmetry generators can be obtained, with
1. g=g0ug1, with g0=±1, g1∈R−{0,−4/3}.
2. g=g0u−4/3, with g0=±1.
3. g=g0eug1, with g0=±1, g1∈R−{0}.
Special functions f and α were considered for each function g presented, yielding extra Lie point symmetries. These results are shown in Tables 4, 5 and 6 for each form of function g, respectively.
Notes:
(1) In this case α,f0 and g1 satisfy
where
(2) Constants c1,c2,c4∈R must verify Eq. (3.8) in relationship to α,f0 and g1.
(3) In this case, α,f0, and g1 satisfy
where H(x) is given by (3.10), and
(4) Analogously to (2) constants c1,c4∈R verify Eq. (3.11) with α,f0, and g1.
(5) In this case α and f0 must verify the following
(6) Parameters α,f0, and g1 must satisfy the condition
where
(7) The constants c1,c5∈R are linked to α,f0, and g1 by condition (3.14).
(8) In this case α,f0, and g1 must satisfy the condition
where H5(x) is given by (3.16), and
(9) The constants c1,c2, and c5 are linked to α,f0 and g1 by condition (3.17).
3.3. Lie symmetry generators for Eq. (1.8)
In our work [39] we presented the corresponding generators for Eq. (1.8), this is,
which are shown in Table 7.
For Cases 2 and 3 from Table 7 we define the following notation.
Case 2. Functions Fi and B are defined as follows:
with F0(t)≠0. As specified in 2, the characteristic form can be written as P=(A1x+A2)ux+B(t,x,u,ut). Here, A1,A2∈R and the function B=B(t,x,u,ut) is written in terms of the derivatives of F:
and
for B0(t)≠0.
Case 3. Functions Gi are defined as follows for f1,f2∈R:
4.
Tumour-related solutions of the studied equations
In this Section we review some biologically meaningful cases for the equations in study, specially those which better represent some cancer features. Equations of this kind have been employed in real data research. For instance, variants of Eqs. (1.4) and (1.6) have been used in problems of interfaces and differential cell motility in the brain [11,40,41] with longitudinal data coming from serial CT scans. Eq. (1.8) was employed in [38] for understanding superexponential growth by means of imaging data from different types of cancer. We focused here on obtaining analytical solutions of Eqs. (1.4), (1.6) and (1.8) by choosing special forms of the general functions included.
4.1. Solutions of a Fisher's Equation whose proliferation term is dependent on density and space
Considering Eq. (1.4) from [29], functions f(u)=f1(u−g2)+f2(g2−u)g1+1 and g(u)=g3(g2−u)g1 from Table 2 and an arbitrary function c(x), we have
This equation has a biological interest in terms of modelling. Specifically, Verhulst's law of growth can be included into the equation to describe cancer cell proliferation dynamics [3,11,40]. The diffusion term is considered as in invasion dynamics for brain cancer [3,13,40]. Focusing on generator from v4a from Table 7, the similarity variable and similarity solution obtained are
as well as the ODE
The changes of variables h(z)=−√v(z), v(z)=eα(z) is made and α′(z) is denoted as α′(z)=w(z). In this work g1 is set as g1=1 as well as
where K4=arctanh(K3g3√4f22K14+K32g32)K1−K2, K2,K3∈R.
This resulted in a family of exact solutions of Eq. (1.4)
depending on paramaters K2, K1≠0, as well as on function c=c(x) from Eq. (4.4).
Parameter g2 can be considered as the carrying capacity, as for the solution (4.5), whenever f1<0, it yields that
this is, eventually in time, the density of cells is limited by g2. This can be observed in Figure 1. Over time, growth function f and diffusion f disappear, as
4.2. Solutions of a Fisher's Equation describing a tumour interface problem
For Eq. (1.6) we consider a special case of applicability to brain cancer, specifically glioma. In a series of papers by Swanson et al. [42,43,44] the Fisher-Kolmogorov equation was adapted to represent the proliferation and diffusion of glioma cells in the brain. In order to investigate its impact on cellularity, hypoxia-induced neoangiogenesis and necrosis, and to account for spacial heterogeneity, the diffusion coefficient was made dependent on space. This represents the distinction between regions of grey and white matter [11,43], which is fundamental to explain macro- and microscopic patterns of growth. Glioma cells tend to migrate along white matter tracts in the brain, and in vitro experiments have shown that white gray matter enhances cell motility [45]. The spatial limitation on cellular proliferation was also included by means of a carrying capacity. In this line, Konukoglu et al. [41] proposed a parameter estimation method for reaction-diffusion models of brain tumours.
Taking into account these biological hypotheses, and following the previously cited works, we can consider a particular case of Eq. (1.6) with proliferation and diffusion terms f and g specified as
where k∈R+ is the proliferation rate, ρ∈R is the diffusion rate, and u∗∈R is the maximum amount of cells that a given volume of tissue can hold (i.e. the carrying capacity of the tissue). Function f(u) is the Verhulst's law of growth, commonly used to model cancer cell proliferation [3,11] and function g(u) also follows usual representations [3,44,46]. With these considerations, Eq. (1.6) would read as follows:
where u(x,t) denotes the density of cells. As explained above, the space-dependent part of the diffusion coefficient is to represent a single interfacial transition region [30]. We thus choose a hyperbolic tangent function to allow for two different levels of migration potential.
This equation can be simplified by setting ρ=1 and making the following change of variables
Eq. (4.10) can then be written as
This equation falls under the second case in Table 4 from Section 3.2, with g0=1, g1=1, f0=k, and f1=−k. In this case, when α(x)=c′(x)c(x) does not satisfy condition (3.11), Eq. (4.11) only admits the additional generator
We then look for a solution of the form
where U(x) is a solution of equation
Setting U=√V and c′(x)/c(x)=c0tanh(x), with c0∈R, then Eq. (4.14) becomes
whose solution is given by associate Legendre functions
Thus, solutions of Eq. (4.10) will have the following form:
Setting c0=2, the transformation
maps Eq. (4.15) into
whose general solution depends on the value of k. For k>0, we have the following solutions:
In order to obtain biologically meaningful solutions, we select Eq. (4.23) with c0=2 and c1=c2=12 and obtain the following family of solutions of Eq. (1.6):
with k being the free parameter.
We now try to provide a biological interpretation of the previous solution. Figures (2) and (3) show the effect of the transition region placed at x0=0, for x∈(−20,20) and t∈(0,100). The family of solutions for different values of k represents a higher cellular density for x<0 and a decreased density for x>0. Both zones represent white and grey matter respectively, as in [11,13]. For a fixed t=t0, the situation on both sides of the interface becomes symmetrical as k→12. When proliferation rate k decreases, the situation changes and the positive region becomes less populated. This effect is more pronounced as k→0. This is more clearly seen when exploring asymptotic behaviour for a fixed x: Density grows faster for higher values of k. Note that, according to expressions (1.4) and (4.9), when u→0 diffusion and proliferation increase. This would be consistent with the fact that, when passing through the interface, diffusion grows [13].
Overall, different values of k yield different behaviours for a tumour crossing an interface. From Figure (3) we get that when k decreases the density recovery rate after the interface is lower, as well as the density minimum value. For higher values, proliferation dominates diffusion, simulating a non-diffusive, proliferative tumour. Also, over time tumour density grows faster. For lower k, the reverse happens, yielding infiltrative but non-proliferative tumours.
4.3. Solutions of Fisher's Equation with a proliferation term involving tumour development
We now focus on Eq. (1.8) in [39], which we recall as
Tumour mass F=F(t) can be considered to behave as a tanh, modelling transition regions [30].
Considering the Case 4 from Table 7, F=F(t) with a tanh form belongs to this case for
This function may represent a tumour mass growing quickly as a tanh, and eventually reaching an upper bound. For Eq. (1.8) and Eq. (4.26) omitting negative signs, v1 and v4∗ are the symmetries obtained, where
The similarity variable and solution obtained are
This yields the following reduction
A particular solution of Eq. (4.29) is
This implies that for F=btanh(bt)+b) (see Eq. (4.26)), a family of exact solutions off Eq. (1.8), is the following
These solutions are in accordance to the ones already found in the previous Sections. The behaviour of solution (4.31) is shown in Figure 5, depending on parameter b, simulated for b∈{1,1.2,1.4,1.6,1.8,2}. This parameter represented the mass influence on the proliferation. In Figure 5 (A), tumour density increases with space x. In general, in Figure 5 (B), stabilisation of the tumour density is observed over time.
5.
Conclusions
In this review, we have examined generalised Fisher's equations that can model biological phenomena. Mathematical biology has faced since its inception the issue of the non-linearity of the systems that it aims to describe. This has been a motivation for the development of new methods in areas like mathematical analysis and partial differential equations, thereby becoming one of the most active areas of mathematical research over the last decades. One particular issue that has attracted the attention of the mathematical community is the derivation of exact solutions, which is challenging in the case of non-linear systems and PDEs. The Lie symmetry method has gained recognition as a tool for the simplification of these systems and has been widely employed for the finding of exact solutions.
We therefore focus on the achievement of biologically meaningful solutions from generalised Fisher's equations applied to cancer modelling. While some of these equations have been related to real data [11,14,38,41,43] we intended here to review specifically a theoretical approach to the mathematical models. We first provided Lie symmetries for a number of equations. Eq. (1.3) described a population with general density-dependent proliferation and diffusion, analysed in [26]. Generalised Fisher's Equation (1.4) included an explicit space dependence in the diffusion term, which was considered as a tool for cancer modelling and cell dynamics in [29]. Moving to cylindrical coordinates, generators of Eq. (1.5) were obtained in [31]. Its generalization (Eq. (1.6)) was recently studied [34]. Finally, we described the effect of temporal, size-dependent variation of the proliferation rate (Eq. (1.8)) [39].
We then applied the classical Lie group method in Section 4 and provided one-parametric families of solutions with biological meaning, especially for Eqs. (1.4), (1.6) and (1.8). This involved the choice of specific forms for functions f(u) and g(u), which was made following known biological processes in tumour dynamics such as uncontrolled proliferation and potential for invasion and metastasis. In these equations we considered a tanh dependence of the diffusion term, which allowed us to study single transition regions and tumour progression at the interface. This is particularly useful for the discussion of the effect of white and grey matter on brain tumour. Finally, we explored the behaviour of tumours when proliferation rate grows in time according to a tanh. These results are an example of the application of Lie symmetries in the field of mathematical oncology, and supports the use of mathematical models as a predictive tool or as a means to understanding tumour growth dynamics.
Acknowledgments
We would like to acknowledge group FQM-201 from Junta de Andalucía. We also would like to acknowledge Profs. Rita Tracinà and Mariano Torrisi from the University of Catania (Italy) and Víctor M. Pérez García from the University of Castilla-La Mancha (Spain) for discussions. This work was partially supported by the Fundación Española para la Ciencia y la Tecnología [UCA PR214], the Asociación Pablo Ugarte (APU, Spain) and Inversión Territorial Integrada de la Provincia de Cádiz [ITI-0038-2019].
Conflict of interest
The authors declare no conflict of interests.