用吸收马尔可夫链估算药物在体内平均转运时间总论

用吸收马尔可夫链估算药物在体内的平均转运时间丁勇南京医科大学数学教研室，江苏南京（210029）摘要：本文提出了导出矩阵的概念，解决了药物转运速率常数转化为转移概率的问题，从而可用吸收马尔可夫链对药物在体内的转运过程进行研究。推导了用于估算药物在体内各部位平均转运时间的基本矩阵。关键词：导出矩阵，吸收马尔可夫链，基本矩阵，药物动力学中图分类号：O211.6Q3320引言近年来，随机过程在生命科学中的研究得到了越来越多的关注[1~2]。本文对吸收马尔可夫链在药动学房室模型中的应用进行了探讨，提出了导出矩阵的概念，解决了药物转运速率常数转化为转移概率的问题，为估算药物在体内转运过程的吸收、分布和消除的平均时间提供了新的方法。1导出矩阵定义设有矩阵[]nmijaA×=和任一常数δ，对整数ml≤≤1构造一个新矩阵[]nmijaA×=~~，其中：jljlaa,,~δ=ljnj≠=;,,2,1L∑∑∑∑≠==≠==−=−=nljjjlnjjlnljjjlnjjlllaaaaa1,1,1,1,,~~δjijiaa,,~=nlli,,1,1,,2,1LL+−=；nj,,2,1L=称A~是由A、δ和l导出的矩阵，记为]|[~lAAδ=；称此变换为导出变换。导出矩阵A~的元素除第l行外，都与A相同；第l行元素除主对角线元素外，与A只差一个常数因子δ；各行元素之和与A相同。对多个iδ和il，ti≤≤1，可对A进行多次导出，记为]|[]|][|[~2211ttlllAAδδδL=。容易证明，导出矩阵有如下的一些简单性质：（1）对于任意矩阵A，有AlA=]|1[；对0≠δ，有AllA=−]|][|[1δδ。（2）对于任意对角阵A，有AlA=]|[δ（3）对于任意n阶方阵A，有]|0[]2|0][1|0[nAL为对角阵，对角线元素为=∑=。（4）对任两个大小相同的矩阵A和B，有]|)[(]|[]|[lBAlBlAδδδ+=+。2吸收马尔可夫链典范式矩阵的性质设吸收马尔可夫链有ts+个状态，其中s个吸收状态，t个转移状态，其转移概率矩阵可表示为如下的典范式[2]：个转移状态个吸收状态tsQROEP⎥⎦⎤⎢⎣⎡=（1）上述各分块矩阵分别为：E为ss×的单位阵，因为从吸收状态到自身的转移概率为1，到其余状态的转移概率为0；O为ts×的零矩阵，因为从吸收状态到转移状态的概率为0；R为st×的矩阵，表示从转移状态到达吸收状态的概率；Q为tt×的矩阵，表示从一个转移状态到另一个（包括自身）转移状态的概率。定理1吸收马尔可夫链的典范式矩阵有如下性质[2]：（1）矩阵QE−可逆，这里的E为tt×的单位阵。（2）在t个转移状态中，从转移状态i到达吸收状态之前，进入转移状态j的平均步长为矩阵1)(−−QE的第i行、第j列元素之值；（3）在t个转移状态中，从转移状态i到达吸收状态之前，在所有转移状态之间传递的平均步长为矩阵1)(−−QE的第i行元素值之和。称1)(−−QE为吸收马尔可夫链的基本矩阵。由于导出矩阵不改变原矩阵各行元素之和，故对转移概率矩阵P和10δ，导出矩阵P~仍为转移概率矩阵，其直观意义可理解为：当系统的环境或某些条件有所改变时，从某个转移状态到其余各状态的转移概率发生了成比例的变化时，马尔可夫链的转移概率矩阵就从P改变为P~。对于矩阵（1）式，当sl≤≤1时，显然有PlP=]|[δ。定理2设由（1）导出的矩阵为⎥⎦⎤⎢⎣⎡=+++=QROEtsssPPt~~]|[]2|][1|[~21δδδL，其中，10iδ，ti≤≤1，则)~(−−QE=1)(−−QE⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡−−−1121100tδδδO（2）从而当δδ=i时，1)~(−−QE=1−δ1)(−−QE（3）证：设⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡=++++++++++++++++++tstsstsststsssssstsssssspppppppppQ,2,1,,22,21,2,12,11,1OML由于11,=∑+=+tsjjlsp，所以∑++≠=+++=−tslsjjjlslslspp1,,1，tl≤≤1，从而⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡−−−−−−=−∑∑∑−+=+++++++++≠=+++++++++≠=+11,2,1,,221,21,2,12,111,1tsjjtsstsststsstssjjjssstsssstssjjjspppppppppQEOML由定义⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡−−−=∑∑∑−+=+++++++++≠=+++++++++≠=+11,2,1,,2221,221,22,112,1111,11111~tsjjtstststststtsstssjjjssstsssstssjjjspppppppppQδδδδδδδδδOML从而⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡−−−−−−=−∑∑∑−+=+++++++++≠=+++++++++≠=+11,2,1,,2221,221,22,112,1111,11~tsjjtstststststtsstssjjjssstsssstssjjjspppppppppQEδδδδδδδδδOML)(0021QEt−⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=δδδO从而定理得证。3药动学房室模型的马尔可夫链基本矩阵在药物代谢动力学中，研究药物在体内转运过程最常用的模型为线性房室模型[3-4]，一般的n室模型如附图所示。该模型有一个中央室（图中的圆1），分别与其余的1−n个周边室（图中的圆2，……，圆n）相连通，任意两个周边室不直接连通，机体吸收的药物直接进入中央室，并且只从中央室消除。图中各药动学参数的意义为：ak为药物吸收速率常数，k为药物消除速率常数，其余jk1和1jk为中央室与周边室j之间的药物转运速率常数（nj,,3,2L=）。1nknk121kak12kk图1药动学房室模型示意图为了能用马尔可夫链对房室模型进行研究，首先将n个房室与n个状态一一对应，中央室对应状态1，其余1−n个周边室分别与状态2、3、……、n对应；另外再12n增加两个状态，状态0表示药物在吸收部位（如胃、肠等），状态n+1表示药物消除后被排出了体循环，故可视为马尔可夫链的吸收态。然后假设这些药动学参数满足条件：nikkkknjjia,,2,1;1,,0211L=+∑=（*）从而将药物吸收速率常数、消除速率常数和转运速率常数分别视为药物从一个状态进入另一个状态的转移概率，根据附图，可得如下的典范式的一步转移概率矩阵：}个转移状态个吸收状态状态状态状态状态111011000000100000100100001000000111313121211131221+⎪⎪⎪⎪⎭⎪⎪⎪⎪⎬⎫+⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡−−−−−−=∑=nnnkkkkkkkkkkkkkkPnnnniiaaMMOMML从而⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡−−−−−−+−=−∑=113131212111312210000000000000nnnniiaakkkkkkkkkkkkkQEOML为求1)(−−QE，对矩阵⎥⎦⎤⎢⎣⎡−EQE做如下一系列的初等列变换。（i）将第1、3、……、n列分别加到第2列上（ii）将第1列除以ak，将第2列除以k，将第j列除以1jk（nj≤≤3）（iii）将第2列分别乘以111−jjkk再加到第j列上（nj≤≤3）可得⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡+++=−−112131132112113113312112113121122111311321121131132112111011031101011)(nnnnnnnnnnakkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkQEOML（4）由定理1的性质（2）、（3），从上述基本矩阵可得药物在体内各部位的平均转运时间。例如由第一行、第一列的元素值可知，药物在吸收状态的平均转运时间为ak1，其值可认为是药物的平均吸收时间；由第二行、第二列的元素值可知，药物在中央室的平均转运时间为k1，由于药物只从中央室消除，故其值可认为是药物的平均消除时间；由第一行、第三列的元素值可知，药物从吸收状态进入周边室2的平均转运时间为2112kkk。由于开始时，药物都在吸收状态，故药物在体内分布的总平均转运时间为第一行元素值之和，即∑=++njjjakkkkk21111。4应用与讨论10名健康志愿者口服750mg国产舒他西林胶囊后，体内氨苄西林的药－时曲线符合二室模型，其药动学参数（均数）为[5]：1min0818.0−=ak，1min0231.0−=k，112min0188.0−=k和121min0143.0−=k，由此可求出基本矩阵（4）式为：⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡84.12629.43091.5629.43091.5629.4322.12由此可知，氨苄西林的平均吸收时间ak1为12.22min，平均消除时间k1为43.29min，在体内分布的总平均转运时间211211kkkkka++为112.43min，将近为2小时，这些数据为临床给药方案拟定提供了参考依据。在前面的讨论中，我们是将药物转运速率常数作为转移概率，从而要求这些参数满足条件（*），但在实际问题中，情况并非都是如此，例如药动学参数大于1的情况。如下方法可解决这一问题：由于药物转运速率常数的单位是时间倒数，故总可以采用较小的时间单位而将它们化为满足所需的条件（*）。如上例，若给出的参数以小时为单位，则1908.4−=hrka，1386.1−=hrk，112128.1−=hrk和121858.0−=hrk，其中有三个参数大于1，这时需要将时间单位进行调整，例如采用以分钟为单位，则这些参数就满足了所需的条件（*），如上例所示。但实际计算时，并不需要这么做，可直接进行计算。若直接以小时为单位的参数进行计算，我们由（4）式可求得氨苄西林的平均吸收时间ak1为0.204hr，平均消除时间k1为0.722hr，在体内总平均滞留时间211211kkkkka++为1.874hr，结果是一致的，仅是时间单位不同而已。不同时间单位的换算处理，实际上是将参数乘以同一常数，相当于前述定义中的导出变换，由定理2可知，其结果是再除以这个常数。如上例，以分钟为单位的矩阵看成P，以小时为单位的矩阵P~其实是在导出变换中取60=iδ，3,2,1=i，由定理2的（3）式可知，两者的结果只相差一个常数6011=−iδ，即时间换算单位，从而结果是一致的。由此可知，利用导出矩阵，基本矩阵（4）式的计算公式适用于任意的药动学参数，从而可用吸收马尔可夫链估算药物在体内各部位的平均转运时间。参考文献[1]孙荣恒编著．随机过程及其应用．北京：清华大学出版社[Ｍ]．2004．38-60[2]徐克学著．生物数学．北京：科学出版社[M]．1999．181-205[3]周怀梧编著．数理医药学．上海：上海科技出版社[M]．1983．98-139[4]Ｍ.Gibaldi，D.Perrier．Pharmacokinetics.Secondedition[M]．NewYork：MarcelDekker．1982．409-416[5]金锐，刘红梅，马霖，等．国产舒他西林胶囊的相对生物利用度．中国药学杂志[J]．1999；34（4）：256-258