基于两阶段随机仿真优化算法的体检顾客预约调度

刘 丹，耿 娜

（1.上海交通大学工业工程与管理系，上海200240；2.上海交通大学中美物流研究院，上海200230）

0 概述

随着社会经济的发展，我国人民健康意识增强，大众的健康观念从看病转向保健，从治病转向防病，有越来越多的人定期进行健康体检。据国家卫生健康委员会数据显示，2018年，我国体检人次数突破5.8 亿次，相较2012年的3.7 亿人次增长了56.8%。体检机构顾客排队现象日益严重，在上海的一家大型医院的体检中心进行调研发现，顾客在排队等待中耗费的时间甚至超过了在体检机构花费的总时长的80%。因此，采取合理的方式减少顾客的等待时间并提高资源利用率己经成为完善体检机构管理工作的当务之急。

医疗服务系统的预约调度可以预先分配医疗资源的服务能力，匹配供求关系，缓解系统拥堵，提高资源利用效率。目前我国大型公立医院已经普及预约就诊服务，然而对于体检机构等其他体疗服务机构的预约调度却没有得到广泛应用。顾客在进行体检时，除了空腹项目由于可用检查时间结束较早需要尽快完成外，多数体检项目没有硬性的顺序约束。目前，大多数体检机构不对顾客体检顺序进行安排，在实际体检服务过程中，由顾客自行选择检查项目的顺序。由于缺乏顾客预约调度及体检项目顺序的优化，因此顾客的等待时间过长，使顾客满意度下降。同时，对于体检机构而言，医学检查仪器是昂贵的成本构成。如何在不增加设备和人员投入的基础上，提升服务能力，减少顾客等待时间，提高医疗资源的利用效率，成为各体检机构面对的重要管理问题。

本文从预约调度及检查项目顺序调度2 个方面进行研究，考虑随机服务时长，提出了一种两阶段随机仿真优化算法，包括粗糙仿真评估阶段和精确仿真评估阶段。运用序优化（Ordinal Optimization，OO）策略，将基于亲和度评估的多种群遗传算法（Multipopulation Genetic Algorithm，MPGA）作为迭代策略，并在粗糙评估阶段中使用最优计算量分配

（Optimal Computational Budget Allocation，OCBA）

方法优化仿真预算分配来减少计算时间。

1 相关研究

目前关于医学检查系统中的顾客调度问题的研究大都仅考虑一个检查科室［1-2］，对多个检查科室串联的系统研究较少。BARON 等［3］研究了加拿大的一家预防保健服务中心的顾客调度问题，该中心的检查服务由10～20 个不同项目组成，使用了几种基本的调度策略来实时调度顾客，其目标有2 个层次：传统宏观层面的目标，如最小化总系统时间或最大限度减少延迟顾客总数；非传统的微观层面的目标，即减少流程中任何工作站的等待时间过长的事件的发生。刘阳等［4］考虑了2 类检查项目和2 类不同紧急程度的门诊患者的预约调度问题，建立了有限时域马尔科夫决策过程模型，利用启发式策略对门诊患者进行调度。HU 等［5］研究了门诊患者的4 项基本检查，为门诊检查制定了一些调度规则，以最大程度地减少总等待时间。CHERN 等［6］提出了一个基于最长处理时间优先调度策略的启发式模型来解决体检顾客项目顺序调度问题，包含20 多个不同的检查项目，目标是在考虑资源限制的情况下尽可能地减少顾客或医生的等待时间和整体检查时间。上述研究大都使用启发式调度策略，其调度的质量取决于所选择的策略。

