IL⁃2驱动的肿瘤⁃免疫动力学模型及数值模拟

IL-2 Driven Tumor⁃Immune Dynamic Model and Numerical Simulation

摘要

患者在接受白细胞介素-2（IL-2）免疫疗法时会经常出现一定的副作用， 这些副作用不仅影响了治疗效果， 而且对患者生命产生了严重的威胁。为了探究IL-2对免疫系统杀伤肿瘤细胞的影响， 建立了一类包含细胞因子IL-2的肿瘤-免疫动力学模型。计算了模型的无肿瘤平衡点和基本再生数， 并从数值上刻画了模型的动力学性态。数值模拟结果表明： IL-2诱导的细胞毒性T淋巴细胞（CD8+T细胞）的增殖分化是免疫系统更有效地杀伤肿瘤细胞的主要原因； IL-2诱导的CD8+T细胞增殖率的增加虽然会促进免疫细胞更快地发挥作用但同时也会增加肿瘤复发的几率； 过多IL-2的存在反而会促进肿瘤的生长。本文通过大量的数值模拟表明， 对癌症患者进行IL-2免疫疗法时， 存在IL-2注射的最佳剂量。然而， 肿瘤细胞的种类纷繁复杂， 针对不同的肿瘤细胞， 选取不同的参数进行模拟， 才会得出更具有现实意义的结论。

Abstract

Patients often experience certain side effects when receiving interleukin-2 （IL-2） immunotherapy， which not only affect the treatment effect， but also pose a serious threat to the patient’s life. In order to investigate the effect of IL-2 on the immune system’s killing of tumor cells， a tumor-immune dynamic model containing the cytokine IL-2 was established. We calculated the tumor-free equilibrium point and basic reproduction number of the model， and numerically characterized the dynamic behavior of the model. The results of numerical simulation show that IL-2-induced proliferation and differentiation of CD8+T cells is the main reason for the immune system to kill tumor cells more effectively. The increase in the proliferation rate of CD8+T cells induced by IL-2 can promote the faster function of immune cells but also increase the probability of tumor recurrence. The presence of too much IL-2 can promote tumor growth. The extensive numerical simulation also demonstrates that there is an optimal dose of IL-2 injection in administering IL-2 immunotherapy to cancer patients. However， the types of tumor cells are diverse and complex， and selecting different parameters for simulation based on different tumor cells will lead to more practical conclusions.

Graphical abstract

关键词

细胞因子IL-2 / 肿瘤 / 免疫 / 动力学模型 / 基本再生数

0 引　言

癌症是导致人类死亡的主要原因之一， 每年有数百万新的癌症病例被确诊。根据2015年世界卫生组织发布的数据显示， 新癌症病例将会在2035年增加70%。根据世界卫生组织国际癌症研究机构发布的全球最新癌症数据显示， 2020年全球癌症新发病例约为1 929万例， 癌症及癌症相关的死亡病例约为1 000万例^［1］。2022年9月， 美国科学家进行的一项研究显示， 早发型癌症的发病率， 自1990年开始在世界范围内急剧上升。随着癌症发病群体的年轻化， 以及发病率和死亡率的不断升高， 人类生活和社会发展受到了严重影响^［2］。

目前， 治疗癌症的手段有手术疗法^［3］、 放射疗法^［4］、 化学疗法^［5］和免疫疗法^［6］等。其中， IL-2免疫疗法作为新兴的治疗方法， 主要通过增强免疫反应， 刺激免疫系统对抗肿瘤， 受到了研究者的广泛关注。Kirschner等^［7］建立了包含肿瘤细胞、 效应细胞和细胞因子IL-2的数学模型， 探究了过继细胞免疫治疗的效果， 结果表明， 在一定情况下， 肿瘤会消失。Smith等^［8］和Yang等^［9］研究了IL-2免疫疗法在转移性黑色素瘤和转移性肾癌治疗中的应用， 研究发现， 持续注射IL-2可以产生有益的效果。Cho等^［10］研究了IL-2免疫疗法在宫颈癌治疗中的应用， 结果表明， 宫颈癌的治疗效果依赖于IL-2的注射剂量。此外， 临床实验发现， IL-2免疫疗法会产生一定的副作用， 这常常与人体内IL-2与免疫细胞之间十分复杂的相互作用关系有关， 而有关此类作用的研究却很少。本文扩展了Heyrim^［10］的模型， 考虑了不同的效应细胞类型， 建立了一个包含细胞因子IL-2的肿瘤-免疫系统相互作用的数学模型， 探究了IL-2对免疫系统杀伤肿瘤细胞的影响。

