跳转至

第 11 章:脑肿瘤的生长与控制(Growth and Control of Brain Tumours)

11.1 医学背景(Medical Background)

脑肿瘤(胶质瘤 / 胶质母细胞瘤)约占所有原发性脑肿瘤的一半,预后极差。胶质瘤高度侵袭性,会浸润周围组织。CT 与 MRI 成像技术的进步使早期检测成为可能,但早期治疗的收益仍然很小;即便手术切除范围超出肉眼可见的肿瘤边界,残余浸润细胞仍会在切除边缘附近复发,最终导致死亡——这与"试图从火线后方扑灭森林大火"的困境相似,火势(肿瘤生长)主要发生在边缘。

大脑主要由灰质与白质两类组织构成:灰质由神经元和胶质细胞体组成,覆盖在大脑表层(即皮层);白质则是髓鞘化的神经轴突束,构成灰质区域之间的通路,其中连接左右大脑半球的厚带状结构称为胼胝体(图 11.2、11.3)。胶质瘤源自可分裂的神经胶质细胞,多发于大脑上部的脑半球,但也可见于脑内各处。最常见的是星形细胞瘤,按侵袭程度(等级)分为低等级星形细胞瘤、中等级的间变性星形细胞瘤和高等级的胶质母细胞瘤(多形性)。肿瘤等级反映恶性程度,基于显微镜下癌细胞的间变(行为与形态的异形程度);同一肿瘤内可含不同等级的细胞,以最高等级的细胞决定整体等级,且临床上对分级尚无统一共识。等级越高,侵袭正常组织的能力越强,但即便如此,胶质瘤也很少转移到脑外。

神经系统肿瘤患者的预后取决于多种因素,关键之一是对肿瘤时空浸润的定量评估,要兼顾解剖位置与各种治疗手段的疗效。由于本章的建模工作对临床治疗有实际意义,需要给出更详细的医学背景,因为这被认为是现实医学建模的重要组成部分。

治疗脑肿瘤的困难(Difficulties in Treating Brain Tumours) 已有大量实验与部分理论工作致力于弄清胶质瘤难以治疗的原因。与多数肿瘤不同,胶质瘤可以是高度弥散的:在大鼠脑内植入肿瘤后 7 天内,胶质瘤细胞即可在整个中枢神经系统中被识别(Silbergeld and Chicoine 1997)。植入处保留局部致密的肿瘤生长,但整个中枢神经系统分布着散在的肿瘤细胞(Silbergeld and Chicoine 1997;Silbergeld 个人通信 1998)。多数胶质瘤治疗聚焦于局部的块状肿瘤,而事实上肿瘤生长与侵袭的作用发生在别处。

胶质瘤的常规治疗包括化疗、放疗与手术干预。切除(resection,即可及肿瘤的外科摘除)历来成功率极差。切除边界处肿瘤再生长是胶质瘤研究中公认的现象(Silbergeld and Chicoine 1997;Woodward et al. 1996;Kreth et al. 1993;Kelley and Hunt 1994 等)。实验学家与理论家都相信,正是远处被侵袭的细胞导致术后肿瘤再生(Chicoine and Silbergeld 1995;Silbergeld and Chicoine 1997)。切除后残余癌细胞的密度在切除边界处最高,因此再生最可能发生于该位置。另一种假说(Silbergeld and Chicoine 1997 提出并正在检验)是:切除位点受损的脑组织释放细胞因子,召集弥散侵袭的肿瘤细胞。两种解释都与"胶质瘤的弥散特性是切除边界附近复发的根本原因"这一论断相符;区别在于前者是物理模型,后者更具生化色彩。本章研究切除治疗的模型(基本上是一个极其简单、甚至是线性的模型),并阐明其为何通常失败——结果表明它对生存期的延长极小,我们会把这些结果与预测与患者数据进行比较。

化疗本质上是用专门化学品毒杀肿瘤细胞。脑通过血脑屏障的精细毛细管结构天然抵御这些以及其他类型化学品——水溶性药物、离子和蛋白不能透过血脑屏障,但脂溶性药剂可以;最近已设计出能暂时打破血脑屏障的药剂。许多化疗是细胞周期依赖性的:药物由细胞周期的某些特定阶段触发。Silbergeld and Chicoine (1997) 观察到远离块状肿瘤的迁移细胞似乎不进入有丝分裂,因此细胞周期特异性药物与标准放疗效果有限。不仅如此,胶质瘤通常是异质性肿瘤;能够抵达癌细胞的药物又受到癌细胞异质性常伴的耐药性的阻碍——当某一类细胞对治疗有响应而死亡时,其他类型细胞正等待占据主导。这种现象需要一个包含细胞突变为耐药细胞(即多克隆模型)的模型。下面将描述并分析这样一个化疗模型,并再次将其结果与患者数据进行比较。

胶质瘤的生物学复杂性使治疗成为艰巨任务。为规划有效的治疗策略,关于肿瘤生长率与侵袭特性的信息至关重要。数学建模有助于量化切除、化疗和放疗(试图以辐射杀死肿瘤细胞)对恶性胶质瘤生长与扩散的作用。本章旨在揭示脑肿瘤治疗的若干方面,以帮助确定更佳乃至最优的疗法方案。主要目标是开发交互式计算机模型,借以考察针对特定肿瘤的各种治疗策略的效果。尽管如此,上述所有治疗的成功率都很低,迫切需要全新的治疗思路——其中数种目前正在研究中。

11.2 胶质瘤生长与侵袭的基础数学模型(Basic Mathematical Model of Glioma Growth and Invasion)

与所有肿瘤一样,胶质瘤的生物学与临床方面都很复杂,其时空生长细节仍未被充分理解。因此在建模时我们不得不做出若干重大假设。这一过程中我们重度依赖华盛顿大学神经病理学教授 E.C. Alvord 博士的医学投入,J.D. Murray 与其自 1990 年代初开始的长期合作囊括了许多学生与博士后,他们的工作都见诸下文所引论文,其中给出更多医学细节。本章将主要讨论他们的工作与相关的脑肿瘤医学研究。

秉承本书一贯的哲学,我们从尽可能简单的模型出发逐步构建。最简单的理论模型只涉及肿瘤细胞总数,肿瘤生长通常被假设为指数、Gompertzian 或 logistic(Swan 1987;Marusic et al. 1994)。这些模型没有考虑特定解剖位置上细胞的空间排布,也没有考虑癌性细胞的空间扩散;这两类空间属性在估计肿瘤生长时至关重要,因为它们决定了侵袭性与肿瘤的可见边界。在多数治疗情境(如估计手术切除的可能收益)中,必须设法确定肿瘤的浸润范围(Alvord 1991)。在出现本章所述工作之前,即便连一个解释人类胶质瘤生长与复发的简单模型都缺失,使得难以解释为何手术切除的结果如此令人失望(Nazzarro and Neuwelt 1990;Kreth et al. 1993)。本章所述工作令人惊讶的一面在于:一个非常简单的(线性)确定性模型就能提供有临床意义的、有助于患者照护的信息。

本节建立肿瘤生长时空动力学的数学模型。重要的是,我们可以从患者连续 CT 扫描与独立实验工作获得的临床数据估计模型参数(包括增殖率 / 生长率与细胞扩散系数)。这将在 11.10 节化疗部分详述,那里将描述 CT 扫描所示内容及我们的使用方法。一旦确立从组织学切片重建侵袭动力学事件的可行性,就可以用同样的数据类型研究其他不同特征生长模式、几何形状的胶质瘤,以及各种疗法(手术、化疗和/或放疗)的效果。生长模式本质上决定了不同恶性程度经典肿瘤的大体与微观特征,也决定了"混合性胶质瘤"和"多中心胶质瘤"(Alvord 1992)的特性。

先前的数学建模(Tracqui et al. 1995;Cruywagen et al. 1995;Cook et al. 1995;Woodward et al. 1996;Burgess et al. 1997)使用理论框架描述胶质瘤的侵袭性(有或没有治疗),方法是把两个特征——增殖与扩散——分离出来;这里的扩散代表胶质瘤细胞的活动迁移性。这些模型表明扩散对决定生存期比增殖率更重要。体内研究显示植入大鼠的恶性胶质瘤细胞经由白质纤维束迅速侵袭对侧脑半球(Kelley and Hunt 1994;Silbergeld and Chicoine 1997)。白质中的胶质瘤细胞扩散与灰质中不同,这一点被纳入更现实的模型。

基础模型认为胶质瘤细胞群的增长主要由增殖与扩散支配,假设肿瘤细胞呈指数增长——这对于我们关心的尺度(即患者至死亡的时间)而言是合理的生物学近似。虽然部分胶质瘤可见坏死核(死细胞区)形成,但只有高增殖、低扩散的肿瘤才显著受坏死影响;本模型暂未纳入,但可以改造以包含此项。典型的 logistic 生长在严格意义上会更准确,但我们最关心的是理解肿瘤行为的扩散性本质以及所考虑的医学尺度。logistic 生长可以平凡地纳入。(在乳腺癌的情形中,Hart et al. 1998 与 Shochat et al. 1999 已证明指数增长并非有效近似。)Silbergeld and Chicoine (1997) 提示扩散是对肿瘤细胞迁移性的良好近似。Giese and Westphal (1996) 给出关于胶质瘤侵袭的优秀综述。我们将在后文展示,扩散合理地模拟了在体外观察到的细胞扩散动力学。

设 c̄(x̄, t̄) 为位置 x̄、时间 t̄ 的细胞数。基础模型取为维数形式的守恒方程

\[ \frac{\partial \bar{c}}{\partial \bar{t}} = \bar{\nabla} \cdot J + \rho \bar{c} \tag{11.1} \]

其中 ρ(时间⁻¹)表示细胞的净增长率(含增殖与死亡 / 损失)。细胞扩散通量 J 取为与细胞密度梯度成正比:

\[ J = \bar{D} \bar{\nabla} \bar{c} \tag{11.2} \]

其中 D̄(距离²/时间)是细胞在脑组织中的扩散系数。上述理论模型把脑组织视为均匀的,因此肿瘤细胞的扩散与增长率在整个脑中都取常数。但考虑肿瘤从灰质侵袭到白质时当然并非如此。在 D̄ 为常数下,方程 (11.1) 与 (11.2) 合并为

\[ \frac{\partial \bar{c}}{\partial \bar{t}} = \bar{D} \bar{\nabla}^2 \bar{c} + \rho \bar{c} \tag{11.3} \]

后文将看到这个模型与作为其基础的 CT 扫描合理吻合,并在预测不同治疗情境下的生存时间方面给出了令人惊讶的好结果(Tracqui et al. 1995;Cruywagen et al. 1995;Cook et al. 1995;Woodward et al. 1996;Burgess et al. 1997)。虽然模型相对其粗略简化假设而言给出了惊人的好结果,但它们包含若干可重新考虑的基本简化(Swanson 1999;Swanson et al. 2000 已这样做)。例如,给定某处胶质瘤细胞源,多数先前模型出于数值简化(并忽略解剖边界)而假定"可检测细胞前沿"从源头对称地向外扩展;它们当然知道临床与实验观察表明实际上肿瘤生长的对称性并不成立。本节讨论的第一个模型也处理这一点以及组织异质性。

白质充当灰质区之间胶质瘤细胞侵袭的通道。胶质瘤细胞在白质中的扩散系数(迁移性)大于灰质。体内研究显示植入大鼠的恶性胶质瘤细胞经白质纤维束迅速侵袭对侧半球(Chicoine and Silbergeld 1995;Silbergeld and Chicoine 1997;Kelley and Hunt 1994)。本节研究的模型纳入组织异质性对细胞扩散与肿瘤增长率的影响,以更准确地模拟临床与实验观察到的可见肿瘤边界的不对称性。

含空间异质性的模型(Model with Spatial Heterogeneity) 我们可以通过令扩散系数 D̄ 为空间变量 x̄ 的函数来在模型中纳入空间异质性,由此区分灰质与白质区域。这给出替代 (11.3) 的

\[ \frac{\partial \bar{c}}{\partial \bar{t}} = \bar{\nabla} \cdot (\bar{D}(\bar{x}) \bar{\nabla} \bar{c}) + \rho \bar{c} \tag{11.4} \]

在脑与脑室的解剖边界上取零通量边界条件。于是,若 B 为 (11.4) 待解的脑区域,则边界条件为

\[ n \cdot \bar{D}(\bar{x}) \bar{\nabla} \bar{c} = 0 \quad \text{for } x \in \partial B \tag{11.5} \]

其中 n 是 ∂B 的单位法向。在解剖学真实的脑(我们实际上会采用)的几何复杂度下,这即便在二维情形也是极难的解析问题与不平凡的数值问题。我们在二维与三维情形下研究了这一问题,所用白质灰质分布由 EMMA(Extensible MATLAB Medical Analysis)从解剖学准确的图像得到;这是 Collins et al. (1998) 在蒙特利尔神经学研究所 McConnell 脑成像中心开发的用以辅助医学图像分析的程序包。图 11.4 显示 EMMA 给出的脑水平切片白质灰质分布一例。

在二维情形下,Swanson (1999) 映射出白质灰质区域以模拟原初二维模型 (11.3) 所用 CT 扫描的临床数据。她的初始目标是在质量上比仅考虑均匀脑组织的原模型更准确地复现这些 CT 图像。如上所述,即便使用这种简化,结果也与多种大体医学量及患者数据相符。我们对灰白质分布如何影响胶质瘤生长与侵袭尤其感兴趣。

我们首先使用简化的大鼠脑解剖区域;该较简单情形的数值模拟与解析结果将与华盛顿大学医学院神经病理学实验室的(合作者)Ellsworth Alvord Jr. 博士(神经病理学主任)与 Daniel Silbergeld 博士(神经外科)所取得与分析的实验数据作比较。随后我们将在解剖学准确的人脑区域上详尽考虑该模型的二维与三维情形。

我们先对含空间异质性的模型作无量纲化,如常所做的那样,这同时减少系统中的有效参数数量,并(在不考虑单位的前提下)让我们对各项的相对重要性有所认识。为对相关数字有所概念,一个小肿瘤中可以有 10¹¹ 个癌细胞,而扩散系数的数量级为 10⁻⁴ cm²/day。参数估计将在后文若干表格中给出。

我们假设在两种组织(白质与灰质)中扩散系数分别为常数(但彼此不同)。于是我们须求解

\[ \frac{\partial \bar{c}}{\partial \bar{t}} = \bar{\nabla} \cdot (\bar{D}(\bar{x}) \bar{\nabla} \bar{c}) + \rho \bar{c} \tag{11.6} \]

其中

\[ \bar{D}(\bar{x}) = \begin{cases} D_w & \text{for } \bar{x} \text{ in white matter} \\ D_g & \text{for } \bar{x} \text{ in grey matter} \end{cases} \tag{11.7} \]

并以零通量边界条件 \(n \cdot \bar{D}(\bar{x}) \bar{\nabla} \bar{c} = 0\) for x ∈ ∂B (11.8) 以及初始条件 c̄(x̄,0) = f̄(x̄) 为约束;f̄(x̄) 的适当形式将在下面给出。

引入无量纲变量 \(x = \sqrt{\rho/D_w} \bar{x}\)\(t = \rho \bar{t}\)\(c(x, t) = (D_w / (\rho N_0)) \bar{c}(\sqrt{\rho/D_w} \bar{x}, \rho \bar{t})\),其中 \(N_0 = \int \bar{f}(\bar{x}) d\bar{x}\) 是 t̄=0 时脑中肿瘤细胞的初始总数。由此 (11.6) 成为

\[ \frac{\partial c}{\partial t} = \nabla \cdot (D(x) \nabla c) + c \tag{11.10} \]

其中

\[ D(x) = \begin{cases} 1 & \text{for } x \text{ in white matter} \\ \gamma = D_g / D_w & \text{for } x \text{ in grey matter} \end{cases} \tag{11.11} \]

并有 \(c(x,0) = f(x) = (D_w/(\rho N_0)) \bar{f}(\sqrt{\rho/D_w} \bar{x})\)\(n \cdot D(x) \nabla c = 0\) for x ∈ ∂B。在该无量纲化下,扩散以白质扩散为空间尺度,时间以肿瘤生长的尺度度量。

现在考察初始条件的形式及 f(x) 的适当形式。理论上肿瘤源自一个癌细胞,但该原始细胞出现的时间、生长类型与早期癌细胞的扩散均未知。我们假设在首次扫描时扩散过程已经打破了此前任何可能的均匀分布。因此,细胞呈正态分布,在肿瘤中心 x₀ 处有最大细胞密度 a,即

\[ c(x, 0) = a \exp\left(-\frac{|x - x_0|^2}{b}\right) \tag{11.12} \]

其中 b 是肿瘤细胞扩散范围的度量。