本文考虑的随机服务时长下的体检顾客调度问题属于离散随机优化问题范畴。在解决离散随机优化问题方面，仿真优化方法［7-8］受到越来越多的关注。离散随机优化问题的解空间往往非常大，并且对解性能的准确评估也会非常耗时。OO 策略和OCBA 方法通过从大量解中抽样，并评估较少数量的解来克服解空间过大的难题［9］。多服务台串联的体检调度问题与多工作站串联的车间调度问题具有类似的特征。在车间调度问题中，已有许多文献将OO、OCBA 和全局搜索方法集成运用。SONG 等［10］使用基于OO 的进化策略来减少调度解评估的计算时间。HORNG 等［11］提出将进化策略嵌入到OO 中以解决作业车间调度问题，缩写为ESOO，其目的是最小化仓储费用和延期罚款。YANG 等［12］将OCBA技术嵌入到ESOO 算法的探索阶段，开发了用于解决作业车间调度问题的ESOO-OCBA 算法，以最小化提前和延迟成本。本文问题解空间更加复杂，对解空间搜索和抽样更困难，在进化策略的选择上，使用改进的MPGA 以改善算法搜索性能，在进化过程中维持种群染色体多样性［13］。车间调度问题的目标大都是最小化最大完工时间［14-15］、总完成时间［16］、总拖期［17］等，而在医疗服务系统中，顾客等待时间是影响服务质量的关键因素［18-19］。最常见的等待时间目标是顾客总等待时间［20］，这属于宏观层面的目标。此外，本文也考虑了微观层面的目标，即顾客在任一项目前的等待时间不超过某阈值。

2 问题描述

本文使用多人时间槽预约调度规则进行顾客预约调度，即对顾客进行分批次、等间隔预约安排，每批次预约人数相同，同时安排每位顾客的体检项目顺序。主要需要做出3 个决策：1）每批次预约人数；2）每个项目上顾客的排队顺序；3）每个顾客的检查项目顺序。

问题详细描述如下：将体检机构允许顾客进入的时间划分为长度相同的预约时间槽，每个预约时间槽内安排n位顾客，nmin≤n≤nmax。所有预约时间槽共可安排N位顾客，每位体检顾客需要进行m项检查，每位顾客必须访问每个检查项目1 次，顾客可以按任何顺序进行检查。1 项检查每次只能服务1 位顾客，在完成当前顾客的服务之前，不允许检查项目接受新顾客。顾客i在项目j上的服务时间为si，j，顾客完成所有检查后离开系统。制定以下假设：

假设1所有顾客是同质的，需要接受相同的检查。检查项目不存在任何顺序约束，也不存在完成时间约束。

假设2体检机构每日08：00-12：00 工作，但在11：00 后便不再允许新的顾客进入系统。假设计划工作时间为240 min，顾客必须在前180 min 进入系统。

假设3每个检查项目只有1 个服务台，服务时长是独立同分布的随机向量。

假设4预约顾客严格按照安排的时间到达系统。

对于顾客调度问题，顾客等待时间指标是优化的重点。据研究，体检机构中顾客在任一项目前等待超过某一时长时，会显著降低其满意度［3］。本文的优化目标同时考虑人均等待时间、系统超时加班时间、所有顾客两项检查之间的等待时间超过给定阈值TW的总时长、服务台平均空闲时间等。随着每批次调度人数的增多，系统服务人数增加，前3 个指标值会增大，服务台空闲时间指标值会减小，问题的本质是优化权衡这2 个部分指标，目标函数为其加权和，表示如下：

其中：i表示顾客索引，i∈｛1，2，…，N｝；j表示检查项目索引，j∈｛1，2，…，m｝；Fi表示顾客i完成服务并离开系统的时间；H表示体检机构每日计划工作时间；Wˉ表示所有顾客的平均等待时间；wi，j表示顾客i在项目j前的等待时间；TW表示顾客两项检查之间的等待时间阈值；Kˉ服务台平均空闲时间；γ1、γ2表示目标权重系数。

3 两阶段仿真优化算法

由于检查顺序不受限制，且服务时间随机，因此该问题的解空间非常大。运用仿真方法对调度解进行评估十分耗时，同时为了消除随机因素对系统的扰动，需要对调度解重复仿真多次。为了有效地缩小解空间并解决随机仿真耗时长的问题，本文给出了基于序优化策略的两阶段仿真优化算法框架，包含粗糙仿真评估阶段和精确仿真评估阶段，如图1所示。粗糙仿真评估阶段是仿真优化算法的关键，将MPGA 作为迭代搜索策略，使用OCBA 方法自适应地在搜索空间分配仿真次数。

