藥物經(jīng)濟(jì)學(xué)評價中Markov模型的周期內(nèi)校正方法探討

沃田 陳磊 席曉宇

摘 要 目的：為減少藥物經(jīng)濟(jì)學(xué)研究中Markov模型的誤差提供參考。方法：參考國外相關(guān)文獻(xiàn)，解釋Markov模型誤差產(chǎn)生的原理，并介紹常用的半周期法、梯形法、Simpsons 1/3 法、Simpsons 3/8法、生命表法及其在Excel、TreeAge軟件中的實現(xiàn)方法。結(jié)果與結(jié)論：Markov模型因離散化過程而產(chǎn)生誤差，使用周期內(nèi)校正方法可校正這一誤差。半周期法是最常用的校正方法，是通過增加第1周期一半的結(jié)果并減去最后1周期一半的結(jié)果進(jìn)行校正;梯形法的校正是以區(qū)間首尾端點值的均值代表該區(qū)間，以直角梯形的面積作為累計的結(jié)果;Simpsons 1/3法和Simpsons 3/8法則在梯形法的基礎(chǔ)上在區(qū)間中插值并以插值點所在的連續(xù)曲線代表整個曲線;壽命表法的校正中，生存人年數(shù)為該年齡組的組距與該年齡組和下一年齡組尚存人數(shù)的均值的乘積，生存人年數(shù)則為各年齡組生存總?cè)四陻?shù)的總和。在Excel中，可根據(jù)方法原理設(shè)置函數(shù)來實現(xiàn)校正;而在TreeAge軟件中，可通過設(shè)置Init Rwd、Incr Rwd、Final Rwd的函數(shù)表達(dá)式來實現(xiàn)校正。在使用Markov模型進(jìn)行藥物經(jīng)濟(jì)學(xué)評價時，對于周期校正方法的選擇，若從建模的易用性、適用情形的廣泛性等角度出發(fā)，建議使用梯形法;從結(jié)果的精確性角度出發(fā)，則建議使用Simpsons 1/3 法，以提高M(jìn)arkov模型的準(zhǔn)確度。

關(guān)鍵詞 藥物經(jīng)濟(jì)學(xué);Markov模型;周期內(nèi)校正;半周期法;梯形法;Simpsons 1/3 法;Simpsons 3/8法;壽命表法

ABSTRACT? ?OBJECTIVE： To provide reference for reducing the error of Markov model in pharmacoeconomic evaluation. METHODS： Referring to foreign literatures， the errors in Markov model were explained. The commonly used correction methods were introduced： half-cycle correction， trapezoidal rule， Simpsons 1/3 rule， Simpsons 3/8 method and life table method and their implementation in Excel and TreeAge software. RESULTS & CONCLUSIONS： The error of Markov model was caused by the discretization process and could be corrected by the within-cycle correction methods. Half-cycle correction， was the most commonly used correction method and could be corrected by adding half of the results of the first cycle and subtracting half of the results of the last cycle. By trapezoidal rule， the interval was represented by the mean value of the first end point value of the interval， and the area of right trapezoid was the cummulative result. Simpsons 1/3 rule and Simpsons 3/8 rule took another point in the interval on the basis of trapezoidal correction and the continuous curve where these three points were represented the whole curve. By life table method， the person-years of survival was the product of the group distance of the age group and the mean value of the number of survivors of the age group and the next age group， and the total person-years of survival was the sum of person-years of survival of each age group. In Excel， the fuction was set according to the method principle to realize the correction; in TreeAge software， the function expression of Init Rwd， Inor Rwd and Final Rwd were set to realize the correction. When using Markov model for pharmacoeconomic evaluation， if the selection of periodic correction method is based on the ease of modeling and the universality of application， the trapezoid rule method is recommended; based on the accuracy of results， it is suggested to use Simpsons 1/3 method， so as to improve the accuracy of Markov model.