生存时间的初始解析估计(Initial Analytical Estimate of Survival Time) 临床上极低密度的肿瘤细胞无法被检出。在 CT 扫描上,肿瘤轮廓由某一非零分辨率水平定义,对应细胞密度 c∗(约为 40,000 cells/cm²)。低于该阈值,癌细胞无法被成像技术检出。典型的连续 CT 扫描本质上是肿瘤所在脑的几个不同层面的一系列图像;据此可重构大致的三维形状。图 11.5 是一个三层肿瘤的理想化扫描。利用这些可估计肿瘤大小与等效球面半径;详见 11.10 节,那里将直接使用患者扫描估计参数以确定肿瘤大小与生长。

在常增长率 ρ 与均匀扩散 D̄(x̄) = D 的情形下,二维空间方程 (11.3) 在 x̄=0 处的初始 δ 函数源 N₀ 个肿瘤细胞下有解

\[ \bar{c}(\bar{x}, \bar{t}) = \frac{N_0}{4\pi D \bar{t}} \exp\left(\rho \bar{t} - \frac{\bar{r}^2}{4D\bar{t}}\right) \tag{11.13} \]

其中 r̄ 为轴对称径向坐标。若可检出阈值密度为 c̄∗,则由 (11.13) 与 c̄ = c̄∗ 给出肿瘤轮廓半径 r̄∗:

\[ \bar{r}_* = 2\sqrt{\frac{D}{\rho} \bar{t}} \left(1 - \frac{1}{\rho \bar{t}} \ln \frac{4\pi D \bar{t} \bar{c}_*}{N_0}\right) \sim 2\sqrt{\frac{D}{\rho} \bar{t}} \quad \text{for } \bar{t} \text{ large} \tag{11.14} \]

最后一个表达式当然就是轴对称 Fisher–Kolmogoroff 方程径向行进波(速度为 \(2\sqrt{D\rho}\),参见第 I 卷第 13 章)的渐近形式。若肿瘤在半径 r̄_detect 时被识别、在半径 r̄_lethal 时致死,则未治疗生存时间可近似为

\[ \text{survival time} = \bar{t}_{\text{lethal}} - \bar{t}_{\text{detect}} \approx \frac{1}{\sqrt{D\rho}} (\bar{r}_{\text{lethal}} - \bar{r}_{\text{detect}}) \]

这表明 D 与 ρ 都是决定生存时间的重要参数:增大 ρ 或 D 中的任一个都会缩短生存时间。肿瘤被识别的平均半径(或更准确说等效半径)为 1.5 cm,致死半径为 3.0 cm(Burgess et al. 1997)。这两个半径只是临床上观察到的平均值(Blankenberg et al. 1995),范围相当大。对高等级肿瘤(参数值见 11.5 节表 11.5),模型所暗示的未治疗生存时间约为 200 天。这与已观察到的胶质母细胞瘤经治疗后的 6–12 个月中位生存期一致(Alvord 1991)。

在 (11.9) 所定义的无量纲变量下,

\[ r_* = 2\sqrt{t} \left(1 - \frac{1}{t} \ln \psi t\right), \quad \psi = \frac{4\pi D \bar{c}_*}{\rho N_0} \tag{11.15} \]

取决于参数 ψ 的值,观察到的肿瘤半径可能随时间单调增加,或在肿瘤起始与可见进展之间出现一段延迟——这段时间有时称为"建立期"。基本上它是肿瘤在被扩散推动前扎根所需的时间。建立期长度 t_e 可由令 r∗=0 与 t=t_e 从 (11.15) 给出:

\[ t_e = \ln \psi t_e \tag{11.16} \]

建立期只在 ψ>1 时存在。立即可由 ψ = e^{t_e}/t_e 的曲线看出——其极小值在 t_e=e 处,因此仅当 ψ>e 时才有解。这定义了所有参数间的关系,可在关联实验观察与模型参数时有用。由 (11.16),建立期长度随 ψ 增大而增长;又因 ψ ∝ D/ρ,故随 D/ρ 增大而增长。

对固定 ψ,无量纲可检肿瘤半径 r∗ 由单条曲线给出。因此可以设想一个高等级肿瘤(高 ρ,高 D)与一个低等级肿瘤(低 ρ,低 D)使比值 D/ρ(因而 ψ)固定。两条肿瘤的可检半径 r∗ 以同样的方式增长,但低等级肿瘤需要更长时间才能达到致死大小。虽然固定 ψ 下空间侵袭特征相同,所发生的时间尺度可以截然不同。

异质环境中的空间扩散主要在生态学背景下被研究。例如 Cantrell and Cosner (1991) 用带空间变化增长率的扩散 logistic 方程建模了异质环境下种群动力学;他们关心有利斑块与不利斑块的空间排布对整体环境适宜度的影响。Cruywagen et al. (1996) 用反应扩散模型在空间异质环境上建模了基因工程生物体的扩散风险;对空间周期扩散率与承载力的情形确定了生物体的侵袭条件;该工作在第 1 章讨论。Shigesada and Kawasaki (1997) 在其生物入侵专著中讨论了带分段常数空间变化扩散系数的扩散;环境由有利与不利斑块排布而成,扩散与 logistic 增长空间周期时,他们证明存在行进周期波解。

关于脑肿瘤在异质区域内的空间扩散,Swanson (1999) 对模型方程 (11.10) 与 (11.11) 进行了 1D 分析,其中空间区域由嵌入无穷灰质域的小 O(ε) 白质域构成。边界条件包括灰白区域边界处细胞密度与细胞通量的连续性。从渐近观点看,ε≪1 时由于白质区域扩散系数大于灰质区域,直觉上预期细胞密度在白质区域近似为常数——这一点已被解析地确认。

11.3 体外肿瘤扩散:参数估计(Tumour Spread In Vitro: Parameter Estimation)

体外实验常被用于帮助表征体内行为。必然地,多个参数估计来自此类体外实验。在基于实验数据估计参数的背后,我们一般假设建模方法与模型适用于这些实验。本节考察两个这样的实验以确定参数估计。我们将较为详细地讨论此事,因为参数估计在预测可能的患者治疗时尤为重要。

正常胶质细胞的迁移率很低(Silbergeld and Chicoine 1997),而胶质瘤细胞可表现出异常高的迁移率(Chicoine and Silbergeld 1995;Giese et al. 1996a,b,c;Giese and Westphal 1996;Pilkington 1997a,b;Silbergeld and Chicoine 1997;Amberger et al. 1998;Giese et al. 1998)。利用延时视频显微术和其他技术,Chicoine and Silbergeld (1995) 量化了脑肿瘤细胞在体内外的迁移与侵袭能力,结果提示:体外人类胶质瘤细胞的平均线速度为 12.5 µm/hr,体内最低线速度为 4.8 µm/hr。

Chicoine and Silbergeld (1995) 的体外实验:细胞迁移性 Chicoine and Silbergeld (1995) 发展了一种称为"径向皿测定"的肿瘤细胞迁移测定法。简单来说,将 2×10⁴ 个细胞铺在一个 8 cm 直径培养皿中心一个 2 cm 直径的小圆盘上。抑制细胞有丝分裂,每日显微观察皿上细胞群的空间扩散。我们可以用上述模型(适用于该实验)估计胶质瘤细胞在体外的扩散系数。

该情形下无生长(ρ=0),且扩散无异质性 D̄(x̄)=D。带径向对称的维数化模型 (11.17) 为

\[ \frac{\partial \bar{c}}{\partial \bar{t}} = D \bar{\nabla}^2 \bar{c} \tag{11.17} \]

以零通量边界条件

\[ n \cdot \bar{\nabla} \bar{c}(\bar{r}, \bar{t}) = 0 \quad \text{for } \bar{r} = R_0 \tag{11.18} \]

其中 R₀ 为培养皿半径;在该测定中 R₀=4 cm。最初细胞均匀分布在半径为 R 的中心圆面内:

\[ \bar{c}(\bar{r}, 0) = \bar{c}_0 H(R - \bar{r}) \tag{11.19} \]

其中 H 为 Heaviside 函数,R=1 cm 在 Chicoine et al. (1995) 的实验中。由于有丝分裂被阻断,无生长,皿中共 N = c̄₀ πR² 个细胞贯穿整个实验。

现通过令 \(x = \bar{x}/R_0\)\(t = D \bar{t}/R_0^2\)\(c(x, t) = \bar{c}(\bar{x}/R_0, D\bar{t}/R_0^2)/\bar{c}_0\) 作无量纲化,其中 N = c̄₀ πR² 是 t=0 时肿瘤细胞的初始数。无量纲模型仅涉及径向坐标 r(在径向对称下),成为

\[ \frac{\partial c}{\partial t} = \frac{\partial^2 c}{\partial r^2} + \frac{1}{r} \frac{\partial c}{\partial r}, \quad 0 < r < 1; \quad c(r, 0) = H(\lambda - r), \quad \frac{\partial c}{\partial r} = 0 \text{ at } r = 1 \tag{11.21} \]

其中 λ = R/R₀;在径向皿测定实验中 λ = 1/4。

细胞密度 c(r, θ, t) 的渐近近似(Asymptotic Approximation for the Cell Density) 若假定培养皿足够大(R₀ ≫ R),可用 (11.13) 令 ρ=0 时的解——无穷域上简单扩散方程的解——作为近似。然而 λ=1/4 不算小,故必须使用更精确的近似。

(11.21) 的解是涉及 Bessel 函数的经典解,可由标准方法得到(例如见 Kevorkian 1999 或 Carslaw and Jaeger 1959)。它可通过叠加半径 r₀ 处的环的基本解并从 0 积分到 λ 得到。

由于我们主要关心大 t 情形,只需解的渐近近似。对位于 r₀ 的细胞环,其解为

\[ c_{\text{ring}}(r, t; r_0) = \frac{1}{4\pi t} \exp\left(-\frac{r^2 + r_0^2}{4t}\right) I_0\left(\frac{r r_0}{2t}\right) \tag{11.22} \]

其中 I₀ 为修正 Bessel 函数。对 r r₀/2t 较小情形,该解近似为

\[ c_{\text{ring}}(r, t; r_0) \sim \frac{1}{4\pi t} \exp\left(-\frac{r^2 + r_0^2}{4t}\right)\left\{1 + \frac{1}{4}\left(\frac{r r_0}{2t}\right)^2\right\} + O\left(\left(\frac{r r_0}{2t}\right)^4\right) \tag{11.23} \]

当 r₀=0 时,这正好给出原点处点源细胞的精确解。

对半径 R 的细胞盘,我们有 \(c(r,0) = H(\lambda - r) = 2\pi \int_0^\lambda c_{\text{ring}}(r,0;r_0) r_0 dr_0\)。结合近似 (11.23),完整问题 (11.21) 带盘状初始条件的解为

\[ c(r, t) = 2\pi \int_0^\lambda c_{\text{ring}}(r, t; r_0) r_0 dr_0 \sim e^{-r^2/4t}\left(1 - e^{-\lambda^2/4t}\right) + \frac{r^2}{4t} e^{-r^2/4t}\left[1 - \left(1 + \frac{\lambda^2}{4t}\right) e^{-\lambda^2/4t}\right] + \cdots \tag{11.24} \]

对 ν = r λ/2t 较小成立。在 Chicoine and Sibergeld (1995) 的实验中 ν = r/(8t)。

现用近似解 (11.24) 得到胶质瘤细胞的体外扩散系数估计。我们预期该估计对小 r 与 t~O(1) 或更大时最精确。图 11.6(a) 是 Chicoine and Silbergeld (1995) 径向皿测定细胞密度的实验测量图,叠加上 (11.24) 的(维数化)渐近解(D = 0.002 cm²/hr)在不同实验时间的曲线。Swanson (1999) 以 D=0.002 cm²/hr 对完整问题 (11.22) 作数值求解并用同一实验数据,见图 11.6(b)。两图比较表明渐近近似 (11.24) 用于确定参数估计相当精确。

细胞距初始接种物的平均距离 ⟨r⟩ 的渐近近似 Chicoine and Silbergeld (1995) 也计算了径向皿测定法中胶质瘤细胞距原点的平均距离 ⟨r⟩。为确定 ⟨r⟩,须对所有 r 积分。渐近解 (11.24) 仅对 r λ/2t 较小严格成立,而对任意大 r 不算小。因此实际上必须完整计算 ⟨r⟩。

若 c(r,t) 是 t 时刻细胞的分布,则细胞距原点的平均径向距离 ⟨r⟩ 为

\[ \langle r \rangle = \frac{\int_0^\infty r^2 c(r,t) dr}{\int_0^\infty r c(r,t) dr} = \frac{2\pi \int_0^\infty r^2 c(r,t) dr}{\int_0^\infty r c(r,t) dr} \tag{11.25} \]

其中 c(r,t) 由 (11.24) 的积分给出。该答案含一个积分,我们用 Laplace 方法(涉及小量与大量积分的渐近评估;例如可参见 Murray 1984 的教学讨论)对小 t 与大 t 分别作近似。对小 t(A=1/(4t))有

\[ \langle r \rangle \sim \frac{3 + 2A\lambda^2}{6A\lambda}\left(1 + \text{erf}(\sqrt{A}\lambda)\right) - \frac{5}{6\lambda^2 \sqrt{A^3 \pi}} + \frac{5 + 2A\lambda^2}{6\sqrt{\pi} A\lambda^2} e^{-A\lambda^2} = \frac{6t + \lambda^2}{3\lambda}\left(1 + \text{erf}\left(\frac{\lambda}{2\sqrt{t}}\right)\right) - \frac{20}{3\lambda^2 \sqrt{t^3/\pi + 10t}} + \frac{\lambda^2}{3\sqrt{\pi}\lambda^2} e^{-\lambda^2/4t} \sim \frac{2\lambda}{3} + \frac{4t}{\lambda} + \cdots \tag{11.26} \]

随 A=1/(4t) → ∞。当 t→0(A→∞),⟨r⟩ 趋于 2λ/3,这对应于初始细胞分布的精确平均径向距离。

对大 t(小 A),平均位移为

\[ \langle r \rangle = \sqrt{\pi t}\left(1 + \frac{\lambda^2}{16t} - \frac{\lambda^4}{768t^2} + \cdots\right) \tag{11.27} \]

对 t 较大,平均半径收敛到原点处点源情形:⟨r⟩ = √(πt)(以维数化形式 ⟨r̄⟩=√(π D t̄))。其含义是长时间后扩散将细胞群扩散得足够开,以至于初始分布不能被识别。然而对长时间,培养皿边界效应变得重要。由于 Chicoine and Silbergeld (1995) 实验中皿内细胞不被允许生长,t→∞ 时细胞密度趋于均匀稳态 c→λ²(以维数化形式 c̄→c̄₀ λ²)。因此对长时间行为的更精确渐近估计应使用该稳态来计算平均位移:

\[ \langle r \rangle = \frac{2\pi \int_0^\infty r^2 c(r,t) dr}{2\pi \int_0^\infty r c(r,t) dr} \to \frac{2\pi \int_0^1 \lambda^2 r^2 dr}{2\pi \int_0^1 \lambda^2 r dr} = \frac{2\pi/3 \lambda^2}{\pi \lambda^2} \tag{11.28} \]

这是培养皿实验的长时间极限。与实验结果比较,这些近似都不算非常好,但小 t 近似更好。

Chicoine and Silbergeld (1995) 实际计算"平均半径"的方式与我们这里做的略有不同。他们的计算忽略了初始半径 R 内细胞的贡献,相当于把 (11.25) 积分的下限从 0 改为 λ。若把这个新的"平均半径"记为 ⟨r∗⟩,则它满足

\[ \langle r_* \rangle = \frac{\int_\lambda^\infty r(r - \lambda) c(r,t) dr}{\int_\lambda^\infty r c(r,t) dr} = \frac{2\pi}{1 - \lambda^2} \int_\lambda^\infty r(r - \lambda) c(r,t) dr \tag{11.29} \]

Swanson (1999) 用与上面 ⟨r⟩ 相似的近似方法证明了

\[ \langle r_* \rangle \sim \frac{2\pi t^2}{1 - \lambda^2} + \frac{40}{3(1 - \lambda^2)} \sqrt{\pi t^3} + \cdots \quad \text{for small } t (\text{large } A) \tag{11.30} \]

注意 R=1 cm 在实验中给出 λ=1/4。

观察到胶质母细胞瘤细胞在 96 小时内即可抵达培养皿边缘,这界定了图 11.7(c) 中对大 t 值渐近结果的高估。

在图 11.7 中我们画出由实验观察得到的 ⟨r∗⟩ 与渐近表达式 (11.30)(严格地只对小 t 成立)的曲线。该图暗示了间变性星形细胞瘤、混合性胶质瘤与胶质母细胞瘤多形性细胞的扩散系数估计分别为 1.6×10⁻⁴ cm²/hr、2×10⁻³ cm²/hr 与 3×10⁻³ cm²/hr。因此,恶性程度越高,细胞迁移性越大。事实上这些结果表明中等级间变性星形细胞瘤与高等级胶质母细胞瘤的迁移性在体外大约有两倍差异。