图1 两阶段仿真优化算法结构Fig.1 Structure of two-stage simulation optimization algorithm

3.1 多种群遗传算法迭代策略

使用MPGA 从每一代各子种群选择优秀的候选解是迭代地进行OO 策略的一种方法。本文在传统MPGA 的基础上加入亲和度评价过程，以在迭代搜索过程中保持各子种群染色体的多样性。

3.1.1 染色体编码和解码

本文采用分段实数编码方式：1 个批次内顾客数量为n，nmin≤n≤nmax，m表示检查项目数，调度批次数目为d，每条染色体长度为d×nmax×m，染色体从左到右每n×m个基因代表1 个调度批次，1 个批次内每位顾客的每个项目用1 到n×m的实数唯一表示，不足的位用0 补齐；假设有2 项检查，2 个调度批次，每批次最多调度3 人，图2 给出了每批次调度2 位顾客的情况下生成的一条染色体。Oi，j表示顾客i在项目j的检查服务。

图2 2×3×2 问题染色体表示Fig.2 Chromosome representation of 2×3×2 problem

染色体编码完成后，需要对其进行解码以得到适应度值。对于同一批次内的顾客，采用半活动解码方法［21］，使每个项目在不延迟其他已解码的项目的情况下尽可能早地被加工。解码过程描述如下：

步骤1从染色体中提取1 个基因，求出其对应的Oi，j，服务时长为si，j。

步骤2按时间先后顺序遍历当前项目j上所有的空闲时间区间［Sj，Ej］。根据式（2）可以得到Oi，j的最早开始时间θi,j。其中ci，j-1表示顾客i上一个被解码的项目的结束时间。如果Oi，j是顾客i被安排的第1 项检查，ci，j-1=0：

步骤3假设Oi，j的开始时间为其最早开始时间θi,j，根据式（3）判断Oi，j是否可以插入到项目j的空闲时间内，若可以插入，则更新项目j的这一空闲区间的开始时刻Sj为(θi,j+si,j)；如果不满足式（3），则根据式（4）重新计算Oi，j的开始时间：

其中：Uj表示当前项目j上最后一位顾客完成检查的时刻。更新项目j的空闲区间。

步骤4若染色体中所有基因都已经安排，则结束程序，否则返回步骤1。

图3 为半活动解码方式的示意图。项目1 的空闲时间为［S1，E1］。若O2，1的开始时间等于S1，当O2，1的结束时间不大于E1时，则可将O2，1插入区间［S1，E1］，如图3（a）所示；相反，当O2，1的结束时间大于E1时，则不能将O2，1插入此时间区间，如图3（b）所示。此时需重新计算O2，1的开始时间。在项目1 上安排的最后一位为3 号顾客，其完成检查的时刻为U1，于是令O2，1的开始时间为U1。

图3 半活动解码示意图Fig.3 Schematic diagram of semi-active decoding

在不同批次间，若某项目上前一名顾客的服务超过当前顾客的预约到达时间，则当前顾客需要等待前一名服务结束后才能开始服务；如果一个项目前一批次所有服务结束时还未到当前顾客预约时间，在当前顾客到来之前服务台空闲。

3.1.2 初始种群

依据实际问题的需要，建立（nmax-nmin+1）个辅助子种群和1 个主种群，每个子种群对应1 种批次人数。对于每个子种群，依据最短队列优先和先到先服务的启发式规则生成2 个初始解。另外2 个初始解通过完全随机的方法生成。将这些初始解交叉、变异产生新的个体，新的个体再交叉、变异，最终生成大小为Z的初始种群。主种群又称为最优保存种群，用来存放子种群的最优解。

3.1.3 交叉和变异

变异算子为随机交换一条染色体上代表某1 批次的基因段上的2 个基因。由于不同初始子种群中的个体代表的批次人数不同，并且顾客不能重复检查同1 个项目，本文设计了基于批次和顾客的交叉方法。对于2 条亲本染色体，随机确定1 个调度批次，在2 条染色体代表这一调度批次内的基因段上进行交换。将1 条染色体基因段上代表同一顾客全部项目的基因，与另1 条染色体基因段上对应同一顾客全部项目的基因进行交换。交叉和变异的概率分别为pc和pm。