1 肿瘤-免疫动力学模型的建立

为了更好地建立数学模型， 首先介绍机体杀伤肿瘤的过程。肿瘤细胞是一种功能异常的畸形细胞， 具有强大的增殖和分化能力。当机体发现肿瘤细胞时， 首先发挥作用的是自然杀伤细胞（NK细胞）， NK细胞作为机体抗肿瘤的第一道防线， 能够迅速识别、 进入并破坏实体组织中的肿瘤细胞^［11］。其次发挥作用的是细胞毒性T淋巴细胞（CD8+T细胞）。CD8+T细胞识别肿瘤抗原后， 会迅速释放大量包含穿孔素的颗粒， 这些颗粒能够破坏肿瘤细胞的细胞膜结构， 使颗粒酶能够进入肿瘤细胞中， 从而诱导肿瘤细胞发生凋亡^［12］。与此同时， CD8+T细胞与肿瘤细胞的相互作用会刺激IL-2的增殖， IL-2作为人体内的一种细胞因子， 它具有激活并促进NK细胞和CD8+T细胞增殖分化的能力。值得注意的是， 当CD8+T细胞和NK细胞与肿瘤细胞之间相互作用时会发生失活， 且被NK细胞裂解的肿瘤细胞碎片会对CD8+T细胞起到招募的作用。人体内肿瘤和免疫系统之间的相互作用关系如图 1 所示， 实线箭头表示促进作用， 虚线箭头表示抑制作用。

本文在文献［10］的基础上， 考虑不同的效应细胞类型， 建立了包含细胞因子IL-2的肿瘤-免疫系统相互作用的微分方程模型