在图 11.7 中看到,对胶质母细胞瘤细胞,渐近展开 (11.26) 在 t̄>96 小时时高估了平均半径。正如预期:Chicoine and Silbergeld (1995) 观察到胶质母细胞瘤细胞移动如此之快以至于 96 小时就能抵达皿边缘,因此边界效应变得重要。间变性星形细胞瘤与混合性胶质瘤细胞直到较晚时间才能抵达培养皿边缘。

由 ⟨r⟩ 的计算,存在皿中肿瘤细胞群的稳态极限——即对大 t 时 c→λ²。我们可以计算该稳态对 ⟨r∗⟩ 极限值的影响:

\[ \langle r_* \rangle \to \frac{\int_\lambda^1 r(r - \lambda) \lambda^2 dr}{\int_\lambda^1 \lambda^2 r dr} = \frac{2 - 3\lambda + \lambda^3}{3(1 - \lambda^2)} \tag{11.31} \]

随 t→∞。该极限值包含在图 11.7 中,界定了 (11.30) 的渐近近似对图 11.7(c) 中胶质母细胞瘤多形性细胞的高估。

Giese et al. (1996a,b,c) 的体外实验:细胞迁移性与增殖 如前述,胶质瘤细胞在白质中的迁移性高于灰质。Giese et al. (1996b) 进行了体外实验以表征胶质瘤细胞在髓鞘(白质的关键组分)上的增强迁移性。这些实验与 Chicoine and Silbergeld (1995) 类似,但时间尺度更短(约 40 至 100 小时),主要区别在于 Giese et al. (1996b) 允许细胞在实验中增殖。他们还以表格形式给出了肿瘤细胞可见前沿半径的增加。此外,细胞的初始接种物覆盖区域小于 Chicoine and Silbergeld (1995) 的实验。因此预期实验可由下式更恰当地建模:

\[ \frac{\partial \bar{c}}{\partial \bar{t}} = \bar{D} \bar{\nabla}^2 \bar{c} + \rho \bar{c} \quad \text{for } \bar{r} \le R_0 \tag{11.32} \]

其中 R₀ 为细胞被允许迁移的皿半径,Laplacian 为轴对称,并有零通量边界条件

\[ n \cdot \bar{\nabla} \bar{c}(\bar{r}, \bar{t}) = 0 \quad \text{for } \bar{r} = R_0 \tag{11.33} \]

n 为皿边外向法向。最初在原点有点源 N 个细胞:

\[ \bar{c}(\bar{r}, 0) = N \delta(\bar{r}) \tag{11.34} \]

对充分大的 R₀,解近似为 (11.13),即

\[ \bar{c}(\bar{r}, \bar{t}) \sim \frac{N}{4\pi D \bar{t}} \exp\left(\rho \bar{t} - \frac{\bar{r}^2}{4D\bar{t}}\right) \tag{11.35} \]

为从实验数据估计扩散系数,我们使用 Fisher–Kolmogoroff 近似。在第 I 卷第 13 章详细研究的一维 Fisher–Kolmogoroff 方程中,行进波速度为 \(v = 2\sqrt{\rho D}\),D 与 ρ 分别为扩散系数与线性增长率。于是仅受增长与扩散支配的种群在大时间后以 \(2\sqrt{\rho D}\) 的速率扩展。在二维情形 (11.32)–(11.34) 下,若边界 R₀ 足够远可近似为无穷域,则肿瘤细胞浓度由 (11.35) 给出。由于实验上我们可以追踪某检测阈值 c∗ 以上的肿瘤细胞轮廓,代入 (11.35) 并对 D 求解给出

\[ D = \frac{\bar{r}^2}{4\bar{t}} \left[\rho \bar{t} - \ln(4 c_* \pi D \bar{t}/N)\right]^{-1} \approx \frac{\bar{r}^2}{4\rho \bar{t}^2} \quad \text{for large } \bar{t} \text{ since } \rho \bar{t} \gg \ln \bar{t} \]

\(v = \bar{r}/\bar{t}\) 为轮廓前沿速度,则

\[ D \approx \frac{v^2}{4\rho} \tag{11.36} \]

因此我们用 \(D_w = v_w^2/(4\rho)\)\(D_g = v_g^2/(4\rho)\)(其中 v_w 与 v_g 分别为实验观察到的线速度)来关联白质与灰质中的扩散系数。为导出扩散系数的估计,需确定肿瘤细胞的增长率。现将此近似 (11.36) 与实验结果一起使用。

这里应插入一个警示。该推导下,Fisher–Kolmogoroff 估计意味着很低的增长率配合很大的线速度会对灰白质中的肿瘤细胞给出很高的扩散率估计。由于增殖率与迁移率一般(正向地)与恶性程度相关(如下文所示),具高线速度的肿瘤细胞群可具有高(或至少不小)的增长率。当然,某些 ρ 小、D 大的中等级肿瘤可能产生 D 更大、ρ 仅微弱增长(甚至不增长)的子群。

Giese et al. (1996b) 测定了髓鞘上种群随时间的范围扩展。即给定某检测阈值 c∗,可检测肿瘤细胞区域的半径被记录为时间的函数。由 (11.14) 可检测肿瘤细胞群的半径 r∗ 在大 t 时为

\[ r_* \sim 2\sqrt{D \rho t} \tag{11.37} \]

图 11.8 是 Giese et al. (1996b) 实验观察(三个胶质母细胞瘤细胞系 G-112、G-140、G-168)的可检测半径对时间的线性最小二乘拟合。直线的斜率被取为等于 \(2\sqrt{\rho D}\)。因此在估计出增长率 ρ 后,可从这些线性最小二乘拟合反推扩散系数 D 的值。

表 11.1 给出图 11.8 中线性拟合的斜率。假设当细胞在对照细胞外基质(ECM)上或髓鞘上迁移时 ρ 不变,可由表 11.1 的斜率导出 ECM 上与髓鞘上扩散系数之间的关系。记 ECM 上的扩散系数与线速度分别为 D_ECM 与 v_ECM;类似地,记髓鞘上的分别为 D_m 与 v_m。线速度则满足

\[ v_{\text{ECM}} = 2\sqrt{\rho D_{\text{ECM}}}, \quad v_m = 2\sqrt{\rho D_m} \implies \frac{D_m}{D_{\text{ECM}}} = \left(\frac{v_m}{v_{\text{ECM}}}\right)^2 \tag{11.38} \]

由表 11.1 给出的速度值可推出髓鞘上扩散系数与 ECM 上扩散系数的比值约为 2。若把 D_ECM 与灰质中的扩散系数关联,则可推测白质中的扩散系数是灰质的两倍。体外实验条件与体内观察到的不同,但从该结果可预期灰白质间扩散系数至少有两倍差异。

与所有实验结果一样,存在一个问题:培养皿中观察到的行为是否与体内行为类比。为处理这一问题,常考虑一系列实验条件。表 11.2 给出 Giese et al. (1996b) 另一系列实验的结果。表 11.2 的每一行对应于供应给细胞的胎牛血清量略有不同的实验条件。随血清浓度增加(实验 1→4),肿瘤细胞的增长率随血清浓度上升,但迁移性达到极大值后下降。由各实验条件,我们现用 Fisher–Kolmogoroff 近似对长时间下扩张率的估计(即 \(2\sqrt{\rho D}\),本质上为 (11.37) 关于 t 的导数)从 \(D = v^2/(4\rho)\) 反推胶质母细胞瘤细胞的扩散系数。

11.4 大鼠脑内的肿瘤侵袭(Tumour Invasion in the Rat Brain)

大鼠通常用于研究肿瘤扩散的体内动力学。对胶质瘤而言,典型做法是将肿瘤植入大鼠脑皮层并让之生长扩散。在某较后时间处死大鼠并近似肿瘤动力学。大鼠胶质瘤的实验研究提供了肿瘤细胞最低线速度的估计。虽然胶质瘤细胞并不沿线性路径行进,但若在处死时某些细胞在距原始植入位置某距离处被识别,这些细胞必定以某最低线速度行进至该距离。我们可将该线速度对应于可检测肿瘤半径的速度,此处检测阈值浓度非常低。描述大鼠体内肿瘤细胞动力学的模型方程是 11.2 节推导并由 (11.4) 给出(维数形式)的含空间异质性模型。

大鼠体内参数估计(In Vivo Parameter Estimation for the Rat Brain) 利用延时视频显微术与其他技术,Chicoine and Silbergeld (1995) 与 Silbergeld and Chicoine (1997) 量化了大鼠体内脑肿瘤细胞的迁移与侵袭能力。结果提示:体内最低线速度为 4.8 µm/hr。

我们对白质与灰质区域中参数值的变化尤其感兴趣。据 D.L. Silbergeld 个人通信 1998,灰质(v_g)与白质(v_w)中的平均线速度分别为 36 µm/hr 与 70 µm/hr。这些值高于上节给出的估计,但代表了评估迁移性的新实验技术的演化。这凸显了确定决定性参数估计的一个持续问题。即使使用最好的理论技术,我们仍须依赖不断变化的最佳实验技术,因此参数估计不可避免地会改变。

胶质瘤细胞不沿线性路径行进,这就是我们以随机游走建模胶质瘤迁移性的原因。要将线速度 v、增殖率 ρ 与随机游走扩散 D 关联,我们再次使用 Fisher–Kolmogoroff 近似 (11.36) \(D \approx v^2/(4\rho)\),Burgess et al. (1997) 也这样用过。我们再次以 \(D_w = v_w^2/(4\rho)\)\(D_g = v_g^2/(4\rho)\)(v_w、v_g 为实验观察到的线速度)关联白质灰质中的扩散系数。为导出扩散系数的估计,须确定肿瘤细胞的增长率。Alvord and Shaw (1991) 引用体内胶质瘤倍增时间为一周至一个月。由此得出的大鼠体内参数估计列于表 11.4。这些数值下,使用 Fisher–Kolmogoroff 估计(低增长率与高线速度配合)对灰白质中的肿瘤细胞给出很高的扩散率估计。由于增殖率与迁移率一般(正向地)与恶性程度相关(自上所见),具高线速度的肿瘤细胞群可具有高(或至少不小)的增长率。因此表 11.4 给出的很大扩散系数可能不那么精确。

大鼠脑内肿瘤侵袭的数值模拟(Numerical Simulations of Tumour Invasion in the Rat Brain) 在先前未治疗肿瘤侵袭的研究(Burgess et al. 1997)中,脑被视作均匀的。Swanson (1999) 研究了因引入白质灰质异质性而增加的侵袭能力。我们首先聚焦于图 11.9 所示简单大鼠脑拓扑。为示意目的,她假设白质中的扩散系数约为灰质中的 10 倍。表 11.4 的参数估计表明实验结果指示白质中的扩散系数为灰质中的 2 至 100 倍。

通常认为胼胝体是肿瘤细胞侵袭对侧半球的常见通道。在图 11.9 较简单大鼠脑拓扑的冠状面中,胼胝体是连接左右脑半球的紧凑白质纤维拱。在模拟中我们将大鼠脑的胼胝体近似为图 11.10 中曲线(一对角状)所围的规则区域。

在初步模拟中,胼胝体通过由胼胝体辐射出的白质纤维与灰质皮层相连,呈现为回旋镖形状,如图 11.10(a)–(d) 所示。注意胼胝体是一大束白质纤维。从胼胝体延伸到皮层的这些纤维束比胼胝体窄得多。它们被表示为从胼胝体延伸的直线。在人脑中,部分连接皮层与胼胝体的(曲线)纤维束可被一条细的癌细胞轨迹所识别。Swanson (1999) 模拟了 (11.10) 与 (11.11) 所定义的模型系统,零通量边界条件与 11.2 节讨论的连续与守恒条件,证明肿瘤侵袭如何被这些细纤维束与胼胝体所促进——这些纤维帮助肿瘤细胞更快地侵袭更远的区域。Giese and Westphal (1996) 的综述展示了若干肿瘤扩散被人脑胼胝体所促进的例子。

图 11.10 的模拟(其中植入一小肿瘤并允许其按模型扩散与生长)显示了白质的存在使侵袭能力显著增强——肿瘤细胞能在相对短的时间内基本填满胼胝体,正如实验和临床上在人脑中所观察到的那样(Giese and Westphal 1996)。

实验中常注意到肿瘤细胞横向越过胼胝体。正是因为这种横向侵袭,故认为胶质瘤细胞优先沿白质胼胝体的纤维行进(例如见 Giese and Westphal 1996)。图 11.10 的模拟显示扩散的各向异性并非造成肿瘤边缘沿胼胝体横向大量移动的必要条件。可以认为灰质与白质中扩散率的不同是一种各向异性。该横向侵袭被脑在胼胝体区域附近的物理结构所强化——胼胝体下方的脑室与上方皮层的裂隙共同促成了横向侵袭。

图 11.11 是大鼠脑上带解剖边界的模拟,以及图 11.9 中空间异质灰白质分布。由于大鼠脑皮层的限制,肿瘤细胞在边界附近累积。此外,由于胼胝体位于脑室与皮层之间,可预期见到穿越胼胝体至对侧半球的细胞侵袭显著增加。

11.5 人脑内的肿瘤侵袭(Tumour Invasion in the Human Brain)

大鼠脑的简单几何是胚胎人类脑拓扑的良好模型(彼时皮层平滑)。成年人脑的复杂几何只会强化我们在大鼠脑几何中已看到的效果。我们预期由于灰白质的异质分布,边界效应会增加、侵袭模式更复杂。本节讨论模型在人脑几何上的应用,并考察异质扩散对引入虚拟肿瘤扩散的影响。

将肿瘤侵袭模型应用于人脑的一项主要进展由 BrainWeb 数据库(Collins et al. 1998)与前述 EMMA 提供。简言之,EMMA(Extensible MATLAB Medical Analysis)是用于处理医学图像的工具。BrainWeb 数据库使用 MRI 模拟器创建,定义了在三维 181×217×181 网格上脑中灰质与白质的位置与分布。它被创建以便在 MATLAB 中可视化冻结图像,从而可被处理与研究。本节描述 Swanson et al. (2000) 给出的模型系统 (11.10) 与 (11.11) 的若干数值模拟。结果取脑的二维切片以便于表示。三维研究也已进行,部分将在下面讨论。

参数估计取自中高等级星形细胞瘤(胶质母细胞瘤多形性)的数据。高等级星形细胞瘤约占所有星形细胞瘤的 50%,而后者又约占所有原发性脑肿瘤的 50%(Alvord and Shaw 1991)。我们使用 Alvord and Shaw (1991) 给出的胶质瘤倍增时间 1 周至 12 个月。因此对高等级星形细胞瘤的细胞增殖,取增长率 ρ=0.012/day,对应 60 天的倍增时间。

上一节我们讨论了胶质瘤细胞在体外与大鼠体内的显著迁移性,并使用 Chicoine and Silbergeld (1995) 与 Silbergeld and Chicoine (1997) 的实验数据估计模型中的参数。

对人脑模型,我们关心白质与灰质区域中参数值的变化。实践中,给定患者的肿瘤生长与侵袭数据取自 CT 扫描与 MRI 等医学图像。这些扫描大致定义了肿瘤可检测部分的轮廓。虽图 11.5 是 CT 扫描的理想化形式,实际扫描远不那么清晰。该研究最初始于一名间变性星形细胞瘤患者终末期年内一系列 CT 扫描的分析;患者经历了将在 11.10 节详细讨论的多种化疗与放疗方案。图 11.12 显示其中一张扫描。

与前节相同,我们使用 Fisher–Kolmogoroff 近似用线速度 v 与增殖率 ρ 近似 D。因此有 \(D_w = v_w^2/(4\rho)\)\(D_g = v_g^2/(4\rho)\),v_w、v_g 分别为实验观察到的线速度。我们现在须为给定患者确定 v_w 与 v_g。关于人肿瘤侵袭最容易获得的信息来自 CT、MRI 或其他成像扫描。我们用发展原始模型所用 CT 扫描(见图 11.12)确定肿瘤前沿在灰质与白质中的速度。在右半球内,可检测肿瘤边缘在大约 180 天内移动约 1.5 cm(Tracqui et al. 1995;Woodward et al. 1996),即平均速度 v=8.0×10⁻³ cm/day。对增长率 ρ=0.012/day,Fisher–Kolmogoroff 近似给出扩散系数 \(D = v^2/(4\rho) = 1.3×10⁻³\) cm²/day。由于侵袭邻近右半球深部核团(灰质),我们把该值与灰质扩散关联,故 v_g ≈ 8.0×10⁻³ cm/day(Tracqui et al. 1995)且 D_g ≈ 1.3×10⁻³ cm²/day。由 CT 扫描,肿瘤边缘穿越胼胝体(白质)的速度比(主要在)灰质中快 2 至 3 倍,故估计 v_w > 2 v_g ≈ 1.6×10⁻² cm/day,D_w > 4 D_g ≈ 4.2×10⁻³ cm²/day。在模拟中我们通常假设灰白质间扩散系数有 5 倍差异 \(D_w = 5 D_g\)