KEYWORDS? ?Pharmacoeconomics; Markov model; Within-cycle correction; Semi-cycle correction; Trapezoidal rule; Simpsons 1/3 rule; Simpsons 3/8 rule; Life table method

藥物經(jīng)濟(jì)學(xué)（Pharmacoeconomics）是為合理配置健康資源而新興的一門交叉學(xué)科，在近幾年由國家醫(yī)保部門組織的藥品談判中發(fā)揮了重要作用。藥物經(jīng)濟(jì)學(xué)可分為基于臨床數(shù)據(jù)的研究和基于模型的研究：其中，基于臨床數(shù)據(jù)的研究通過收集臨床數(shù)據(jù)進(jìn)行分析以判斷藥物的已知療效和經(jīng)濟(jì)性;基于模型的研究則通過二手?jǐn)?shù)據(jù)、專家訪談等構(gòu)建模型分析預(yù)測藥物的療效和經(jīng)濟(jì)性。目前藥物經(jīng)濟(jì)學(xué)應(yīng)用的主流模型包括：模擬群體層面病程短且復(fù)雜程度低的靜態(tài)決策樹模型、模擬群體層面病程長且復(fù)雜程度高的動態(tài)Markov模型、模擬個體層面疾病進(jìn)展的離散事件模型、基于生存分析模擬群體層面疾病進(jìn)展的分區(qū)生存模型以及模擬傳染病傳播的動態(tài)傳染模型[1]。其中，Markov模型因其能模擬復(fù)雜慢性疾病（如癌癥、慢性腎病、糖尿病等）的多階段多狀態(tài)進(jìn)程，得到了學(xué)界的認(rèn)可。該模型在20世紀(jì)80年代就已應(yīng)用于藥物經(jīng)濟(jì)學(xué)評價中，現(xiàn)已成為復(fù)雜慢性病藥物經(jīng)濟(jì)學(xué)研究建模最常用的方法之一。Markov模型的基本組成包括Markov健康狀態(tài)、循環(huán)周期、轉(zhuǎn)移概率、健康產(chǎn)出及成本，其基本原理是將疾病發(fā)生發(fā)展的進(jìn)程劃分為多個互斥的Markov健康狀態(tài)，每個循環(huán)周期內(nèi)各狀態(tài)間按照轉(zhuǎn)移概率相互轉(zhuǎn)移，當(dāng)絕大多數(shù)研究對象死亡或完成預(yù)設(shè)循環(huán)次數(shù)時模型停止運(yùn)行，從而得出不同健康狀態(tài)對應(yīng)的健康產(chǎn)出和成本[1]。

Markov模型本質(zhì)上是將連續(xù)的復(fù)雜的疾病進(jìn)程離散化，從而實現(xiàn)對慢性疾病進(jìn)程的模擬，然而這一離散化的過程也同時造成了健康產(chǎn)出和成本計算的誤差。周期內(nèi)校正（Within-cycle correction）可用于校正該誤差，且已被廣泛應(yīng)用于國外相關(guān)研究中[2]，但在國內(nèi)相關(guān)文獻(xiàn)中尚未見明確提及。本研究對Markov模型中健康產(chǎn)出和成本誤差的產(chǎn)生和校正方法進(jìn)行介紹和分析，并給出誤差校正方法在Excel軟件和TreeAge軟件中的具體實現(xiàn)方法，旨在為減少M(fèi)arkov模型的誤差提供參考。

1 Markov模型周期內(nèi)誤差產(chǎn)生原理