\left\{\begin{array}{l}{M}^{\text{'}}\left(t\right)=a_{\mathrm{1}}M\left(\mathrm{1}-b_{\mathrm{1}}M\right)-c_{\mathrm{1}}MC-d_{\mathrm{1}}MN,\\ C^{\text{'}}\left(t\right)=-b_{\mathrm{2}}C-c_{\mathrm{2}}CM+d_{\mathrm{2}}NM+\frac{e_{\mathrm{2}}IC}{g_{\mathrm{2}}+I},\\ N^{\text{'}}\left(t\right)=a_{\mathrm{3}}-b_{\mathrm{3}}N-c_{\mathrm{3}}NM+\frac{e_{\mathrm{3}}IN}{g_{\mathrm{3}}+I},\\ I^{\text{'}}\left(t\right)=\frac{e_{\mathrm{4}}CM}{g_{\mathrm{4}}+M}-\mu I\mathrm{。}\end{array}\right.（1）

模型（1）中： M表示肿瘤细胞群， C表示CD8+T细胞群， N表示NK细胞群， I表示细胞因子IL-2。假设肿瘤细胞的自然增殖遵循logistic增长， 则a_{\mathrm{1}}为肿瘤生长速率， \mathrm{1}/b_{\mathrm{1}}为肿瘤承载能力。肿瘤细胞与CD8+T细胞之间的相互作用用-c_{\mathrm{1}}MC表示。同样地， NK细胞也能诱导肿瘤细胞死亡， 其速率为d_{\mathrm{1}}。假设CD8+T细胞只有当肿瘤细胞存在时才存在， 则CD8+T细胞的自然增殖只包含自然死亡率。与CD8+T细胞不同的是， NK细胞在人体内有一个恒定的产生速率a_{\mathrm{3}}。假设CD8+T细胞和NK细胞的死亡率与其在人体内的数量成正比， 分别用-b_{\mathrm{2}}C和-b_{\mathrm{3}}N表示。此外， 参数c_{\mathrm{2}}和c_{\mathrm{3}}分别表示CD8+T细胞和NK细胞与肿瘤细胞之间相互作用时的平均失活率。NK细胞裂解的肿瘤细胞碎片对CD8+T细胞起到的招募作用， 用d_{\mathrm{2}}NM表示。\frac{e_{\mathrm{2}}IC}{g_{\mathrm{2}}+I}和\frac{e_{\mathrm{3}}IN}{g_{\mathrm{3}}+I}分别表示由于IL-2的刺激， CD8+T细胞和NK细胞的克隆扩增。e_{\mathrm{2}}和e_{\mathrm{3}}分别表示由于IL-2刺激， CD8+T细胞和NK细胞的增殖率。g_{\mathrm{2}}和g_{\mathrm{3}}分别表示由于IL-2刺激， CD8+T细胞和NK细胞的生长动力学达到饱和。关于细胞因子IL-2， \frac{e_{\mathrm{4}}CM}{g_{\mathrm{4}}+M}表示CD8+T细胞与肿瘤细胞的相互作用刺激了IL-2的产生。-\mu I表示人体血液中IL-2浓度呈指数衰减。模型中的\frac{e_{\mathrm{2}}IC}{g_{\mathrm{2}}+I},\frac{e_{\mathrm{3}}IN}{g_{\mathrm{3}}+I}和\frac{e_{\mathrm{4}}CM}{g_{\mathrm{4}}+M}均为功能反应项， 表示有限的免疫应答。表 1 给出了参数详细的生物学描述及其取值， 所有参数值的选取均来源于Kirschner^［7］、 Heyrim^［10］和de Pillis^［13-14］的研究中肿瘤-免疫动力学模型中的参数取值。

2 动力学性态分析

通过计算可知， 系统（1）始终存在唯一的无肿瘤平衡点

E\left(\mathrm{0,0},\frac{a_{\mathrm{3}}}{b_{\mathrm{3}}},\mathrm{0}\right)\mathrm{。}

为了研究无肿瘤平衡点E的稳定性， 首先给出系统（1）的基本再生数R_{\mathrm{0}}。根据van den Driessche等^［15］给出的下一代矩阵的方法计算基本再生数， 系统（1）中保留的仓室有M,C和I,构造矩阵\mathit{F}和\mathit{V}

\mathit{F}=\left[\begin{array}{c}{a}_{\mathrm{1}}M\left(\mathrm{1}-b_{\mathrm{1}}M\right)\\ d_{\mathrm{2}}NM\\ \mathrm{0}\end{array}\right],

\mathit{V}==\left[\begin{array}{c}{c}_{\mathrm{1}}MC+d_{\mathrm{1}}MN\\ b_{\mathrm{2}}C+c_{\mathrm{2}}CM-e_{\mathrm{2}}IC/(g_{\mathrm{2}}+I)\\ -e_{\mathrm{4}}CM/(g_{\mathrm{4}}+M)+\mu I\end{array}\right]\mathrm{。}==

矩阵\mathit{F}和\mathit{V}在E处的雅可比矩阵分别表示为\mathit{f}和\mathit{v}

\mathit{f}=\left[\begin{array}{ccc}{a}_{\mathrm{1}}& \mathrm{0}& \mathrm{0}\\ d_{\mathrm{2}}{a}_{\mathrm{3}}/b_{\mathrm{3}}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mathrm{0}\end{array}\right],\mathit{v}=\left[\begin{array}{ccc}{d}_{\mathrm{1}}{a}_{\mathrm{3}}/b_{\mathrm{3}}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& b_{\mathrm{2}}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mu \end{array}\right]\mathrm{。}

则\mathit{f}是一个非负矩阵， \mathit{v}是非奇异M-矩阵， 且

\mathit{f}{\mathit{v}}^{-\mathrm{1}}=\left(\begin{array}{ccc}{a}_{\mathrm{1}}{b}_{\mathrm{3}}/d_{\mathrm{1}}{a}_{\mathrm{3}}& \mathrm{0}& \mathrm{0}\\ d_{\mathrm{2}}/d_{\mathrm{1}}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mathrm{0}\end{array}\right)\mathrm{。}

从而

R_{\mathrm{0}}=\rho \left(\mathit{f}{\mathit{v}}^{-\mathrm{1}}\right)=\frac{a_{\mathrm{1}}{b}_{\mathrm{3}}}{d_{\mathrm{1}}{a}_{\mathrm{3}}},

式中： \rho \left(\mathit{f}{\mathit{v}}^{-\mathrm{1}}\right)表示非负矩阵\mathit{f}{\mathit{v}}^{-\mathrm{1}}的谱半径。

定理 1　当时， 系统（1）的无肿瘤平衡点E\left(\mathrm{0,0},\frac{a_{\mathrm{3}}}{b_{\mathrm{3}}},\mathrm{0}\right)局部渐近稳定， 当R_{\mathrm{0}}>\mathrm{1}时， 无肿瘤平衡点E不稳定。

根据Wang等^［16］研究中的定理4.1可知，

定理 2　当R_{\mathrm{0}}>\mathrm{1}时， 肿瘤一致持续生存。

由于系统（1）的理论分析比较困难， 本文选取适当的参数对系统（1）进行数值模拟， 观察系统可能出现的复杂性态。

1） 在系统（1）中取a_{\mathrm{1}}=0.18， b_{\mathrm{1}}=1×10^-9， c_{\mathrm{1}}=3×10^-7， d_{\mathrm{1}}=3.23×10^-7， b_{\mathrm{2}}=0.02， c_{\mathrm{2}}=3.24×10^-10， d_{\mathrm{2}}=1.1×10^-7，e_{\mathrm{2}}=0.124 5， g_{\mathrm{2}}=2×10⁷，a_{\mathrm{3}}=1.3×10⁴， b_{\mathrm{3}}=4.12×10^-2， c_{\mathrm{3}}=1×10^-7， e_{\mathrm{3}}=0.124 5， g_{\mathrm{3}}=2×10⁷， e_{\mathrm{4}}=5， g_{\mathrm{4}}=1×10³， \mu=10， 此时系统（1）的正平衡点是渐近稳定的， 仿真结果如图 2 所示。

2） 在系统（1）中取a_{\mathrm{1}}=0.18， b_{\mathrm{1}}=1×10^-9， c_{\mathrm{1}}=3×10^-7， d_{\mathrm{1}}=3.23×10^-7， b_{\mathrm{2}}=0.02， c_{\mathrm{2}}=3.42×10^-10， d_{\mathrm{2}}=1.1×10^-7， e_{\mathrm{2}}=0.124 5， g_{\mathrm{2}}=2×10⁷， a_{\mathrm{3}}=1.3×10⁴， b_{\mathrm{3}}=4.12×10^-2， c_{\mathrm{3}}=1×10^-7， e_{\mathrm{3}}=0.124 5， g_{\mathrm{3}}=2×10⁷， e_{\mathrm{4}}=30， g_{\mathrm{4}}=1×10³， \mu=10， 此时系统（1）的正平衡点不稳定， 仿真结果如图 3 所示。由图 3 可以看出， 系统（1）的解都趋向于周期解。

3 IL-2对免疫系统杀伤肿瘤细胞的影响

相关参数值或范围见表 1。首先对肿瘤数量进行敏感性分析， 结果如图 4 所示。由图 4 可以看出， 肿瘤生长速率、 CD8+T细胞死亡率和IL-2的死亡率是影响肿瘤数量的主要因素。然后通过数值模拟的方式来探讨IL-2对免疫系统杀伤肿瘤细胞的影响。

1） 由于肿瘤细胞多种多样， 不同肿瘤细胞的生长速率也不同。肿瘤细胞的生长速率越大， 恶性程度往往越高， 免疫细胞杀伤肿瘤细胞的效果也就越不明显。对于不同种类的肿瘤细胞， IL-2的存在是否对免疫细胞杀伤肿瘤细胞的效果都有影响？针对这一问题， 分别选取两组具有不同生长速率的肿瘤细胞进行数值模拟， 结果如图 5 所示。

图 5 中， a_{\mathrm{1}}=0.12为小肿瘤状态， a_{\mathrm{1}}=0.18为大肿瘤状态， 实线表示免疫系统中不存在IL-2， 虚线表示免疫系统中存在IL-2。由图 5 可以看出， 对于具有不同生长速率的肿瘤细胞， IL-2的存在都能够促进免疫系统更快地杀伤肿瘤细胞且使肿瘤细胞最终稳定在一个更低的水平。

2） CD8+T细胞和NK细胞是杀死肿瘤细胞的免疫细胞， 而IL-2存在的主要作用是诱导CD8+T细胞和NK细胞的增殖分化。当IL-2对CD8+T细胞不存在诱导作用时或者IL-2对NK细胞不存在诱导作用时， 肿瘤数量有何变化？针对这一问题进行数值模拟， 结果如图 6 所示。

图 6（a） 中， 实线表示IL-2对CD8+T细胞或NK细胞不存在诱导作用， 虚线表示IL-2对CD8+T细胞存在诱导作用。图 6（b） 中， 实线表示IL-2对NK细胞不存在诱导作用， 虚线表示IL-2对NK细胞存在诱导作用。由图 6 可以看出， 当IL-2对CD8+T细胞不存在诱导作用时， 肿瘤数量变化有明显的不同， 最终稳定在一个更高的水平。当IL-2对NK细胞不存在诱导作用时， 肿瘤数量变化没有明显的不同。这表明IL-2诱导的CD8+T细胞的增殖分化是IL-2能够促进免疫细胞杀伤肿瘤细胞的主要原因， 这是因为CD8+T细胞在对抗肿瘤细胞过程中发挥着主要作用。

3） IL-2诱导CD8+T细胞的增殖分化是免疫细胞更快更有效地杀伤肿瘤细胞的主要原因， 而IL-2诱导的CD8+T细胞增殖率大小是影响CD8+T细胞增殖分化能力强弱的主要因素。CD8+T细胞的增殖分化能力越强， 免疫系统杀伤肿瘤细胞的效果是否越明显？针对这一问题， 分别选取3组不同的IL-2诱导的CD8+T细胞增殖率进行数值模拟， 结果如图 7（a） 所示。

同样地， 模拟IL-2诱导的NK细胞增殖率对免疫系统杀伤肿瘤细胞效果的影响， 结果如图 7（b） 所示。可以看出， IL-2诱导的CD8+T细胞增殖率的增加比IL-2诱导的NK细胞增殖率的增加更能够促进免疫系统杀伤肿瘤细胞， 这合理地解释了之前的模拟结果。此外， 还可以看出IL-2诱导的CD8+T细胞增殖率的增加也会引起肿瘤细胞的复发， 这是因为当肿瘤细胞变小时， CD8+T细胞的数量也会随之变小， 从而给了肿瘤细胞生存的空间。

4） IL-2的主要来源是肿瘤与免疫细胞之间的相互作用， 当人体内IL-2的数量越多， IL-2诱导产生的CD8+T细胞和NK细胞的数量也会越多。人体内IL-2的数量越多， 肿瘤细胞的数量是否就会越少？针对这一问题， 选取3组不同的肿瘤与免疫细胞之间的相互作用对IL-2的招募率进行数值模拟， 结果如图 8 所示。

由图 8 可以看出， 当人体内的IL-2达到一定数量时， 肿瘤甚至会接近消失， 但与此同时也会出现肿瘤复发的现象， 且肿瘤数量会越来越多， 最终不受控制。这也意味着， 当对患者进行IL-2免疫疗法时， 存在IL-2注射的最佳剂量。

4 结论和讨论

在本文中， 根据IL-2在免疫系统杀伤肿瘤细胞过程中的作用， 建立了一类包含细胞因子IL-2的肿瘤-免疫系统相互作用的动力学模型。对于建立的动力学模型， 给出了无肿瘤平衡点E和基本再生数R_{\mathrm{0}}， 并从数值上刻画了模型的动力学性态。通过数值模拟的方式研究了细胞因子IL-2对免疫系统杀伤肿瘤细胞效果的影响。数值模拟的结果表明， IL-2能够促进免疫细胞更快地发挥作用， 使得肿瘤数量最终处于一个更低的水平， 且IL-2诱导的CD8+T细胞的增殖分化是免疫系统更有效地杀伤肿瘤细胞的主要原因。此外， IL-2诱导的CD8+T细胞增殖率的增加虽然能够促进免疫细胞更快地发挥作用， 但同时也给了肿瘤细胞生存的空间， 增加了肿瘤复发的几率。当人体内的IL-2达到一定数量时， 也会出现复发的现象， 且肿瘤数量会越来越多， 因此在对癌症患者进行IL-2免疫疗法时， 存在IL-2注射的最佳剂量。

本文模型仅仅考虑了细胞因子IL-2， 然而人体内的细胞因子是多种多样的， 肿瘤细胞的种类也是错综复杂的， 因此， 后续工作中， 需要考虑更多的细胞因子和免疫细胞， 建立一个更符合实际的动力学模型， 且针对不同的肿瘤细胞， 选取不同的参数进行模拟， 以得出更具有现实意义的结论。

参考文献

[1] SUNG H， FERLAY J， SIEGEL R L， et al. Global cancer statistics 2020： GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries［J］. CA-A Cancer Journal for Clinicians， 2021， 71： 209-249.

[2] HALKOLA A S， PARVINEN K， KASANEN H. Modelling of killer T-cell and cancer cell subpopulation dynamics under immuno-and chemotherapies［J］. Journal of Theoretical biology， 2020， 488： 110136.

[3] NAKAYAMA K， ORIHATA K， YAMAGUCHI K. Surgical treatment combined with preoperative concentrated irradiation for esophageal cancer［J］. Cancer， 2015， 20： 778-788.

[4] YU D， HU W， CHEN L， et al. Effect of radiotherapy interruption due to COVID-19 outbreak［J］. Radiotherapy and Oncology， 2020， 155： 1-2.

[5] KANAVOS P. The rising burden of cancer in the developing world［J］. Annals of Oncology， 2006， 17： 15-23.

[6] ALLEN-VERCOE E， COBURN B. A microbiota-derived metabolite augments cancer immunotherapy responses in mice［J］. Cancer Cell， 2020， 38： 452-453.

[7] KIRSCHNER D， PANETTA J C. Modeling immunotherapy of the tumor-immune interaction［J］. Journal of Mathematical Biology， 1998， 37： 235-252.

[8] SMITH F O， DOWNEY S G， KLAPPER J A， et al. Treatment of metastatic melanoma using interleukin-2 alone or in conjunction with vaccines［J］. Clinical Cancer Research， 2008， 14（17）： 5610-5618.

[9] YANG J C， SHERRY R M， STEINBERG S M， et al. Randomized study of high-dose and low-dose interleukin-2 in patients with metastatic renal cancer［J］. Journal of Clinical Oncology， 2003， 21（16）： 3127-3132.

[10] CHO H， WANG Z P， LEVY D. Study of dose-dependent combination immunotherapy using engineered T cells and IL-2 in cervical cancer［J］. Journal of Theoretical Biology， 2020， 505： 110403.

[11] MESSERSCHMIDT J L， PRENDERGAST G C， MESSERSCHMIDT G L. How cancers escape immune destruction and mechanisms of action for the new significantly active immune therapies： Helping nonimmunologists decipher recent advances［J］. The Oncologist， 2016， 21（2）： 233-243.

[12] GUILLEREY C， HUNTINGTON N D， SMYTH M J. Targeting natural killer cells in cancer immunotherapy［J］. Nature Immunology， 2016， 17（9）： 1025-1036.

[13] DE PILLIS L G， FISTER K R， GU W， et al. Optimal control of mixed immunotherapy and chemotherapy of tumors［J］. Journal of Biological Systems， 2008， 16（1）： 51-80.

[14] DE PILLIS L G， RADUNSKAYA A E， WISEMAN C L. A validated mathematical model of cell-mediated immune response to tumor growth［J］. Cancer Research， 2005， 65： 7950-7958.

[15] VAN DEN DRIESSCHE P， WATMOUGH J. Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission［J］. Mathematical Biosciences， 2002， 180（1-2）： 29-48.

[16] WANG W D， ZHAO X Q. Threshold dynamics for compartmental epidemic models in periodic environments［J］. Journal of Dynamics and Differential Equations， 2008， 20（3）： 699-717.

基金资助

国家自然科学基金资助项目（12071217）