图4所示为包含2 项检查的交叉算子。亲本染色体P1 代表每批次预约3 人，P2 代表每批次2 人，交换第2 个批次内代表编号1 顾客项目的全部基因，得到子代C1 和C2。交叉和变异的概率分别为pc和pm。

图4 基于批次和顾客的交叉方法Fig.4 Crossover method based on batch and customer

3.1.4 选择机制

本文采用α+β选择机制。在每个子种群中，根据每个个体仿真获得适应度值均值的排名，选择最优的α个个体进行精确仿真评估。据观测，对一个体进行LS=105次仿真，可以得到十分稳定的观测值，将LS=105作为精确评估的仿真次数。其他β个个体通过轮盘赌选择获得。

3.1.5 基于亲和度评估的子种群交流方法

基于亲和度评估的子种群交流方法如图5所示。从每个子种群中选择那些可以提高算法性能的个体，被称为移民算子，在子种群之间进行交流，从而打破不同子种群之间的隔离机制，提高种群染色体多样性。在交流时采用单向环状拓扑，即1 个种群的交流对象为其后继种群，最后1 个种群的后继种群为第1 个种群。

图5 基于亲和度的子种群交流方法Fig.5 Population communication method based on affinity

在算法的每一代，子种群选出的α个个体精确仿真后，其最优解为该子种群这一代的最优解，作为移民算子。计算移民算子与目标种群这一代选择的（α+β）个个体的亲和度值。依据亲和度，评估目标种群中的（α+β）个个体与该移民算子的亲缘关系，选择关系最远且适应度值不如该移民算子的个体，用移民算子替换掉它。这样做使得每次交换为目标种群引入外来优势基因的同时，保留一部分自身的基因，外来优秀个体不会很快成为主导，能够保持种群的多样性。

假设有某2 条长为L的染色体：1 条染色体是X（x1，x2,…,xL）；另外1 条染色体是Y（y1，y2,…,yL），xi、yi分别代表X、Y染色体上的第i个位置对应的基因，则它们的亲和度为：

其中：maxi代表所使用的编码方式在这一位上可取到的最大值；mini代表在这一位上可取到的最小值。在本文中，maxi为m×nmax，mini为0。AAcev（X，Y）越大说明染色体X、Y之间的亲缘关系越近。

3.1.6 算法迭代终止条件

在每一代中，各子种群的最优解贡献给主种群。主种群的最优解为算法这一代的最优解。重复迭代优化过程，直到一定数量的连续代，genmax、最优适应度值均没有改进，算法停止。

3.2 仿真资源分配策略

在算法的粗糙评估阶段中，对于每一子种群中的个体，使用OCBA 方法根据个体的均值和方差，渐进地优化分配给定仿真计算预算。在使用OCBA 方法时，可能存在一些个体会吸收大部分仿真预算，而其他个体则很难获得任何仿真预算，这些个体被称为“超级个体”。“超级个体”的存在将降低正确选择的可能性，因此，需要对OCBA 方法进行改进。

将LS设置为仿真次数的阈值来检测超级个体。“超级个体”的集合为Θ，其规模为nnum。从原始种群中删除所有“超级个体”，从给定的仿真总次数中减去“超级个体”所占用的仿真次数。最后，将现有的仿真次数分配给剩余个体。将Lη，η=1,2,…,κ定义为分配给调度解η的仿真次数，并将T定义为每次迭代分配的总仿真次数。由以上改进的OCBA 方法得到的调度解η的仿真次数如式（6）所示：

其中：Ta表示已经消耗的仿真次数；是调度解Sη的样本方差；Sb是目前的最佳调度解；δb,η是Sη和Sb适应度值的均值之差。

3.3 仿真优化算法的实施步骤

本文仿真优化算法的实现步骤如下：

步骤1生成初始种群，包含（nmax-nmin+1）个子种群和1 个用于保存最优解的主种群，种群规模为Z。

步骤2对每个初始子种群的每个个体都分配相同的基本仿真次数T0，计算每个个体的样本均值和方差，获得初始信息。