慢性疾病的進(jìn)程是連續(xù)的，即每時每刻都有可能發(fā)生不同狀態(tài)間的轉(zhuǎn)移。以標(biāo)準(zhǔn)“健康-患病-死亡”三狀態(tài)的Markov模型為例，患者由健康人群按照患病概率轉(zhuǎn)移而來，同時患者還會按照死亡概率轉(zhuǎn)移至死亡狀態(tài)。隨著疾病的進(jìn)展，患者人數(shù)將會先增加后減少。若患病概率和死亡概率均已知，則可列出微分方程求解當(dāng)前患者數(shù)量占總?cè)藬?shù)的比例關(guān)于時間（x）的函數(shù)f（x），該函數(shù)以時間為橫坐標(biāo)、患病人數(shù)比例為縱坐標(biāo)，其曲線下面積即為患者在該狀態(tài)下的平均存活時長。這一連續(xù)時間連續(xù)狀態(tài)的建模方法即為隨機(jī)微分方程模型（SDE）。SDE的優(yōu)勢在于其是“金標(biāo)準(zhǔn)”，結(jié)果最為準(zhǔn)確，但是其易用性較差，難于解釋，不夠直觀[3]。

Markov模型將時間劃分為多個周期（t），假定患者人數(shù)在周期內(nèi)保持不變，下一周期的人數(shù)取決于上一周期所有狀態(tài)的人數(shù)與對應(yīng)的轉(zhuǎn)移概率。以圖1為例，第1周期開始時所有研究對象均處于健康狀態(tài)，到第1周期結(jié)束時，部分研究對象患病進(jìn)入疾病狀態(tài)，則第1周期的期初患病人數(shù)為0，第1周期的期末患病人數(shù)為f（t）。SDE計算第10周期平均存活時長的方法是以積分求出f（0）與f（10t）圍成的曲線的面積;而Markov模型計算第1周期平均存活時長時，若設(shè)定狀態(tài)轉(zhuǎn)移發(fā)生在周期結(jié)束時，則如圖1所示按照期初值f（0）作為第1周期下的患病人數(shù)，以f（0）為矩形的長、周期長度為矩形的寬計算矩形的面積;若設(shè)定狀態(tài)轉(zhuǎn)移發(fā)生在周期開始時，則如圖2所示按照期末值f（t）作為第1周期下的患病人數(shù)，以f（t）為矩形的長、周期長度為矩形的寬計算矩形的面積。由圖1、圖2可以看出，Markov模型在第1周期即出現(xiàn)了明顯的偏差，在后續(xù)的周期中誤差依然存在，多周期下f（x）的曲線與矩形相差的面積即為Markov模型離散化過程中產(chǎn)生的誤差，其誤差的方向如表1所示[3]。

2 常見的周期內(nèi)校正方法

目前，國際上已有多種成熟的周期內(nèi)校正方法，包括半周期法（Half-cycle correction）、梯形法（Trapezoidal rule）、Simpsons 1/3 法（Simpsons 1/3 rule）、Simpsons 3/8法（Simpsons 3/8 rule）、壽命表法（Life table method），并已有廣泛應(yīng)用[4-6]。

2.1 半周期法

半周期法是目前最常用的周期內(nèi)校正方法，國際藥物經(jīng)濟(jì)學(xué)與結(jié)果研究學(xué)會醫(yī)學(xué)決策制定學(xué)會下屬的建模研究規(guī)范工作組（ISPOR-SMDM Modeling Good Research Practices Task Force）建議使用半周期法校正Markov模型的臨床產(chǎn)出和成本[7]。半周期法假設(shè)狀態(tài)的轉(zhuǎn)變發(fā)生在每一循環(huán)周期的中點（如圖3所示）[8]。在計算累計結(jié)果時，通過增加第1周期一半的結(jié)果并減去最后1周期一半的結(jié)果，即可得出校正后的累計結(jié)果。例如，未校正的總質(zhì)量調(diào)整生命年（QALYs）為9.90年，其中第1周期和最后1周期中的QALYs分別為1和0.32年。半周期法的校正運(yùn)算過程即為：9.90+1/2-0.32/2≈10.24年。

盡管半周期法操作簡單、應(yīng)用廣泛，但是其原理和計算過程使得其僅能基于首末周期校正，若周期長度不一致、f（t）分段不連續(xù)，則無法校正。