我们也可以通过使用距原点(二维)的均方根距离估计扩散系数

\[ \langle \bar{r}^2 \rangle = 4D \bar{t} \]

由于生长未包含在该定义中,用该公式结合肿瘤轮廓进展的实验数据计算出的 D 将是实际扩散系数的过高估计。仅在大 t̄ 时成立,Fisher–Kolmogoroff 估计可能是胶质瘤扩散系数的较好近似。表 11.5 列出由患者数据产生的高等级胶质瘤参数估计。

模型系统的数值模拟:虚拟肿瘤的侵袭(Numerical Simulations of the Model System: Invasion of Virtual Tumours) 在解剖学准确的区域上求解模型方程是相当复杂的数值任务。在我们的模拟(Swanson 1999;Swanson et al. 2000)中允许 ρ 与 D 有 10 倍变化以模拟不同肿瘤等级:高等级(高 ρ、高 D)、中等级(高 ρ 低 D 或低 ρ 高 D)和低等级(低 ρ、低 D)。Woodward et al. (1996) 发现该范围与各种胶质瘤生长与侵袭非常吻合。数值模拟过程允许我们追踪任意初始大小、分布的虚拟肿瘤从任意位点的侵袭。

图 11.13 显示由 EMMA 获得的人脑冠状切片作为我们模型模拟的区域。图中白质灰质分别呈白与灰。具体考察的三个肿瘤位点在图 11.13 中标出:位置 1 为额顶下肿瘤,位置 2 为额顶上肿瘤,位置 3 为颞叶肿瘤。在每个位置考虑代表 ρ 与 D 的 10 倍变化的 4 个肿瘤等级:高等级(高 ρ、高 D)、中等级(高 ρ 低 D 或低 ρ 高 D)、低等级(低 ρ、低 D)。

肿瘤侵袭作为初始位置的函数(Tumour Invasion as a Function of Initial Tumour Location) 由于脑组织组成的异质性与皮层的几何(这在图 11.13 中非常清楚),侵袭动力学在不同的块状肿瘤位置下可以非常不同。对给定肿瘤等级(即固定 ρ 与 D),肿瘤侵袭范围在每个肿瘤位置下不同(如图 11.14–11.16 所见),取决于灰白质的局部分布以及与解剖边界的邻近程度。

图 11.15 表示高等级(高 ρ、高 D)肿瘤在脑内位置 2(图 11.13 中定义)处的模拟。图 11.15(a) 代表增强 CT 诊断时的可检测肿瘤部分。这里假设诊断发生在可检测肿瘤(即密度高于检测阈值部分的肿瘤)覆盖相当于直径 3 cm 圆的面积时。虽然肿瘤看似局部化,但若把检测能力提高 20 倍(数学上我们可把阈值设为任意正值),对应 500 cells/cm²,在图 11.15(c) 中可见肿瘤已显著侵袭整个右脑叶并穿越胼胝体进入对侧半球。140 天后,图 11.15(b) 表示死亡时增强 CT 可检测的肿瘤部分。此处假设当可检测肿瘤(仍为密度高于检测阈值的部分)覆盖相当于直径 6 cm 圆的面积时发生死亡。把检测能力提高到对应 500 tumour cells/cm² 的分辨率,图 11.15(b) 与 (d) 清楚显示了当前增强 CT 技术无法看到的附加侵袭。注意图 11.15(b) 与 (d) 并不代表肿瘤侵袭的全部范围,而是代表与理论成像技术相关的检测能力的提升。

Swanson (1999) 还把图 11.13 中位置 3 处高等级肿瘤的模拟侵袭与 Burger et al. (1988) 通过临床显微镜尸检整个脑所观察到的侵袭作了比较,发现形状与高等级肿瘤的亚阈值 CT 侵袭(图 11.16(b) 与 (d))相似。

未治疗生存时间作为肿瘤位置的函数(Untreated Survival Time as a Function of Tumour Location) 表 11.6 量化了肿瘤位置对生存时间的影响,针对代表 ρ 增长率与 D_g、D_w 扩散系数 10 倍变化的 4 个假想肿瘤等级。

在表 11.6 中,诊断时间定义为可检测肿瘤覆盖 π(1.5)² cm² 面积时,即相当于半径 1.5 cm 的圆形肿瘤。死亡时间假设为可检测肿瘤覆盖 π(3)² cm² 面积时,即相当于半径 3.0 cm 的圆形肿瘤。可临床定义的生存时间 t_S 是死亡与诊断之间的时间差。注意每个肿瘤等级下白质中的扩散系数仅是灰质中的 5 倍。还需注意位于颞叶的肿瘤对应于诊断后未治疗预期寿命增加。

在均匀脑组织中均匀扩散系数的假设下,11.2 节得到一个近似的生存估计,即

\[ \text{survival time} \approx \frac{1}{\sqrt{D\rho}} (\bar{r}_{\text{lethal}} - \bar{r}_{\text{detect}}) \tag{11.39} \]

若取 \(\bar{r}_{\text{lethal}} - \bar{r}_{\text{detect}} = 1.5\) cm,对高等级肿瘤 ρ=1.2×10⁻²/day 与扩散系数 \(1.3×10⁻³ \le D \le 6.5×10⁻³\) cm²/day(即灰白质之间的迁移性),则得到 170 天至 380 天的估计。对低等级肿瘤 ρ=1.2×10⁻³/day 与 \(1.3×10⁻⁴ \le D \le 6.5×10⁻⁴\) cm²/day,得到 1698 天至 3798 天。这些生存时间的下限合理接近表 11.6 中位置 3 高等级与低等级肿瘤的取值。一般而言,若使用两个(灰质与白质)扩散系数的较小者,这一相当粗略的估计不算太离谱。Burgess et al. (1997) 求解了完整球对称三维模型,对 ρ=1.2×10⁻²/day 与 D=1.3×10⁻³ cm²/day 的高等级肿瘤得到生存期 179 天,对 ρ=1.2×10⁻³/day 与 D=1.3×10⁻⁴ cm²/day 的低等级肿瘤得到 1796 天。把估计与均匀情形比较时吻合得很好。位置 3 有最好的生存期预后,故对其他位置的肿瘤,近似生存期估计是相当大的过高估计。对位于位置 1 的高等级情形,过高估计生存期约 50%;对低等级情形约 85%。

由图 11.13 可见,颞叶白质的几何倾向于将肿瘤细胞导向深部核团(灰质),从而限制肿瘤从位置 3 向外侵袭并提供最好预后;见表 11.6。类似地比较图 11.13 与 11.15,可见位置 2 的肿瘤被上下与内侧所限,且可直接通向白质 corona radiata 胼胝体。位置 1 的肿瘤预后最差,因其邻接 corona radiata 白质且与解剖边界的相互作用很少,仅在外侧受限。

表 11.6 的结果(可预期)表明高等级胶质瘤一般预后最差,低等级肿瘤最好,中等级肿瘤介于其间。然而,高增殖的中等级胶质瘤预后优于高扩散的中等级肿瘤,位置 1 处高扩散的中等级肿瘤预后甚至比高等级胶质瘤更差。这暗示肿瘤的扩散能力是比增长率更强的预后指标。没有一项预后是令人振奋的。

现外推基于肿瘤等级与位置的理论上生存曲线。对每个肿瘤等级,让增长 ρ 与迁移性 D 在表 11.6 的值上下变化 50%。这为 4 个 ρ 与 D 组合中的每一个定义 9 个假想患者。各位置肿瘤的生存时间示于图 11.17。给定肿瘤等级的 3 条曲线之每一条定义了未治疗胶质瘤的 9 个假想患者的生存时间。图 11.17(a) 表示高等级肿瘤(HH)——高增长率 ρ 与高扩散系数 D。图 11.17(d) 显示对立的极端情形(LL)。高–低(HL)与低–高(LH)组合的中等级种群示于图 11.17(b) 与 (c)。注意颞叶肿瘤(位置 3)始终预后最好。虽然曲线看起来很相似,每幅图中的时间尺度却很不同。据我们所知,由于我们的模拟肿瘤是未治疗的,目前尚无可用的临床数据用于比较。(下文讨论考虑各种治疗的模型修正。)然而如上所述,即使在手术切除后也能找到一些与不同位置胶质瘤相关的临床结果。

增强 CT 上可检测肿瘤的比例随扩散系数增加与增长率减小而下降,符合预期。对固定增长率,扩散增加时肿瘤细胞迁移更大距离,因此扩散倾向于平滑胶质瘤细胞的空间分布。对固定扩散系数,增长率减小时肿瘤细胞较少机会累积到可检测水平。有趣的是,若固定增长率与扩散系数的比值而改变 ρ 与 D,则可检测肿瘤比例固定。这类关系对个体患者病例研究可能非常重要,因为无法独立估计扩散与增长率参数。

我们使用对应 40,000 cells/cm² 细胞密度的检测阈值。超过此阈值,肿瘤在增强 CT 上可被观察到;低于此阈值则不能。由于胶质瘤的弥散本质,我们只能用现有成像技术"检测"实际肿瘤的一小部分。这一事实已由 Chicoine and Silbergeld (1995) 与 Silbergeld and Chicoine (1997) 通过他们对远离块状肿瘤位置的组织学正常脑中培养的肿瘤细胞所证实。

数学建模的一个优点是能够在理论上改变任意参数并分析其影响,从而模拟假想实验。通过在数学上降低模型中的检测阈值,我们可以可视化增加医学成像检测能力的用处。图 11.18 展示降低检测阈值对实际检测肿瘤比例的影响。增强检测能力必然增加可识别肿瘤量,但增长速率因肿瘤位置与等级而异。注意对高等级与低等级肿瘤的检测一致。这暗示增长率 ρ 与扩散系数 D 的比值对决定可检测肿瘤的比例至关重要。再次说明肿瘤的侵袭性由 D/ρ 决定。因此该模型对患者的另一潜在实际应用得以简化,因为不需要对 D 与 ρ 单独估计。这里的问题当然是如何获得 D/ρ 的估计。

注意图 11.18 中缓慢增长、快速扩散的中等级肿瘤对应于诊断时可检测肿瘤比例的最小值。这显然是由于这些肿瘤广泛的侵袭能力。高等级与低等级肿瘤对应于相同的可检测肿瘤比例(作为检测阈值的函数)。高等级与低等级肿瘤之间的差异在于生长与侵袭发生的时间尺度。高等级肿瘤在很短的时间内侵袭;低等级肿瘤可能沿高等级肿瘤的侵袭路径行进,但慢得多。

本节展示的所有结果都可以对脑的任何矢状、冠状或水平切片获得。

我们现在必须问的关键问题是模型的三维模拟会如何影响这些结果。数值模拟显然更为复杂。Burgess et al. (1997) 研究了三维常扩散球对称情形,并研究了在有或没有切除情况下增长率与扩散系数的相互作用。他们特别关心量化切除体积大小对增加预期寿命的影响。结合此处更精确的模型,三维模拟也已进行。图 11.19 是在位置 3(图 11.13)处引入未治疗虚拟肿瘤的此类模拟示例,显示死亡时的细胞分布。

模型假设与局限性的回顾(Recap of the Model Assumptions and Limitations) 在讨论各种治疗情境之前,或许应回顾建模中的各种假设。模型假设胶质瘤细胞以常数指数率 ρ 增长。我们未纳入任何坏死核形成(可以通过对细胞使用 logistic 增长动力学来纳入):这在某些胶质瘤中可见,但只有高增殖、低扩散的肿瘤才显著受坏死影响。事实上我们的实验数据取自无坏死的肿瘤(间变性星形细胞瘤)。Tracqui et al. (1995)、Woodward et al. (1996) 与 Burgess et al. (1997) 已证明我们的模型用 logistic 增长所得结果与指数增长差别不大。Tracqui et al. (1995) 表明 logistic 增长减少了肿瘤被识别前的时间,并随坏死细胞密度减少而降低预期寿命。

扩散约在 23 代后开始——此时肿瘤包含约 4×10³ 个细胞;我们认为模型方程 (11.6) 从此开始支配肿瘤生长与侵袭。没有这一限制,局部肿瘤团块似乎无法形成,对应某些大脑胶质瘤病的情形。

我们把扩散系数取为白质中约比灰质中大 5 倍(举例而言)。沿单条白质纤维的扩散可能有一定各向异性,但纤维网络在全脑尺度上非常复杂。存在若干促进横向扩散的通道,如胼胝体,以及在有限程度上的额顶方向。整体扩散系数的改变也许足以模拟胶质瘤在白质中的增强迁移性。鉴于 Burgess et al. (1997) 与患者数据的良好比较,这种可能性有一定根据。

肿瘤在高于临界密度 40,000 cells/cm² 时可被检测——对应于增强 CT 的可检测阈值。该值通过对尸检所得组织学切片与终末期 CT 扫描的患者数据分析确定(Tracqui et al. 1995)。肿瘤平均在可检测面积覆盖 π(1.5)² cm²(相当于半径 1.5 cm 圆形肿瘤的面积)时被诊断。当可检测面积为诊断时的两倍直径(即覆盖 9π cm² 面积)时致死。这些选择基于大量扫描或尸检患者中胶质瘤的平均大小,但范围很宽。

虽然脑边界与组织异质性已纳入考虑,我们并未定义任何正常组织因肿瘤细胞侵袭而被破坏。在实验上目前仍有争议肿瘤细胞是否能破坏远离肿瘤位点的正常脑实质并影响至死亡的时间。

我们的模拟可在二维与三维下进行。三维当然更好,相对于二维 CT 或 MRI 扫描也更易可视化。实践中,我们预期把模型拟合到患者二维 CT 扫描(图 11.5 示意)所得数据上。使用二维计算当然并非不合适。

尽管涉及诸多假设,模型已让我们获得许多临床上潜在有用的结果。尤其是,我们已识别颞叶肿瘤可能比额顶叶肿瘤预后更好。虽然颞叶与额顶叶肿瘤生存时间的差异在临床上可能统计上不显著,但它至少可以解释一些临床变异性。模型还让我们在脑组织中纳入显著的异质性度量,以解释肿瘤几何的不对称性。我们也已能模拟临床上观察到的肿瘤几何并预测患者生存时间。当然这需要进一步的临床分析以观察是否与临床数据一致。鉴于可获得灰白质在脑中局部组成的详细描述,模型仍有大量未探索的潜力,下文将讨论其中一些。

11.6 治疗情境建模:一般性评述(Modelling Treatment Scenarios: General Comments)

如前所述,胶质瘤患者的预后不容乐观,取决于多种因素,其中之一是肿瘤类型与恶性等级。我们提出的相对简单的数学模型考虑胶质瘤的生长与扩散率以及组织异质性。如前述章节所见,结果与临床数据合理一致(Cruywagen et al. 1995;Tracqui et al. 1995;Cook et al. 1995;Woodward et al. 1996;Burgess et al. 1997;Swanson et al. 2000),因此我们相信模型可用于预测各种治疗(包括手术切除、放疗与化疗)后患者的生存时间。

当前治疗策略的一个主要问题是治疗的局部聚焦,但肿瘤生长与侵袭的作用发生在别处(例如见 Gaspar et al. 1992;Liang and Weil 1998)。胶质瘤治疗方案的失败占据胶质瘤临床与实验文献的很大部分(例如 Yount et al. 1998 等及其中文献)。局部治疗当然是理想的,可减少主要引起压迫症状的块状肿瘤。然而它们无法控制造成复发的迁移侵袭细胞(Silbergeld and Chicoine 1997)。

对给定胶质瘤肿瘤的最合适治疗往往完全不清楚,无论是否考虑随后剩余生存的质量。由于它们的侵袭特性,恶性胶质瘤很少能仅通过手术或放射切除治愈。关于完全不同治疗(如肿瘤攻击病毒)有许多新研究。不同程度的切除已被证明对胶质母细胞瘤多形性仅轻微增加生存时间,对较低等级(间变性)星形细胞瘤的增加一般更显著。在单独手术切除通常失败的情况下,为增加患者生存时间,已发展出结合切除、放疗与化疗及其他疗法的多模态治疗。Giese and Westphal (1996) 的综述讨论了抗侵袭治疗的前景。大量临床研究试图展示各种治疗与组合的疗效(例如见 Ramina et al. 1999 等及其中文献)。反直觉的是,由于各种原因,某些治疗组合已被证明比每种治疗单独使用时更无效。最近发现电离辐射能抑制某些胶质母细胞瘤细胞中化疗诱导的细胞死亡(凋亡)(Yount et al. 1998)。(也有人提出抗氧化剂可能被癌细胞利用以帮助防止其被破坏。)这种多模态治疗失败可归因于暴露于强烈化学品或辐射的细胞所诱发的突变。癌细胞本身就是突变过的"正常"细胞,因此随着突变的累积,癌细胞逐渐变得更加恶性与耐药。我们将在后文讨论多克隆模型时讨论化疗治疗时的基本模型修正。