步骤3在每个子种群中，根据个体均值和方差的统计指标，迭代地使用OCBA 方法将仿真次数预算分配给每个个体。每次分配的仿真次数增量为T。计算2 次分配后仿真获得的最优值之间的差ggap，与给定阈值DOCBA比较。如果ggap＞DOCBA，则继续执行迭代OCBA；如果ggap≤DOCBA，则停止OCBA 迭代并执行步骤4。

步骤4从总体中选择最优的α个个体。首先对α个个体进行仿真次数为LS=105的精确评估。然后选择最优个体作为该子种群这一代的最优解，贡献给主种群。当主种群数量达到Z时，需要对子种群贡献个体进行评估，若优于目前主种群中的最差调度解，则替代该调度解。每一代主种群的最优解为算法这一代的最优解。通过轮盘赌选择方法获得其他β个个体。

步骤5将每个子种群的最优个体作为移民算子，与目标子种群的（α+β）个个体之间基于亲和度进行交流。

步骤6在每个子种群中，交流后的（α+β）个个体使用交叉和变异操作来产生新种群，种群大小是Z。新种群的每个个体都被分配了相同的基本仿真次数T0，计算每个个体的样本均值和方差。

步骤7重复步骤3～步骤6，直到在一定数量的连续代，genmax、最优值没有改进，停止迭代输出当前最优解。

4 数值实验与结果分析

在MATLAB 中实现上述仿真优化算法，在配置为（Intel®CoreTMi7-9700 3.00 GHz CPU 16 GB RAM）的计算机上运行。使用来自某大型医院体检中心的真实数据，该体检中心常规项目中包含8 项检查。分析并拟合实际的测量数据，每项检查的服务时间遵循指数分布，服务时间的数学期望分别为｛3，4，5，2，3，2，6，10｝min。假设每个检查科室只有1 个服务台。顾客可以进入系统的时间为系统开放的前180 min。将调度批次时间间隔设置为30 min，故可以安排7 个批次顾客到达。优化一个批次内顾客到达人数n及每位顾客的检查项目顺序。

4.1 参数设置

对于每个批次内顾客到达人数n的取值范围，考虑到服务时间最长的项目的时长均值为10 min，调度批次时间间隔30 min 内可安排3 位顾客。为了避免系统超时加班时长过多及服务台空闲时间过长，将n的取值范围定为2≤n≤5，子种群个数为4。在基础算例中，等待时间阈值TW设置为15 min，之后会针对此阈值进行算法敏感度分析。OCBA 初始仿真次数T0以及每次迭代中增加的仿真次数预算T的设置为T0=30，T= 66 000［11］。

设计对比实验以得到以下参数值：1）种群大小Z=100，150 或200；2）OCBA 阈值DOCBA=0.5 或1；3）交叉概率pc=0.7，0.8，0.9 或1.0，以及变异概率pm=0.2，0.3，0.4，0.5；4）每代选择的最好的调度解数量α=3，5，以及通过轮盘赌选择的解数量β=Z/2-α或Z/5-α；5）停止规则中没有改进的代数genmax=30，40，50；6）2 组不同权重系数：γ1=1,γ2=5 或γ1=1,γ2=10。同时，考虑解的质量和运行时间，所获得的最适当参数值如表1所示。

表1 参数设置Table1 Parameter settings

4.2 数值实验结果

根据上述参数设置，使用仿真优化算法计算问题实例得到：在单服务台情况下，每个批次的最优人数为3，工作时间内服务21 人，顾客平均等待时间为29.1 min，顾客超时等待时间之和为120.6 min，科室平均空闲时间为99.8 min。图6所示为适应度的最优值随着迭代次数的变化情况。最优值在100 代左右收敛。

图6 适应度函数最优值收敛情况Fig.6 Convergence of fitness function optimal value