2.2 梯形法

為了解決半周期法僅能通過首末周期校正的缺陷，可應(yīng)用數(shù)值積分中的梯形法。其基本原理是將f（t）每一周期的端點值連線，由此將標(biāo)準(zhǔn)Markov模型中未經(jīng)校正的矩形變成梯形，將折線下的面積作為累計的結(jié)果（如圖4所示）[3]。若f（t）為凸函數(shù)，則梯形法校正的結(jié)果將略低于曲線下面積;若f（t）為凹函數(shù)，則梯形法計算的結(jié)果將略高于曲線下面積。在計算健康產(chǎn)出和成本時，梯形法對Markov模型的每一周期都進(jìn)行了調(diào)整，即使出現(xiàn)周期長度不同、f（t）分段不連續(xù)的情況也能予以校正，克服了半周期法僅能基于首末周期校正的缺陷，能得到相較半周期法更精準(zhǔn)的結(jié)果。

2.3 Simpsons 1/3 法和Simpsons 3/8法

數(shù)值積分中，Simpsons 1/3法和Simpsons 3/8法是相較梯形法精度更高的計算方法，也可用于周期內(nèi)校正。數(shù)值積分是將連續(xù)的定積分離散為不同特定的形狀，諸如矩形、梯形、曲邊梯形等，通過計算這些形狀的面積從而獲得近似的積分值，因其能計算難以獲得積分函數(shù)或難以計算原函數(shù)的積分，故在計算機(jī)領(lǐng)域被廣泛應(yīng)用。數(shù)值積分的基本方法是使用插值法，在特定區(qū)間內(nèi)選擇1組點以代表該區(qū)間的函數(shù)值取值情況，其構(gòu)造的基本插值型求積公式包括Simpsons 1/3公式和Simpsons 3/8公式[3]。

Simpsons 1/3公式基于兩次插值推導(dǎo)，能擬合較為平滑的曲線的定積分，設(shè)Markov模型單位同期為n，k為任一正整數(shù)，當(dāng)k取1，2…，T時，tk=n，2n，…nT，即總周期為nT。函數(shù)f（x）在區(qū)間[tk-1，tk+1]中的定積分則約等于通過（tk-1，f（tk-1））、（tk，f（tk））、（tk+1，f（tk+1））的曲線與x軸圍成的曲邊梯形面積，函數(shù)f（x）在區(qū)間[tk-1，tk+1]圍成的面積Ck為：

Simpsons 3/8公式基于3次插值推導(dǎo)，能更精確地擬合較為平滑的曲線的定積分，函數(shù)f（x）在區(qū)間[tk，tk+3]（k=1，2，3，…，T）中的定積分則約等于通過（tk，f（tk））、（tk+1，f（tk+1））、（tk+2，f（tk+2））、（tk+3，f（tk+3））的曲線與x軸圍成的曲邊梯形面積，規(guī)定Markov模型中單位周期為n，函數(shù)f（x）在區(qū)間[tk，tk+3]圍成的面積Ck為：

2.4 壽命表法

壽命表法借鑒了保險精算和生存分析中常用的壽命表。壽命表由年齡組、死亡概率（qx）、尚存人數(shù)（lx）、死亡人數(shù)（dx）、生存人年數(shù)（Lx）和總生存人年數(shù)（Tx）構(gòu)成。某一年齡組的尚存人數(shù)與該年齡組的死亡概率的乘積即為這一年齡組的死亡人數(shù)（dx=lx×qx），下一年齡組的尚存人數(shù)即為上一年齡組的尚存人數(shù)減去死亡人數(shù)（lx+1= lx-dx），生存人年數(shù)取該年齡組的組距t與該年齡組和下一年齡組尚存人數(shù)的均值的乘積[Lx=t×（lx+ lx+1）/2]，生存總?cè)四陻?shù)即為所有年齡組生存人年數(shù)的總和。

