局部脑血流的测定

局部脑血流的测定

一. 问题简介

脑血流量是诊断和治疗脑梗塞,脑出血,动脉瘤和先天性动,静脉血管畸形等脑血管疾病的主要依据。测定脑血流量可为研究人脑在不同的病理和生理条件下的功能提供客观指标，它对研究脑循环药物的药理作用也很有帮助。所以人们长期致力于寻找有效地测定脑血流量的方法。

近年来出现了以放射性同位素作示踪剂测定人脑局部血流量的方法。这种方法大致可描述如下:由受试者吸入某种放射性同位素的气体，然后将探测器置于受试者头部某固定处，定时测量该处放射性同位素的计数率(简称计数率)，同时测量他呼出气的计数率。

由于动脉血将肺部的放射性同位素输送至大脑，使脑部同位素增加，而脑血流又将同位素带离，使同位素减少，实验证明由脑血流引起局部地区计数率下降的速率与当时该处的记数率成正比，其比例系数反映了该处的脑血流量，被称为脑血流量系数，只要确定该系数即可推算出脑血流量。动脉血从肺输送同位素至大脑引起脑部计数率上升的速率与当时呼出气的计数率成正比。

试建立确定脑血流系数的数学模型并计算上述受试者的脑血流系数。

备注:该题目是市(1990 年)大学生数学建模竞赛A 题。二.

模型的假定

= − 1. 脑部计数率(记为 h (t ) )的上升只与肺部的放射性同位素有关，上 升速度与呼出气的记数率(记为 p (t ) )成正比，比例系数记为 k ；

2. 脑部记数率 h (t ) 的下降只与该处脑血流量有关，其下降速度正比 于 h (t ) ，比例系数为脑血流系数，记为 K ，这里忽略了放射性元素的衰 变和其它因素；

3. 脑血流量在测定期间恒定，心脏博动，被测试者大脑活动，情感 波动等带来的变化忽略不予考虑；

4. 每次仪器测量为相互独立事件，各测量值无记忆相关；

5. 放射性同位素在人体传递是从吸入气体(含有放射物)开始的， 并假定一次吸入，因此认为同位素在肺中瞬时达到最大浓度；

6. 在吸入气体瞬时，脑中放射物记数率为零；

7. 脑血流量与脑血流系数 K 成单值函数关系，求得后者即可确定前 者。

三. 模型的建立与分析 由于已知脑局部同位素的增减与已测定的头部

记数率 h (t ) 和呼出气 记数率 p (t ) 成正比关系,于是很自然地确定以脑部同位素量,即脑部记数 率作为讨论对象.。

1.原始模型的建立 设某时刻 t ≥ 0 时,头部记数率为 h (t ) ,在 Δt 时段后记数率 h (t + Δt ) ， 由假定可知, 头部记数率的增量 Δh = h (t + Δt ) − h (t ) 仅与三个因素有关： (i) 肺动脉血将肺部的放射性同位素送至大脑,使脑部记数率增量 为 Δh 1 ；

(ii) 脑血流将同位素带离,脑记数率下降为 Δh 2 ；

(iii)放射性同位素自身有衰减引起记数率下降量为 Δh 3 ,设其半衰期 为τ .

因此,由医学试验及假定有

dh 1 = kp (t ), dh 2

= Kh (t ),

dh 3 ln 2

h (t ),

τ

dt

dt

dt

而

Δh (t ) = Δh 1 (t ) − Δh 2 (t ) + Δh 3 (t ) ，

于是

dh = dh 1 − dh 2 + dh 3 , dt dt dt dt

即

dh = − K h (t ) + kp (t ) − ln 2 h (t ) . (1)

dt τ

由于在测试时放射性同位素(如133 Xe )的半衰期τ 一般很大,而测试

时间又很短(大约十几分钟左右),由此假定τ → +∞ ,于是(1)式变为:

dh

= − K h (t ) + kp (t ) . (2) dt

2. 算法模型的建立与改进

在建立算法模型之前,首先必须对 p (t ) 进行预测。作 p (t ) ~ t 和

ln p (t ) ~ t 的离散图(图 1 和图 2),由此发现 p (t ) 与 t 有近似于 Ae − λt 的函 数关系。

通过对 ln p (t ) ~ t 的离散图 2 的观察,去掉时刻 6 及 6.5 以后的样本点 (这样作的原因见文章后面的评注),再利用最小二乘法进行拟合得

图

1 图2

ln p (t ) = 9.15158 − 1.47577t

其相关系数 r = 0.999887 ,由此得知 p (t ) = 9429.33e −1.47577t .

图3

我们作出p(t ) = 9429.33e −1.47577t 的图形,并将此图和图1 放在一起见图3,由图3 及相关系数r = 0.999887 可以认为p(t )确实是负指数曲线p(t ) = Ae −λt , (A = 9429.33, λ= 1.47577) .

由(2)及假设f，即h(0) = 0 ，解得

h(t ) =

kA

K −λ

(e −λt −e −Kt )，(3)

此式从数学上来看并不复杂，但要利用此式求出参数K 和k 却并非易事，而参数K 则需要在测试中使用，因此我们的问题归结为：如何利用实际测量值和(2)及h(0)=0 去决定参数k 和K ．这类数学问题称为参数辨识问题。

下面建立几个算法模型:

算法模型Ⅰ.一般差分拟合法：将方程(2)离散化，记时间步长为T，利用前插公式得：

h −h

n+1 n = −K h

n + kp n

T

即

h n+1 = (1 −KT )h n + kTp n ，(4)其中

h n = h(t0 + nT ), p n = p(t0 + nT ).

用差分法求解，其截断误差为o(T 2 ) ，显然大了些，为了提高精度和准确度，最直接的方法是由插值方法获得更多的结点，缩短步长，使截断误差减少。如用三次样条插值法在每两个结点的中点进行插值，可使截断误差减少到原来的1/4，但仍然为o(T 2 ) ，且继续缩短步长，计算