显然存在某给定患者可持续的治疗阈值,确定最优最小策略至关重要。不仅如此,即便肿瘤大小、组织学恶性程度与解剖位置相似,某些治疗对一些患者比对另一些更有效。没有一种通用治疗,但在为患者决定最佳治疗过程时,必须考虑关于该肿瘤的所有可用信息再继续。我们相信现实的建模有助于通过量化某些有或没有治疗的虚拟胶质瘤行为来应对这一复杂问题。

如下文所见,对模型进行相对微小的修正即可考虑化疗、放疗与切除对肿瘤时空行为的影响。这一能力让我们能比较不同治疗情境下肿瘤的预期生长与侵袭,从而对最优治疗过程给出一些洞见。给定某特定肿瘤的位置、大小、形状、扩散系数与增长率的认知,我们的模型可帮助建议最佳治疗类型以最大化生存时间——即诊断与死亡之间的时间差。

11.7 均匀组织中肿瘤切除建模(Modelling Tumour Resection in Homogeneous Tissue)

Cook(J. Cook 博士个人通信 1994)、Woodward et al. (1996) 与 Burgess et al. (1997) 在二维情形下使用我们的空间均匀模型 (11.3) 较早研究了切除对患者生存时间的作用。他们的分析提示手术仅比不手术增加约两个月的生存时间。肿瘤在 CT 扫描上的复发是由于远离大体肿瘤位点的低密度浸润肿瘤细胞。

数学上,手术切除通过去除某确定位置内的所有肿瘤细胞来模拟。假设死亡发生在可见肿瘤达到某确定半径时,或总癌细胞数达到某临界值时。给定初始肿瘤大小、等级与位置,这些模型可模拟各种切除大小与几何对生存时间的作用。Woodward et al. (1996) 取 (11.12) 为初始条件,然后当肿瘤具有相当于半径 1.5 cm 圆的可(CT)检测面积时切除一个区域。Woodward et al. (1996) 以不同于前述方法估计扩散系数;该方法在 Cruywagen et al. (1995) 与 Tracqui et al. (1995) 中详细描述。他们使用穿过肿瘤的三层患者 CT 扫描(图 11.5 中理想化)计算每层的肿瘤面积。在给定增殖参数 ρ 后,他们数值求解模型方程并确定在每层细胞密度处于扫描上观察到的阈值时所对应的扩散系数。这对每层给出扩散系数的三个估计,他们取平均值:有趣的是(事实上相当关键)每层 D 的估计差异不大。我们将在 11.10 节更详细讨论该程序。

以模型参数的估计值,Woodward et al. (1996) 对肿瘤达 CT 可检测大小(直径 3 cm)时进行不同直径的假想切除。以该切除后的状态作为初始条件,他们数值求解模型方程,直至复发肿瘤达到等同于致死大小(即直径 6 cm 的圆)。肿瘤形状当然不再是圆盘,而是内低密度、外远场的环。模拟在包含脑室(Tracqui et al. 1995 与 Cruywagen et al. 1995 在其化疗治疗建模中也这样做)的解剖约束下进行;如前节所示并阐明,这些约束对结果有重要影响。限制区域内细胞的累积实际上减缓了肿瘤的空间扩散。Cook(J. Cook 博士个人通信 1994)发展了一种解析程序以研究该问题,并在前期分布与之后生长特征相同时获得了参数值的极限。我们将在下节讨论他的结果。

Woodward et al. (1996) 考虑了三种切除直径,记为 S,即 S=3、4、5 cm。每种情形下切除区域外的细胞密度当然低于可检测水平。医学上 S=3 cm 模拟诊断与治疗时的"大体全切除"。因此更大的切除要广泛得多。假设死亡发生在视觉可检测面积相当于半径 3 cm 的圆时,对应约 113 cm³ 的体积。

图 11.20(a) 显示由 9 个具所有可能三 ρ 值与三 D 值组合的理想化患者构成的均匀种群的累积生存曲线。各种群经历四种程度的手术:S=0 代表"仅活检"(无切除),S=3 代表"大体全切除",S=4 与 S=5 代表对肿瘤边缘递增的切除;我们不能考虑 S=6,因为那被视为致死。图 11.20(b) 把计算的生存曲线与 Kreth et al. (1993) 报告的 115 名实际患者的曲线相比较——两组胶质母细胞瘤患者,一组接受活检加放疗(58 名患者),另一组接受手术切除加放疗(57 名患者)。所有人在 100 周内死亡。活检组的中位生存时间为 32 周,切除组为 39.5 周。模拟中 S=0 cm 与 S=4 cm 切除之间的中位生存时间增加 7.7 周,与 Kreth et al. (1993) 发现的中位值 7.5 周相符。但我们的模拟要求多 1 cm 的更广泛切除,这可以由我们对阈值细胞数的估计或 Kreth et al. (1993) 切除实际多广来解释。令人惊讶的是,即便使用这样一个相对简单的模型——均匀组织但包含脑边界——结果与患者数据的相关性也这么好。

Woodward et al. (1996) 进一步探究了 ρ 与 D 在他们图 11.20 所用高低值上下变化 10 倍的效果,以模拟前述四个肿瘤等级——从高 ρ 与高 D(HH)到低等级肿瘤(LL)。对四个组合中的每一个,通过让 ρ 与 D 各有 50% 变化形成 9 个假想患者的种群,并计算生存时间。图 11.21 显示结果。图 11.21(a) 为最高等级肿瘤,与图 11.20(a) 相同,仅用于比较;图 11.21(d) 给出最低等级结果。

虽然图 11.21 的曲线看起来相当相似,时间尺度却极为不同。曲线斜率从高 ρ–高 D 等级肿瘤到低 ρ–低 D 等级肿瘤逐步减小。当所有曲线放在同一时间尺度上,它们就可以与 Alvord (1992) 总结的实际患者结果直接比较。其结果展示了胶质瘤行为的范围很宽,模型结果也发现了这一点——这并不令人意外。还应补充,临床上对如何对胶质瘤进行医学分级并无普遍共识。理论或许可以提供此类分级的基础。

假设可检测肿瘤大小的致死半径为 3 cm,图 11.22 显示各种切除半径对后续生存时间的影响。该图提示只有缓慢扩散的胶质瘤才被切除显著影响,正如我们所预期。同样,我们发现扩散是胶质瘤治疗困难与复发的根本原因。

上述详细分析与患者数据比较得出的一个主要结论是:没有明确证据表明切除手术能增加生存——至少对中等级与高等级星形细胞瘤患者如此。这也是 Nazzarro and Neuwalt (1990)(回顾了过去 50 年 33 篇主要报告)的结论。Kreth et al. (1993) 注意到 7.5 周在统计上不显著的增长。Kreth et al. (1993) 指出即使是低等级肿瘤,他们也发现同样难以定义切除的益处。

虽然我们改变了肿瘤细胞的增长率,对给定肿瘤我们仍假设其常数。正如 Woodward et al. (1996) 评论的,一般认为肿瘤边缘包含看起来恶性较低的细胞,但它们仍然是癌性的,无论切除(CT)可见肿瘤外多少组织都会生长。

11.8 切除后肿瘤复发的解析解(Analytical Solution for Tumour Recurrence After Resection)

Cook(J. Cook 博士个人通信 1994)以解析方式考虑了均匀组织中的切除效应。如前所述,数学上手术切除通过去除以肿瘤为中心的某确定位置内的所有肿瘤细胞来模拟。这里我们展示他的分析。如前所述,假设癌症在可见大小达到等价于半径 3 cm 的圆时或总癌细胞数达到某临界值时致死。我们预期切除的主要效应可被均匀情形合理量化,但与肿瘤位置和局部组织异质性相关的细节可能调整这些结果;下节讨论这种情形。

无量纲模型为 11.2 节推导并由 (11.10) 与 (11.11) 定义,以边界条件 (11.8) 为约束。无量纲化由 (11.9) 给出。在扩散系数无空间异质性时(即 D(x) ≡ 1),模型方程为

\[ \frac{\partial c}{\partial t} = \nabla^2 c + c \tag{11.40} \]

切除前后肿瘤细胞群满足该方程,但初始条件不同。

为简化起见,并鉴于上述二维研究与患者数据相比给出的令人惊讶的好结果,我们考虑极坐标下的切除后二维问题。假设肿瘤最初由原点处 N 个细胞的点源引发,故切除前问题对 \(0 < t < t_r\) 满足 (11.38),初始条件 \(c_{\text{preresect}}(r,\theta,0) = N\delta(r)\)。在无穷域上解由 (11.13) 给出,其在 t_r 时刻的无量纲形式为

\[ c_{\text{preresect}}(r, \theta, t_r) = \frac{N}{4\pi t_r} \exp\left(t_r - \frac{r^2}{4t_r}\right) \tag{11.41} \]

切除时去除半径 R_r 的中心核,故切除后问题对 t>t_r 满足 (11.38),初始条件(t=t_r)为

\[ c_{\text{postresect}}(r, \theta, t_r) = F(r, \theta) = N H(r - R_r) c_{\text{preresect}}(r, \theta, t_r) = H(r - R_r) \frac{N}{4\pi t_r} \exp\left(t_r - \frac{r^2}{4t_r}\right) \tag{11.42} \]

其中 H 为 Heaviside 函数。

由叠加原理,带这些初始条件的切除后问题的解可表示为积分

\[ c_{\text{postresect}}(r, \theta, t) = \int_0^{2\pi} \int_0^\infty K(r, \theta, t; \xi, \alpha, t_r) F(\xi, \alpha) \xi d\xi d\alpha \quad \text{for } t > t_r \]

其中 K 为在 t=t_r 引入 \((r,\theta)=(\xi,\alpha)\) 处点源的 (11.38) 基本解:

\[ K(r, \theta, t; \xi, \alpha, t_r) = \frac{1}{4\pi(t - t_r)} \exp\left((t - t_r) - \frac{r^2 + \xi^2 - 2r\xi\cos(\theta - \alpha)}{4(t - t_r)}\right) \tag{11.43} \]

\[ c_{\text{postresect}}(r, \theta, t) = N \exp\left(t - \frac{r^2}{4(t - t_r)}\right) \frac{1}{(4\pi)^2 t_r (t - t_r)} \int_{R_r}^\infty \xi \exp\left(-\frac{\xi^2}{4t_r} - \frac{\xi^2}{4(t - t_r)}\right) \left[\int_0^{2\pi} \exp\left(\frac{2r\xi}{4(t - t_r)} \cos(\theta - \alpha)\right) d\alpha\right] d\xi \]

考虑方括号中的第二个积分,记为 I,引入变量替换 ν=θ−α,则

\[ I = \int_0^{2\pi} e^{A\cos(\theta - \alpha)} d\alpha = \int_{\theta-2\pi}^\theta e^{A\cos\nu} d\nu \tag{11.44} \]

其中 \(A = 2r\xi/(4(t - t_r))\)。由于我们最终要把 I 从 ξ=R_r 积分到 ∞,且 A 与 ξ 成正比,A 大是我们关心的极限。此外从数值模拟可知,对我们考虑的增长率 ρ 范围,无量纲时间 t=ρ t̄ 至多为 O(1)。这两个观察支持使用 A 大的极限。对 A 大,我们可以用 Laplace 方法(见例如 Murray 1984)近似 I;该方法利用积分由指数最大值附近的贡献支配:cos v=1 的极大值在 v=0 或 2π。在 v=0 处重标度,对 A 大有

\[ I = 2\int_0^\epsilon e^{A(1 - v^2/2 + \cdots)} dv \sim e^A \int_{-\infty}^\infty e^{-Av^2/2} dv \sim e^A \sqrt{\frac{2\pi}{A}} \]

由此

\[ c_{\text{postresect}}(r, \theta, t) \sim N \sqrt{\frac{2\pi}{e}} \exp\left(t - \frac{r^2}{4(t - t_r)}\right) \frac{1}{16\pi^2 t_r (t - t_r)} \int_{R_r}^\infty \exp\left(-\frac{\xi^2}{4t_r} - \frac{\xi^2}{4(t - t_r)} + \frac{r\xi}{2(t - t_r)}\right) \sqrt{\frac{r\xi}{2(t - t_r)}} \xi d\xi \sim N \frac{e^t}{8t_r} \frac{r\sqrt{\pi}^3 (t - t_r)}{} \int_{R_r}^\infty \xi \exp\left(-\frac{\xi^2}{4t_r} - \frac{(r - \xi)^2}{4(t - t_r)}\right) d\xi \]

随 A→∞。设

\[ J = \int_{R_r}^\infty \xi \exp\left(-\frac{\xi^2}{4t_r} - \frac{(r - \xi)^2}{4(t - t_r)}\right) d\xi = \int_{R_r}^\infty g(\xi) e^{-x h(\xi)} d\xi \tag{11.45} \]

其中 \(g(\xi) = \sqrt{\xi}\)\(h(\xi) = \xi^2/r + r t_r/t - 2\xi t_r/t\)\(x = r t/(4 t_r (t - t_r))\),我们可以对 x 大再次用 Laplace 方法近似。这等价于假设 \(t \approx t_r\)\(r \sim O(1)\)\(t \gg t_r\)\(r \gg 4t_r\)。渐近贡献来自 h(ξ) 的极小值附近,其在 \(\xi = r t_r/t\) 处。有两种可能。

(i) 极小值位于积分范围内:\((r t_r/t) > R_r\)。引入新变量 \(w = \xi - r t_r/t\) 并对 x 大展开,给出

\[ J \sim \sqrt{\frac{4\pi r(t - t_r) t_r}{t}} e^{-r^2/(4t_r)} + O\left(\sqrt{\frac{4 t_r (t - t_r)}{rt}}\right)^3 e^{-r^2/(4t_r)} \tag{11.46} \]

随 x→∞。

(ii) 极小值不在积分范围内:\(r t_r/t < R_r\)。引入新变量 \(w = \xi - R_r\) 并对 x 大展开,给出

\[ J \sim \frac{2 t_r (t - t_r) \sqrt{R_r}}{R_r t - r t_r} \exp\left(-\frac{R_r^2}{4t_r} - \frac{(r - R_r)^2}{4(t - t_r)}\right) + O\left(\exp\left(-\frac{R_r(R_r t - t_r r)}{4 t_r (t - t_r)}\right) \frac{4 t_r (t - t_r)}{rt}\right)^2 \tag{11.47} \]

随 x→∞。

现汇总这些结果得到切除后的解

\[ c_{\text{postresect}}(r, \theta, t) \sim \begin{cases} \dfrac{N}{4(R_r t - r t_r)} \sqrt{\dfrac{(t - t_r) R_r}{r \pi^3}} e^{t - R_r^2/(4t_r) - (r - R_r)^2/(4(t - t_r))} & \text{for } r < (t/t_r) R_r \\ \dfrac{N}{4\pi t} e^{t - r^2/(4t_r)} & \text{for } r > (t/t_r) R_r \end{cases} \tag{11.48} \]

仅在 \((t/2 t_r)(r R_r/(2(t - t_r)))\)\((r R_r/(2(t - t_r)))\) 大时严格成立。从该解可以清楚看到,无论切除区域如何,解最终都将(指数地)增长到临界大小。

我们主要关心切除与不切除下肿瘤细胞群生长与侵袭的差异。因此考察超出块状肿瘤边缘侵袭组织的肿瘤细胞前沿。由渐近解 (11.43) 可推出切除后侵袭前沿如何被减缓。对 \(r > (t/t_r) R_r\),切除后的渐近解可写为

\[ c_{\text{postresect}}(r, \theta, t) \sim \frac{N}{4\pi t} \exp\left(t - \frac{r^2}{4t_r}\right) = c_{\text{noresect}}(r, \theta, t) \exp\left(-\frac{r^2}{4t_r} + \frac{r^2}{4t}\right) \tag{11.49} \]

其中(对照 (11.41))

\[ c_{\text{noresect}}(r, \theta, t) = \frac{N}{4\pi t} \exp\left(t - \frac{r^2}{4t}\right) \tag{11.50} \]

为若未进行切除时的解析解。因此切除后切除区域之外的肿瘤细胞浓度被因子 \(\exp(-r^2/(4t_r) + r^2/(4t))\) 所抑制;这只在上述渐近参数充分大时成立。假设细胞密度在阈值 c∗ 处(由 CT)可检测,则可检测肿瘤的外半径满足

\[ r_{\text{postresect}}^2 \sim 4 t t_r \left(1 - \frac{1}{t} \ln \frac{4\pi t c_*}{N}\right) < 4 t^2 \left(1 - \frac{1}{t} \ln \frac{4\pi t c_*}{N}\right) = r_{\text{noresect}}^2 \tag{11.51} \]