表2為以上述原理按照1 000人規(guī)模且下年死亡率為0.52時編制的壽命表，依照壽命表計算的總生存人年數(shù)與依照半周期法計算的總生存人年數(shù)基本一致，表明其精度尚可。實際應(yīng)用時，年齡組即為Markov模型中的循環(huán)周期，死亡概率可類比為轉(zhuǎn)移概率;另外，在實質(zhì)上，壽命表法和梯形法原理完全一致，可以認(rèn)為二者等價[9]。

3 Markov模型周期內(nèi)校正的實現(xiàn)方法

Markov模型一般通過Excel或者TreeAge軟件建立，為了方便敘述，本文僅展示考慮貼現(xiàn)情形下的成本校正。成本關(guān)于周期和貼現(xiàn)率的公式為：

式中，t為周期，St為t周期后處于疾病狀態(tài)的人數(shù)，c表示單位周期內(nèi)的治療費(fèi)用，rD表示年折現(xiàn)率?？紤]到壽命表法和梯形法在數(shù)學(xué)運(yùn)算上完全一致[2]，本部分不予以重復(fù)介紹。

3.1 基于Excel軟件的周期內(nèi)校正實現(xiàn)方法

半周期法：若構(gòu)建的模型假設(shè)狀態(tài)的轉(zhuǎn)變發(fā)生在每一循環(huán)周期開始時，則校正方式為在計算總成本加上（C0-CT）/2;若發(fā)生在每一循環(huán)周期的結(jié)束時，則減去（C0-CT））/2。式中，T為總周期數(shù)，C0為第1周期的費(fèi)用，CT為最后1周期的費(fèi)用。

梯形法和壽命表法：梯形法在每一周期期初、期末兩端點下的均值代表這一周期，即建立函數(shù)ft=（Ct-1+Ct）/2，最后對所有的ft求和即可得到總費(fèi)用。

3.2 基于TreeAge軟件的周期內(nèi)校正實現(xiàn)方法

TreeAge軟件中，Markov狀態(tài)的值統(tǒng)稱為Rwd，通過設(shè)置Init Rwd、Incr Rwd、Final Rwd的函數(shù)表達(dá)式即可實現(xiàn)對成本、效果的校正。TreeAge軟件內(nèi)置的“半周期法”校正其原理同梯形法完全一致，因此本部分將二者合并介紹。

半周期法和梯形法：初始值Init Rwd=C0/2，期末值Final Rwd=CT/2，每一期增加的值Incr Rwd=discount（cost;rD;_stage）;式中，discount是TreeAge軟件內(nèi)置的貼現(xiàn)函數(shù)，cost是給定狀態(tài)下的每周期成本，_stage是TreeAge軟件內(nèi)置的模型當(dāng)前運(yùn)行周期t的計數(shù)變量（下同）。

值得注意的是，上述介紹的實現(xiàn)方法中，還有一些技術(shù)細(xì)節(jié)值得討論：第一，周期數(shù)若為奇數(shù)時使用Simpsons 1/3法校正會造成不能被2整除的周期未納入校正，而周期數(shù)若不能被3整除時Simpsons 3/8法也會出現(xiàn)類似的情況。解決方法為對不能整除的周期采用梯形法校正，或者將周期數(shù)延長至可整除的周期。第二，上述校正方法實現(xiàn)的前提是每個周期的長度相同;若不相同，則每周期的結(jié)果要考慮周期長度的差異，因此需要將周期長度也作為模型的參數(shù)。

4 討論