为了验证所提出的仿真优化算法的有效性，模拟实际中不对体检顾客进行调度的情况建立离散事件仿真模型，将仿真优化算法（SOA）的调度结果与离散事件仿真结果对比。在仿真模型1 中，顾客自由到达，到达过程服从泊松分布，顾客完成一项检查后，在未完成检查中选择队列长度最短的加入队列；在仿真模型2 中，顾客自由到达，对于顾客体检项目顺序的调度使用几种启发式调度策略。在每个项目上，顾客优先规则为：最长剩余服务时间优先（A1）；最短剩余服务时间优先（A2）。当顾客完成一项检查时，下一项检查的选择规则为：最长服务时间优先（B1）；最短服务时间优先（B2）。考虑到实际中顾客在系统开放的初期集中到达的特点，设置自由到达的顾客在系统开放前30 min 遵循参数为2λ的泊松分布，随后在150 min 内遵循参数为λ的泊松分布。仿真实验考虑与仿真优化算法最优调度解相同的顾客总数，利用与仿真优化算法相同的服务时间随机种子，运行LS=105次。表2所示为仿真优化算法（SOA）、不进行任何调度的仿真（DES）以及顾客自由到达情况下使用最优顺序调度策略（A2B2）的结果对比。

表2 仿真优化算法与离散事件仿真结果对比Table 2 Comparison of simulation optimization algorithm and discrete event simulation results

由表2 可以看出，考虑随机服务时间，与DES得到的结果相比，SOA 得到的结果虽然加班时间有所增加，但顾客平均等待时间减少25.8%，顾客在一个科室前等待超过15 min 的总超时等待时长减少71.2%。与A2B2 调度策略的结果对比，SOA算法的顾客平均等待时间减少19.4%，使顾客在一个科室前等待超过15 min 的总超时等待时长减少70.9%。可以看出，本文所提的调度策略及仿真优化算法显著缓解了体检系统的顾客排队现象。

4.3 敏感度分析

为了进一步验证算法的有效性，针对不同的最长等待时间阈值进行敏感度分析。保持其他参数不变，取最长等待时间阈值TW为5 min、10 min、15 min、20 min。仿真优化算法（SOA）、不进行任何调度的仿真（DES）以及使用最优顺序调度策略（A2B2）的结果对比如表3所示。

表3 等待时间阈值敏感度分析Table 3 Sensitivity analysis of waiting time threshold

由表3 可以看出，在4 组不同的最长等待时间阈值实验中，SOA 的调度结果均优于DES 和A2B2 算法的结果。在2 个等待时间指标上，无论等待时间阈值设置为多少，SOA 的结果均优于另外2 种算法的结果。尤其是超时等待指标，SOA 算法的结果较启发式调度策略结果至少减少55.5%。这说明本文所提出的调度策略及仿真优化算法，可以有效地解决体检机构顾客排队时间过长的问题。最长等待时间阈值为5 min、10 min 和15 min 时，仿真优化算法所求出的最优批次人数为3；当等待时间阈值为20 min 时，最优批次人数为4。这说明随着最长等待时间阈值的增加，顾客在一个项目前等待时间超过阈值时间减少，每批次可以调度的顾客数量增加。同时也可以看出，当每批次调度人数增加时，系统加班时间、顾客平均等待时间、顾客超时等待时间均增加，科室平均空闲时间减少。在实际应用中，管理者可以根据需要适当放宽顾客等待时间限制，以服务更多的顾客。

5 结束语

本文考虑了服务时间随机的体检顾客调度问题，使用多人时间槽预约调度规则，从预约调度批次人数和顾客体检项目顺序2 个方面进行调度优化。提出了一种两阶段仿真优化算法，采用序优化策略提高仿真效率，多种群遗传算法作为迭代优化策略，最优计算量分配算法被嵌入到算法的粗糙仿真评估阶段，形成了仿真资源的全局和自适应优化分配机制。与不进行任何调度的真实情况以及启发式调度规则对顾客体检顺序的结果进行比较，实验结果验证了所提出的调度策略及仿真优化算法的有效性。本文研究不仅可以为类似体检机构这样多科室串联的服务机构管理者提供有效的顾客调度方法，还可以辅助体检机构管理者进行决策，增加关键科室的设备或人员投入，平衡拥堵科室和非拥堵科室的服务能力。下一步研究将多服务台、检查项目的部分顺序约束或某项检查完成时间的约束以及顾客的不按时到达或爽约行为，以此提高体检系统的整体服务水平。