\(t_r < t < (r/R_r) t_r\)。图 11.23 给出由上式渐近形式得到的可检测肿瘤半径(有无切除)。

11.9 含脑组织异质性的手术切除建模(Modelling Surgical Resection with Brain Tissue Heterogeneity)

Swanson (1999) 在解剖学准确的脑上重新考虑了切除(如上节所述)。原理与上节相同,但数学问题由于白质灰质分布带来的扩散系数空间异质性而显著更复杂。与之前一样,切除后肿瘤细胞群继续由 11.2 节发展的模型系统支配。含切除的数学问题由无量纲系统

\[ \frac{\partial c}{\partial t} = \nabla \cdot (D(x) \nabla c) + c \tag{11.52} \]

给出,对 t>t_r,其中

\[ D(x) = \begin{cases} 1 & \text{for } x \in \text{White Matter} \\ \gamma = D_g / D_w & \text{for } x \in \text{Grey Matter} \end{cases} \tag{11.53} \]

\(c(x, t_r) = F(x)\)(切除后细胞的初始分布)与 \(n \cdot D(x) \nabla c = 0\)(x 在 ∂B 即脑的边界上)为约束。

我们关心切除对肿瘤生长与侵袭的作用。使用上文对未治疗肿瘤生长所描述的程序,我们可以考察切除对脑内不同位置肿瘤的作用。解必然是数值的。作为详细示例,我们专门关注位置 1(额顶下)肿瘤且仅考虑高等级肿瘤——即高增长率与高扩散;Swanson (1999) 给出了所有等级肿瘤的进一步结果。增长率 ρ 并不出现在无量纲模型中,但当然出现在无量纲化中(见 (11.9))。

图 11.24 是位置 1(见图 11.13)处高等级胶质瘤在大体全切除(即去除肿瘤的可(CT)检测部分)后的模拟响应。图 11.25 显示同一肿瘤在广泛切除(即去除相当于可检测部分两倍的区域)后的演化。回想大体全切除对应去除相当于半径 1.5 cm 圆的肿瘤面积;也就是说总切除面积(在 CT 扫描上可见)为 π(1.5)² cm²。如上所述,"大体全切除"是神经外科医生常用的术语,描述手术刚结束时增强 CT 扫描上无可检测肿瘤的肿瘤成功切除。最大广泛切除对应去除相当于半径 3 cm 圆的肿瘤面积;也就是说总切除面积(在 CT 扫描上可见)为 π·3² cm²。Woodward et al. (1996) 工作中所考虑的最大切除相当于半径 2.5 cm 的圆。这一情形定义了最近的趋势——不仅切除可检测肿瘤,还切除环绕大体可见肿瘤的某些半径的正常外观组织。对最终结果影响不大。

图 11.24(a) 表示诊断时(即 t=t_r)肿瘤的可(CT 扫描)检测部分。这之后立即进行大体全切除,图 11.24(b) 显示使用 11.7 节阈值的肿瘤实际范围。可见切除仅移除了肿瘤总量的很小一部分——事实上 36.9%(见表 11.7)。60 天后,可检测肿瘤部分已从面积 0 增至 π(1.5)² cm²;换言之肿瘤已复发。复发肿瘤表示为环绕切除床的可检测肿瘤环——这是人类胶质瘤的常见特征。可检测复发性肿瘤与实际肿瘤在图 11.24(c) 与 (d) 中分别显示。假设不考虑其他治疗,肿瘤在切除后 130 天致死。

若进行广泛切除(见图 11.25(b)),肿瘤总量的 86.7% 被去除(见表 11.7)。肿瘤在 225 天后复发。复发肿瘤不在肿瘤床周围形成环,而是在切除床边缘附近形成孤立的肿瘤岛(图 11.25(c))。显然广泛切除后的复发显著受到脑边界的影响。三个复发性肿瘤岛位于边界处或附近。由于切除后肿瘤细胞浓度相当低,单一块状肿瘤并不形成。切除后残余的肿瘤细胞继续扩散,并在遇到边界时累积。图 11.26 显示临床观察到的多灶复发一例。模拟肿瘤在切除后 253 天致死,假设不尝试其他治疗方案——即增加生存时间 123 天或约 100%。广泛切除后的这种复发常被视为与所切除肿瘤分离的"新"肿瘤。模型提示事实上这种再生长是同一肿瘤的一部分——即已成为可检测的弥散侵袭部分。

异质介质中肿瘤切除后的生存时间(Survival Time Following Tumour Resection in Heterogeneous Medium) 虽然显然弥散侵袭的肿瘤细胞直接贡献于患者的生存时间,但考虑我们全文使用的理想化致死肿瘤大小的模型含义是有教益的。我们对图 11.13 中其他标记位置(即位置 2 与 3)处的肿瘤执行了与上述类似的模拟;Swanson (1999) 给出了这些模拟结果。图 11.27 显示三位置处高等级与两个中等级肿瘤的生存时间。未给出低等级胶质瘤,因为它们的行为本质上与高等级胶质瘤相同但时间尺度慢得多。表 11.7 与 11.8 总结生存时间与切除的结果。我们以表格形式给出大体全切除或广泛切除后去除的肿瘤总体积比例。在我们的模型假设下,大体全切除可少至去除 12% 的总肿瘤——因诊断时它已如此弥散。广泛切除可去除超过一半的肿瘤。切除的有效性当然取决于扩散系数 D 与增长率 ρ 的比值。亦如上所述,肿瘤位置相对于边界的情况显著影响模型预测的生存时间。

切除量与生存时间的关系显然取决于许多因素。特别重要的是了解肿瘤位置的物理环境——如与边界的邻近程度以及白质灰质分布。有了这些信息,我们可以用这种方法更准确地推断切除治疗的有效性。还有其他因素,下文讨论其中一个主要因素——化疗。

11.10 化疗对肿瘤生长作用的建模(Modelling the Effect of Chemotherapy on Tumour Growth)

本节描述 Cruywagen et al. (1995) 与 Tracqui et al. (1995) 发展的用以量化化疗对肿瘤生长作用的模型。它部分基于人脑胶质瘤组织学切片的定量图像分析,特别是化疗期间连续 CT 图像的横截面积/体积测量。我们用优化技术估计模型参数以使模拟肿瘤面积与 CT 扫描数据最佳拟合。我们在二维区域上进行了数值模拟,考虑了脑(仅脑室与颅骨)的几何及其对扩散的自然屏障。结果被用于确定化疗对肿瘤时空生长的作用。(Shochat et al. (1999) 用基于计算机的基本常微分方程模型模拟评估了乳腺癌化疗方案的有效性。)

讨论此工作的原因之一是它是一种稍有不同的方法,且展示了如何从一名间变性星形细胞瘤患者(接受放疗与化疗)的连续 CT 图像估计模型参数。当其他方法或独立数据不可用时,这是一种可行的方法。如前所述,治疗的时间进程被纳入模型中,通过用优化技术拟合模型响应到实验数据,并通过识别刻画肿瘤细胞死亡率的参数来获得化疗有效性的评估。

与上述讨论一样,一旦从组织学切片重建侵袭动力学事件的可行性已建立,希望是能用这样的建模研究其他不同特征生长模式、几何形状的胶质瘤以及使用同一患者的其他数据类型的各种疗法的作用。

实验数据(Experimental Data) 由于许多脑肿瘤患者在检测后很快死亡而没有足够的后续扫描,我们幸运地能够研究一名间变性星形细胞瘤患者——该患者 3 年前被诊断并接受 X 射线治疗。肿瘤复发并以适度成功再次接受化疗,同时在患者死亡前 12 个月内反复进行 CT 扫描。从这些扫描中观察到一块或多块肿瘤。最大块(块 a)的肿瘤面积在脑的三个不同层面(称为层面 1、2、3;回顾图 11.3)上使用数字化平板与计算机化面积测量技术测量。由于记录面积的层面在所有扫描中并非总是完全相同,某些数据点通过相邻上下层面之间的线性插值获得。

在同一终末期年内,患者接受两种不同的化疗(图 11.28)。第一种是 6 种药物(6-thioguanine、procarbazine、dibromodulcitol、CCNN、5-氟尿嘧啶与 hydroxyurea)组成的方案,给予 15 天并每 6 至 8 周重复一次以允许骨髓恢复,共实施 5 次。第二种包括以月间隔给予的两程顺铂。此外,患者在最后三周接受中子束照射。

肿瘤细胞密度由首次扫描时活检的图像分析确定。肿瘤组织中每 mm² 检测到的细胞核数在 5 张图像上取平均值。该细胞密度被用作模型的起始条件。对尸检材料进行另一次肿瘤细胞密度评估。在 5 张图像上对正常组织中检测到的每 mm² 细胞核数取平均值。该细胞密度被用作模型中确定某区域为肿瘤性的阈值。

数学模型(Mathematical Model) 此处我们把组织视作均匀的——即扩散系数与增长率参数取常数。不含化疗作用的基础数学模型由 (11.3) 给出,以脑边界上的零通量边界条件为约束。

现在须量化化疗的作用。我们假设由化疗引起的细胞死亡可以由线性去除率 \(K(t) c(x,t)\) 建模,其中函数 K(t) 描述治疗的时间轮廓。我们进一步假设对所考虑的时间尺度,化疗作用可作为一阶近似用具有常数振幅的阶跃函数建模,使得

\[ K(t) = \begin{cases} k_1 & \text{if } t_{1,i} \le t \le t_{1,i+1}, i = 1, 3, \ldots, 9 \\ k_2 & \text{if } t_{2,j} \le t \le t_{2,j+1}, j = 1, 3 \\ 0 & \text{otherwise} \end{cases} \tag{11.54} \]

其中 \(t_{1,i}\)\(t_{1,i+1}\) 分别为第一次六药化疗每一治疗的开始与结束时间,而 \(t_{2,j}\)\(t_{2,j+1}\) 对应于第二次顺铂化疗的开始与结束时间,如图 11.29 所示。结合此与方程 (11.3),基础(单细胞类型)方程为

\[ \frac{\partial c}{\partial t} = D \nabla^2 c + \rho c - K(t) c \tag{11.55} \]

现在需要初始与边界条件。理论上肿瘤最初源自一个癌细胞,但该原始细胞出现的时间、生长类型与早期癌细胞扩散均未知。由于须给出近似初始条件,我们如上假设首次扫描时扩散过程已打破了任何先前的可能均匀分布。因此假设细胞以正态分布,在所考虑层面的肿瘤中心 \(x_0\) 处有最大细胞密度 a,即

\[ c(x, 0) = a \exp\left(-\frac{|x - x_0|^2}{b}\right) \tag{11.56} \]

其中 b 是肿瘤细胞扩散范围的度量。

作为边界条件,要求细胞对脑外或对脑室内的通量为零,故在区域边界上仍为

\[ n \cdot \nabla c(x, t) = 0 \tag{11.57} \]

其中 n 为边界法向单位向量。

最初以单细胞群假设加上上述模型重现肿瘤时空变化的尝试未获成功,因为它无法解释临床数据中看到的定性行为变化。特别是它无法捕捉在图 11.28 中看到的第一次化疗末肿瘤面积的急剧增加——这表明该治疗失效。这提示需要一个含细胞突变的多克隆模型——一个合理的第一修正。

因此必须修正模型以包含两个细胞群。第一类细胞记为 \(c_1(x,t)\),对第一次化疗(UW 六药方案)与第二次(顺铂)均敏感,构成首次扫描时肿瘤的大部分(实际计算超过 90%)。第二类细胞记为 \(c_2(x,t)\),假设对第一次化疗耐药,但可能对第二次治疗敏感,构成其余癌细胞群。为建模化疗的作用,我们考虑细胞死亡与化疗期间对特定治疗敏感的现存细胞数成正比。第二类肿瘤细胞可能源自早期基因突变——或是自发的或由 3 年前的放疗诱导的,或是化疗开始时诱导的类似突变。在此假设下,第二细胞群的增殖由该细胞类型的初始密度决定。另一种模型可考虑包含 Goldie and Coldman (1979) 提出的每代突变概率的假设。我们进一步假设两个细胞群关于扩散具有相同性质,但可能具有不同增长率。

双细胞群模型(Two-Cell Population Model) 此处模型中我们不把第二细胞群视为突变,而视为独立的细胞系。在下一节关于多克隆性的内容中,我们描述这样的模型。这里我们考虑两个细胞群 \(c_1(x,t)\)\(c_2(x,t)\),具相同扩散系数 D 但可能具有不同增长率 \(r_1\)\(r_2\)。以同一模型 (11.55) 为基础,我们对两细胞群取的模型为

\[ \frac{\partial c_1}{\partial t} = D \nabla^2 c_1 + r_1 c_1 - K_1(t) c_1 - K_2(t) c_1, \quad \frac{\partial c_2}{\partial t} = D \nabla^2 c_2 + r_2 c_2 - K_2(t) c_2 \tag{11.58} \]

这里假设第一类细胞 \(c_1\) 对两种治疗都敏感,第二类 \(c_2\) 对第一次化疗耐药但对第二次治疗敏感,即

\[ K_1(t) = \begin{cases} k_1 & \text{if } t_{1,i} \le t \le t_{1,i+1}, i = 1, 3, 5, 7, 9 \\ 0 & \text{otherwise} \end{cases}, \quad K_2(t) = \begin{cases} k_2 & \text{if } t_{2,1} + 4 \le t \le t_{2,2} + 6 \\ k_2 & \text{if } t_{2,3} \le t \le t_{2,4} + 2 \\ 0 & \text{otherwise} \end{cases} \tag{11.59} \]

时间 \(t_{1,i}\)\(i=1,\ldots,10\))与 \(t_{2,j}\)\(j=1,\ldots,4\))与图 11.29 中所给时间直接相关。这里(在数值模拟中)我们考虑第一次化疗从治疗开始到结束起作用,但第二次化疗在治疗日期后 0 到 4 天才开始起作用,并在最后一次给药后 2 天仍有效。确切时序未知,但这些假设对临床数据给出了更好的拟合。

相应边界条件(对照 (11.57))为

\[ n \cdot \nabla c_1 = 0, \quad n \cdot \nabla c_2 = 0 \tag{11.60} \]

初始空间条件(对照 (11.56))为

\[ c_1(x, 0) = a_1 \exp\left(-\frac{|x - x_0|^2}{b}\right), \quad c_2(x, 0) = a_2 \exp\left(-\frac{|x - x_0|^2}{b}\right) \tag{11.61} \]

参数 \(a = a_1 + a_2\) 为肿瘤中心 \(x_0\) 处细胞的最大初始密度,参数 b 衡量肿瘤细胞扩散程度。a 与 b 的初始值由首次扫描时活检测定的细胞密度与首次扫描面积估计,以便首次扫描时肿瘤面积的模拟值与实验测定值接近。参数 \(a_1\)(第一细胞群的最大初始密度)仍须估计。假设从边界开始的初始三维均匀扩散,对应于首次扫描测量的不同初始面积,为不同层面选取不同的初始细胞浓度。因此层面 2 中心处肿瘤的初始细胞密度 a 大于层面 1 与 3;相应地刻画细胞初始扩散的参数 b 对层面 2 略大于其他两层。表 11.9 展示这些参数估计在各层面的比较。

数值方法与结果(Numerical Methods and Results) 两个细胞群 \(c_1(x,t)\)\(c_2(x,t)\) 由模型系统 (11.58)–(11.61) 的数值积分得到。一个直接问题是为参数赋值。在从 CT 扫描评估肿瘤面积时,模拟面积的数值确定涉及一个细胞密度阈值 \(c_{th}\)。数学上我们可以任意选择。但既然我们在与数据结合的情况下使用该模型,必须用不同的方法。根据尸检所得组织学切片与终末期 CT 扫描比较的实验数据,归一化阈值的取值取为活检密度的 40%。这是一个合理的估计,但显然最好能从脑切片扫描独立估计。然而我们可以定量说明阈值变化对模型参数的影响。例如,检测阈值下降 25% 之多对其他参数的影响小于 10%,表现为 D 与 \(k_1\) 略增,增长率 \(r_1\)\(r_2\) 略减;下节回到此点。数值模拟细节见 Tracqui et al. (1995) 与 Cruywagen et al. (1995)。

模型参数的初始值(与观察数据良好吻合的)由启发式方法确定;从前述章节的文献(如 Alvord and Shaw 1991 关于细胞倍增时间的实验数据与 Chicoine and Silbergeld 1995 关于体外迁移性的数据)中我们对这些值有一些概念,但此处它们最初是由这些实验数据估计的。表 11.9 给出的模型参数是优化程序所得。随后通过全局程序通过对 6 个未知模型参数的优化来改善模型解对数据的拟合。

