
Citation: Michael Sebastiano, Xiaolei Chu, Fikret Aydin, Leebyn Chong, Meenakshi Dutt. Interactions of Bio-Inspired Membranes with Peptides and Peptide-Mimetic Nanoparticles[J]. AIMS Materials Science, 2015, 2(3): 303-318. doi: 10.3934/matersci.2015.3.303
[1] | Chibuike Chiedozie Ibebuchi . Can synoptic patterns influence the track and formation of tropical cyclones in the Mozambique Channel?. AIMS Geosciences, 2022, 8(1): 33-51. doi: 10.3934/geosci.2022003 |
[2] | Miyuru B Gunathilake, Thamashi Senerath, Upaka Rathnayake . Artificial neural network based PERSIANN data sets in evaluation of hydrologic utility of precipitation estimations in a tropical watershed of Sri Lanka. AIMS Geosciences, 2021, 7(3): 478-489. doi: 10.3934/geosci.2021027 |
[3] | Manogaran Madhiarasan . Long-term wind speed prediction using artificial neural network-based approaches. AIMS Geosciences, 2021, 7(4): 542-552. doi: 10.3934/geosci.2021031 |
[4] | Shuo Yang, Dong Wang, Zeguang Dong, Yingge Li, Dongxing Du . ANN prediction of the CO2 solubility in water and brine under reservoir conditions. AIMS Geosciences, 2025, 11(1): 201-227. doi: 10.3934/geosci.2025009 |
[5] | Brian E. Bunker, Jason A. Tullis, Jackson D. Cothren, Jesse Casana, Mohamed H. Aly . Object-based Dimensionality Reduction in Land Surface Phenology Classification. AIMS Geosciences, 2016, 2(4): 302-328. doi: 10.3934/geosci.2016.4.302 |
[6] | Paolo Dell’Aversana, Gianluca Gabbriellini, Alfonso Iunio Marini, Alfonso Amendola . Application of Musical Information Retrieval (MIR) Techniques to Seismic Facies Classification. Examples in Hydrocarbon Exploration. AIMS Geosciences, 2016, 2(4): 413-425. doi: 10.3934/geosci.2016.4.413 |
[7] | Leszek J. Kaszubowski . Seismic profiling of the sea-bottom in recognition of geotechnical conditions. AIMS Geosciences, 2020, 6(2): 199-230. doi: 10.3934/geosci.2020013 |
[8] | Efthymios Panagiotis, Alena Zhelezova, Irene Rocchi . Analysis of CPTu data from the North Sea region to estimate the frequency of inaccurate pore pressure measurements. AIMS Geosciences, 2025, 11(2): 517-527. doi: 10.3934/geosci.2025021 |
[9] | Pingping Chen, Mingyang Qi, Long Chen . Distributed sensors and neural network driven building earthquake resistance mechanism. AIMS Geosciences, 2022, 8(4): 718-730. doi: 10.3934/geosci.2022040 |
[10] | Konstantinos X Soulis, Evangelos E Nikitakis, Aikaterini N Katsogiannou, Dionissios P Kalivas . Examination of empirical and Machine Learning methods for regression of missing or invalid solar radiation data using routine meteorological data as predictors. AIMS Geosciences, 2024, 10(4): 939-964. doi: 10.3934/geosci.2024044 |
The elaboration of mathematical models that describe interactions between cancer cells and the immune system is an actively developing field in mathematical medicine. One of the important goals of this modeling is to analyze the ultimate dynamics which is necessary for monitoring the long-term tumor dynamics and, in the most optimistic case, for the creation of tumor eradication algorithms as well. Oncolytic virotherapy of cancer attracts attention of many researchers, see e.g. papers [1,2,3,4,5,6,8,14] and so on. Eftimie et al. in [5] described and studied the seven-dimensional model characterizing behavior of the cancer population in the case when the dual immunization protocol, by vaccine viruses and oncolytic viruses, is applied. They obtained local stability conditions for equilibrium points of their model and numerically investigated multi-stability phenomenon which may exist in the case of four types of equilibrium points. Later, in [6], Eftimie et al. presented the five-dimensional system which is obtained by a simplification of the model from [5]. This system was derived from an experimental protocol of the paper [3] and can be written as
$
⋅x1=rx1(1−x1+x2κ)−a1x1x5b1+x1−a2x1x4b2+x4⋅x2=a1x1x5b1+x1−δx2−a2x2x4b2+x4⋅x3=p1x3x5b3+x5(1−x3m)⋅x4=p2x3x5b3+x5−a3x4−a4x1x4+σ⋅x5=δβx2−γx5.
$
|
(1.1) |
Below the system (1.1) is considered on the nonnegative orthant ${\bf R}_{+, 0}^{5}$ and it is denoted by $x_{1}$/$x_{2}$/$x_{3}$ /$x_{4}$/ $x_{5}$ the density of uninfected tumor cells/ infected tumor cells/ memory cells/ effector immune cells/ oncolytic (Vesicular Stomatitis) virus particles (OVPs) correspondingly. Parameters are supposed to be positive excepting $\sigma $, $\sigma \geq 0\text{, }\; $ and are defined as follows: $r$ is proliferation rate for uninfected tumor cells; $\kappa $ is carrying capacity for tumor cells; $a_{1}$ is infection rate of tumor cells with the oncolytic virus; $b_{1}$ is half-saturation constant for the tumor cells infected with the oncolytic virus; $\delta $ is rate at which the oncolytic virus kills the tumor cells; $a_{2}$ is lysis rate of tumor cells (infected and uninfected) by the immune cells; $b_{2}$ is half-saturation constant for the effector cells that support half the maximum killing rate; $p_{1}$ is proliferation rate of memory cells following secondary encounter with tumor antigens carried by virus particles; $b_{3}$ is half-saturation constant of viral antigens that induce half the maximum proliferation rate of immune cells; $m$ is carrying capacity for memory cells; $p_{2}$ is rate at which memory cells become effector cells following secondary encounter with tumor antigens carried by virus particles; effector cells die at rate $a_{3}$; $a_{4}$ is inactivation rate of effector immune cells by tumor cells; $\beta $ is number of OVPs released from an infected cell, capable of forming plaques; $\gamma $ is decay rate for the concentration of OVPs in the blood; $\sigma $ is the influx rate into the lymphoid tissue at which the naive cells become effector cells, see [5]. It is pointed out that σ has sufficiently small values and may be neglected, [6]. In this paper it is assumed that the rate σ reflects the presence of the influx of immune cells from the external source and the dependence of global behavior of (1.1) is explored with a change in $\sigma $.
Authors of [6] studied only the case $\sigma = 0$. They found equilibrium points and performed local stability analysis. As well, authors of [6] fulfilled numerical studies transient and asymptotic behavior of the system which are based on using bifurcation diagrams. In particular, they showed the existence of periodic orbits and a chaotic attractor for some values of parameters. Nevertheless, the rigorous exploration of long time-dynamics is of interest for better understanding potentialities of cancer virotherapy described by equations (1.1).
The purpose of this paper and its novelty consist in studies of the case of convergence dynamics of (1.1), with $\sigma \geq 0$. In particular, conditions are found under which every trajectory in the nonnegative orthant tends to one of equilibrium points. It is worth to remind that systems with convergence dynamics have been investigated in articles of Hirsch; R.A. Smith; Li and Muldowney, see [7,15,16], in the broader perspective of studies of competitive or cooperative systems; higher dimensional generalizations of Bendixon's criterion and Dulac's criterion. Our results are based on the LaSalle theorem and the localization method of compact invariant sets (LMCIS), see [10,11]. It is worth to mention that the LMCIS has been efficiently utilized in papers [12,13,17,18,19,20] for analysis of ultimate dynamics of various cancer models.
The dynamics of (1.1) has some distinctive features: the unbounded set of equilibrium points, positive Lagrange stability of trajectories and several invariant planes that are useful in this analysis. In the case of a sufficiently large $\sigma >0$ it is found that (1.1) acquires a new global qualitative property consisting in the fact that all trajectories are attracted to the tumor-free plane $x_{1} = x_{2} = 0$. This property may be characterized quantitatively as well. Namely, it is computed a bound $\sigma _{\min }^{(1)}$ called the first eradication bound such that if $\sigma \geq \sigma _{\min }^{(1)}$ then each $\omega $-limit set is located in $x_{1} = x_{2} = 0$ plane. This signifies that the population of tumor cells is under control after sufficiently large time of observations. From the biological point of view, this case may correspond to the immunotherapy injection of some agent like IL-2 which is combined with virotherapy. Although the tumor may be cleared for $\sigma \geq \sigma _{\min }^{(1)}$, but the harmful side effects from an over-activated immune system may eliminate the benefits of tumor clearance. We propose an approach for reducing this risk. We show that there is the sequence of eradication bounds $\{\sigma _{\min }^{(k)}\}, k = 1, 2, .., $with the property that for $\sigma \geq \sigma _{\min }^{(k)}$ each $\omega $-limit set is located in $x_{1} = x_{2} = 0$ plane. The value $\sigma _{\min }^{(k)}$ appears as the minimal positive root of some algebraic equation of degree $k.$ The existence of $\sigma _{\min }^{(k)}$ is established in Theorem 5.1. The formula for $\sigma _{\min }^{(2)}$ is provided; here we have that $\sigma _{\min }^{(1)}>\sigma _{\min }^{(2)}$. It is demonstrated numerically that $\sigma _{\min }^{(2)}>\sigma _{\min }^{(3)}.$ We expect that $\{\sigma _{\min }^{(k)}\}$ is a decreasing sequence. Figs. 1-3 show the tumor eradication dynamics with $\sigma = \sigma _{\min }^{(i)}, i = 1, 2, 3.$
In the present paper ultimate upper bounds are computed for variables of (1.1). The domain defined by these bounds contains all compact invariant sets including a chaotic attractor found in [6] for some values of parameters. Further, it is noted that in the case $\sigma = 0$ tumor-free TF-equilibrium points are always locally unstable [6]. Main results of this paper are
1) Theorem 4.2 containing conditions under which the $\omega $-limit set of any positive half trajectory in ${\bf R}_{+}^{5}$ is either one of TF-equilibrium points or one of tumor-only TO-equilibrium points and
2) Propositions 3;4 and Theorem 5.1 in which we give minimal bounds for $\sigma $ guaranteeing convergence dynamics of any positive half trajectory in ${\bf R}_{+}^{5}$ to one of tumor-only TO-equilibrium points.
Therefore conditions of this paper may be utilized for a prediction of the possible health scenario for a sufficiently long time of observations.
The results of this work are stated in terms of inequalities imposed on parameters of the model. This may be important in applications, because in medical models, as a rule, biological parameters are known in a certain range.
The present paper is organized as follows. Helpful auxiliary results are provided in Section 2 and in Appendix. In Section 3 ultimate dynamics of the system (1.1) is studied via computing ultimate bounds in case $\sigma \geq 0$. These bounds form a domain in which the chaotic attractor is detected for some values of parameters given in [6]. Further, Section 4 contains conditions under which (1.1) possesses convergent dynamics to one of TO- or TF-equilibrium points in case $\sigma = 0$. In Section 5 we describe the case when the immunotherapy is utilized $\sigma >0$. Here we show that the tumor clearance phenomenon may be achieved with a sufficiently large value of $\sigma $. The tumor eradication bounds for $\sigma $ are derived. Section 6 contains concluding remarks. Appendix contains the derivation of the cubic equation for finding $\sigma _{\min }^{(3)}$.
With the goal of describing helpful results the following objects are introduced:
1) a nonlinear system
$ \dot{x} = v (x) $ | (2.1) |
where $v$ is a $C^{1} -$differentiable vector field; $x \in {\bf R}^{n}$ is the state vector;
2) a $C^{1} -$differentiable function $h (x)$ such that $h$ is not the first integral of (2.1);
3) the set $S (h): = \{x \in {\bf R}^{n} \mid L_{v} h (x) = 0\}$ where $L_{v} h$ is the Lie derivative;
4) values $h_{\inf }: = \inf \{h (x) \mid x \in S (h)\}$; $h_{\sup }: = \sup \{h (x) \mid x \in S (h)\}\text{.}$
If it is set up the localization problem of compact invariant sets contained in the set $M \subset {\bf R}^{n}$ then the sets
$
S(h;M):=S(h)∩M={x∈M∣Lvh(x)=0};hinf(M):=inf{h(x)∣x∈S(h;M)};hsup(M):=sup{h(x)∣x∈S(h;M)}.
$
|
are defined.
Assertion 1 [9,10,11] For any $h(x)\in C^{1}({\bf R}^{n})$ all compact invariant sets of the system (2.1) located in $M$ are contained in the localization set $K(h; M)$ defined by the formula $K(h; M) = \{x\in M\mid h_{\inf }(M)\leq h(x)\leq h_{\sup }(M)\}$ as well. If $M\cap S(h) = \varnothing $ then there are no compact invariant sets located in $M$.
The function $h$ mentioned in this assertion is called localizing.
Below the following notations are used: if the system is examined on the positive orthant $M = {\bf R}_{+, 0}^{n}$ then $K(h): = K(h; M)\text{.}\; $ If the system is considered on the domain $M: = {\bf R}_{+, 0}^{n}\cap \{x_{5}>0\}$ then this $M$ is indicated.
Further, it is utilized
Assertion 2 [12]. Let $h$ be a differential function which is bounded below on the set $D\subset {\bf R}^{n}$ and $L_{v}h(x)\mid _{D}\leq -\varepsilon $ for some $\varepsilon >0$. Then the solution $\varphi (x, t)$ of the system (2.1), with initial condition $\varphi (x, 0) = x\in D\text{, }$ leaves the set $D$ in a finite time.
To formulate the following useful result, the nonautonomous system
$ \dot{y} = G (y , t) , y \in {\bf R}^{n} , $ | (2.2) |
is introduced; here $G$ is a $C^{1}-$differentiable vector function. The system (2.2) is called asymptotically autonomous with the limit system (2.1) if $G(x, t)\rightarrow v(x)$, with $t\rightarrow\infty$, for any compact subset in ${\bf R}^{n}$, see e.g. [21].
Theorem 2.1. (Markus), see e.g. [21]. Let us consider systems (2.1) and (2.2) in case $n = 2.$ Now if $\Omega $ is the $\omega $ -limit set of a forward bounded solution $x$ of (2.2) then $\Omega $ either contains equilibria of (2.1) or $\Omega $ is the union of periodic orbits of (2.1).
The Markus theorem is applicable to the nonnegative orthant ${\bf R}_{+, 0}^{2}$ in case when ${\bf R}_{+, 0}^{2}$ is a positively invariant domain.
The planes $x_{1} = 0$; $x_{1} = x_{2} = 0$; $x_{3} = 0$; $x_{3} = m$; $x_{2} = x_{5} = 0$ are invariant.
Further, if $\sigma = 0$ then there are four types of equilibrium points, [6]:
a) TF-equilibrium points are given by ${\bf E}_{1}: = (0, 0, e_{1 3}, 0, 0)^{T}; e_{1 3} \in {\bf R}_{ +, 0}^{1};$
b) TO-equilibrium points are given by ${\bf E}_{2}: = (\kappa, 0, e_{2 3}, 0, 0)^{T}; e_{2 3} \in {\bf R}_{ +, 0}^{1};$
c) a tumor with virus but no immune response TV-equilibrium point is given by ${\bf E}_{3}: = (e_{3 1}, e_{3 2}, 0, 0, e_{3 5})^{T}$, with
$
e31=b1γa1β−γ;e32=rδβ(κ−e31)(b1+e31)a1δβκ+b1rγ+e31rγ;e35=e35γβδ;
$
|
d) an internal tumor virus ITV-equilibrium point is given by ${\bf E}_{4}: = (e_{41}, e_{42}, m, e_{44}, e_{45})^{T}$, with components $e_{4j}, j = 1, 2, 4, 5\text{, }\; $defined by very cumbersome formulas and omitted here.
Since TF-equilibrium points and TO-equilibrium points form an unbounded set the system (1.1) does not possess the bounded positively invariant domain.
It is noted that if $\sigma = 0$ then each TF-equilibrium point $(0, 0, e_{1 3}, 0, 0)^{T}$ is connected with TO-equilibrium point $(\kappa, 0, e_{1 3}, 0, 0)^{T}$ by a heteroclinic orbit which is located in the invariant plane $x_{i} = 0, i = 2, 4, 5;x_{3} = e_{13}$ and is an arc directed from $(0, 0, e_{13}, 0, 0)^{T}$ to $(\kappa, 0, e_{13}, 0, 0)^{T}$. So heteroclinic orbits cover half strip $x_{3} \geq 0, x_{2} = x_{4} = x_{5} = 0.$
Next, upper ultimate bounds are got for all state variables of the system (1.1) with help of computing the localization domain containing all compact invariant sets.
1. Bound $x_{1\max }$. The localizing function $h_{1} = x_{1}$ is applied. Then
$
S(h1)∩{x1>0}={r(1−x1+x2κ)=a1x5b1+x1+a2x4b2+x4≥0}⊂{r(1−x1+x2κ)≥0}.
$
|
Therefore the localization set $K_{1\max }: = \left \{0 \leq x_{1} \leq \kappa : = x_{1\max }\right \}$ is derived.
2) Bound $x_{2\max }\; \text{.}\; $ Here the function $h_{2} = x_{1} +x_{2}$ is utilized. Then
$ L_{f} h_{2} = r x_{1} \left (1 -\frac{x_{1} +x_{2}}{\kappa }\right ) -a_{2} \frac{(x_{1} +x_{2}) x_{4}}{b_{2} +x_{4}} -\delta x_{2} $ |
Therefore the inequality
$ \delta (x_{1} +x_{2}) \leq (\delta +r) x_{1} -\frac{r}{\kappa } x_{1}^{2} $ |
holds on $S (h_{2})\text{.}$ Considering the last inequality within $K_{1\max }$ the following inequality holds
$
\delta (x_{1} +x_{2}) \leq \zeta : = \left \{δκ,δ>rκ(δ+r)24r,δ≤r \right \}\;\text{.}\;
$
|
Thus the localization set $K_{21}: = \left \{h_{2} \leq h_{2\max }: = \zeta \delta ^{ -1}\right \}$ is derived. So the upper bound $x_{2\max }$ given in
$ K_{2\max }: = \left \{x_{2} \leq \zeta \delta ^{ -1}: = x_{2\max }\right \} $ | (3.2) |
is obtained.
3) Bound $x_{5\max }$. Here using $h_{3} = x_{5}$ it is noted that
$ S \left (h_{3}\right ) \cap K_{2\max } \subset \left \{x_{5} \leq \frac{\delta \beta }{\gamma } x_{2\max }\right \} $ |
and so the localization set with the bound $x_{5\max }$ is obtained:
$ K_{5\max }: = \left \{x_{5} \leq x_{5\max }: = \frac{\delta \beta }{\gamma } x_{2\max }\right \} . $ | (3.3) |
4) The bound $x_{3\max }$ is computed with respect to the domain ${\bf R}_{ +, 0}^{5} \cap \{x_{5} >0\}\; \text{.}\; $ The function $h_{4} = x_{3}$ provides the localization set:
$ K_{3\max } (h_{4} , {\bf R}_{ + , 0}^{5} \cap \{x_{5} \gt 0\}): = \{0 \leq x_{3} \leq m ;x_{5} \gt 0\}\;\text{.} $ |
5) The bound $x_{4\max }$ is computed with respect to the domain ${\bf R}_{ +, 0}^{5} \cap \{x_{5} >0\}\; \text{.}\; $ Here the function $h_{5} = x_{4}$ is applied. Then
$ L_{f} h_{5} = p_{2} \frac{x_{3} x_{5}}{b_{3} +x_{5}} -a_{3} x_{4} -a_{4} x_{1} x_{4} +\sigma $ |
and therefore the set $S (h_{5}) \cap K_{5\max } \cap \{x_{5} >0\}$ is contained in the set defined by the inequality
$ x_{4} \leq \frac{p_{2} x_{3} x_{5}}{a_{3} (b_{3} +x_{5})} +\frac{\sigma }{a_{3}} \leq x_{4\max }: = \frac{p_{2} m x_{5\max }}{a_{3} (b_{3} +x_{5\max })} +\frac{\sigma }{a_{3}}\;\text{.}\; $ |
Therefore the localization set:
$ K_{4\max }(h_{5} , {\bf R}_{ + , 0}^{5} \cap \{x_{5} \gt 0\}): = \left \{0 \leq x_{4} \leq x_{4\max } ;x_{5} \gt 0\right\} $ |
is got.
As a result, it is established
Proposition 1. All compact invariant sets in${\bf R}_{ +, 0}^{5} \cap \{x_{5} >0\}$ are located in the domain
$ \Pi : = \cap _{i = 1 , 2 , 5}K_{i\max } \cap K_{3\max } (h_{4} , {\bf R}_{ + , 0}^{5} \cap \{x_{5} \gt 0\}) \cap K_{4\max } (h_{5} , {\bf R}_{ + , 0}^{5} \cap \{x_{5} \gt 0\})\;\text{.}\; $ |
Let $\varphi (x, t)$ be a solution of (1.1), with $\varphi (x, 0) = x$ for all $x$. Further, the following fact is used in this paper:
Lemma 3.1. Let $\sigma \geq 0$. The system (1.1) is positively Lagrange stable in ${\bf R}_{+, 0}^{5}$, i.e. for any point $x\in {\bf R}_{+, 0}^{5}$ the time functions $\varphi _{i}(x, t), t\geq 0$, $i = 1, ..., 5\text{, }\; $ are upper bounded. So there are no escaping to infinity trajectories in positive time.
Proof. 1. Taking the function $h_{1} = x_{1}$ and the domain $D_{1} = \{x_{1}\geq \kappa +\varepsilon _{1}\}\cap {\bf R}_{+, 0}^{5}$ for $\varepsilon _{1}>0$ it is easy to see that
$ L_{f}h_{1}\mid _{D_{1}}\leq -r\varepsilon _{1}-\frac{r\varepsilon _{1}^{2}}{\kappa } \lt 0\text{.}\; $ |
So applying Assertion 2 it is established that each trajectory leaves $D_{1} $ in a finite time. So $\varphi _{1}(x, t)$ is upper bounded as a time function.
2. Since $x_{3} = m$ is an invariant plane any trajectory with an initial condition satisfying the inequality $x_{3}\leq m$ has the property that $\varphi _{3}(x, t)$ is upper bounded as a time function. From the other side, if a trajectory has an initial condition satisfying $x_{3}>m$ then since $\dot{x}_{3}(t)\leq 0$ on the domain $x_{3}>m$ then $\varphi _{3}(x, t)$ is upper bounded as a time function.
3. We have in ${\bf R}_{+, 0}^{5}$ that
$ \dot{x}_{4}\leq -a_{3}x_{4}+p_{2}\varphi _{3}(x, t)+\sigma . $ |
Integrating this inequality we obtain that the function $\varphi _{4}(x, t)$ is upper bounded.
4. Taking the function $h_{2} = x_{1}+x_{2}$ and the domain
$ D_{2} = \{x_{1}+x_{2}\geq h_{2\max }+\varepsilon _{2};x_{1} \lt \kappa +\varepsilon _{1}\}\cap {\bf R}_{+, 0}^{5}, $ |
for $\varepsilon _{2}>0$ it is easy to see that
$ L_{f}h_{2}\mid _{D_{2}}\leq -r\varepsilon _{1}-\frac{r\varepsilon _{1}^{2}}{\kappa }-a_{2}(h_{2\max }+\varepsilon _{2})\inf\limits_{t \gt 0}\frac{\varphi _{4}(x, t)}{a_{2}+\varphi _{4}(x, t)} \lt 0\text{.}\; $ |
So applying Assertion 2 it is established that each trajectory leaves $D_{2}$ in a finite time. So both of $\varphi _{i}(x, t), i = 1, 2\; \text{, }\; $ are upper bounded time functions.
5. At last, it follows from integration of $\dot{x}_{5}$ with the upper bounded time function $\varphi _{2}(x, t)$ that $\varphi _{5}(x, t)$ is the upper bounded time function as well.
The main purpose of this section is to present attraction conditions to either TF- equilibrium points or TO- equilibrium points in case $\sigma = 0$. All statements are established in the form of simple inequalities imposed on the decay rate $\gamma $ of VSV particles in the blood. Here $\gamma $ can be considered as a control parameter. Firstly, it is introduced a system obtained from (1.1) by a restriction on the invariant plane $x_{3} = x_{4} = 0$:
$
⋅x1=rx1(1−x1+x2κ)−a1x1x5b1+x1;⋅x2=a1x1x5b1+x1−δx2;⋅x5=δβx2−γx5.
$
|
(4.1) |
By $g$ it is denoted its vector field. This system has the following equilibrium points
$
e1=(0,0,0)T;e2=(κ,0,0)T;e3=(e31,e32,e33)T,
$
|
with
$
e31=b1γa1β−γ;e32=rb1γ(κ−e31)(e31+γ)rb1γ(e31+γ)+e31κa1βδ;e33=δβγe32.
$
|
Before considering the system (1.1) it is formulated
Lemma 4.1. The system (4.1) has the localization domain$\Pi _{1}$ in ${\bf R}_{+, 0}^{3}$ defined by
$ x_{1}\leq x_{1\max }: = \kappa ;x_{2}\leq x_{2\max };x_{5}\leq x_{5\max }\text{, }\; $ |
where $x_{2\max }$ is defined by (3.2). This polytope contains all compact invariant sets and the attracting set of (4.1).
These bounds are derived as in case of Proposition 1.
Proposition 2. Assume that
$ \gamma \gt \gamma_{min}: = \frac{a_{1}\beta \kappa }{\kappa +b_{1}}. $ | (4.2) |
Then each trajectory of the system (4.1) in ${\bf R}_{+, 0}^{3}$ tends to ${\bf e}_{1}$ or ${\bf e}_{2}\text{.}\; $
Proof. Let us apply the function $h_{6} = x_{2}+\beta ^{-1}x_{5}$. Then
$ L_{g}h_{6} = (\frac{a_{1}x_{1}}{b_{1}+x_{1}}-\frac{\gamma }{\beta })x_{5}. $ | (4.3) |
It is noted that all $\omega $-limit sets of the system (4.1) are located in $\Pi _{1}$. So one may examine the equality (4.3) on $\Pi _{1}$ and proceed as follows
$ L_{g}h_{6}\mid _{\Pi _{1}}\leq (\frac{a_{1}\kappa }{b_{1}+\kappa }-\frac{\gamma }{\beta })x_{5}\mid _{\Pi _{1}}\leq 0. $ |
By using Lemma 4.1 and the LaSalle theorem it is established that the $\omega $-limit set of any point of the system (4.1) is not empty and belongs to the set $\{x_{5} = 0\}\text{.}\; $Therefore, in order to analyze the ultimate dynamics of (4.1) it should be examined the ultimate dynamics of the system
$
⋅x1=rx1(1−x1+x2κ)−a1x1φ5(x,t)b1+x1⋅x2=a1x1φ5(x,t)b1+x1−δx2,
$
|
(4.4) |
with the input $\varphi _{5}(x, t)$ satisfying the condition
$ \underset{t\rightarrow \infty }{\lim }\varphi _{5}(x, t) = 0. $ |
The system (4.4) is asymptotically autonomous with the limit system
$
⋅x1=rx1(1−x1+x2κ)⋅x2=−δx2,
$
|
(4.5) |
All solutions of the system (4.5) in ${\bf R}_{+, 0}^{2}$ are convergent to the equilibrium point $(\kappa, 0)^{T}$. In virtue of the Markus theorem, all solutions of the system (4.4) are convergent to the equilibrium point $(\kappa, 0)^{T}\text{.}\; $ Therefore any solution of the system (4.1) tends to ${\bf e}_{1}$ or ${\bf e}_{2}.$
Remark 1. If the inequality
$ 1 \lt \frac{a_{1}\beta }{\gamma } \lt 1+\frac{b_{1}}{\kappa } $ |
then there is the third equilibrium ${\bf e}_{3}$ which is not contained in ${\bf R}_{ +, 0}^{3}$ because $\kappa -e_{3 1} < 0$ and, therefore, $e_{3 2} < 0$. It follows from formulas for equilibrium points.
Now it is established
Theorem 4.2. Suppose that (4.2) is fulfilled. Then the$\omega $ -limit set in ${\bf R}_{+, 0}^{5}$ is either one of TF-equilibrium points or one of TO-equilibrium points.
It should be noted that the condition (4.2) is the condition of local stability of TO-equilibrium points, see [6].
Proof. Let us apply the function $h_{6} = x_{2}+\beta ^{-1}x_{5}$. In this case
$ L_{f}h_{6} = a_{1}\frac{x_{1}x_{5}}{b_{1}+x_{1}}-\frac{\gamma }{\beta }x_{5}-a_{2}\frac{x_{2}x_{4}}{b_{2}+x_{4}}. $ | (4.6) |
Since all $\omega $-limit sets of the system (1.1) in ${\bf R}_{+, 0}^{5}\cap \{x_{5}>0\}$ are nonempty and are located in $K_{1\max }$ one may restrict the equality (4.6) on $K_{1\max }$ and proceed as follows
$ L_{f}h_{6}\mid _{K_{1\max }}\leq \{(\frac{a_{1}\kappa }{b_{1}+\kappa }-\frac{\gamma }{\beta })x_{5}-a_{2}\frac{x_{2}x_{4}}{b_{2}+x_{4}}\}\mid _{K_{1\max }}\leq 0. $ |
Using the condition (4.2) and the LaSalle theorem it is possible to conclude that any $\omega $-limit set belongs to the set $\{x_{2} = x_{5} = 0\}\cup \{x_{4} = x_{5} = 0\}$. Hence, any $\omega $-limit set is one of TF-equilibrium points or TO-equilibrium points.
Below everwhere we take the case $a_{2}>r$. We shall utilize the notation $x_{1\max }^{(1)}: = x_{1\max }$.
Proposition 3. Then all compact invariant sets are located in the domain $K_{4\min }^{(1)}\cap K_{1\max }^{(2)}$, with
$
K(1)4min:={x4≥x(1)4min:=σa3+a4κ};K(2)1max:={x1≤x(2)1max:=κ−a2κx(1)4minr(b2+x(1)4min)}.{}
$
|
Here these formulas are valid under condition
$ \sigma \lt \sigma _{\min }^{(1)}: = \frac{b_{2}r(a_{3}+a_{4}\kappa )}{a_{2}-r}. $ | (5.1) |
2) Further, if
$ \sigma \geq \sigma _{\min }^{(1)} $ | (5.2) |
then all $\omega $ -limit sets are located in the plane $x_{1} = 0.$ Moreover, all $\omega $ -limit sets are located in the plane $x_{1} = 0;x_{2} = 0$ as well.
The value of $\sigma _{\min }^{(1)}$ is called the first eradication bound.
Proof. 1). Using the function $h_{5}$ one may obtain inequalities
$ \sigma \leq p_{2}\frac{x_{3}x_{5}}{b_{3}+x_{5}}+\sigma = a_{3}x_{4}+a_{4}x_{1}x_{4}\leq (a_{3}+a_{4}\kappa )x_{4} $ |
within the set $K_{1\max }\cap S(h_{5})$ and, as a result, the localization set $K_{4\min }^{(1)}$ is got. Next,
$ \dot{x}_{1}\mid _{x_{1} \gt 0}\leq (r-\frac{r}{\kappa }x_{1}-\frac{a_{2}x_{4}}{b_{2}+x_{4}})\mid _{x_{1} \gt 0} $ | (5.3) |
and the inequality (5.3) on the domain $K_{4\min }^{(1)}$ leads to the inequality
$ 0 \lt \frac{r}{\kappa }x_{1}\leq r-\frac{a_{2}x_{4\min }^{(1)}}{b_{2}+x_{4\min }^{(1)}}\;\text{.}\; $ |
This implies the desirable assertion.
2) It follows from (5.3), the LaSalle theorem and the fact that the second equation of (1.1) with $x_{1}\equiv 0$ is globally asymptotically stable.
Now we discuss how we can decrease the value for $\sigma $ for which the tumor clearance can be reached.
One may iterate this argument and obtain sequences of localization sets
$
K(1)4min⊂K(2)4min⊂…K(2)1max⊃K(3)1max⊃…,
$
|
with
$ x_{1\max }^{(k)} (\sigma ): = \kappa -\frac{a_{2} \kappa x_{4\min }^{(k -1)} (\sigma )}{r (b_{2} +x_{4\min }^{(k -1)}) (\sigma )} , k = 3 , 4 , \ldots $ |
and
$ x_{1\max }^{(k)}(\sigma ) \gt x_{1\max }^{(k+1)}(\sigma ), k = 1, 2, ... $ |
So similarly as above, we establish
Proposition 4. If for given value of the influx rate $\sigma $ there is $k$ such that
$ x_{1\max }^{(1)}(\sigma ) \gt x_{1\max }^{(2)}(\sigma ) \gt ... \gt x_{1\max }^{(k-1)}(\sigma ) \gt 0 $ |
but $x_{1\max }^{(k)}(\sigma)\leq 0$ then we have the tumor clearance for this value of $\sigma $.
The minimal positive root of the equation (5.4)
$ x_{1\max }^{(k)}(\sigma ) = 0 $ | (5.4) |
is called the $k-$th eradication bound which will be denoted by $\sigma _{\min }^{(k)}$. In what follows, we examine the issue concerning the existence of a positive root of (5.4) for any $k$. Finding bounds for the value $\sigma _{\min }^{(k)}$ encounters computational difficulties and is not the subject of this study. We notice that though the function $x_{1\max }^{(k)}(\sigma)$ may possess more than one positive root, this fact does not contradict the biological meaning of the use of the treatment $\sigma $. This is because in Propositions 3 and 4 only sufficient conditions of the tumor clearance are given.
Below we shall utilize notations
$ A_{k}(\sigma ) = \frac{P_{k}(\sigma )}{Q_{k}(\sigma )} = x_{1\max }^{(k)}(\sigma );B_{k}(\sigma ) = x_{4\min }^{(k)}(\sigma ) $ | (5.5) |
and for simplicity we shall omit $\sigma $ for functions introduced in (5.5). Further, our goal is to establish
Theorem 5.1. The equation (5.4) for $k\geq 2$ has at least one positive root.
Proof. Firstly, we prove
Lemma 5.2. We have for any $k\geq 2$ that 1) $\deg P_{k} = \deg Q_{k} = k-1;$ 2) coefficients of $P_{k}$ and $Q_{k}$ satisfy the following condition:
$
coefσ0(Pk)>0;coefσ0(Qk)>0;coefσk−1(Pk)<0;coefσk−1(Qk)>0.
$
|
Proof. We shall prove this statement by induction on $k.$ For $k = 2$ this result is valid since
$
A2=P2Q2=κ−κa2B1r(b2+B1)=rb2κ(a3+a4κ)+σκ(r−a2)r[b2(a3+a4κ)+σ].
$
|
Suppose that this assertion is proven for $k.$ Let us consider this assertion for $k+1.$ Then
$ A_{k+1} = \frac{P_{k+1}}{Q_{k+1}} = \kappa -\frac{\kappa a_{2}B_{k}}{r(b_{2}+B_{k})} = \frac{\kappa rb_{2}+B_{k}\kappa (r-a_{2})}{r(b_{2}+B_{k})}. $ | (5.6) |
Using the formula
$ B_{k} = \frac{\sigma }{a_{3}+a_{4}A_{k}} = \frac{\sigma Q_{k}}{a_{3}Q_{k}+a_{4}P_{k}} $ |
we proceed (5.6) as follows:
$ A_{k+1} = \frac{\kappa rb_{2}(a_{3}Q_{k}+a_{4}P_{k})+\sigma Q_{k}\kappa (r-a_{2})}{r[b_{2}(a_{3}Q_{k}+a_{4}P_{k})+\sigma Q_{k}]}. $ | (5.7) |
Now we get from the inductive hypothesis and the formula (5.7) that the first assertion of this lemma is valid and, besides, we notice that
$
coefσ0(Pk+1)=κrb2[a3coefσ0(Qk)+a4coefσ0(Pk)]>0;coefσ0(Qk+1)=rb2[a3coefσ0(Qk)+a4coefσ0(Pk)]>0;coefσk(Pk+1)=κ(r−a2)coefσk−1(Qk)<0;coefσk(Qk+1)=rcoefσk−1(Qk)>0.
$
|
Therefore Lemma 5.2 is proven.
As a result, in virtue of Lemma 5.2 and the formula (5.7) we obtain the desirable assertion.
Now we compute the value $\sigma _{\min }^{(2)}$. Then we have that
$ B_{2} = \frac{\sigma r(b_{2}+B_{1})}{r(a_{3}+a_{4}\kappa )(b_{2}+B_{1})-a_{2}a_{4}\kappa B_{1}}. $ |
Examining the equation $A_{3}(\sigma) = 0$ we get the equation
$ rb_{2} = B_{2}(a_{2}-r). $ | (5.8) |
For our convenience we introduce notations
$
α1=b2a2−r;α2=rb2(a3+a4κ)2;α3=r(a3+a4κ)−a2a4κ
$
|
in (5.8) and after this we come to the equation
$ \sigma ^{2}+\sigma (b_{2}(a_{3}+a_{4}\kappa )-\alpha _{1}\alpha _{3})-\alpha _{1}\alpha _{2} = 0. $ |
Its unique positive root has the form
$ \sigma _{+} = \frac{\alpha _{1}\alpha _{3}-b_{2}(a_{3}+a_{4}\kappa )}{2}+\sqrt{\frac{(\alpha _{1}\alpha _{3}-b_{2}(a_{3}+a_{4}\kappa ))^{2}}{4}+\alpha _{1}\alpha _{2}} $ |
We notice that $\sigma _{\min }^{(2)} = \sigma _{+}.$ Finally, by the routine computations omitted here it is easy to establish that
$ \sigma _{\min }^{(1)}+\frac{b_{2}(a_{3}+a_{4}\kappa )-\alpha _{1}\alpha _{3}}{2} \gt 0 $ |
and, its square exceeds
$ \frac{(\alpha _{1}\alpha _{3}-b_{2}(a_{3}+a_{4}\kappa ))^{2}}{4}+\alpha _{1}\alpha _{2}. $ |
So, we have that $\sigma _{\min }^{(1)}> \sigma _{\min }^{(2)}$ and in case of using $\sigma \geq \sigma _{\min }^{(2)}$ we still achieve the tumor eradication.
Now taking parameters from Table 2 [6], we demonstrate by calculations that the tumor eradication bounds $\sigma _{\min }^{(1)}; \sigma _{\min }^{(2)}; \sigma _{\min }^{(3)}$ are decreasing: $\sigma _{\min }^{(1)} = 1118.3;\sigma _{\min }^{(2)} = 530.3;\sigma _{\min }^{(3)} = 419.9$.
On Figures 1; 2; 3 we show the tumor eradication dynamics corresponding these values of $\sigma$ at $\gamma = 1.26; 1.46$ and $1.66:$
It should be said that one can get improved bounds for $x_{2}$ and $x_{5}$ with help of using formulas (3.2)-(3.3). In this case the function $h_{2} = x_{1}+x_{2}$ it should be employed in the same way as above in the formula (3.1). So considering $S(h_{2})$ within $K_{1\max }^{(1)}\cap K_{4\min }^{(1)}$ it is easy to see that this set is located in the set defined by the inequality
$ (\delta +\frac{a_{2} x_{4\min }}{b_{2} +x_{4\min }}) (x_{1} +x_{2}) \leq \zeta \;\text{.}\; $ |
Thus the localization set $K_{2 2}: = \left \{h_{2} \leq h_{2\max }: = \zeta \delta _{1}^{ -1}\right \}$ is obtained, with
$ \delta _{1} = \delta +\frac{a_{2} x_{4\min }}{b_{2} +x_{4\min }}\;\text{.}\; $ |
As a result, one can get upper bounds indicated in the formula
$
K(1)2max:={x2≤ζδ−11:=x(1)2max}, K(1)5max:={x5≤x(1)5max:=σβγx(1)2max}.
$
|
Hence one may conclude:
Proposition 5. All compact invariant sets are located in the domain $K_{1\max }^{(1)}\cap K_{2\max }^{(1)}\cap K_{5\max }^{(1)}\cap K_{4\min }^{(1)}$.
Convergence conditions for trajectories of the model (1.1) are presented in this paper. Such kind of ultimate dynamics may be important for predicting the patient's health scenario, its monitoring, correction and a proper application of treatments. For this purpose one may vary the parameters $\gamma $ (the decay rate for the concentration of OVPs in the blood) and $\sigma $ (the rate of influx of immune cells from the external source).
1) Ultimate upper bounds for $x_{1}$-, $x_{2}$-cells populations and for concentration of OVPs ($x_{5}$) are obtained. Upper bounds for $x_{3}$-, $x_{4}$-cells populations are given provided we examine compact invariant sets outside the plane $x_{5} = 0.$ As a result, it is shown that ultimate dynamics of the system (1.1) outside the plane $x_{5} = 0$ is contained in the domain defined by these bounds. Then it is established that in case $\sigma >0$ upper bounds for $x_{1}$-, $x_{2}$-, $x_{5}$- variables and the lower bound for $x_{4}$- variable can be improved by some number of iterations.
2) In the case $\sigma = 0$ we examine one specific type of the convergence dynamics (1.1), which excludes the existence of TV-, ITV-equilibrium points, periodic orbits and the chaotic attractor. All results are expressed in terms of $\gamma $. It is computed the lower bound $\gamma _{\min }$ such that for $\gamma >\gamma _{\min }$ any trajectory is attracted to either one of TO-equilibrium points or one of TF-equilibrium points. This problem is solved in Theorem 4.2.
3) In Proposition 3 conditions for the global tumor elimination are presented for sufficiently large values of the influx rate $\sigma $ of immune cells from the external source $\sigma, \sigma \geq \sigma _{\min }^{(1)}$. The formula (5.2) describes how the eradication bound $\sigma _{\min }^{(1)}$ depends on parameters of the model. It is shown that the bound $\sigma _{\min }^{(1)}$ may be decreased if we apply the iterational procedure. This leads to a reduced risk of over-immunization. As a result, we get the sequence of bounds $\sigma _{\min }^{(k)}$ which are characterized in Theorem 5.1. It is demonstrated that $\sigma _{\min }^{(2)} < \sigma _{\min }^{(1)}$. The decrease of the sequence of eradication bounds $\{\sigma _{\min }^{(k)}\}$ seems clear from the point of view of common sense and is confirmed by calculations. However, the additional analysis is still required and deserves a separate research. On Figures 1; 2; 3 it is shown numerically that there is the tumor eradication dynamics at different values of $\gamma $ and $\sigma = \sigma _{\min }^{(1)}; \sigma _{\min }^{(2)}; \sigma _{\min }^{(3)}\text{.}\; $Our statements about the tumor eradication are presented in a convenient form for applications and may be helpful in administering immunotherapy.
Further, let us derive an equation for finding $\sigma _{\min }^{(3)}$. With this goal we notice that $_{.}$
$ A_{2} = \frac{\theta _{1}+\sigma \theta_{2}}{\theta _{3}+\sigma \theta_{4}}, $ | (6.1) |
with
$
θ1=κrb2(a3+a4κ);abcdθ2=κ(r−a2);θ3=rb2(a3+a4κ);abcdeθ4=r.
$
|
Utilizing these formulas for $B_{2}$ we have that
$ B_{2} = \frac{\sigma }{a_{3}+a_{4}A_{2}} = \frac{\theta _{3}\sigma +\theta _{4}\sigma ^{2}}{\theta _{5}+\theta _{6}\sigma }, $ |
with
$ \theta _{5} = a_{3}\theta _{3}+a_{4}\theta _{1}; \phantom{abcd}\theta _{6} = a_{3}\theta _{4}+a_{4}\theta _{2}. $ |
Next, we compute
$ A_{3} = \kappa -\frac{\kappa a_{2}B_{2}}{rb_{2}+rB_{2}} = \kappa \frac{\theta _{7}+\sigma \theta_{8}+\sigma ^{2}\theta_{9}}{\theta _{10}+\sigma \theta _{11}+\sigma ^{2}\theta _{12}}, $ |
with
$
θ7=rb2θ5;abcdθ8=rb2θ6+(r−a2)θ3;abcdθ9=(r−a2)θ4;abcdθ10=rb2θ6+rθ3;abcdθ11=rθ4.
$
|
Further,
$ B_{3} = \frac{\sigma }{a_{3}+a_{4}A_{3}} = \frac{\sigma \theta _{7}+\sigma ^{2}\theta _{10}+\sigma ^{3}\theta _{11}}{\theta_{12}+\sigma \theta _{13}+\sigma ^{2}\theta _{14}}, $ |
with
$ \theta_{12} = (a_{3}+a_{4}\kappa )\theta _{7};\theta _{13} = a_{3}\theta _{10}+a_{4}\kappa \theta _{8};\theta _{14} = a_{3}\theta _{11}+a_{4}\theta _{9} $ |
Finally, we notice that $\sigma _{\min }^{(3)}$ is a minimal positive root of the cubic equation $rb_{2}+B_{3}(r-a_{2}) = 0$ which can be written in the form
$ \sigma ^{3}\theta _{11}(r-a_{2})+\sigma ^{2}(rb_{2}\theta _{14}+(r-a_{2})\theta _{10})+\sigma (rb_{2}\theta _{13}+(r-a_{2})\theta _{7})+rb_{2}\theta _{12} = 0. $ |
All authors declare no conflicts of interest in this paper.
[1] | Alberts B, Johnson A, Lewis J, et al. (2007) Molecular Biology of the Cell, Garland Science: New York. |
[2] | Brannigan G, Brown FLH (2005) Composition Dependence of Bilayer Elasticity. J Chem Phys 122: 07490. |
[3] | Shillcock JC, Lipowsky R (2002) Equilibrium structure and lateral stress distribution of amphiphilic bilayers from dissipative particle dynamics simulations. J Chem Phys 117: 5048-5061. |
[4] | Lipowsky R, Sackmann E (1995) Structure and dynamics of membranes, Handbook of biological physics, Elsevier, Amsterdam. |
[5] | Petelska AD, Figaszewski ZA (2002) Effect of pH on the Interfacial Tension of Lipid Bilayer Membrane. Biophys J 1561:135-146. |
[6] |
Cooke IR, Kremer K, Deserno M (2005) Tunable Generic Model for Fluid Bilayer Membranes. Phys Rev E 72: 011506. doi: 10.1103/PhysRevE.72.011506
![]() |
[7] |
Laradji M, Kumar PBS (2004) Dynamics of Domain Growth in Self-assembled Fluid Vesicles. Phys Rev Lett 93: 198105. doi: 10.1103/PhysRevLett.93.198105
![]() |
[8] |
Laradji M, Kumar PBS (2005) Domain Growth, Budding, and Fission in Phase Separating Self-assembled Fluid Bilayers. J Chem Phys 123: 224902. doi: 10.1063/1.2102894
![]() |
[9] |
Ramachandran S, Laradji M, Kumar PBS (2009) Lateral Organization of Lipids in Multi-component Liposomes. J Phys Soc Jpn 78: 041006. doi: 10.1143/JPSJ.78.041006
![]() |
[10] |
Taniguchi T (1996) Shape Deformation and Phase Separation Dynamics of Two-component Vesicles. Phys Rev Lett 76: 4444-4447. doi: 10.1103/PhysRevLett.76.4444
![]() |
[11] | Fan J, Han T, Haataja M (2010) Hydrodynamic Effects on Spinodal Decomposition Kinetics in Planar Lipid Bilayer Membranes. J Chem Phys 133: 235101. |
[12] |
Stanich CA, Honerkamp-Smith AR, Putzel GG, et al. (2013) Coarsening Dynamics of Domains in Lipid Membranes. Biophys J 105: 444-454. doi: 10.1016/j.bpj.2013.06.013
![]() |
[13] |
Veatch SL, Keller SL (2003) Separation of Liquid Phases in Giant Vesicles of Ternary Mixtures of Phospholipids and Cholesterol. Biophys J 85:3074-3083. doi: 10.1016/S0006-3495(03)74726-2
![]() |
[14] |
Esposito C, Tian A, Melamed S, et al. (2007) Flicker Spectroscopy of Thermal Lipid Bilayer Domain Boundary Fluctuations. Biophys J 93: 3169-3181. doi: 10.1529/biophysj.107.111922
![]() |
[15] | Lipowsky R (1992) Budding of Membranes Induced by Intramembrane Domains. J Phys II 2: 1825. |
[16] |
Bagatolli LA, Gratton E (2001) Direct Observation of Lipid Domains in Free Standing Bilayers Using Two-photon Excitation Fluorescence Microscopy. J Fluorescence 11: 141-160. doi: 10.1023/A:1012228631693
![]() |
[17] |
Ramachandran S, Komura S, Gommper G (2010) Effects of an Embedding Bulk Fluid on Phase Separation Dynamics in a Thin Liquid Film. EPL 89: 56001. doi: 10.1209/0295-5075/89/56001
![]() |
[18] |
Ursell TS, Klug WS, Phillips R (2009) Morphology and Interaction between Lipid Domains. Proc Natl Acad Sci U S A 106: 13301. doi: 10.1073/pnas.0903825106
![]() |
[19] |
Bagatolli L, Kumar PBS (2009) Phase Behavior of Multicomponent Membranes: Experimental and Compuatational Techniques. Soft Matter 5: 3234-3248. doi: 10.1039/b901866b
![]() |
[20] | Marrink SJ, de Vries AH, Tieleman DP (2009) Lipids on the Move: Simulations of Membrane Pores, Domains, Stalks and Curves. Biochim Biophys Acta Biomembr 1788: 149-168. |
[21] |
Lipowsky R (2002) Domains and Rafts in Membranes—Hidden Dimensions of Selforganization. J Biol Phys 28: 195-210. doi: 10.1023/A:1019994628793
![]() |
[22] | Simons K, Vaz WLC (2004) Model Systems, Lipid Rafts, and Cell Membranes. Annu Rev Biophys Biomol Struct 3: 269. |
[23] |
Barberousse A, Franceschelli S, Imbert C (2009) Computer Simulations as Experiments. Synthese 169: 557-574. doi: 10.1007/s11229-008-9430-7
![]() |
[24] | Farago O (2003) “Water-free” Computer Model for Fluid Bilayer Membranes. O J Chem Phys 119: 596-605. |
[25] |
Brannigan G, Brown FLH (2004) Solvent-free Simulations of Fluid Membrane Bilayers. J Chem Phys 120: 1059. doi: 10.1063/1.1625913
![]() |
[26] | Shillcock JC (2012) Spontaneous Vesicle Self-Assembly: A Mesoscopic View of Membrane Dynamics. Langmuir 28: 541-547. |
[27] | Tieleman DP, Leontiadau H, Mark AE, et al. (2003) Simulation of Pore Formation in Lipid Bilayers by Mechanical Stress and Electric Fields. J Am Chem Soc125: 6382-6383. |
[28] |
Damodaran KV, Merz KM (1994) A Comparison of DMPC- and DLPE-based Lipid Bilayers. Biophys J 66: 1076-1087. doi: 10.1016/S0006-3495(94)80889-6
![]() |
[29] |
Moore PB, Lopez CF, Klein ML (2001) Dynamical Properties of a Hydrated Lipid Bilayer from a Multinanosecond Molecular Dynamics Simulation. Biophys J 81: 2484-2494. doi: 10.1016/S0006-3495(01)75894-8
![]() |
[30] |
Essmann U, Perera L, Berkowitz ML (1995) The Origin of the Hydration Interaction of Lipid Bilayers from MD Simulation of Dipalmitoylphosphatidylcholine Membranes in Gel and Liquid Crystalline Phases. Langmuir 11: 4519-4531. doi: 10.1021/la00011a056
![]() |
[31] |
Cooke IR, Deserno M (2005) Solvent-free Model for Self-assembling Fluid Bilayer Membranes: Stabilization of the Fluid Phase based on Broad Attractive Tail Potentials. J Chem Phys 123: 224710. doi: 10.1063/1.2135785
![]() |
[32] |
West B, Schmid F (2010) Fluctuations and Elastic Properties of Lipid Membranes in the Gel L-beta State: A Coarse-grained Monte Carlo Study. Soft Matter 6: 1275. doi: 10.1039/b920978f
![]() |
[33] |
Farago O (2008) Mode Excitation Monte Carlo Simulations of Mesoscopically Large Membranes. J Chem Phys 128: 184105. doi: 10.1063/1.2918736
![]() |
[34] |
Farago O (2010) Fluctuation-induced Attraction between Adhesion Sites of Supported Membranes. Phys Rev E 81: 050902. doi: 10.1103/PhysRevE.81.050902
![]() |
[35] |
Farago O (2008) Membrane Fluctuations near a Plane Rigid Surface. Phys Rev E 78:051919. doi: 10.1103/PhysRevE.78.051919
![]() |
[36] | Dutt M, Nayhouse MJ, Kuksenok O, et al. (2011) Interactions of End-Functionalized Nanotubes with Lipid Vesicles: Spontaneous Insertion and Nanotube Self-organization. Current Nanoscience 7: 699-715. |
[37] |
Dutt M, Kuksenok O, Nayhouse MJ, et al. (2011) Modeling the Self-Assembly of Lipids and Nanotubes in Solution: Forming Vesicles and Bicelles with Transmembrane Nanotube Channels. ACS Nano 5: 4769-4782. doi: 10.1021/nn201260r
![]() |
[38] |
Dutt M, Kuksenok O, Little SR, et al. (2011) Forming Transmembrane Channels Using End-Functionalized Nanotubes. Nanoscale. 3: 240-250. doi: 10.1039/C0NR00578A
![]() |
[39] | Dutt M, Kuksenok O, Little SR, et al. (2012) Designing Tunable Bio-nanostructured Materials via Self-assembly of Amphiphilic Lipids and Functionalized Nanotubes. MRS Spring 2012 Conference Proceedings; 1464. |
[40] | Ludford P, Aydin F, Dutt M (2013) Design and Characterization of Nanostructured Biomaterials via the Self-assembly of Lipids. MRS Fall 2013 Conference Proceedings; 1498. |
[41] | Koufos E, Dutt M (2013) Design of Nanostructured Hybrid Inorganic-biological Materials via Self-assembly. MRS Spring 2013 Conference Proceedings; 1569. |
[42] |
Smith KA, Jasnow D, Balazs AC (2007) Designing Synthetic Vesicles that Engulf Nanoscopic Particles. J Chem Phys 127: 084703. doi: 10.1063/1.2766953
![]() |
[43] | Goetz R, Lipowsky R (1998) Computer Simulations of Bilayer Membranes: Self-assembly and Interfacial Tension. J Chem Phys 108: 7397-7409. |
[44] | Kranenburg M, Venturoli M, Smit B. (2003) Phase Behavior and Induced Interdigitation in Bilayers Studied with Dissipative Particle Dynamics. J Phys Chem 41: 11491. |
[45] |
Kranenburg M, Laforge C, Smit B (2004) Mesoscopic Simulations of Phase Transitions in Lipid Bilayers. Phys Chem Chem Phys 6: 4531-4534. doi: 10.1039/b410914g
![]() |
[46] | Yamamoto S, Maruyama Y, Hyodo S (2002) Dissipative Particle Dynamics Study of Spontaneous Vesicle Formation of Amphiphilic Molecules. J Chem Phys 116: 5842. |
[47] |
Yamamoto S, Hyodo S (2003) Budding and Fission Dynamics of Two-Component Vesicles. J Chem Phys 118: 7937-7943. doi: 10.1063/1.1563613
![]() |
[48] |
Stevens MJ, Hoh JH, Woolf TB (2003) Insights into the Molecular Mechanism of Membrane Fusion from Simulation: Evidence for the Association of Splayed Tails. Phys Rev Lett 91: 188102. doi: 10.1103/PhysRevLett.91.188102
![]() |
[49] | Stevens MJ (2004) Coarse-grained Simulations of Lipid Bilayers. Chem Phys 121: 11942-11948. |
[50] | Arkhipov A, Yin Y, Schulten K (2009) Membrane-bending Mechanism of Amphiphysin N-BAR Domains. Biophys J 97: 2727-2735. |
[51] |
Shih AY, Arkhipov A, Freddolino PL, et al. (2006) A Coarse-grained Protein-lipid Model with Application to Lipoprotein Particles. J Phys Chem 110: 3674-3684. doi: 10.1021/jp0550816
![]() |
[52] |
Marrink SJ, Risselada HJ, Yefimov S, et al. (2007) The MARTINI Forcefield: Coarse-grained Model for Biomolecular Simulations. J Phys Chem B 111: 7812-7824. doi: 10.1021/jp071097f
![]() |
[53] | Wang Z, Frenkel DJ (2005) Modeling Flexible Amphiphilic Bilayers: A Solvent-free Off-lattice Monte Carlo Study. Chem Phys 122: 234711. |
[54] |
Brannigan G, Philips PF, Brown FLH (2005) Flexible Lipid Bilayers in Implicit Solvent. Phys Rev E 72: 011915. doi: 10.1103/PhysRevE.72.011915
![]() |
[55] |
Noguchi H, Takasu M (2001) Self-assembly of Amphiphiles into Vesicles: A Brownian Dynamics Simulation. Phys Rev E 64: 041913. doi: 10.1103/PhysRevE.64.041913
![]() |
[56] |
Noguchi H (2002) Fusion and Toroidal Formation of Vesicles by Mechanical Forces: A Brownian Dynamics Simulation. J Chem Phys 117: 8130-8137. doi: 10.1063/1.1510114
![]() |
[57] |
Katsov K, Mueller M, Schick M (2004) Field Theoretic Study of Bilayer Membrane Fusion I Hemifusion Mechanism. Biophys J 87: 3277. doi: 10.1529/biophysj.103.038943
![]() |
[58] | Schick M (2012) Membranes: A Field-theoretic Description. Encyclopedia of Biophysics. Roberts, G.C.K., Ed., Springer-Verlag: Berlin Heidelberg. |
[59] |
May S, Kozlovsky Y, Ben-Shaul A, et al. (2004) Tilt Modulus of a Lipid Monolayer. Eur Phys J E 14: 299-308. doi: 10.1140/epje/i2004-10019-y
![]() |
[60] |
May S (2000) A Molecular Model for the Line Tension of Lipid Membranes. Eur Phys J E 3: 37-44. doi: 10.1007/s101890070039
![]() |
[61] |
Lee WB, Mezzenga R, Fredrickson GH (2008) Self-consistent Field Theory for Lipid-based Liquid Crystals: Hydrogen Bonding Effect. J Chem Phys 128: 074504-074510. doi: 10.1063/1.2838624
![]() |
[62] | Ginzburg VV, Balijepalli S (2007) Modelling the Thermodynamics of the Interaction of Nanoparticles with Cell Membranes. Nano Lett 7: 3716-3722. |
[63] |
Ayton G, Voth GA (2002) Bridging Microscopic and Mesoscopic Simulations of Lipid Bilayers. Biophys J 83: 3357-3370. doi: 10.1016/S0006-3495(02)75336-8
![]() |
[64] |
Wang ZJ, Deserno MA (2010) Systematically Coarse-grained Solvent-free Model for Quantitative Phospholipid Bilayer Simulations. J Phys Chem B 114: 11207. doi: 10.1021/jp102543j
![]() |
[65] |
Wang ZJ, Deserno M (2010) Systematic Implicit Solvent Coarse-graining of Bilayer Membranes: Lipid and Phase Transferability of the Force Field. New J Phys 12: 095004. doi: 10.1088/1367-2630/12/9/095004
![]() |
[66] |
Ge Z, Li Q, Wang Y (2014) Free energy Calculation of Nanodiamond-Membrane Association—The Effect of Shape and Surface Functionalization. J Chem Theory Comput 10: 2751-2758. doi: 10.1021/ct500194s
![]() |
[67] | Reid CVL, Ricci M, Silva PHJ, et al. (2014) Lipid tail protrusions mediate the insertion of nanoparticles into model cell membranes. Nat Commun 5: 4482. |
[68] |
Wong-Ekkabut J, Baoukina S, Triampo W, et al. (2008) Computer simulation study of fullerene translocation through lipid membranes. Nat Nanotechnol 3: 363-368. doi: 10.1038/nnano.2008.130
![]() |
[69] |
Li Y, Chen X, Gu N (2008) Computational Investigation of Interaction between Nanoparticles and Membranes: Hydrophobic/Hydrophilic Effect. J Phys Chem B 112: 16647-16653. doi: 10.1021/jp8051906
![]() |
[70] |
Huang C, Zhang Y, Yuan H, et al. (2013) Role of Nanoparticle Geometry in Endocytosis: Laying Down to Stand Up. Nano Lett 13: 4546-4550. doi: 10.1021/nl402628n
![]() |
[71] |
Shi X, Bussche AVD, Hurt RH, et al. (2011) Cell entry of one-dimensional nanomaterials occurs by tip recognition and rotation. Nat Nanotechnol 6: 714-719. doi: 10.1038/nnano.2011.151
![]() |
[72] |
Illya G, Lipowsky R, Shillcock JC (2006) Two-component membrane material properties and domain formation from dissipative particle dynamics. J Chem Phys 125: 114710. doi: 10.1063/1.2353114
![]() |
[73] |
Groot RD, Warren PB (1997) Dissipative Particle Dynamics: Bridging the gap between atomistic and mesoscopic simulation. J Chem Phys 107: 4423-4435. doi: 10.1063/1.474784
![]() |
[74] |
Chou H, Tsao HK, Sheng YJ (2006) Morphologies of multicompartment micelles formed by triblock copolymers. J Chem Phys 125: 194903. doi: 10.1063/1.2390716
![]() |
[75] |
Ortiz V, Nielsen SO, Discher DE, et al. (2005) Disipative Particle Dyanmics simulations of polymerosome. J Phys Chem B 109: 17708-17714. doi: 10.1021/jp0512762
![]() |
[76] | Boek ES, Coveney PV, Lekkerkerker HNW, et al. (1997) Simulating rheology of dense colloidal suspensions using dissipative particle dynamics. Phys Rev E 55: 3124-3131. |
[77] | Spenley NA (200) Scaling laws for polymers in dissipative particle dynamics, Europhys Lett 49: 534-540. |
[78] |
Fan XJ, Phan-Thien N, Chen S, et al. (2006) Simulating flow of DNA suspension using dissipative particle dynamics. Phys Fluids 18: 063102. doi: 10.1063/1.2206595
![]() |
[79] |
Chem S, Phan-Thien N, Fan XJ, et al. (2004) Dissipative particle dynamics of polymer drops in periodic shear flow. J Non-Newtonian Fluid Mech 118: 65-81. doi: 10.1016/j.jnnfm.2004.02.005
![]() |
[80] |
Arai N, Yasuoka K, Zeng XC (2013) A vesicle cell under collision with a Janus or homogeneous nanoparticle: translocation dynamics and late-stage morphology. Nanoscale 5: 9089-9100. doi: 10.1039/c3nr02024j
![]() |
[81] |
Yang K, Ma Y (2010) Computer simulation of the translocation of nanoparticles with different shapes across a lipid bilayer. Nat Nanotechnol 5: 579-583. doi: 10.1038/nnano.2010.141
![]() |
[82] |
Ding H, Tian W, Ma Y (2012) Designing Nanoparticle Translocation through Membranes by Computer Simulations. ACS Nano 6: 1230-1238. doi: 10.1021/nn2038862
![]() |
[83] |
Chen X, Tian F, Zhang X, et al. (2013) Internalization pathways of nanoparticles and their interaction with a vesicle. Soft Matter 9: 7592-7600. doi: 10.1039/c3sm50931a
![]() |
[84] |
Arnarez C, Uusitalo JJ, Masman MF, et al. (2015) Dry Martini, a Coarse-Grained Force Field for Lipid Membrane Simulations with Implicit Solvent. J Chem Theory Comput 11: 260-275. doi: 10.1021/ct500477k
![]() |
[85] | Hall BA, Chetwynd AP, Sansom MSP (2011) Exploring Peptide-Membrane Interactions with Coarse-Grained MD Simulations. Biophys J 100: 1940-1948. |
[86] |
Gkeka P, Sarkisov L (2010). Interactions of Phospholipid Bilayers with Several Classes of Amphiphilic α-Helical Peptides: Insights from Coarse-Grained Molecular Dynamics Simulations. J Phys Chem B 114: 826-839. doi: 10.1021/jp908320b
![]() |
[87] | Allen MP, Tildesley DJ (2001) Computer simulations of liquids, Clarendon Press, Oxford. |
[88] |
Koufos E, Muralidharan B, Dutt M (2014) Computational Design of Multi-component Bio-Inspired Bilayer. AIMS Materials Science 1: 103-120. doi: 10.3934/matersci.2014.2.103
![]() |
[89] |
Milletti F (2012) Cell-Penetrating Peptides: Classes, Origin, and Current Landscape. Drug Discov Today 17: 850-860. doi: 10.1016/j.drudis.2012.03.002
![]() |
[90] |
Aydin F, Ludford P, Dutt M (2014) Phase segregation in bio-inspired multi-component vesicles encompassing double tail phospholipid species. Soft Matter 10: 6096-6108. doi: 10.1039/C4SM00998C
![]() |
[91] |
Aydin F, Uppaladadium G, Dutt M (2015) The Design of Shape-Tunable Hairy Vesicles. Colloids Surf B Biointerfaces 128: 268-275. doi: 10.1016/j.colsurfb.2015.01.049
![]() |
[92] |
Monticelli L, Kandasamy SK, Periole X, et al. (2008) The MARTINI Coarse-Grained Force Field: Extension to Proteins. J Chem Theory Comput 4: 819-834. doi: 10.1021/ct700324x
![]() |
[93] |
Feller SE, Pastor RW (1999) Constant surface tension simulations of lipid bilayers: The sensitivity of surface areas and compressibilities. J Chem Phys 111: 1281-1287. doi: 10.1063/1.479313
![]() |