本文介紹了Markov模型常見的周期內(nèi)校正方法，周期內(nèi)校正的目的是解決Markov模型離散化連續(xù)疾病進(jìn)程中存在的誤差問題。標(biāo)準(zhǔn)Markov模型運(yùn)算過程實則為誤差較大的以端點值代表單位區(qū)間取值的數(shù)值積分過程，而數(shù)值積分領(lǐng)域已有成熟的誤差較小的積分方法，應(yīng)用這些積分方法即可更為精確地求解Markov模型，即所謂的周期內(nèi)校正。上述幾種校正方法中，半周期法的中點假設(shè)實則為數(shù)值積分中的中點法，即以區(qū)間中的中點代表整個區(qū)間;梯形法則以區(qū)間的兩端點代表整個區(qū)間;Simpsons 1/3法則在梯形法的基礎(chǔ)上在區(qū)間中再取一點并以這3點所在的連續(xù)曲線代表整個曲線。因而上述方法的精度是逐漸遞增的，如圖5所示（圖中的曲邊為函數(shù)圖像，圓圈表示取值點，取值點的連線指該校正法所代表的整個區(qū)間）[2]。

5 結(jié)語

在這些周期內(nèi)校正方法中，半周期法應(yīng)用最早但誤差也最大;Simpsons 1/3 法和Simpsons 3/8法最為精確但計算略為復(fù)雜，且無法應(yīng)用于出現(xiàn)轉(zhuǎn)移概率不連續(xù)的情況[10]。因此，在方法選擇時，若從建模的易用性、適用情形的廣泛性（特別是模型存在不連續(xù)性時）等角度出發(fā)，建議使用梯形法;從結(jié)果的精確性角度出發(fā)，則建議使用Simpsons 1/3法[2];其他周期內(nèi)校正方法可在敏感性分析中呈現(xiàn)結(jié)果并予以討論。

參考文獻(xiàn)

[ 1 ] BECK JR，PAUKER SG. The Markov process in medical prognosis[J]. Med Decis Making，1983，3（4）：419-458.

[ 2 ] ELBASHA E，CHHATWAL J. Myths and misconceptions of within-cycle correction：a guide for modelers and decision makers[J]. Pharmacoeconomics，2016，34（1）：13-22.

[ 3 ] ELBASHA E，CHHATWAL J. Theoretical foundations and practical applications of within-cycle correction methods[J]. Med Decis Making，2016，36（1）：115-131.

[ 4 ] ELBASHA E，GREAVES W，ROTH D，et al. Cost-effectiveness of elbasvir/grazoprevir use in treatment-naive and treatment-experienced patients with hepatitis C virus genotype 1 infection and chronic kidney disease in the United States[J]. J Viral Hepatitis，2017，24（4）：268-279.

[ 5 ] ELBASHA E，ROBERTSON MN，NWANKWO C. The cost-effectiveness of testing for NS5a resistance-associated polymorphisms at baseline in genotype 1a-infected （treatment-naive and treatment-experienced）subjects tr- eated with all-oral elbasvir/grazoprevir regimens in the United States[J]. Aliment Pharmacol Ther，2017，45（3）：455-467.

[ 6 ] SAMUR S，KLEBANOFF M，BANKEN R，et al. Long- term clinical impact and cost-effectiveness of obeticholic acid for the treatment of primary biliary cholangitis[J]. Hepatology，2017，65（3）：920-928.

[ 7 ] WEINSTEIN M，OBRIEN B，HORNBERGER J，et al. Principles of good practice for decision analytic modeling in health-care evaluation：report of the ISPOR task force on good research practicesmodeling studies[J]. Value Health，2003，6（1）：9-17.

[ 8 ] SONNENBERG FA，BECK JR. Markov models in medical decision making：a practical guide[J]. Med Decis Making，1993，13（4）：322-338.

[ 9 ] BARENDREGT J. The half-cycle correction：banish rather than explain it[J]. Med Decis Making，2009，29（4）：500-502.

[10] PREZ-MART?N J，BERMEJO I，DEZ FJ. Evaluation of Markov models with discontinuities[J]. Med Decis Making，2019，39（4）：414-420.

（收稿日期：2019-10-12 修回日期：2020-03-06）

（編輯：孫 冰）