该模型下未知参数为扩散系数 D、两细胞群的两个增长率 \(r_1\)\(r_2\) 与刻画两细胞群在两种不同化疗下死亡率的 \(k_1\)\(k_2\)。起始时间(即首次扫描时)第一类细胞与第二类细胞的百分比也是未知参数。

癌细胞总密度 \(c(x,t) = c_1(x,t) + c_2(x,t)\) 由数值解得到。每层相关的颅骨与脑室的解剖边界被纳入考虑。

图 11.29 显示层面 1、2、3 处肿瘤面积的模拟时间演化。在模拟中,肿瘤面积在第一次化疗开始前增大。如预期,第一次化疗作用 \(k_1\) 的值越高,肿瘤面积减小越多;而增大扩散系数 D 或增长率 \(r_1\) 的值有相反效果。为在第一次化疗末之前产生肿瘤面积的增长,\(k_1\) 的作用必须由第二细胞群增长率 \(r_2\) 平衡。第二次化疗的作用由 \(k_2\) 度量。这些参数的值以及肿瘤的初始组成因此对确定肿瘤的时空变化至关重要。在模型中纳入物理脑边界与脑室位置对结果有相当影响:边界不仅由于限制区域内细胞的累积而对肿瘤细胞密度的提升有贡献,显然也存在肿瘤面积的最大值——即脑除脑室外的面积。

由模型模拟的肿瘤面积演化与实验数据良好吻合。对三个层面都得到良好拟合——考虑到实验曲线后半段(即第二次化疗期间)观察到的数据点大幅变化(图 11.28),这对模型有效性是令人鼓舞的。同样令人鼓舞的是从表 11.9 看到对所考虑的三个层面模型参数获得相似值——如果假设不同层面间细胞行为与性质无特定差异,这是预期的。化疗有效性度量参数 \(k_1\)\(k_2\) 显示出较大的层面间变化。然而图 11.30 中由连续治疗导致的锯齿状曲线形状对模拟化疗每阶段的阶跃函数的开始与结束时间的任何小变化都很敏感。因此在缺乏关于治疗动力学的更详细信息(这些信息将允许对化疗作用随时间进行更精确建模)的情况下,可以预期这样的变化。

Tracqui et al. (1995) 用检测阈值降低 25% 与基本相同的优化准则值重新运行了优化,以观察其对参数估计的影响。他们发现这导致对其他估计参数的变化小于 10%,D 与 \(k_1\) 略增,增长率 \(r_1\)\(r_2\) 略减。他们还进行了模拟,以关于总癌细胞以及两亚群的肿瘤时空变化。他们发现,正如从模型基础所预期的那样,尽管第一类细胞减少,第二类细胞却增加——显然患者的死亡是由于第二类细胞亚群的出现。

图 11.30 显示模拟肿瘤的二维图像,并与第 113 与 300 天扫描所得图像作比较。均匀介质的假设导致模拟肿瘤形状比观察到的更规则。扫描层面的变化伴随着脑室位置的修改,使得不可能精确追踪肿瘤形状的演化——这一问题在 Swanson (1999) 的方法中不会遇到。

如上所述,模型参数由临床数据估计。基本参数——细胞增殖与扩散参数——是稳健的,对它们的小幅改变给出在量上相似的结果。另一方面,度量化疗治疗强度与持续时间的参数对小变化更敏感,如上所述。比各层面参数之间的一致性更重要的是模拟值与观察值之间的一致性——后者由其他实验以及表 11.5 与 11.6 给出的各种其他生物学数据确定。这进一步支持了用于估计参数的方法。

初始细胞密度(第 1 天)对应于首次扫描后活检测得的平均值 1026 cells/mm²。模型在模拟时间末(第 356 天)给出的细胞密度平均值为 710 cells/mm²,与尸检测量的平均值 750 cells/mm² 非常接近。

由模型估计的第二突变细胞群的增长率与第一细胞群没有显著差异。相应的倍增时间为第一细胞群 66 天与第二细胞群 60 天(表 11.9),与该类型肿瘤的报道值一致(Alvord and Shaw 1991)。

表 11.9 中的扩散系数值(以 cm²/day 给出)由扫描与真实脑大小之间 \((5/2.2)^2\) 的比例因子确定。模型识别的平均值为 \(1.25×10⁻⁷\) cm²/sec。该值可与体外癌细胞迁移率的估计(Chicoine and Silbergeld 1995)比较。使用实验速度值 v=15 µm/hr 与模型识别的平均增长率 ρ=1.2×10⁻²/day,得到 D 的估计 \(v^2/(4\rho) = 3.13×10⁻⁷\) cm²/sec;亦见表 11.5。

显然其他病例必须被分析以确定各种参数的范围。要分析其他病例,至少须有跨越一段时间(最好至少一年)的几次 MRI 或 CT 扫描,以及来自活检或尸检的肿瘤与边缘的组织学样本。由于星形细胞瘤的估计倍增时间从一周到一年或以上不等,一般与组织学恶性程度相关,我们预期在其他病例中可以用模型定义 \(r_1\)\(r_2\) 的对应差异。当然星形细胞瘤对放疗和/或化疗的响应程度难以预料,因此我们预期参数 \(k_1\)\(k_2\) 也有对应差异。

或许最关键的问题之一是扩散系数 D 在不同胶质瘤间有多大变化。另一问题是细胞浓度或细胞类型是否各自独立贡献于 D。

当前胶质瘤的治疗受限于正常脑对 X 射线放疗以及造血组织对化疗的耐受,但合理预期随着额外治疗的发展,模型中参数的定义将有助于规划这些治疗的剂量与时机。

预测肿瘤对治疗的响应是数学模型可帮助达到的目标。这一益处对化疗应特别明显——化疗是一种因对药物药代动力学以及大细胞群中耐药亚群存在的知识严重缺乏而受困的治疗(两个因素都因患者而异)。对给定类型脑肿瘤的模型参数典型值的了解将提供对常规治疗下肿瘤生长响应的第一估计。基于每位被治疗患者的每次新 CT 或 MRI 扫描所得特定数据对参数的空间在线适应将提供接近最佳药物递送模式所需的必要修正。

我们应该再次就模型的使用提出警示。在能达到最优临床策略之前,须对模型预测进行大量谨慎的探索与验证。但这种方法允许癌症治疗的改进——这不意味着或依赖于完全新药的发现,而是建议对已有药物的最优使用与原位评估。

当然模型存在一些明显局限。肿瘤异质性可通过考虑具有特定动力学与可能竞争性相互作用的多个种群纳入,这会影响肿瘤生长动力学(Michelson et al. 1987)。不同治疗时间安排的作用可通过移动建模化疗的阶跃函数的时间边界值容易地考察,关于药物药代动力学的额外知识可进一步引入以调制这些函数的形状。不同程度手术切除的作用已由 Woodward et al. (1996) 与最近的 Swanson (1999) 与 Swanson et al. (2000) 用这类模型探究,后者关键地纳入了沿白质纤维束方向的增强扩散。虽然我们用扩散建模迁移性,肿瘤细胞不太可能被动地扩散穿过正常脑组织——脑组织太过黏稠而不允许太多被动扩散。相反,如实验所示(Chicoine and Silbergeld 1995),肿瘤细胞可能通过主动变形虫式过程迁移。本章所发展的模型中数学框架足够一般化以形式化这些特性——举例来说,本书中若干种群动力学中的不同模型即如此。Tracqui (1995) 通过对比本节讨论的方法与关于细胞迁移与细胞牵引力的最小假说讨论了这一方面。Tracqui et al. (1995) 给出了其他局限与改进模型的建议。

在参数估计的这种方法下,还须额外工作以确定把肿瘤建模为三维对象在何种程度上显著地修改模型参数的估计,特别是扩散系数。Swanson (1999) 与 Swanson et al. (2000) 的工作对此有重要贡献。她也已用她的方法(Swanson 1999)应用于化疗治疗问题。然而更准确描述现象的益处将需要脑与脑室边界的三维分析与灰白质的异质分布。

11.11 肿瘤多克隆性与细胞突变建模(Modelling Tumour Polyclonality and Cell Mutation)

除上一节外,我们假设肿瘤细胞群是均匀的,只涉及一类细胞。尽管与数据(如关于切除后预期寿命)有良好的定量比较,上节中胶质瘤对化疗药物可变敏感性所展示的现象需要一个含肿瘤细胞异质性的模型。已知胶质瘤是异质性的(多克隆的),且异质性一般随等级增加。恶性程度更高的细胞被认为具有更强的突变倾向,从而增加异质性。因此预期在肿瘤内看到不同类型的细胞(例如 Pilkington 1992)。基础模型可扩展为考虑多克隆肿瘤的情形,如上一节所见,只需在肿瘤内创建两(或更多)个可能具不同扩散性与增长率的细胞群。在模型方程 (11.3)(均匀组织情形)中,细胞密度 c̄ 成为细胞密度向量 c̄,扩散系数 D 与增长率 ρ 现在成为扩散系数与增长率的对角矩阵。在某给定的时空点 \((x̄, t̄)\),总肿瘤种群为 c̄ 向量各分量之和。在上一节化疗模型中我们取两个细胞群相互独立。然而亚群常非独立,而是通过把癌细胞从细胞群 i 转移到细胞群 j 的突变事件相连。一个相当一般但基本的考虑种群多克隆性的模型可写为维数形式

\[ \frac{\partial \bar{c}}{\partial \bar{t}} = \bar{\nabla} \cdot (\bar{D} \bar{\nabla} \bar{c}) + \rho \bar{c} + T \bar{c} \tag{11.62} \]

其中 T 是表示亚群间转移的矩阵。为完成数学问题,我们使用初始与边界条件

\[ n \cdot \bar{D} \bar{\nabla} \bar{c} = 0 \quad \text{for } \bar{x} \text{ on } \partial B \text{ (brain boundary)}, \quad \bar{c}(\bar{x}, 0) = \bar{f}(\bar{x}) \quad \text{for } \bar{x} \text{ in } B \text{ (brain domain)} \]

我们预期在肿瘤中引入多个细胞群会在模拟肿瘤的生长模式中引入更多异质性。临床与实验结果已显示可见肿瘤的指状分支。

本节作为示例,考虑肿瘤中存在两个带突变转移的克隆亚群。具体而言假设一个种群具高增长率与低扩散系数,而第二个种群具中等增长率与高扩散系数;我们可以取许多其他情形。设

\[ \bar{c} = \begin{pmatrix} \bar{u} \\ \bar{v} \end{pmatrix}, \quad \bar{D} = \begin{pmatrix} D_1 & 0 \\ 0 & D_2 \end{pmatrix}, \quad \rho = \begin{pmatrix} \rho_1 & 0 \\ 0 & \rho_2 \end{pmatrix}, \quad T = \begin{pmatrix} -k & 0 \\ k & 0 \end{pmatrix} \]

其中 D、ρ 与 k 为常数参数,(11.62) 成为

\[ \frac{\partial \bar{u}}{\partial \bar{t}} = \bar{\nabla} \cdot (D_1 \bar{\nabla} \bar{u}) + \rho_1 \bar{u} - k \bar{u}, \quad \frac{\partial \bar{v}}{\partial \bar{t}} = \bar{\nabla} \cdot (D_2 \bar{\nabla} \bar{v}) + \rho_2 \bar{v} + k \bar{u} \tag{11.63} \]

以给定初始条件 \(\bar{u}(\bar{x}, 0) = \bar{f}(\bar{x})\)\(\bar{v}(\bar{x}, 0) = \bar{g}(\bar{x})\)。以 u 为更快速增殖的群、v 为更快速扩散的群,\(D_2 > D_1\)\(\rho_1 > \rho_2\)。我们进一步假设 u 细胞是最初唯一存在的肿瘤细胞(\(\bar{f}(\bar{x}) > 0\)\(\bar{g}(\bar{x}) = 0\))。以小概率 \(k \ll \rho_1\),u 细胞突变为形成 v 细胞。虽然此处未纳入,但每个细胞群都可保留在白质区域中更快扩散的能力。

最初假设存在由健康细胞突变而来的 u 肿瘤细胞源,能比邻近正常细胞更快增殖从而开始形成肿瘤。可把 k 视为 u 肿瘤细胞突变为快速扩散肿瘤细胞群 v 的概率度量。

引入无量纲变量

\[ x = \sqrt{\frac{\rho_1}{D_1}} \bar{x}, \quad t = \rho_1 \bar{t}, \quad \beta = \frac{\rho_2}{\rho_1} < 1, \quad \alpha = \frac{k}{\rho_1}, \quad \nu = \frac{D_2}{D_1} \tag{11.64} \]
\[ u(x, t) = \frac{D_1}{\rho_1 u_0} \bar{u}\left(\sqrt{\frac{\rho_1}{D_1}} \bar{x}, \rho_1 \bar{t}\right), \quad v(x, t) = \frac{D_1}{\rho_1 u_0} \bar{v}\left(\sqrt{\frac{\rho_1}{D_1}} \bar{x}, \rho_1 \bar{t}\right) \tag{11.65} \]

其中 \(u_0 = \int \bar{f}(\bar{x}) d\bar{x}\) 为原始癌细胞总群。增长以 u 种群增殖的时间尺度度量,扩散以 u 细胞扩散的空间尺度度量。

方程 (11.63) 现变为

\[ \frac{\partial u}{\partial t} = \nabla^2 u + u - \alpha u \tag{11.66} \]
\[ \frac{\partial v}{\partial t} = \nu \nabla^2 v + \beta v + \alpha u \tag{11.67} \]

参数 \(\alpha = k/\rho_1 \ll 1\) 是第一亚群增长中因突变损失的比例。

在一维无穷域上我们可以写下解析解,从中可推导出一些有趣且高度相关的结论。在一维上

\[ \frac{\partial u}{\partial t} = \frac{\partial^2 u}{\partial x^2} + u - \alpha u \tag{11.68} \]
\[ \frac{\partial v}{\partial t} = \nu \frac{\partial^2 v}{\partial x^2} + \beta v + \alpha u \tag{11.69} \]

取 u 肿瘤细胞的初始源为 \(u(x,0) = \delta(x)\)\(v(x,0) = 0\)。v 种群扩散比 u 快,故 ν>1。u 种群的增长率大于 v,故 β<1。u 的增长率远高于突变概率,故 α≪1。

(11.68) 的解可独立于 v 方程求解,得到

\[ u(x, t) = \frac{1}{\sqrt{4\pi t}} \exp\left((1 - \alpha) t - \frac{x^2}{4t}\right) \tag{11.70} \]

代入 v 方程 (11.69) 给出

\[ \frac{\partial v}{\partial t} = \nu \frac{\partial^2 v}{\partial x^2} + \beta v + \alpha \frac{1}{\sqrt{4\pi t}} \exp\left((1 - \alpha) t - \frac{x^2}{4t}\right) \tag{11.71} \]

我们用关于空间变量 x 的 Fourier 变换求解该方程。变换及其逆定义为

\[ F[f(x,t)](\omega, t) = \int_{-\infty}^\infty f(x,t) e^{-i\omega x} dx, \quad F^{-1}[F(\omega, t)](x, t) = \frac{1}{2\pi} \int_{-\infty}^\infty F(\omega, t) e^{i\omega x} d\omega \]

变换后 v 的方程为

\[ \frac{\partial V}{\partial t} = -\nu \omega^2 V + \beta V + \alpha e^{(1 - \alpha - \omega^2) t} \]

其中 \(V(t=0;\omega) = 0\)。V 的解为

\[ V(t; \omega) = \alpha \left[\frac{e^{(1 - \alpha - \omega^2) t} - e^{(\beta - \nu \omega^2) t}}{1 - \alpha - \beta + \omega^2 (\nu - 1)}\right] \tag{11.72} \]

取逆变换并经若干代数运算后,\(v(x,t)\)

\[ v(x, t) = F^{-1}[V(t; \omega)] = \alpha F^{-1}\left[e^{(1 - \alpha - \omega^2) t} - e^{(\beta - \nu \omega^2) t}\right] * F^{-1}\left[\frac{1}{1 - \alpha - \beta + \omega^2(\nu - 1)}\right] = \alpha e^{(1 - \alpha) t} \int_{-\infty}^\infty \frac{\exp\left(-\frac{(x - \xi)^2}{4t}\right)}{\sqrt{4\pi t}} \frac{e^{-A|\xi|}}{A(\nu - 1)} d\xi - \alpha e^{\beta t} \int_{-\infty}^\infty \frac{\exp\left(-\frac{(x - \xi)^2}{4\nu t}\right)}{\sqrt{4\pi \nu t}} \frac{e^{-A|\xi|}}{A(\nu - 1)} d\xi \tag{11.73} \]

现假设 \(1 - \alpha - \beta > 0\),对卷积积分进行积分给出

