第 4 章：生长域上的图样形成：鳄鱼与蛇（Pattern Formation on Growing Domains: Alligators and Snakes）
4.1 鳄鱼条纹图样形成的实验（Stripe Pattern Formation in the Alligator: Experiments）
本章研究两个具体的空间图样问题，都与生长过程紧密相关——这与前几章研究的图样形成机制有别，本章的两个问题均涉及胚胎生长的不同方面及其对实际图样过程的影响。第一个是鳄鱼（Alligator mississippiensis）皮肤上的条纹图样，第二个是牙齿原基（teeth primordia，即牙齿前体）的空间排列。可以见到生长在发育过程中起到关键作用。第 1 章卷 I 中曾讨论鳄鱼（鳄形目整体）性别决定的机制，及其对其惊人长寿的可能关键作用；鳄形目还有代谢、生理、捕食技能等其他显著特征，均与其生存相关。
鳄鱼胚胎在所有研究对象中特别便于研究：从鳄鱼的角度看除外，其胚胎不像哺乳动物胚胎那样在发育过程中不可接近；发育在蛋内（外部于母体）进行；蛋壳和胚胎膜的广泛研究（Ferguson 1981a,b,c, 1985）已发展出在整个孵化期内半裸壳培养和操控胚胎的技术。由于蛋的成分特别是钙动力学，鳄鱼胚胎在裸壳培养下发育正常。鳄鱼（鳄形目）的另一特征是呈现某些哺乳类特征——这在研究人类颅面发育时很重要（Ferguson 1981c,d）；本章后续讨论牙齿原基空间图样时会回到这一点。本章描述的实验数据全部来自对胚胎不同发育阶段的直接观察，以及在牙齿研究中通过蛋壳和蛋壳膜开窗手术操控的效应。正因这种可达性，理论与实验得以如此有效地对接。在第一个主题——鳄鱼皮肤条纹——理论（Murray 1989）暗示了特定实验（Deeming–Ferguson 1989a），以解决发育中反复出现的问题（即遗传在图样决定中的作用）在鳄鱼条纹上的应用。在第二个主题——牙齿原基高度规则的顺序空间图样——基于现有发育生物学知识的理论对一系列实验结果作预测，可能有助于进一步理解潜在的生物学图样形成机制。
与其他空间图样形成研究一样，本书反复强调：我们并不知道发育过程中实际的图样形成机制是什么。我们虽然研究了若干可能的图样生成机制（反应扩散、化学趋化-扩散、机械系统），但特定机制的实验证据仍然缺乏。重大障碍之一是：在许多情形下我们不知道发育过程中图样生成机制何时起作用，只能观察其结果。在第 3 章讨论的斑马条纹案例中，我们通过数条纹、测其尺寸和数量相对胚胎当时大小，推断机制可能作用在发育特定时间。然而由于斑马胚胎发育数据匮乏，缺少严格实验验证。鳄鱼胚胎便于操控且生长数据可靠，所以 Murray 等（1990）决定研究鳄鱼条纹图样，尝试确定图样形成机制的启动时间并量化胚胎大小对条纹色素图样的影响。这里描述的部分工作是他们完成的。我们用其结果说明理论（Murray 1989）如何暗示特定实验（Deeming–Ferguson 1989a）以解决发育中的反复问题，即遗传在图样决定中的作用在鳄鱼条纹上的应用。一个有趣的结果是：我们将表明遗传并不像常被声称的那样在条纹细节形成中起作用。
Murray 等（1990）的实验数据来源如下：从颈部背侧到尾尖逐条计数白色条纹（头部无条纹）；记录身体（颈到臀部）和尾部（臀部到尾尖）的条纹数及尾尖颜色；测量总长、颈臀长和臀尾长，精度 0.1 mm。样本来自两种孵化温度——30°C（100% 雌）和 33°C（100% 雄）（与 Deeming–Ferguson 1989a 中同一批动物）。为了研究性别对色素图样（特别是条纹数）的影响，分析了"shift twice"实验（Deeming–Ferguson 1988）的孵化卵——这些实验把卵在 33°C 孵化，但第 7–14 天改用 30°C。该处理产生 23 雄 5 雌，尽管其余孵化时间是雄性诱导温度 33°C（Deeming–Ferguson 1988）；出现雌性是因为此时间段与该温度下性别被决定为雌。第二组实验中，30°C 与 33°C 胚胎在孵化第 10–50 天测量总长、颈臀长、臀尾长（Deeming–Ferguson 1989b）；胚胎被分配发育期（Ferguson 1985）；对两种温度下的胚胎生长作回归估计。第三组胚胎孵化 32、36、40、44、48、52 天后作形态测量并分配发育期，用解剖显微镜寻找最早的宏观色素迹象。28 期时从尾部取含黑白条纹的皮肤样本，对神经和黑色素细胞做染色。
主要结果如下。孵化温度显著影响孵化幼体的色素图样（见表 1 Murray 等 1990）。33°C 个体条纹数高于 30°C；白尾尖个体平均比黑尾尖个体多约一条；典型身体约 8 条、尾部约 12 条。但温度对孵化体长无显著影响，也不存在条纹数与孵化大小的直接关系。条纹数与性别无关联：33°C 孵化（7–14 天用 30°C）的雄性孵化体平均条纹 19.96±1.15 条，同处理雌性为 20.00±0.71 条。30°C 与 33°C 胚胎尾长与颈臀长随时间变化的回归估计见图 4.3、4.4；两种温度下的生长比较见图 4.5。33°C 胚胎长得更快。颈臀/臀尾比例达到 8/12 的时间点：30°C 约 46.5 天，33°C 约 36.5 天；这些时间与 Ferguson（1985）报道的 23 期相近。胚胎总长也与温度相关（图 4.6(a)）；回归分析显示尽管 33°C 胚胎在任何给定时间都更长，但 33°C 与 30°C 胚胎在任何给定发育期是等价的（图 4.6(b)）。色素最早可见时间：33°C 约 36 天，30°C 约 44 天；身体部位色素出现早于尾部。组织学上，黑色素细胞（melanocytes）存在于基底层；黑条纹中黑色素细胞密集、黑色素浓度高；白条纹中黑色素细胞稀少、虽产黑但有限并局限在细胞及其紧邻环境；黑白条纹分界锐利。
4.2 建模思路：确定条纹形成时间（Modelling Concepts: Determining the Time of Stripe Formation）
任何机制——反应扩散、化学趋化、机械模型、神经网络——原则上都能解释鳄鱼条纹。但 Murray 等（1990）的皮肤组织学数据更倾向支持一种细胞-化学趋化-扩散机制：细胞产生自身化学趋化因子，聚集效应（chemotaxis）强于扩散时即形成图样。该机制的第 2 章 2.6 节"行波"解尤其契合：条纹从颈部分别向前后方向按波状铺展。模式机制本身能否在某生物过程中产生图样是一回事，与实际生物学是否相关是另一回事；但不同模型往往会引出不同实验。在鳄鱼情形下，条纹图案约在孵化 40 天可见（孵化期约 70 天），所以图样生成机制几乎可以肯定要在此之前就开始作用。模型机制能与所有上述候选机制兼容，关键是它在量上能给出与生物学观察吻合的预测，并由此暗示可检验的实验设计。
对白条纹有三种可能解释：（i）黑色素细胞缺失；（ii）所有细胞同等产黑，但局部细胞密度过低而颜色不够深；（iii）产黑本身对细胞数有阈值依赖（例如某面积内细胞数需达到阈值才开始产黑色素）。鳄鱼白条纹似乎符合（i）这一类，但更细致地说是黑色素细胞数目低且产黑少——Murray 等（1990）所观察到的。
发育从头侧向尾侧推进；头部分化程度始终大于尾部。黑色素沉积也呈头先尾后的顺序；白条纹先在身体出现、逐渐向尾部推进。这暗示：若模型机制在头端被启动，由于胚胎此时长且细（即近一维），条纹将依波状机制顺序铺下。Murray 等（1990）基于 Myerscough–Murray（1992）的一维解（无量纲化细胞-化学趋化-扩散系统，参数 D=0.25, αb/µ=2，从 n=1, c=0.5 的均匀态在 x=0 处施加小扰动），数值模拟得到了从颈向尾传播的、有规律的细胞密度波。n 扩散并产趋化因子 c；c 扩散且一阶降解；α 项刻画细胞朝 c 浓度梯度方向的趋化运动；b 项为细胞产 c 的速率；n=1 为均匀稳态细胞分布加 x=0 端的小幅增加——波状稳态图样由颈向尾推进，留下规则的细胞密度异质条纹图样。结合上述组织学证据可论证：所涉细胞都是黑色素细胞，仅因细胞密度高才显深色。然而这仍主要是推测。
沿用上述假设——机制生成图样的时间（约数小时）远小于胚胎显著生长的时间——则最终会在胚胎背部形成等间距的条纹（图 4.1）。
要在生物学意义上确定图样形成的真实时间，唯一可行的是利用生长曲线数据。所需信息是发育过程中颈臀/臀尾长之比为 8/12（成体条纹数比）的时间。鳄鱼身上恰好有 Deeming 与 Ferguson（1990）的完整数据。Murray 等（1990）的生长曲线（图 4.3、4.4）显示：30°C 约 46.5 天达到 8:12，33°C 约 35.5 天。这一时间给出了机制最可能作用的时间窗口，把实验搜索范围大幅收窄。发育是顺序进行的，不同机制在不同时间起作用；有必要大致知道何时去观察。在鳄鱼身上，恰好在条纹可见之前才去定位机制几乎肯定太迟。本节主要使用数学建模概念，但这些概念导致了真正可验证（且已被验证）的生物学含义。Murray 等（1990）的整个方法学展示了：即便没有用到精细的数学，只要建模概念本身给出可检验预测、且结果可被实验证实，就能产生有价值的生物学结论。第 2 章 2.6 节中讨论的旅行波生成条纹的模型机制可与 Myerscough–Murray 的实际模拟过程结合——这就是本章后两节将用到的扩散-化学趋化机制的来源。
4.3 条纹与阴影条纹（Stripes and Shadow Stripes on the Alligator）
色素从头侧开始、再向尾部推进，这一观察对机制的最低要求是能从头端生成传播行波。Murray 等（1990）的细胞-化学趋化模型在合适参数下确实能生成这样的顺序条纹铺展。图 4.2(a) 给出了 x=0 对应头侧的典型模拟，细胞密度的峰对应所观察到的黑条纹——细胞密度高的区域能产足够黑色素，低密度则不能。这与 Murray 等（1990）表 1 的观察一致：可容纳更多条纹时尾尖常为白。一个推论是：胚胎的"默认色"本应是均匀深色；若模式机制未被激活，全身本应都成深色（密度足够产黑色素，但略浅于正常条纹间色）。这种 melanistic 形态确实罕见（Ferguson 1985）。
假设机制生成图样相对胚胎显著生长很快，则机制会在背部产生等间距图样（图 4.2(b) 和 4.7(a)）。如果机制作用跨越数日，且细胞响应能力逐渐衰减，则条纹会不规则。这些图样被解释为色素的 prepattern。波长由参数决定；图 4.2(a) 的波形峰十分锐利，正好符合 A. mississippiensis 实测条纹。其他机制也能产生类似顺序图样，但峰通常不那么锐利。
若假设机制在某特定发育期被激活且参数固定，则条纹数主要取决于该期胚胎长度。给定长度的胚胎只能容纳特定数量的波长；胚胎越长、条纹越多；参数和胚胎大小共同决定条纹数。结合鳄鱼数据：高温加速生长和达到某期的时间；在高温下机制激活时胚胎更大；3°C 温差对模型参数影响不大。在细胞图样形成之后还需时间让色素过程显现条纹。所以模式的关键时间应在条纹可见之前。33°C 与 30°C 胚胎达到可辨条纹的时间相差约 8 天；33°C 胚胎在条纹可见前的较早时间点已显著大于 30°C 胚胎。所以雄性（33°C）比雌性（30°C）多条纹完全是长度而不是性别的差异。基因数目有限，让基因规定雄雌鳄鱼条纹精确数和位置既无必要也不高效。
数值上，平均每条条纹（白条纹前缘到下一个白条纹前缘）所对应的胚胎长度 L = （尾长+颈臀长（激活时））/孵化时条纹数。Murray 等（1990）表 2 数据给出：33°C、32 天 L ≈ 4.71 mm；30°C、40 天 L ≈ 4.89 mm。33°C 比 30°C 多出 7.75 mm 长度可容纳约 1.64（7.75/4.71）条额外条纹，与 Murray 等（1990）的实测条纹数差异定性一致。雄性激活时条纹略小，可能源自温度对模型参数的影响。
Nijhout（1980a）在 Precis coenia 蝶翅眼斑上的类似温度效应与第 3 章反应扩散分析一致：参数变化与扩散系数和反应率随温度改变的数量级相符。Harrison 等（1981）对 Acetabularia 毛间距的温度效应也支持反应扩散假设。
阴影条纹。放松"机制作用时间远小于生长时间"假设，意味着图样在生长域上铺展。细胞保持产趋化因子的能力（即便效率下降），作用时间延长；皮肤在长度和宽度（尤其躯干）上都增大。当长细域在主条纹启动后稍变宽时，数学模型预测出现较不明显的、夹在主条纹之间的"阴影条纹"（图 4.7(b)）——这与图 4.1 实际观察一致。因为模式从头部启动，条纹向尾部方向密度变淡，这也常被观察到。
任何能较快形成图样的机制——反应扩散、化学趋化、机械模型——都可给出类似结论。区别机制的是各自暗示的实验。Murray–Oster 机械化学模型和本节细胞-化学趋化模型都直接涉及细胞；操纵细胞密度比操纵未知反应扩散化学物质容易，且更容易证伪，但在证伪过程中通常加深对图样过程的理解。Murray 等（1990）的全部结果都与具体机制无关——关键结论：30°C 鳄鱼孵化得到较少条纹、33°C 较多；机制起作用的时间 33°C 早于 30°C；激活时 33°C 胚胎更长因而条纹更多；条纹只通过温度与性别相关，本身不是性连锁。
幼鱼条纹。Pomacanthus（幼刺尻鱼，刺盖鱼科）初始 3 条、随生长在原条纹之间插入新条纹——这是 Kondo–Asai（1995）提出的解释；简单域增长不能解释插入条纹的细窄特征。Painter 等（1999）从已有实验获取生长率，假设所有皮下细胞都进行有丝分裂（包括色素细胞），提出扩散-化学趋化模型，趋化由独立反应扩散系统的形态发生素调节。Painter 等的二维权衡结果在数量上模拟出刺尻鱼条纹顺序形成。Aragón 等（1998）研究了反应扩散模型中边界条件、域生长与耦合对图样的影响。Denton–Rowe（1998）认为鲭鱼（Scomber scombrus）的条纹主要用作个体间精确运动信号而非伪装：脊背两侧亮暗条纹中央覆盖反射血小板，使鱼体在滚动、偏航、俯仰上的细微变化立即在亮度图样上显现——这可能是鱼群协调运动速度如此之快的机制。
4.4 鳄鱼牙齿原基的空间图样：背景与意义（Spatial Patterning of Teeth Primordia in the Alligator: Background and Relevance）
鳄鱼（以及鳄形目整体）牙齿原基在颌骨上的发育是图样在动态方式下形成的又一例子。细胞协作构建每个原基（最终牙齿的前体），形成精确的时空序列。该过程在颌骨动态生长过程中发生，生长与图样形成的交互对最终间距和出现顺序至关重要。这些牙齿原基的最终图样为功能性齿列奠定基础。
在演化过程中牙齿发育并非没有缺陷。例如先天性畸形如唇裂与腭裂影响全球儿童，可引起喂养、呼吸和言语困难。即使到今日，关于颅面发育的生物学事件细节仍知之甚少（Ferguson 1981c, 1994），部分原因是哺乳动物胚胎发育体内研究的不可接近性。虽然鳄鱼属爬行类，但其牙齿与颌具有某些哺乳类特征（Ferguson 1981c,d），便于作非侵入性与侵入性在体研究，特别是腭发育与牙齿启动同时在胚胎颌骨上发生。这让研究者能够开始调查正常和异常发育过程的生物学机制。牙齿原基发育模式研究与腭研究互补，构成颅面发育的重要方面。鳄鱼牙齿发育和腭闭合研究可能有助于确定人类发育过程，或至少为如何预防此类先天畸形和预测人类产前腭裂修复效果提供线索（Ferguson 1981c, 1994）。结合胎儿伤口愈合（后文专章讨论）相比成人伤口修复的优势，潜在收益可能很大。所以在基础层面上，这些深入研究完全合理——尽管本章主要关心图样形成机制的一般原理和该具体生物学案例。
鳄鱼与人类颌骨及齿列的相对相似性见图 4.9。还有其他相似点，当然也有许多差异。两者都有次生腭、单列牙、相似腭结构；但人类为二套齿（biphodont，一生两套牙）、三型牙（切牙、尖牙、磨牙），鳄鱼为多套齿（polyphydont，许多套牙）、一型牙但尺寸不同。在所有鳄鱼中，A. mississippiensis 具有最像哺乳类的吻部和次生腭（Ferguson 1981b,c,d）——这让它成为与人类牙齿发育比较的特别有用的模型。图 4.9 给出 13 英尺（4 米）鳄鱼与人类颌骨齿列的视觉对比。
Westergaard–Ferguson（1986, 1987, 1990）实验性记录了牙齿原基的启动与空间图样，构成将实验观察与假设纳入理论框架的数据库。本章回答两个基本问题：单个原基启动的机制是什么？精确的空间分布如何确定？通过构造一个反应扩散模型机制（但与之前章节研究的有本质区别），并在第 4.8 节展示广泛模拟与实验的比较结果；然后用模型预测可能的实验结果，以指导进一步的实验。借助对牙齿原基启动可能机制及其如何实现的更深入理解，希望为上述人类颅面先天缺陷的形成提供线索。除了教学价值外，研究鳄鱼牙齿发育还有众多生物学理由。体积 I 已给出研究鳄形目的一般理由（讨论其显著长寿时）；这里鳄鱼条纹研究说明了研究如何增进生物学理解；牙齿研究进一步延伸这一理解。本节主要使用建模概念但不使用具体数学——建模概念与具体生物数据的融合才是工具。这也是本章结构安排的基本思路：先呈现生物学数据，再讨论可能的机制，最后给出数学建模。Westergaard–Ferguson 的工作提供了从第 1 天到第 75 天的准确启动数据，是首批足够详细以揭示胚胎发育期牙齿精确图样的实验。
4.5 牙齿启动生物学（Biology of Tooth Initiation）
脊椎动物牙齿形态大小各异，但都经历相似的发育阶段。颌骨中有两层基本细胞层：上皮层（排列成片）和下方间充质（mesenchyme，含可运动细胞、结缔组织与胶原）。图 4.10 示意启动早期事件：最早可见的牙齿结构征象是原基；原基首先表现为口腔上皮局部增厚形成基板（placode）；通过一系列复杂的上皮-间充质（真皮）相互作用（颌骨生长同时进行），这些上皮细胞团内陷进入间充质，造成局部间充质细胞聚集（牙乳头，papilla），形成牙蕾（tooth bud）。在某些脊椎动物中，早期原基退化进入间充质、被吸收或脱落；另一些则早期原基就发育为功能性牙齿。在鳄鱼身上，后续原基以高度稳定精确的时空序列持续形成，类似过程亦见于羽毛胚芽启动（第 6 章）。本章关心基板（placode）的空间图样。
组织相互作用建模在第 6 章讨论，模型直接与牙齿形成相关——每种组织都有图样生成器，两个耦合系统必须作为整体研究。Cruywagen（1992）、Cruywagen–Murray（1992）给出了组织相互作用模型与简化 caricature 模型；通过牙乳头阶段研究牙齿形成将需要非常复杂的组织相互作用模型。
牙齿启动与腭形成都是胚胎事件，体内观察的不可接近性限制了研究。鳄鱼胚胎在蛋内发育，便于通过蛋壳开窗进行手术操控和全程观察，所以同时能研究腭形成与牙齿启动（Ferguson 1981c,d, 1988）。鳄鱼具备多种其他爬行类不具备的形态特征，是与人类牙齿发育比较的有用模型；所有鳄鱼中，A. mississippiensis 具有最像哺乳类的吻部和次生腭（Ferguson 1981b,c,d）。
牙齿原基的时空序列。Westergaard–Ferguson（1986, 1987, 1990）实验性研究了鳄鱼颌骨牙齿原基精确的时空启动序列。第一个原基（dental determinant）形成于下颌前部，但并非最前端的牙齿。牙齿启动从 dental determinant 同时向前后扩散；间隙原基在空间允许处且靠近较成熟邻居的位置形成。精确间距和出现顺序示意在图 4.11（图 4.13(a) 为真实胚胎；图 4.18 为成体）。Westergaard–Ferguson 的主要结论是：牙齿原基启动直接依赖于颌骨生长，因此必须在生长域上建模动态图样形成。
早期原基在颌前后形成，间隙牙在生长空间中早期原基之间出现。Kulesa 等（Kulesa 1995, Kulesa et al. 1993, 1996a,b, Kulesa–Murray 1995, Murray–Kulesa 1996）研究了这种关系。每个发育天所形成的牙齿数在早期呈明显指数关系（图 4.12(a)，下颌 N(t)=0.0066 exp(0.3077 t)、上颌 N(t)=0.0047 exp(0.3442 t)，t 单位天），整个孵化期则呈 Gompertz 增长（图 4.12(b)，下颌 N1=71.8, N2=8.9, r=0.068/day；上颌 N1=69.6, N2=12.0, r=0.082/day）。Murray–Frenzen（1986）已对指数到 Gompertz 的过渡作过模型解释。早期颌域在该发育期必以恒定速率指数增长，这一假设在下一节模型中明确引入。
颌骨四分之一象限内早期牙齿原基启动：前 7 个原基交替形成；第 8 个开始非交替序列，到第 11 个恢复交替。后续原基形成受前一原基重吸收干扰；第 9a 形成于第一原基重吸收位置。整个 65–70 天孵化期，约 19 个早期原基被重吸收或脱落而未功能性；7 个为短期功能性或被重吸收；36 个构成功能性齿列。下颌数据（第 17 原基启动后）显示：后续基板形成在旧原基重吸收的位置附近，颌骨生长减缓、对新原基启动的重要性降低；早期原基似乎为后续原基的位置作标记。
上皮中标记牙齿基板的初始细胞局部凝聚的位置与精确启动信号机制未知。研究集中于：表皮生长因子（EGF）及其受体（Thesleff–Partanen 1987, Partanen et al. 1985, Kronmiller et al. 1991）、组织相互作用（Mina–Kollar 1987）和同源框基因的局部表达（MacKenzie et al. 1991, 1992）。小鼠实验显示 EGF 引起牙上皮增殖；Kronmiller 等（1991）化学阻断 EGF 后启动停止，证明其必要性。Ferguson（1994）综述了牙齿原基启动和形成的分子机制研究。同源框基因 Msx1 与 Msx2 在牙齿形成局部区域表达（Jowett et al. 1993, Vainio et al. 1993），这是上皮-间充质相互作用的结果。整个齿列的完整形成是信号与组织物理相互作用协调的一系列过程，Kulesa（1995）给出了详尽综述。
4.6 牙齿原基启动建模背景（Modelling Tooth Primordium Initiation: Background）
早期爬行类牙齿发育研究（Röse 1894, Woerdeman 1919, 1921）构成牙齿形成描述性模型的基础。一般分为两类：预图样（prepattern）模型和动态模型。预图样模型基于 Wolpert（1969）的位置信息概念：形态发生素在组织中反应扩散形成异质浓度场，细胞据此反应并分化。对牙齿启动，预图样模型假设原基由细胞对形态发生素浓度差异的响应所形成，核心问题是确定合适的形态发生素空间图样机制。常规形态发生素预图样模型是非动态的，即图样一旦铺设，便不会被系统动态行为（如域生长）改变。
Woerdeman（1919, 1921）研究爬行类 Gongylus ocellatus 与 Crocodylus porosus，观察到第一牙齿原基在颌前部形成。Edmund（1960a,b）据此提出颌骨从前到后行波式刺激的牙齿启动过程：波经过预定或预图样化的牙齿位置即启动基板。Edmund 借用 Woerdeman 的 "Zahnreihe"（德语"齿行"）概念。Osborn（1970, 1971）指出胚胎中没有符合 Zahnreihe 理论的牙齿启动图样，且 Edmund（1960）实际重排了 Woerdeman（1921）的数据。Zahnreihe 理论与实验数据不符，但"化学行波式冲动启动牙齿原基"的想法仍被广泛接受。
相对传统预图样模型，动态模型将图样视为自生成过程，依赖系统动态生长。鳄鱼阴影条纹章节即说明生长起了关键作用。牙齿原基仍可能由形态发生素梯度启动，但图样形成机制对系统生长作出响应。两种动态模型：clone 模型（Osborn 1978）和机械化学模型（Sneyd et al. 1993）。
Osborn（1978）首次尝试把颌骨生长动态纳入牙齿启动描述模型。clone 模型假设牙齿原基从位于上皮下间充质中一个或多个神经嵴细胞克隆启动；克隆源有边界（即 progress zone），可向前或后沿颌扩张；随克隆扩张，细胞分裂产生有龄梯度的胜任组织；在克隆边缘的成熟组织+可用空间满足时新原基启动。基于 Osborn（1971）观察（新原基靠近较老的邻居），Osborn（1978）假设每个新原基产生一个抑制后续原基在局部启动的区域。在该模型中颌内原基总数取决于克隆最终尺寸和抑制区大小。
Sneyd 等（1993）的机械化学模型描述细胞机械运动和相关组织如何产生牙齿原基结构与形态。模型结合机械化学与反应扩散机制。牙决定子（dental determinant，第一牙齿原基）的启动被建模为受颌骨生长、上皮细胞龄结构、以及细胞黏附分子（CAM，如 Chuong–Edelman 1985, Obrink 1986）的年龄依赖性产生所控制。其数值模拟预测：仅当 CAM 浓度在形成基板的邻域且紧邻发育之前迅速上升，dental determinant 才能在正确空间位置形成；除非颌前部 CAM 产生时间尺度远快于颌骨生长，否则会形成错误图样。这一来自模拟的结论表明：把固定域模型结果外推到生长域是危险的。
以上理论均有严重缺陷。Westergaard–Ferguson（1986, 1987, 1990）的系列实验从第 1 到 75 天提供了精确的启动数据，是首批足够细致以揭示胚胎发育期牙齿精确图样、并凸显上述模型在定性和定量上不足的工作。建模视角看，主要发现指出牙齿启动密切依赖于：（1）颌骨生长；（2）已有牙齿之间的距离；（3）已有牙齿的尺寸和发育成熟度。
实验结果印证了 Zahnreihe 理论的不足，并导致 Westergaard–Ferguson（1986）否定 clone 模型——其依据是按生长克隆预测的序列中牙齿不能形成。Osborn（1993）提出鳄鱼颌象限内牙齿可能由多个克隆发育，但初始位置和生长动态未指明。Sneyd 等（1993）的机械化学模型虽然能在某一约束下于正确位置产生牙决定子，但其主要贡献在于证明必须把颌域生长纳入模型，且该模型未涉及后续原基序列。
Westergaard–Ferguson（1986, 1987, 1990）明确指出：牙齿原基的空间图样并非一次铺设，而是在胚胎颌骨生长过程中动态发育。无论何种图样形成机制，至少必须能再现鳄鱼前几个原基的时空序列，这对理论模型并非平凡挑战。依据生物学数据，下四节将构建并分析这样一个动态模型机制，关键在于包含颌骨生长，满足 A. mississippiensis 原基序列出现的要求。
4.7 鳄鱼牙齿图样的模型机制（Model Mechanism for Alligator Teeth Patterning）
鳄鱼胚胎牙齿发育提供了研究牙齿形成多个方面的模型系统。本节聚焦牙齿原基启动过程。以已知生物学事实为指导构造数学模型机制，描述牙齿原基的启动和空间图样。初步目标是再现 A. mississippiensis 颌骨上前 7 个原基的空间图样——七是因为此后这些牙齿开始被系统重吸收。
Westergaard–Ferguson（1986, 1987, 1990）的实验性工作是模型机制构造基础。小鼠方面的近期实验研究补充了生物学数据库。在鳄鱼尚未显示但可合理假设某些相似特征的事实也予以纳入。
模型机制的关键是纳入颌骨物理生长。前 7 个牙齿原基的启动发生在鳄鱼孵化期前三分之一。牙齿原基数对早期呈指数（图 4.12(a)），对孵化全程呈 Gompertz（图 4.12(b)）。必须首先定量刻画颌骨生长再嵌入模型方程。
模型假设。基于 Westergaard–Ferguson（1986）左右颌对称的事实，假设左右颌启动过程对称。原基形成域相对长度极薄（图 4.13(a)），视为沿前后轴的一维行（图 4.14）。Westergaard–Ferguson（1986, 1987, 1990）实验显示牙齿启动位置相对想象的颌前后轴线几乎没有横向偏移，进一步支持该简化。生物学上启动信号如何被"打开"未知；推测受上皮下间充质神经嵴细胞控制。上皮中细胞凝聚标记牙齿原基启动位点，但启动源亦未知；模型假设在颌后端有化学源，第一原基形成后该源作用迅速减弱。
实验显示牙齿启动和形成涉及 EGF、骨形态发生蛋白（BMP-4）、同源框基因 Msx1、Msx2，提示化学机制的存在——特定化学浓度刺激上皮某区域形成基板。模型采用反应扩散系统，但有别于前几章。
模型方程。沿 Kulesa（1995）、Kulesa–Murray（1995）、Murray–Kulesa（1996）、Kulesa 等（1996a,b）。构造动态图样形成系统必须将颌骨物理生长与能形成图样的机制结合：把静态控制化学物质介导的图样形成机制与颌域物理生长结合，得到动态图样机制。
实验证据要求图样是颌骨生长所致动态结果而非牙齿位点的预图样，因此生长是关键。基于指数颌骨生长（图 4.12(a)），假设颌长 L(t) 以恒定应变率 r 指数增长：
dL/dt = rL ⇒ L(t) = L0 exp(rt), (4.1)
无量纲化后 L0=1；实验数据可给出参数估计。生长域稀释化学浓度。
考虑标量反应扩散方程在生长域上：
c_t = D c_{ξξ} + γ f(c), (4.2)
D 为扩散系数，γ 为尺度因子。设 s 为长度 l 上的反应物量，c = s/l。在 t 到 t+Δt 区间，l → l+Δl，s → s+Δs。计算极限得到在生长域上的稀释项：
c_t = D c_{ξξ} + γ f(c) - r c. (4.3)
为数值模拟方便，作坐标变换到固定域：
x = ξ exp(-rt) ⇒ c_t = D exp(-2rt) c_{xx} + γ f(c) - r c, (4.4)
扩散系数随时间指数衰减。
化学图样机制取第 2 章（2.32）反应扩散系统，结合生物学增加抑制控制化学物质 c 对底物 u 的控制方程：
∂u/∂t = γ(hc - u + u² v) + ∂²u/∂x² ∂v/∂t = γ(b - u² v) + d ∂²v/∂x² ∂c/∂t = -δ c + p ∂²c/∂x². (4.5)
u, v 分别为底物和激活子浓度，γ 为尺度因子，b, h 为常数，δ 为 c 的拟一阶降解常数，d 为激活子与底物扩散系数之比，p 为 c 与底物扩散系数之比。假设 c 为与 EGF 相关的抑制剂，浓度随 EGF 增加而降低、随 EGF 降低而升高；c 来源于颌后端。
经 (4.3) 比例变换和无量纲化，得到在固定域 0≤x≤1 上的方程 (4.6)–(4.8)，每条都含稀释项 -r u, -r v, -r c 和时间依赖扩散系数 D exp(-2rt)。边界条件 (4.9) 和 (4.10)：u, v 在两端零通量；c 在 x=1 零通量，x=0 有源 c0(t)（随时间减小，对应后端源）。x=1 对应颌前端的对称条件。
按第 2 章 2.4、2.5 节分析，反应扩散系统在域长大于最小值时能形成稳定空间图样；参数变化可选择异质稳态或规则空间图样。当 c 高于阈值，u, v 的图样被抑制；c 低于阈值，机制通过扩散驱动失稳启动，子阈值域足够大时形成空间图样。模式与子阈值域大小及参数共同决定具体图样；具体参数空间详见第 4.21 图，对 hump-like 图样作了分析（h c 为常数 a、固定域大小）。本节模型中 c 在时空上变化、域可变，等价分析不易进行，但直觉上能看清行为。
生物学强烈支持已形成的牙齿原基成为 c 的源（抑制因子）。Westergaard–Ferguson（1986, 1987, 1990）、Osborn（1971）实验显示：新原基周围有抑制后续局部原基的"区域"。Westergaard–Ferguson（1986）还注意到若新原基在两个较老原基间形成，会偏向较老邻居侧。数学上每个新原基成为生长因子化学源——抑制局部原基形成。新牙齿位置定义为 u(x,t) 在 0≤x≤1 子域上超过阈值的位置。一旦 u 超过阈值即开启新牙齿源 c̄_i（模拟抑制区），按 logistic 增长：
d/dt c̄_i = k1 c̄_i (1 - c̄_i/k2). (4.11)
这些新牙齿源以 delta 函数形式加入 (4.8) 右端：
Σ_i δ(x - x_i) H(u - u_threshold) F(t - t_i). (4.12)
δ(x-x_i) 为 Dirac delta，H(u-u_th) 为阶跃函数（一旦开启即使 u 降回阈值以下仍保持 1），F(t) 为 (4.11) 右端函数。
机制运作。基本系统 (4.5) 中 c = a 为常数时，参数在合适范围内存在最小域长 L_c，小于则无图样、大于则依模式（mode）形成图样（图 4.15(a) 示意参数空间）。图 4.15(b) 是 c 方程在 x=0 有源、x=1 零通量时的典型解。c 随 -δ c 项衰减，存在临界 L_c 使得 L_c<x<1 范围内平均 c 足够低、让 (4.14) 形成 mode 2 图样（图 4.15(b)(c)）。c = a 情形可由第 2 章方法解析求临界 L_c；c(x,t) 情形需全系统数值求解。但直觉上当 L_c 足够大时，生长域 (4.6)–(4.8) 中会出现 mode 2 模式 (4.15(b))，u 增大；当 u 跨过图案化阈值即触发 placode 启动、固定 tooth 1 位置（图 4.16(a)）。Westergaard–Ferguson（1986）显示牙决定子和每个后续原基成为抑制源——当 u 越过阈值，u 峰位置成为 c 的源，相当于每个牙齿处的内部边界条件 c(x_i, t) = c̄_i(t)（(4.15)）。dental determinant 出现后抑制源增至两个（颌后端与牙决定子位置 x1）。颌骨继续生长，c 在两源之间降到临界阈值以下，u 在颌后端再出现 hump-like 模式；第二个原基在 u 再次跨阈值处形成并成为新 c 源（图 4.16(b)）。如此反复：c(x,t) 跌破阈值 → 域足够大时形成局部图样 → u 跨阈值 → 产生 c 源 → 新原基形成。基于图 4.10 牙齿形成场景，基板（由底物 u 启动）引发间充质中细胞聚集（papilla），二者位置相同；但 placode 与 papilla 哪个先发生的细节尚未定论——第 6 章再议。
4.8 结果与实验数据比较（Results and Comparison with Experimental Data）
按上述启发式思路数值求解前 7 个牙齿原基的时空序列，详见 Kulesa（1995）、Kulesa–Murray（1995）、Murray–Kulesa（1996）、Kulesa 等（1995, 1996a,b）。先估计模型参数：r 由 Westergaard–Ferguson（1986, 1990）实验给出；p（与 EGF 相关）有真实生物学估计；δ（与 EGF 降解有关）也有生物学估计；其余参数通过线性分析简化形式 (4.6)(4.7)（即 (4.14)）按所需 u 模式确定。具体参数集选取使 mode 2 模式在扩散驱动失稳下被启动。当不能解析求时，Bentil–Murray（1991）的逻辑数值程序可通用地确定参数集——这种逻辑数值程序对一般模式生成器机制都适用。整个过程本质上是：先由实验给出关键参数 r、p、δ；再由线性分析确定能给出正确模式的余下参数。这是参数估计中常见的策略。
边界条件由 (4.9)、(4.10) 给出，初始条件 (4.16) 指定初始 u, v, c 分布：
u(x,0) = u0 = h + b v(x,0) = v0 = b / (u0)² c(x,0) = a exp(-k x), (4.16)
其中 a, k 为正参数。颌后端 c 源取为时间连续函数：
c(0, t) = c0(t) = -m tanh[(t-f)/g] + j, (4.17)
m, f, g, j 为常数；这是初始图案化过程开关的平滑阶跃函数形式。对 u0, v0 施加小随机空间扰动后，用 Crank-Nicholson 有限差分格式（x=0.01, Δt=(Δx)²）数值求解非线性模型系统 (4.6)–(4.8)，得到下颌前 7 个牙齿原基的时空序列。所有参数值见图 4.19 图例；该图例给出的参数集是大量参数空间搜索的结果，能使模型再现实验观察到的牙齿启动序列。
图 4.17 显示 u 和 c 数值模拟的浓度剖面（直至前 5 颗牙齿形成时间 t_i）。图 4.17(a) 中 t1 是 u 达到阈值 uth 启动 tooth 1 位置 x1 的时间，并在该位置打开 c 源、构成零通量壁垒；随后以新初始场在 0<x<x1 和 x1<x<1 两个区域重启模拟，直至 c 在足够大的域内再次下降、u 在某处达到 uth（位置 x2，记为 tooth 2 时间 t2）；再在三个区域 0<x<x2, x2<x<x1, x1<x<1 内继续，求 u 在某位置达到 uth（t3, x3，即 tooth 3）；如此递推至所有原基位置与时间确定。模拟结果画在 [0,1] 上，但实际域大小为 [0, exp(rt_i)]（r 为颌骨生长率参数）。每个新模拟在牙齿形成的左右两侧区域重启过程持续进行。该方法模拟上下颌完整序列的能力依赖于：将每次重启的初始条件取自前一次模拟结果，并把新增牙齿源 (4.11)(4.12) 加入反应扩散方程右端。
图 4.18(a) 给出 u 随时间（天）和沿颌前后轴距离的三维演化；图 4.18(b) 为孵化天数与模拟颌长关系，叠加原基顺序。Westergaard–Ferguson（1990）的实验数据给出上下颌牙齿出现序列，导出图 4.12。从图 4.12(a) 可见前几个原基中上颌第 6、7 与下颌在空间位置略有差异；上颌生长速率略高。模拟中除上颌 r=0.34/day 较下颌 r=0.31/day 略高外，使用相同参数集；实验数据可把仿真时间 T 与真实时间 t 直接联系（T = a1 t + a2, a1=27.06, a2=-286.6）。无量纲仿真时间通过该线性关系映射到真实孵化天数；该映射是经验性的，由模拟数据与实验数据拟合得到。
图 4.19 显示带参数集的上下颌牙齿原基序列模拟与实验数据对比：上颌数值数据（∗, 实线）N(t)=N3 exp(r3 t), N3=0.0042, r3=0.35/day；上颌实验数据（点划线 +）N(t)=N4 exp(r4 t), N4=0.0047, r4=0.34/day；下颌数值数据（×, 虚线）N(t)=N1 exp(r1 t), N1=0.012, r1=0.28/day；下颌实验数据（点线 ○）N(t)=N2 exp(r2 t), N2=0.0066, r2=0.31/day。模型机制正确再现了下颌前 8 个牙齿的时空序列、上颌前 6 个牙齿的正确序列。模型参数为 h=1, b=1, d=150, r=0.01, p=0.5, γ=40, σ=0.2, k1=0.3, k2=1.0, a=2.21, k=0.9, c0(t)=-m tanh[(t-f)/g]+j, m=0.65, f=200, g=34, j=1.5。Kulesa 等（1996a）还把数值结果与实际颌骨上原基位置做了直接对比，吻合良好。实验空间位置由附着于图 4.13(a) 的 xy 坐标系获取；一条沿牙齿原基的抛物线 f(x)=b1 x²+b2 x+b3（b1=-0.256, b2=0, b3=7.28）拟合实验数据；数值空间位置用 f(x) 从 x=0 到 x=x（f(x)=0）的长度 l 作为前后轴无量纲长度。结果令人鼓舞，因为模型机制成功再现了实验观察到的牙齿启动序列的时空特征，而非仅是粗略定性相似。完整非线性模型方程组的完整数值分析和结果讨论详见 Kulesa（1995）和 Kulesa 等后续论文。
4.9 预测实验（Prediction Experiments）
构建直接关联真实生物学问题的理论模型机制的主要收益之一是将其用作预测工具。数值模拟可确定模型机制能否再现观察到的实验数据；更重要的是可数值模拟可能的实验场景——虚拟实验能覆盖更广参数空间且更快。基于模型预测可识别关键生物学过程、指导进一步实验；也可检验模型假设。
模型能再现主要实验数据后，可设计若干干预实验作数值模拟，预测牙齿启动序列的可能结果。Murray–Kulesa（1996）描述三类：移除牙齿基板并用口腔上皮替代、牙齿基板移植、在颌上皮插入物理屏障影响扩散过程。图 4.20(a) 给出上下颌正常序列（从 x=0 后端计数为 7-3-6-2-8-5-1-4；第 16 天颌约 0.6 mm，tooth 4 形成时约 3.8 mm）。
基板移除并用口腔组织替代。Murray–Kulesa（1996）考虑移除第一牙齿原基后的情形。实验假设细胞凝聚标记第一原基位点后，placode 被移除并替换为另一段口腔上皮。数值模拟中当 u(x,t) 跨阈值标记第一牙齿位置时，不增加新牙齿源，继续模拟相当于 dental determinant 被移除。用图 4.19 同参数集模拟。结果：模型预测新牙齿原基在原移除位置形成（相当于原牙齿再生）。c(x,t) 仍低于阈值，移除与插入口腔组织只延迟该位置后续原基形成。tooth 2 移除替换的情形类似；后续牙齿原基启动序列不变。任何牙齿原基移除后都预测在同一位置再生，符合直觉。
基板移植。Murray–Kulesa（1996）考虑将牙齿原基移植到口腔上皮不同区段。图 4.20(b) 显示将 dental determinant 移植到 x=0.9 的情形。第一牙齿原基形成后，在颌域 0≤x≤1 内前后各位置移植；模拟了 x=0.25, 0.5, 0.9 的情形。x=0.9 情形下，后续 7 个原基全部在移植 placode 的后端形成——移植 placode 与颌前端之间（0.9≤x≤1）的域增长不足以容纳至少 7 个后续原基——结果 tooth 4 在 tooth 1 后端形成；dental determinant 显然重要。移植到 x=0.25 时新序列为 1-3-5-2-4；移植到 x=0.5 时为 5-3-6-1-4-2-7（图 4.20(b) 图例）。其他牙齿原基的移植也可模拟，颌两端移植引起最显著变化。
颌中屏障。可能最简单的干预方式是在颌上皮放置物理屏障。屏障虽然物理侵入性强，但对化学扩散构成阻碍——这在模型机制中起关键作用。在 dental determinant 形成（作为参考时间）之后，把屏障放在 0≤x≤1 内的若干位置，模拟位置 x=0.25, 0.5, 0.9。屏障在数学上表示为内部零通量边界条件 (4.18)：u_x(x_b, t) = v_x(x_b, t) = c_x(x_b, t) = 0, t > t_b。最显著预测仍来自对颌后端附近的干预：x=0.9 处屏障阻挡 tooth 1 源对 0.9≤x≤1 区域的影响，使该处 EGF c(x,t) 更快跨过低阈值，从而让某原基（图中 tooth 3）提前启动（图 4.20(c)）。x=0.25 处屏障的序列为 4-2-3-1-5；x=0.5 处为 5-2-3-1-4（图 4.20(c) 图例）。这些屏障尽管物理侵入，但放在颌上皮中对实验干扰可最小化；屏障产生化学扩散阻碍，在模型机制中起关键作用。屏障位置和时间选择同样重要。
稳健性。在结束牙齿形成讨论前，必须处理图样形成与时序稳健性问题。Kulesa（1995）深入研究，对完整模型机制 (4.6)–(4.8) 在大量参数值上计算了参数子空间的选择。模型方程中参数为 h, b, δ, p, d, γ，空间图样参数域为 6 维。第 2 章 2.5 节已详细讨论了 (4.5) 前两方程（c = a 常数）的 3 维参数空间（Turing 空间）。Kulesa（1995）的目标是找到能给出前 7 个牙齿正确启动序列的参数集。图 4.21(a) 显示 γ-d 平面上主导空间模式 n 的二维截面（n=2 对应单 hump），其他参数见图 4.19（除 b=2.0 外）。模型机制对前 7 个牙齿稳健——参数域不算小，每个点是该参数值下的一次模拟。稳健性对应真实发育中常见的参数小随机变化。图 4.21(b) 中若参数集在 7 范围内（p-δ 平面，r=0.01 在研究 r 范围中部，其他参数见图 4.19），p 和 δ 的小变化仍给出正确间距，除非靠近分岔边界，数字表示按正确顺序出现的牙齿原基数。若图样形成机制参数空间很小，几乎不可能在现实参数变化下保持同一类型图样。
4.10 鳄鱼牙齿空间图样结论性评述（Concluding Remarks on Alligator Tooth Spatial Patterning）
本章讨论的两个 A. mississippiensis 问题在建模视角下形成鲜明对比，但都凸显了现实生物学建模的一些基本原理。条纹形成时机（包括阴影条纹）讨论除一个示意性的扩散-化学趋化模型外，并未采用任何特定模型；真正对实验者有用的是建模概念而非具体模型。数学分析本身相当简单。如 Murray 所言：如果能提出具有合理生物学依据的具体机制、或能在其上进行生物学上可行的虚拟实验以检验其合理性，那会是有趣且有用的。牙齿形成部分则基于牙齿原基形成的已知生物学事实，是一般所需的、明显更复杂的机制。
总结牙齿部分更复杂建模练习的成就。基于已知生物学，构造了一个基于动态反应扩散系统的定量模型机制，关键在于纳入颌骨生长。模型机制能再现鳄鱼颌骨前 7 个牙齿原基的时空序列，结果与实验吻合良好。Kulesa（1995）检验了模型对参数变化的稳健性：在参数空间中有相当大区域能保持前 7 个牙齿的正确序列；机制对 γ, b, d 等图样参数变化稳健，能维持 u, v 单 hump 空间图样并给出正确原基顺序。该部分基于已分析稳健系统，结果并不意外。
数值模拟验证了颌骨生长对牙齿原基精确时空序列至关重要的实验假设。颌骨生长率 r 是取得正确牙齿原基启动序列最敏感的参数；上颌模拟使用相同参数集，只需把 r 提高（由实验数据得来）即可得到与下颌不同的上颌正确序列。这一发现说明模型对生长率变化是敏感的，与生物学观察一致——上下颌生长率差异是决定两颌牙齿图样差异的关键因素之一。
模型构建涉及把颌骨域物理生长纳入底物-激活子图样机制。Kulesa（1995）证明：仅由底物、激活子、抑制剂组成的图样机制（无颌骨生长）即使有抑制剂也不能得到正确序列；底物-激活子系统（含颌骨生长）也无正确序列。因此颌骨生长与底物-激活子-抑制剂图样机制都是牙齿启动序列所必需的；只有全部纳入才能得到前 7 个牙齿原基的精确时空序列。这是模型设计中的一个具体可检验结论——任何简化都会破坏图样。
预测实验显示原基移植可显著改变启动序列。移植到颌域两端的原基对空间顺序影响最大；颌上皮零通量物理屏障也会改变序列，两端屏障造成最显著变化。预测结果说明颌域两端是模式形成的关键区域——屏障和移植引起的变化模式具有空间系统性。
更多实验将使牙齿原基启动机制更精确，可纳入实际化学成分和更精确的生长数据。分子水平实验结果将特定基因表达与牙齿启动生物学事件联系，提供进一步建模信息。新实验有助于构建更详细的牙齿启动模型机制，纳入组织相互作用，使用本章建模的定性结果；这样的模型分析将非常复杂。如第 6 章所讨论的组织相互作用建模将不可回避地需要把本章的反应扩散机制与机械-化学机制耦合。
当然我们不能说牙齿原基的实际生物图样形成机制（甚至逼真 caricature）是这种反应扩散类。但我们可以说：生长在所形成的图样中起着连续且本质的作用；每个牙齿位置几乎肯定不是预先指定的。最终生物学相关性的仲裁须依靠进一步细致的分子水平实验。本章两个主题中实验研究与理论建模的交互都加深了对基础生物学的理解。这两个案例都展示了建模如何引导实验假设、聚焦搜索方向、并产生可被实验证伪的具体定量预测。
4.11 蛇的色素图样形成（Pigmentation Pattern Formation on Snakes）
生物学背景。蛇（squamata 目的蛇亚目）——爬行类和两栖类总体——在形态和生理上极其多样，呈现非常丰富的图样。Greene（2000）的书涵盖演化、多样性、保护、生物学、毒素、社会行为等。Klauber（1998）则更专门地以响尾蛇为对象。同一物种内常出现极端图样多态性；加州王蛇 Lampropeltis getulus californiae 即是好例子（Zweifel 1981），个体上的图样异常也时有发生（图 4.22(e)(f)）。蛇的图样除简单元素（侧条纹、纵条纹、点）外，还有在其他动物身上找不到的复杂图样元素。许多常见蛇图样似乎不属于反应扩散模型通常能形成的图样类别，除非修改边界条件、生长或变化域、空间变化参数等。空间变化参数情形下，模型结果难以与生物学关联。
爬行类皮肤是身体最大器官，有许多有趣的发育问题（Maderson 1985 综述）。皮肤由外层表皮和下层真皮组成；图样固定于真皮，表皮周期性蜕皮后基本图样保持。色素细胞前体（chromatoblasts）由神经嵴迁移而来、在真皮中较均匀分布。如鳄鱼条纹讨论所述，是否产生色素取决于前体细胞是否活化；这些前体细胞的相互作用和可能的导向运动造成色素与非色素细胞在不同区域聚集形成条纹或点（Bagnara–Hadley 1973）。色素细胞何时获得产色素能力未知；从鳄鱼发育看细胞可能在很早就能产色素。已承诺产色素的细胞仍可分裂一段时间。
色素发育和色素细胞前体迁移的实验研究主要集中于两栖类、哺乳类和鸟类；爬行类除前述鳄鱼工作外研究很少。对蛇和爬行类基本假设是：迁移、分裂、分化的基本过程与其他动物类似。蛇的胚胎学、生态学、演化生物学研究都较少（Hubert 1985 部分文献；Greene 2000 是最权威的）。鳄鱼胚胎学研究已较完善（Ferguson 1985, Deeming–Ferguson 1989a）。
Hubert–Dufaure（1968）描绘了 asp viper（Vipera aspis）的发育。约 106 mm 时色素首先出现在体鳞上，向头部扩展。图样几乎可以肯定在可见之前更早铺设。Zehr（1962）对 Thamnophis sirtalis sirtalis 的观察也提示类似过程：色素图样初现时形态不完整，随发育进行变得更清晰（蝴蝶翅也有类似渐变发育）。Treadwell（1962）注意到 Pituophis melanoleucus sayi 胚胎 29 天在体侧出现 3 行点、31 天背中线出现斑。蛇发育事件的时间要谨慎对待——温度对卵生或胎生种类影响很大；asp viper 妊娠期从 90 到 110 天不等。
鳄鱼情形下能估计条纹铺设时间——远早于可见。细胞-化学趋化-扩散机制同样适用于此。鳄鱼白条纹区有黑色素细胞但产黑有限，缺乏 melanogenesis 可能需要黑色素细胞密度达到阈值，这与鳄鱼身上阴影条纹（夹在深条纹间、位于腹侧）观察一致。类似交错条纹也见于鱼类（如 angelfish Pomacanthus；Aragón et al. 1998, Painter et al. 1999, Painter 2000）。鳄鱼研究给出对蛇皮肤图样两点重要启示：（1）模式生成期间胚胎生长显著时，可产生差异颇大的图样；（2）色素细胞前体图样可能远早于它们开始产色素而被铺设。
本章与第 4.1–4.4 节相同兴趣——生长域上的图样形成机制产生的图样。空间异质解可以与固定域机制得到的完全不同、复杂得多。许多胚胎图样形成机制的作用时间尺度与胚胎生长相当。Maini 等（Aragón et al. 1998, Painter et al. 1999, Crampin et al. 1999）专门研究了反应扩散系统中生长与图样形成的相互作用。基础系统在生长域上产生的新颖而复杂的图样很可能是所有图样形成机制的共性。
由于蛇胚胎形态学数据稀缺，尚无法确定蛇胚胎发育中机制的作用时间。Murray–Myerscough（1991）的目的是说明蛇的一些复杂皮肤图样如何被生成——这是任何潜在机制首先必须满足的必要条件。
4.12 细胞-化学趋化模型机制（Cell-Chemotaxis Model Mechanism）
本节模型涉及实际细胞运动。直接涉及细胞的图样形成模型更易于与实验联系——前文所述组织学观察到的色素细胞密度变化提供了一些实验支持。Le Douarin（1982）推测化学趋化可能是色素细胞迁移入皮肤的因素。直观上看化学趋化可促成点、条纹的聚集和锐化。模型假设 chromatoblasts 既响应又自身产生趋化因子。该机制能促进分化细胞在皮肤某些区域聚集——关联所观察到的蛇皮肤图样。细胞除趋化响应外还扩散；细胞有丝分裂、扩散、化学趋化的相互作用可产生空间异质性。
具体机制为 Oster–Murray（1989）提出的发育约束模型：
∂n/∂t = D_n ∇²n - α ∇·(n ∇c) - r n (N - n), (4.19) ∂c/∂t = D_c ∇²c + S n/(β + n) - γ c. (4.20)
n, c 为细胞和趋化因子密度；D_n, D_c 为扩散系数；r 为线性有丝分裂率、N 为初始均匀细胞密度（logistic 增长形式）；α 为化学趋化强度；S 为细胞最大分泌率；γ 为趋化因子自然降解率；β 为趋化因子产生的 Michaelis 常数。
无量纲化 (4.21) 后方程 (4.22) 化为（去掉星号）：
∂n/∂t = D ∇²n - α ∇·(n ∇c) + s r n (N - n) ∂c/∂t = ∇²c + s (n/(1 + n) - c). (4.22)
Winters 等（1990）、Myerscough 等（1990）、Maini 等（1991）报告了在长宽比大的矩形域（细胞和趋化因子边界零通量）上的详细数值模拟和复杂的分岔图样序列。选择长矩形域是因为皮肤图样很可能在胚胎已经长且类蛇形时铺设——即便蜷曲状态已是长柱形。Hubert–Defaure（1968）、Hubert（1985）给出 Vipera aspis 胚胎的细节。理论上在蜷曲柱面上研究更真实，但平面域数值模拟已困难重重；这里关心机制可产生的图样种类，把柱面铺成平面，主要特征类似；当然也可采用周期边界条件。方程 (4.22) 有正均匀稳态 n0 = N, c0 = N/(1+N)。(4.23) 围绕稳态线性化（设 n=N+u, c=c0+v, |u|,|v|<<1）得 (4.24)：∂u/∂t = D ∇²u - α N ∇²v - s r N u, ∂v/∂t = ∇²v + s (u/(1+N)² - v)。边界条件 (4.25) n·∇u = n·∇v = 0。
按第 2 章 2.3 节方法，设 ∝ exp(ik·x + λt) 得到色散关系 λ(k²) 满足特征多项式 (4.26)：
λ² + [(D+1)k² + rsN + s]λ + [D k⁴ + {rsN + Ds - sN α/(1+N)²} k² + r N s²] = 0。
在矩形域（边长 L_x, L_y）上，由零通量边界条件得 k = (k_x, k_y)，k_x = m π/L, k_y = l π/L；k² = π² (m²/L_x² + l²/L_y²)（4.27）。可选择 D, α, s, r, N 隔离唯一不稳定波矢；模式选择即强制初始增长某特定图样。隔离模式的波矢满足 λ(k²) = 0，即 (4.28) 双根条件给出临界波矢模量 k²_c = (s² r N/D)^{1/2} (4.29)；由 (4.27) 找到合适 k²，从 (4.28) 解 α（取较大正根）。这确定了 (N, D, r, s, α) 空间中隔离该模式的点，故合适参数下方程的解是空间异质的。该模式选择过程是标准的——隔离唯一不稳定模式即保证初始增长为特定图样。
若 r 或 N 减小，临界波数减小、图样间距增大；减够则图样消失。这是一条可实验检验的预测——例如基细胞密度的实验操控。肢体发育中通过细胞密度实验变化对空间结构解分岔的影响详见第 7 章；该处理论与实验预测的吻合程度惊人。
现提出该细胞-化学趋化机制 (4.19)(4.20) 作为蛇皮肤图样形成的候选机制，观察图样反映底层细胞密度空间图样。如前所述，发育中何时铺设图样、铺设时间与胚胎显著生长时间的关系（即方程相关空间域大小）未知。域大小是显著参数（第 3 章分析）。若干参数跨过分岔值后图样开始演化。本节展示方程系统能产生的图样并与具体蛇图样关联。
模拟方法与动物毛皮和蝶翅不同。Murray–Myerscough（1991）假设细胞分化速率和胚胎发育使得化学趋化系统在色素图样固定时已达稳态或近稳态，故求解 (4.24) 的稳态形式 (4.30)。在长矩形二维域上用 ENTWIFE 软件包作数值模拟（Winters 等 1990 详述），沿参数变化跟踪分岔解；这里主要关心解的生物学含义和相关性。该软件包允许跟踪不同分岔路径（从不同稳态解分叉）；Grindrod 等（1989）对 (4.30) 在一维 s=0 情形的有限幅稳态解作了部分解析工作。
4.13 蛇图样的简单与复杂元素（Simple and Complex Snake Pattern Elements）
Murray–Myerscough（1991）通过调整有丝分裂率 r 和化学趋化参数 α 生成了多种条纹图样和规则点图样。基本侧条纹（图 4.23(a)）是常见蛇图样元素——bandy-bandy（Vermicella annulata，图 4.23(b)）、珊瑚蛇（Micrurus）、环蛇（Bungarus fasciatus）、加州王蛇的环纹变种（图 4.23(d)）。侧条纹正是第 3 章所预期。
更有趣的是模型能产生沿长边方向的纵条纹（图 4.24(a)），见于 ribbon snake（Thamnophis sauritus sauritus，图 4.24(b)）、garter snake（Thamnophis sirtalis sirtalis）、four-lined snake（Elaphe quatuorlineata）。起初令人惊讶的是对某些参数集而言，域大小或参数值的较小变化足以使侧条纹和纵条纹互换——这些图样并不稳健。这具有重要生物学含义：同一物种的不同个体可呈现侧条纹或纵条纹。例如加州王蛇既有侧条纹（图 4.23(d)）也有纵条纹（图 4.24(c)），图 4.22(a) 亦是一例。
选择合适参数值还可生成规则点与斑图（图 4.25）。规则点见于多种蛇。例如 Cape mountain adder（Bitis atropos atropos）显示与模型结果类似的交替半圆形图样。
在本化学趋化模型中点是条纹取决于初始条件、域形状大小、以及 α, D, r, N 的参数值。Murray–Myerscough（1991）只改变化学趋化参数 α，但其他参数同样可论证。对给定初始条件和域大小，α 较小时最易形成条纹（参见 Maini 等 1991 分岔图）。从无量纲化 (4.39) 看，这对应维数问题中趋化因子的慢产生或快速扩散/降解。直观理解：α 弱时细胞必须聚集在足够大的高密度区域才能产生足够陡的趋化因子梯度以招募细胞平衡 logistic 损失；要存在陡梯度则需邻近低细胞密度区（趋化因子产生低），所以"许多细胞+邻近少细胞区"对应的条纹比点更易在 α 小时出现。α 快时趋化因子梯度不必很陡；少数细胞即可产生足够趋化因子招募其他细胞形成集群——更可能是点。
考虑改变化学趋化因子 α 有实验支持——第 5 章细菌图样工作中讨论过。蛇图样形成中若化学趋化确有作用，细胞的化学趋化响应完全可能在图样形成过程中改变。
一旦稳态图样形成，增大 α 可使现有图样更锐利或在某些情况下引起定性变化（见下）。未发现条纹分裂为点，但有些情形条纹分裂为二（Maini et al. 1991, Myerscough et al. 1990）：这发生在无有丝分裂情况下。对已建立的条纹，α 增加给出很清晰锐利的带。Pseudonaja modesta 环纹褐蛇（图 4.23(c)）即见这种孤立窄色素带。反应扩散机制要产生如此清晰带需极精细调节色素产生阈值；而化学趋化机制只需 α/D 大即可。对点图样增大 α 通常使集群更聚焦聚焦，但有一种情形出现定性变化。
该看似简单的非线性模型系统的图样潜力不限于条纹与规则点。例如取规则点图样，逐渐增大 α 而不改变域形状，图样转变——α 足够大时出现成对点（图 4.26(a)）。该图样见于 leopard snake（Elaphe situla，图 4.26(b)）；该物种也有单点而非成对点的阶段。
图样发展的关键侧面可能是图样形成过程中生长引起的皮肤域变化。这正是本章讨论蛇图样的原因。正如所预期，域形状变化也产生更复杂图样。侧条纹域的单纯纵向生长会形成额外条纹插入既有条纹之间（图 4.27, 4.28）。第一种情形（图 4.27(a)）：侧生长使非对称点图样变对称——色素细胞聚集移至域中央。该类居中点图样非常常见，例如 corn snake（Elaphe guttata，图 4.27(b)）和多种 Vipera。从稍微不同的点图样出发，侧生长可产生菱形图样（图 4.28(a)(i)）；域稍变窄则产生波状条纹（图 4.28(a)(ii)）。菱形图样是多种响尾蛇的特征，如 eastern diamondback rattlesnake（Crotalus adamanteus，图 4.28(b)）；horseshoe snake（Coluber hippocrepis）也见此图样。尾部附近身体变窄，菱形图样可能变为波状条纹——第 3 章讨论的发育约束示例。这些结果说明域生长在色素图样铺设中可能起关键作用。
化学趋化作为合理初始机制的选定基于其在其他发育系统中的重要性——第 5 章细菌图样中其直接生物学相关性毫无疑问。Murray 等（1990）鳄鱼工作也提供模型涉及细胞的其他证据。
虽然前几章经验预期多种复杂图样，Murray–Myerscough（1991）还发现其他意外图样，如成对点和蛇身上常见的波状条纹、菱形图样。许多蛇和蜥蜴体表图样延伸至尾尖，几乎不变，即便尾部明显变细。图 4.27(b) 即一例，对比许多哺乳类（尾部点图样在域变细时常变为侧条纹）。可能强化学趋化的聚集效应使点即便在变细域也能形成，但需进一步数值计算确认。最终若域在图样形成过程中足够细，所有点都将变为条纹。蛇的渐缩度通常比猎豹、豹（尾部尖变条纹）更平缓。
实际蛇图样与化学趋化模型产生的复杂图样之间的现象学相似性令人鼓舞，也为进一步理论研究和实验研究（化学趋化、细胞密度、域大小在色素图样形成中的可能作用）提供动力。当然这一细胞-化学趋化模型机制对色素图样形成如同前述模型一样是推测性的；与纯反应扩散模型不同，它显式包含了细胞运动和细胞间通过化学介质（已有相当生物学证据）的相互作用。该相对简单的模型能再现多种蛇身上观察到的简单和复杂色素图样。
4.14 细胞-化学趋化系统的传播型图样生成（Propagating Pattern Generation with the Cell-Chemotaxis System）
鳄鱼条纹情形（第 4.2 节）观察到条纹图样从头部开始向背部和尾部渐进传播。第 4.12 节给出在 (4.30) 合适参数空间中细胞-化学趋化系统能生成多种空间图样，由特定初始条件产生。直观上若在域一端施加稳态扰动，该区域会启动空间图样解并沿域渐进传播。Myerscough–Murray（1992）在无有丝分裂（(4.22) 中 r=0）情形下研究该问题：
∂n/∂t = D ∇²n - α ∇·(n ∇c) ∂c/∂t = ∇²c + s (n/(1+n) - c). (4.31)
在能产生规则条纹的参数空间中，初始一端空间扰动引发传播行波，行波过后留下规则条纹稳态空间图样。图 4.2 显示该类图样的典型演化。可用第 2 章 2.6 节方法（渐近技术）对图样波长和传播速度作解析估计。
考虑一维问题、关注扰动前沿附近图样刚开始形成的行为。前沿处扰动幅度小，由线性化方程 (4.32)（(4.24) 中 r=0）控制，色散关系 (4.34) 给 λ² + [(D+1)k² + s]λ + [D k²(k² + s) - s N α/(1+N)²]k² = 0。前沿附近线性近似下解为 Fourier 模式积分 (4.35)；n(Vt, t) = ∫ A(k) exp[t g(k)] dk, (4.36) g(k) = i k V + λ(k²)，V 为有限"传播速度"、即扰动前沿行进速度。用最陡下降法 (Murray 1984) 在大 t 即大 x (V=O(1)) 下渐近求积分得 (4.37) n ≈ F(k) / √t · exp[t (i k V + λ(k²))]。k 为 g(k) 在复 k 平面上使实部最大的鞍点，由 dg/dk = i V + dλ/dk = 0 (4.38) 确定。
由图 4.2 等数值解看，远离初始扰动处图样包络形状基本不变，故对 g(k) 加约束 Re[g(k)] = Re[i k V + λ(k²)] = 0 (4.39) 即边际稳定性假设，与 (4.38) 联合（原则上）可解出波数 k 和图样行波速度 V。在行波包络运动坐标系中，前沿之后从包络前沿开始出现振荡图样，振荡频率 ω = Im[g(k)] (4.40)。若峰保守（不合并），ω 也是远离前沿处的振荡频率。设 k' 为远在包络前沿之后图样的波数，则 ω = k' V，即 k' = Re(k) + Im[λ(k²)] / V (4.41)，给出稳态图样波数。
该方法依赖前沿附近线性行为。弱非线性假设下解写成 Fourier 模式积分 (4.35)。但对远离均匀稳态的非线性系统 (4.31)，弱非线性假设未必成立；隐含假设解的行为由前沿事件控制，亦未必严格成立。因此该方法不太可能与非线性系统 (4.31) 数值解所得波长 w = 2π/k' 和 V 在数量上良好吻合；可期望在定性上一致，用于预测 V 和 w。
解析行波图样。将该方法应用到 (4.31) 配上色散关系 (4.34) 不易。将 λ(k²) 代入 (4.38)(4.39) 所得方程非常复杂，无法解析提取所需信息。Myerscough–Murray（1992）转而考虑 caricature 色散关系 λ(k²) = ε k²(k²_0 - k²), ε>0，对精确色散关系 (4.34) 在 λ 正区域作最佳拟合。ε, k²_0 刻画抛物近似的陡度和在 k² 轴上的截距。Myerscough–Murray（1992）完成了详细最陡下降分析，并把近似结果与 (4.31) 全非线性系统数值解比较。比较并不很接近，但这种解析程序提供了对参数的定性依赖。只有两个参数，难期望极好拟合，部分误差源自拟合本身——但有趣的是误差并非来自 caricature 色散关系的拟合而是线性化过程。
该近似方法隐含假设传播速度和最终图样波长由前沿处（线性理论适用）事件决定。比较结果显示：前沿事件虽对 w, V 重要，前沿之后非线性区事件（化学趋化情形）也有显著影响。Myerscough–Murray（1992）给出可能解释。
前沿之后正在发育的细胞密度峰产生趋化因子浓度峰；该浓度梯度作用于该峰前沿侧的细胞，抑制它们向前进入在包络前沿正在发育的下一个峰。线性情形下可能加入新峰的细胞，要么加入之前形成的细胞集群，要么被较慢地招募到前沿峰。第一种解释了为什么解波长比解析方法得到的更长；第二种解释了为什么解的传播速度比解析更慢。这些化学趋化细胞招募效应在反应扩散系统中不存在。尽管如此，对强非线性 (4.31) 该解析方法仍可作为传播速度和最终图样波长（前提是峰不合并）的有用估计。所需信息仅为线性系统色散关系的合理近似。
本章个人批注
本章最大特征是它把"生长"作为图样形成机制的核心变量——既不简化为预图样，也不简化为固定域反应扩散，而是把生长作为模型方程中的一个显式项（-rc 项稀释、时间依赖扩散系数）。Murray 在第 3 章把固定域做到了精细，本章则示范如何把这种精细在生长情形下重新建立。牙原基模型 (4.6)–(4.8) 清楚地说明：把生长剥离后系统不再给出正确序列——这是模型的有机成分而非装饰。
鳄鱼条纹案例给我最深的启示是"机制无关结论"。Murray 一再强调：决定条纹数的是激活时胚胎长度而非性别或基因；具体机制（反应扩散、化学趋化、机械模型）都能给出相同定性预测，关键是参数选择和域几何。换句话说，本书反复出现的"扩散驱动失稳"+"最小域长"框架是普适的——它对机制的具体生物化学身份保持中性。这一点对生物学家是有用的提醒：不要过早投入某一具体机制。
鳄鱼与蛇两套案例的对比也耐人寻味。鳄鱼实验数据极丰富（Westergaard–Ferguson 系列）、胚胎可操控，所以可以构造细节很多的牙原基模型 (4.6)–(4.8)，参数估计半数来自实验、半数从解析模式反推；而蛇几乎完全缺胚胎学数据，所以模型只能停留在定性预测（哪种参数给出哪种图样）。这导致本章两半风格截然不同：鳄鱼部分更像参数估计+数值验证；蛇部分更像计算实验+图样画廊。Murray 在第 4.13 节反复使用"computed ... pattern"对照"observed ... snake"——这是一种"先看模型能生成什么再去找生物学对应"的逆向研究范式。
第 4.14 节对行波传播的渐近分析（Myerscough–Murray 1992）的坦白陈述值得注意："这种方法对非线性系统 (4.31) 不会给出与数值在数量上良好吻合的结果……误差并非来自 caricature 色散关系的拟合而是来自线性化过程"——这是少见的对方法局限性的诚实评价。化学趋化峰值处的非线性反馈（峰前侧被抑制、招募放缓）让前沿事件和后沿事件的耦合变得不可忽略。这种坦白在科学写作中难得——多数教科书只展示拟合结果。
关于预测实验（第 4.9 节）：dental determinant 移植到 x=0.9 导致后续 7 个原基全部在移植位置之后形成——模型给出了一个非常具体且反直觉的预测。这种预测只有在模型与生物学认真耦合时才会出现。Murray 强调这些预测"原则上可在实验室验证"——尽管他从未做这种实验（鳄鱼保护法严格、牙原基手术窗口窄），但模型作为假说生成器的价值已经体现。
读完本章我有两个疑问需要在 questions.md 中记录：（a）Murray 第 4.4 节第 4 段提到"the model mechanism is clearly quite robust in that the parameter domain for the first seven teeth is not small"——具体说，6 维参数空间里 mode 2 模式区域的体积（或者任意子空间投影的面积）是多少？书中只给了 γ-d、p-δ 截面，未给出体积估计。（b）Painter 等（1999）的 angelfish 模型（被本章引用为同时包含扩散、化学趋化、且由独立反应扩散系统提供形态发生素）是否被 Murray 推广到蛇身上？第 4.12 节仅使用 Oster–Murray (1989) 的 logistic-趋化模型；没看到 Murray 把 Painter 模型（多尺度耦合）应用到蛇。如果 Painter 模型的"形态发生素调节趋化"在蛇身上也成立，可能解释更多图样细节。
与上下章的衔接（一段话）
第 3 章（Pattern Formation with Reaction Diffusion Systems）建立了反应扩散机制的全套数学工具——线性化、色散关系、Turing 空间、分岔、域长依赖——并在动物毛皮、蝶翅等固定域问题上展示了该机制如何解释图样。本章则把视角转向生长域：条纹、牙原基、蛇皮等图样都在胚胎生长同时被铺设，所以"机制作用时间 << 生长时间"这一假设（适用于第 3 章多数固定域情形）不再成立。Murray 通过两组具体生物学问题（鳄鱼条纹、鳄鱼牙原基、蛇皮）展示三种不同的"加入生长"方式：鳄鱼条纹用"行波在长细域上传播"（机制快）；鳄鱼牙原基用"-rc 稀释项+时间依赖扩散系数"（机制慢、与生长耦合）；蛇皮用"参数选择+域形状变化驱动分岔"（机制不一定慢但生长改变域几何）。第 5 章（细菌图样与趋化性）会把细胞-化学趋化机制放在实验数据最完备的细菌系统里深入分析，并把蛇皮模型的具体参数选择置于生物化学实验下检验。第 6、7 章则回到组织相互作用、肢体发育等更复杂的形态发生问题。本章在固定域（ch3）与生长域（ch5+）之间起到承上启下的作用——既把第 3 章的数学工具在生长情形下重塑，又为后续章节的具体生物学应用埋下伏笔。