\[ v(x, t) = \frac{\alpha}{2A(\nu - 1)} \exp\left(\frac{\nu(1 - \alpha - \beta)}{\nu - 1} t\right) \left[e^{-Ax} \left(\text{erf}\left(\frac{x - 2At}{2\sqrt{t}}\right) - \text{erf}\left(\frac{x - 2A\nu t}{2\sqrt{\nu t}}\right)\right) - e^{Ax} \left(\text{erf}\left(\frac{x + 2At}{2\sqrt{t}}\right) - \text{erf}\left(\frac{x + 2A\nu t}{2\sqrt{\nu t}}\right)\right)\right] \tag{11.74} \]

其中 \(A = \sqrt{1 - \alpha - \beta}/(\nu - 1)\)。该解的一个重要含义是:对 α 的某些值,可以看出原本不存在的亚群 v 最终能主导肿瘤生长。Swanson (1999) 用该解计算了两个种群,图 11.31 展示了从快速增殖的 u 群到快速扩散的 v 群的支配转移。虽然初始肿瘤仅由 u 细胞组成,但对长时间 v 亚群可以主导——取决于突变概率参数 α 的值。我们再次发现扩散比增殖对胶质瘤生长与侵袭更重要。

看待支配转移有两种相关方式:(i) 总肿瘤细胞群体积在后期由第二、更具侵略性的肿瘤亚群 v 主导;或 (ii) v 在某点或某邻域内主导——例如肿瘤中心附近——但不一定填满整个体积。

在对肿瘤进行活检时,组织样本取自肿瘤内某个相当随机的位置。预期对该组织片的分析能给出肿瘤组成的准确描述。由于不同等级的肿瘤细胞在组织学上(物理上)常可区分,病理学家可以从组织活检中量化某等级肿瘤的比例。特别地,我们预期活检能揭示总肿瘤体积的克隆分布信息,其结果与分析整个肿瘤(如尸检可能的那样)相比不会显著不同。从 Swanson (1999) 的分析与上文关于某点处与体积上支配转移的讨论,活检实际上不一定能定义总肿瘤细胞组成(肿瘤体积)的准确表述,这一点并不清楚。

给定位置处的支配转移(Transition of Dominance at a Given Location) 为确定亚群 v 是否最终在某一位置(如活检可能进行的肿瘤中心 x=0 处)主导肿瘤组成,我们考察两亚群比 \(v(x,t)/u(x,t)\)。显然对小 t,该比小于 1,因 u 最初是唯一存在的群。若 v 亚群比值 \(v/u > 1\) 则 v 成为主导。

对大 t,两亚群的比近似为

\[ \frac{v(x,t)}{u(x,t)} \sim \frac{2\alpha}{1 - \alpha - \beta} - \frac{\alpha(\nu - 1)}{(1 - \alpha - \beta)^2} \left(\frac{1}{t}\right) - \frac{\alpha x^2 [x^2 (1 - \alpha - \beta) - 4(\nu - 1)]}{16 (1 - \alpha - \beta)} \left(\frac{1}{t^2}\right) + O\left(\frac{1}{t^3}\right) \tag{11.75} \]

故对大时间,比趋于常数值 \(2\alpha/(1 - \alpha - \beta)\)。若 \(2\alpha/(1 - \alpha - \beta) > 1\),则第二亚群 v 对大时间主导。某点(x=0)处扩散群 v 主导的条件为

\[ \frac{2\alpha}{1 - \alpha - \beta} > 1 \implies \beta > 1 - 3\alpha \tag{11.76} \]

注意该参数条件不依赖于 ν(即 u 与 v 扩散系数之比)。我们也可以用渐近表达式 (11.75) 估计支配转移发生的时间。令 (11.75) 等于 1 并对 t 求解,发现某点处(x=0)支配转移发生的时间 \(t_{\text{dominance}}\) 近似为

\[ t_{\text{dominance}} \sim -\frac{\alpha(\nu - 1)}{(2\alpha + \beta - 1)^2 - \alpha} \tag{11.77} \]

该时间依赖于 ν,且若某点处支配转移的条件 (11.76) 满足则为正。对大 ν,v 细胞的扩散系数远高于 u 细胞,至支配转移的时间较大。然而,若参数 α 与 β 满足参数条件 (11.76),则若扩散系数近乎相等(ν≈1),至该转移的时间可相当短。

体积上的支配转移(Transition of Dominance in Volume) 现在考虑可以分析整个肿瘤并确定 u 与 v 细胞各自占据的肿瘤体积比例的情形。这在生物学上类似于对整个脑的细致尸检分析。在这种情况下,支配转移发生于若 v 细胞占据的体积超过 u 细胞占据的体积时。

模型方程 (11.68) 与 (11.69) 定义 u 与 v 亚群的时空动力学。对全空间积分这些方程可以确定肿瘤细胞亚群体积的时间行为。设 \(V_u(t)\)\(V_v(t)\) 分别为 u 与 v 细胞占据的肿瘤体积,则

\[ V_u(t) = \int_{-\infty}^\infty u(x,t) dx, \quad V_v(t) = \int_{-\infty}^\infty v(x,t) dx \tag{11.78} \]

对(无穷域上的)所有空间积分方程 (11.68) 与 (11.69) 得到

\[ \frac{dV_u}{dt} = (1 - \alpha) V_u \quad \text{with } V_u(0) = 1, \quad \frac{dV_v}{dt} = \beta V_v + \alpha V_u \quad \text{with } V_v(0) = 0 \]

解第一个方程给出

\[ V_u(t) = e^{(1 - \alpha) t} \tag{11.79} \]

\(V_v\) 满足

\[ \frac{dV_v}{dt} = \beta V_v + \alpha e^{(1 - \alpha) t} \quad \text{with } V_v(0) = 0 \tag{11.80} \]

其解为

\[ V_v(t) = \alpha \frac{e^{(1 - \alpha) t} - e^{\beta t}}{1 - \alpha - \beta} \tag{11.81} \]

体积上的支配转移发生于若比值 \(V_v/V_u = 1\) 在某有限时间时。所以

\[ \frac{V_v}{V_u} = 1 \implies e^{(1 - \alpha - \beta) t} = \frac{\alpha}{2\alpha + \beta - 1} \implies t_{\text{dominance}} = \frac{1}{1 - \alpha - \beta} \ln\left(\frac{\alpha}{2\alpha + \beta - 1}\right) \tag{11.82} \]

其中支配转移时间存在(即为正)当且仅当

\[ 1 - \alpha - \beta > 0 \quad \text{且} \quad 2\alpha + \beta - 1 > 0 \]

等价地,体积上支配转移的参数条件为

\[ 1 - 2\alpha < \beta < 1 - \alpha \tag{11.83} \]

图 11.32 显示(α–β 平面上)某点处、体积上以及两者皆有的支配转移发生的参数域。图 11.32(a) 的阴影区域代表某点处(取 x=0)支配转移发生的参数域,图 11.32(b) 代表体积上支配转移可能发生的参数值。显然存在使某点处支配转移发生而体积上不发生的 α 与 β 值。这指示了肿瘤活检准确性的失效。从肿瘤中心提取的组织片段的分析不一定能揭示实际总肿瘤组成。

含肿瘤多克隆性的多细胞模型(Multi-Cell Model with Tumour Polyclonality) 现简要考虑一个含均匀扩散的简单模型,其中存在一个细胞系到另一个的突变,但有 n 个系;其表述模仿上述做法,但存在可能复杂得多的时空动力学。一般问题则为适当的无量纲形式下

\[ \frac{\partial c}{\partial t} = \nabla \cdot (D \nabla c) + Pc + Tc \tag{11.84} \]

c 满足

\[ n \cdot D \nabla c = 0 \quad \text{for } x \text{ on } \partial B, \quad c(x, 0) = f(x) \tag{11.85} \]

其中 D、P 与 T 分别为扩散系数矩阵、增长率矩阵与种群间的转移矩阵。

定义本征值问题为

\[ \nabla \cdot (D \nabla W) + K^2 W(x) = 0, \quad n \cdot D \nabla W = 0 \tag{11.86} \]

其中 K 表示本征值。现寻找(对照第 2 章)形如

\[ c(x, t; K_i) \propto e^{(P + T - K_i^2 I) t} \tag{11.87} \]

的解,其中 \(K_i\) 为离散本征值,I 为单位矩阵。注意指数的指数是一个矩阵。在线性基上的时空行为由本征值 \(\sigma(\rho_i, K_i, D_i)\) 决定,其中 \(\rho_i\)\(D_i\) 分别为增长率矩阵与扩散矩阵的分量。本征值由

\[ \det |P + T - (K_i^2 + \sigma) I| = 0 \tag{11.88} \]

给出。

这里解的时空行为分析是本书多章中详细讨论的两物种模型系统的一般化。这类模型所支配的肿瘤生长与侵袭的含义仍是未探索的。

一些结语(Some Concluding Comments) 癌症不是单一疾病。本章我们只考虑了脑肿瘤生长与控制的一些基础模型。我们已看到,即便用本章中讨论的如此简单的线性模型,也有令人惊讶的潜在适用结果的丰富性。显然用这种建模方法还有许多工作可做。关于其他肿瘤的模型有大量文献,不可能全部提及。起点之一是 Chaplain et al. (1999) 编辑的书中关于癌症及相关主题的若干文章。其他肿瘤可涉及广泛的生物学现象,如血管生成与毛细血管网络(Chaplain and Anderson 1999)、癌症巨噬细胞动力学的图样形成(Owen and Sherratt 1997)、抑制剂与细胞黏附的作用(Byrne and Chaplain 1995, 1996)以及细胞牵引(Holmes and Sleeman 2000)。非常清楚的是,使用现实的模型有助于阐明恶性肿瘤生长中发生的许多复杂过程(或至少提出有信息量的问题);与医学情境的联系越紧密,结果越有用。

实验与理论如此紧密结合如何富有成效的一个非常好例子是 Jackson 及其同事发展的治疗癌症(特别是黑色素瘤)的新方法(Jackson 1998;Jackson et al. 1999a,b,c)。她的工作旨在提高化疗药物对癌症治疗的疗效。它涉及一个两步过程,旨在使对身体毒性最小化同时使对癌症毒性最大化。治疗包括酶结合抗体(ECA)联合前药。简单来说,与肿瘤相关抗原结合的抗体被偶联到宿主中不存在的一种酶上,然后被注入血流。ECA 在全身组织中分布,过一段时间后由于高结合亲和力被定位于肿瘤内。然后前药被注入体内,酶将前药转化为有毒形式。如此药物以最小正常组织与血流穿透进入肿瘤。建模涉及把此情境与基础生化结合起来,构建一个可在体外实验测试的现实模型。模型机制涉及空间与时间,并导致不同类型的非线性耦合偏微分方程系统。由于实验跨学科合作,Jackson (1998) 能为其模型赋予现实的参数值,使与实验的比较如此相关,并使她的方法重要。她能够与实验达到的密切定量一致(若干实验由理论推动)对此类化疗方案对各种癌症(但不幸的是不包括脑肿瘤)的方法特别令人鼓舞。

本章个人批注

第 11 章对脑肿瘤(胶质瘤)的处理在 Murray 全书中属于我读到的、与临床医学结合最紧密的一章。它不是纯生物数学——模型参数直接从 CT 扫描、活检细胞密度、大鼠体内最低线速度等具体临床与实验数据估计出来,结果也与 Kreth et al. (1993) 的 115 名实际患者生存曲线、Burger et al. (1988) 的尸检显微镜分析做了对比。这种"模型参数来自真实临床数据、模型预测又返回与真实患者对比"的双向闭环,使本章在书的格局里有一种独特的实用分量。Murray 自己反复强调:即便这样极其简单的(甚至是线性的)反应–扩散方程,竟能给临床决策提供有用的方向。这与他对数学建模作为"一种理解 + 一种预测工具"的定位是一致的。

我在本章中得到的最深印象是扩散(迁移性)对胶质瘤预后的决定性意义。这一论断在 11.2 节关于 v∝2√(Dρ) 的行进波速度中初现,在 11.5 节表 11.6 关于位置 1 高扩散中等级肿瘤预后比高等级还差的结论中再次出现,并在 11.11 节关于多克隆模型的解析分析中以"v 群最终可以支配 u 群"的方式被精确地推导出来。三种独立证据指向同一结论——这是数学建模难得出现的"硬"结果之一,提示临床上"高扩散型胶质瘤"应当作为独立的预后分层维度。

第二个让我印象深刻的是 11.8 节 Cook 的解析切除后解。作者通过 Laplace 方法两次(一次对角度积分 A→∞,一次对径向积分 x→∞)剥离出切除后肿瘤浓度与无切除浓度的比 \(\exp(-r^2/4t_r + r^2/4t)\)。这个简洁的比值直接告诉读者:切除后远处的肿瘤浓度被抑制了多少、与时间的关系如何——而且抑制项与 \(t_r\)(切除时间)有显式关系。配合 11.7 节 Woodward 切除模拟(不同切除半径 S=0、3、4、5)与 Kreth 临床数据的中位生存期 7.7 周 vs. 7.5 周的吻合,作者把"线性模型可以量化切除"这一命题从理论推到了临床验证层面。

第三是 11.11 节关于活检代表性失效的洞见。当 v(高扩散)在某点已主导但体积上未主导时,从肿瘤中心活检一块组织将系统性高估 v 群、可能给出错误的细胞组成判读——并由此误导后续治疗(特别是化疗,因为化疗反应由细胞敏感性决定)。这是数学模型反过来给临床操作提出警诫的少数例子,与 11.1 节"异质性随等级上升"的临床观察一起,揭示了"看局部 ≠ 看整体"在胶质瘤中的具体机制。

从建模史的角度,本章是 1990 年代初 Murray 与 Alvord(华盛顿大学神经病理学)合作开始的跨学科项目的总结。论文序列——Tracqui et al. 1995、Cruywagen et al. 1995、Cook 1994 (个人通信)、Woodward et al. 1996、Burgess et al. 1997、Swanson 1999、Swanson et al. 2000、Tracqui et al. 1995——构成一个递进的链条:先均匀组织 (11.3)→加入异质性 (11.6)–(11.7)→解析切除 (11.40)–(11.51)→数值切除 (11.52)→多克隆 (11.62)–(11.66)。作者把这个故事讲得很清楚,读者能跟随模型的演化追踪到参数估计、临床验证、再回到新模型的反馈过程。

需要留意的几处不严谨:(i) 11.1 节提及肿瘤边缘复发被比喻为"从火线后方扑灭森林大火",但章节未给出量化证据切除边缘附近密度确实最高——这是假说。(ii) 11.10 节提到"用 25% 检测阈值下降、其他参数变化 <10%"的稳健性声明,没有给出详细推导表。(iii) 11.11 节多细胞模型最后只到本征值方程 (11.88),未给出具体算例。这些不影响主线,但若写综述需要补充原始论文细节。

本章与第 10 章(真皮创伤愈合)的区别是显著的:10 章关心的是"愈合过程本身的生物力学建模"(收缩、应力、应变),参数以胶原性质、生长因子浓度为输入;本章关心的是"侵袭性生长 + 治疗干预下肿瘤动力学建模",参数以 CT 扫描轮廓、活检细胞密度、化疗时间表为输入。两者都是医学 + 数学,但建模对象决定了建模风格截然不同。第 10 章的力学色彩在第 11 章中几乎不可见——胶质瘤模型不需要应力–应变平衡,因为它假定肿瘤是被动的、可由扩散–反应动力学描述的细胞群。这与神经组织(颅骨包封、压力敏感)与皮肤组织(开放、可压缩)的物理差异相符。

与上下章的衔接(一段话)

第 11 章是 Murray 全书中把数学建模"拉回到临床决策"的最强一击。在它之前的第 9 章(表皮创伤愈合)与第 10 章(真皮创伤愈合)里,生物学过程是"愈合"——即组织朝恢复平衡的稳态演进;数学模型主要刻画组织如何响应伤口(特别是收缩、生长因子作用、纤维化)。第 11 章彻底转入另一个领域——生长与失控:胶质瘤是指数增长的细胞群、以扩散向外侵袭、对治疗具有耐受性。这一转换既反映了 Murray 与 Alvord 自 1990 年代初的长期临床合作,也对应于 1990 年代医学成像(CT、MRI)的成熟,使基于扫描轮廓估计肿瘤参数成为可行。第 11 章之后,第 12 章是"动态疾病中的网络与模式"——从神经生物学中的癫痫传播、动物体表色素花纹形成等案例切入。回头看,第 11 章是全书在"组织尺度上的具体器官(脑)+ 具体疾病(胶质瘤)"层次上的顶峰案例——它把反应–扩散、空间异质性、参数估计、临床验证、与治疗模型组装在一起。第 12 章则跳出来,把视野扩展到一般网络动力学,回归到 Murray 擅长的"图样形成"主题。第 11 章的位置因此既是第 9–10 章临床建模路线的延续,也是第 12 章抽象网络主题的必要过渡。