githubに公開されているCPLEXのサンプルnurse_schedulingを読み解きます。
問題
このモデルはナース・スケジューリングを扱います。ナーススケジューリングでは一般的に様々なスキルとシフト割り当ての制約を満たしたシフトにナースをアサインしなければなりません。
このモデルのゴールは以下の2つの目的の間でバランスが取れた解を見つけることです。
- 全体のコストの最小化
- 可能な限り公平なシフト割り当て
利用データ
サンプルでは以下のデータが与えられています。
部門
割り当て対象の部門です。ただし、シフト・データではここで定義されている部門以外の部門も記載されており、そちらに従ってシフトのアサインを行うため、このデータはサンプルコードでは使用されません。
| name |
|---|
| Emergency |
| Consultation |
スキル
ナースが保有するスキルの種類です。
サンプルコードではCardiac Careしか使われていません。
| name |
|---|
| Anaesthesiology |
| Cardiac Care |
| Geriatrics |
| Oncology |
| Pediatrics |
シフト
どの部門で何曜日の何時から何時までの最低必要人数と最大必要人数が記載されたシフトの要件の抜粋です。
サンプルでは部門データに記載されているEmergency, Consultation以外の部門のシフトも記載されており、実際のアサイン対象の部門はEmergency, Consultation, Cardiac Care, Geriatricsの4種類になります。
| department | day | start_time | end_time | min_req | max_req |
|---|---|---|---|---|---|
| Emergency | Monday | 2 | 8 | 3 | 5 |
| Emergency | Monday | 8 | 12 | 4 | 7 |
| Emergency | Monday | 12 | 18 | 2 | 5 |
| Emergency | Monday | 18 | 2 | 3 | 7 |
| Consultation | Monday | 8 | 12 | 10 | 13 |
| Consultation | Monday | 12 | 18 | 8 | 12 |
各部門のシフト要件を図にすると以下のようになります。
ここで塗り潰された部分までが最低必要人数(min_req)で枠線のみの部分までが最大必要人数(max_req)です。
ナース
ナースの一覧の抜粋です。seniority, qualificationの記載がありますがサンプルコードでは使用されません。
また、サンプルコードでは時給(pay_rate)が15以上25未満のナースに絞って問題を解いています。
| name | seniority | qualification | pay_rate |
|---|---|---|---|
| Anne | 11 | 1 | 25 |
| Bethanie | 4 | 5 | 28 |
| Betsy | 2 | 2 | 17 |
| Cathy | 2 | 2 | 17 |
| Cecilia | 9 | 5 | 38 |
| Chris | 11 | 4 | 38 |
ナース・スキル
各ナースの持つスキルの一覧の抜粋です。一人のナースが複数のスキルを持つことができます。
サンプルでは様々なスキルが記載されていますがサンプルコードではCardiac Careのみが使用されます。
| nurse | skill |
|---|---|
| Anne | Anaesthesiology |
| Anne | Oncology |
| Anne | Pediatrics |
| Betsy | Cardiac Care |
| Cathy | Anaesthesiology |
| Cecilia | Anaesthesiology |
| Cecilia | Oncology |
| Cecilia | Pediatrics |
必要スキル
部門に必要なスキルの一覧です。
| department | skill | required |
|---|---|---|
| Emergency | Cardiac Care | 1 |
休暇
ナースの休暇希望の一覧の抜粋です。サンプルコードではナースの休暇希望にはシフトを入れることはできません。
| nurse | day |
|---|---|
| Anne | Friday |
| Anne | Sunday |
| Cathy | Thursday |
| Cathy | Tuesday |
| Joan | Thursday |
| Joan | Saturday |
一緒に勤務させたいペア
一緒に勤務させたいペアの一覧です。サンプルコードではこれらのペアが同じシフトで勤務させる必要があります。
| nurse1 | nurse2 |
|---|---|
| Isabelle | Dee |
| Anne | Patrick |
一緒に勤務させたくないペア
一緒に勤務させたくないペアの一覧の抜粋です。サンプルコードではこれらのペアは同じシフトで勤務させることができません。
| nurse1 | nurse2 |
|---|---|
| Patricia | Patrick |
| Janice | Wendie |
| Suzanne | Betsy |
| Janelle | Jane |
| Gloria | David |
| Dee | Jemma |
| Bethanie | Dee |
決定変数
- 各ナースのシフト割り当て
- 各ナースの勤務時間合計
- 各ナースの超過勤務時間
- 各ナースの過小勤務時間
- ナースの勤務時間平均
目的関数
"給与コスト + フェアネス + アサイン総数"を最小化します。
- 給与コストは全ナースの給与の合計です
- フェアネスは全ナースの超過勤務時間+全ナースの過小勤務時間の合計です。これを最小化することで各ナースの勤務時間を平均に近づけることができ、ナース間の勤務時間の差が小さくなるようにします
- アサイン総数は全ナースのシフト総数です。ナースの給与を最小化するのに加えて、シフト総数を最小化することで過剰なシフト割り当てを抑制します
制約
制約条件は以下の項目となります。
サンプルコードでは目的関数の解を得る際に下記の制約条件全てを満たすことはできず、制約条件を緩和して最適解を求めています。
制約条件の緩和優先度はmedium, highの順となっており、以下の制約条件の中でサンプルコードにて緩和可能な制約条件として指定されたものには緩和優先度を記載をしています。
- ナースの勤務時間平均*全ナースの数は全ナースの勤務時間合計と等しくなること
- 各ナースの勤務時間合計は各ナースのシフト割り当て*シフトの勤務時間の合計と等しくなること
- 各ナースの勤務時間合計はナースの勤務時間平均に各ナースの超過勤務時間を足したものまたは各ナースの過小勤務時間を引いたものと等しくなること
- 各ナースの勤務時間合計はナースの最大勤務可能時間を下回ること。サンプルではナースの最大勤務時間は40時間としている
- 各ナースの休暇にシフトを入れないこと（medium）
- 各ナースは同じ時間帯に複数のシフトがアサインされないこと（high）
- 各シフトに対して最小ナース数を上回るナースがアサインされること（high）
- 各シフトに対して最大ナース数を上回るナースがアサインされないこと（medium）
- 各シフトに最低一人のナースがアサインされること
- 各シフトに対して必須スキルを持つナースが必要人数アサインされること（high）
- 一緒に勤務させたいペアは同じシフトにアサインすること（medium）
- 一緒に勤務させたくないペアは同じシフトにアサインしないこと（medium）
問題の種類
混合整数計画
解
この問題は元データから特定の給与レンジのナースに絞り込んでいるため、上記の制約条件の元では解なしになります。
サンプルコードではCPLEXのRelaxerを使用して制約条件を緩和して解を求めています。
データ準備
曜日の定義
月曜日から日曜日までの曜日を定義し、曜日のID(0〜6)を返すようにday_to_day_week関数を定義します。
# utility to conevrt a weekday string to an index in 0..6
_all_days = ["monday", "tuesday", "wednesday", "thursday", "friday", "saturday", "sunday"]
def day_to_day_week(day):
day_map = {day: d for d, day in enumerate(_all_days)}
return day_map[day.lower()]
NamedTupleの定義
後のデータ処理のためにNamedTupleを定義します。
TWorkRules = namedtuple("TWorkRules", ["work_time_max"])
TVacation = namedtuple("TVacation", ["nurse", "day"])
TNursePair = namedtuple("TNursePair", ["firstNurse", "secondNurse"])
TNurseSkill = namedtuple("TNurseSkill", ["nurse", "skill"])
TSkillRequirement = namedtuple("TSkillRequirement", ["department", "skill", "required"])
さらにNamedTupleのクラスを定義します。
# subclass the namedtuple to refine the str() method as the nurse's name
class TNurse(namedtuple("TNurse1", ["name", "pay_rate"])):
""" A subclass to redefine the default str() of namedtuple.
This class is used in variable naming, so we need to redefine the str() method
used by variable naming.
"""
def __str__(self):
return self.name
class TShift(
namedtuple("TShift1", ["department", "day", "start_time", "end_time", "min_requirement", "max_requirement"])):
""" specialize namedtuple to redefine its str() method
"""
def __str__(self):
# keep first two characters in departement, uppercased
dept2 = self.department[0:4].upper()
# keep 3 days of weekday
dayname = self.day[0:3]
return '%s_%s_%02d' % (dept2, dayname, self.start_time)
class ShiftActivity(object):
@staticmethod
def to_abstime(day_index, time_of_day):
""" Convert a pair (day_index, time) into a number of hours since Monday 00:00
:param day_index: The index of the day from 1 to 7 (Monday is 1).
:param time_of_day: An integer number of hours.
:return:
"""
time = 24 * (day_index - 1)
time += time_of_day
return time
def __init__(self, weekday, start_time_of_day, end_time_of_day):
assert (start_time_of_day >= 0)
assert (start_time_of_day <= 24)
assert (end_time_of_day >= 0)
assert (end_time_of_day <= 24)
self._weekday = weekday
self._start_time_of_day = start_time_of_day
self._end_time_of_day = end_time_of_day
# conversion to absolute time.
start_day_index = day_to_day_week(self._weekday)
self.start_time = self.to_abstime(start_day_index, start_time_of_day)
self.day_start_time = self.to_abstime(start_day_index, 0)
end_day_index = start_day_index if end_time_of_day > start_time_of_day else start_day_index + 1
self.end_time = self.to_abstime(end_day_index, end_time_of_day)
assert self.end_time > self.start_time
@property
def duration(self):
return self.end_time - self.start_time
def overlaps(self, other_shift):
if not isinstance(other_shift, ShiftActivity):
return False
else:
return other_shift.end_time > self.start_time and other_shift.start_time < self.end_time
データ読み込み
元データとなるExcelファイルをDataFrameとして読み込んで各シートのデータを定義します。
data_url = "https://github.com/IBMDecisionOptimization/docplex-examples/blob/master/examples/mp/jupyter/nurses_data.xls?raw=true"
nurse_xls_file = pd.ExcelFile(data_url)
SkillTable = nurse_xls_file.parse('Skills')
DeptTable = nurse_xls_file.parse('Departments')
ShiftTable = nurse_xls_file.parse('Shifts')
SkillRequirementTable = nurse_xls_file.parse('SkillRequirements')
NurseTable = nurse_xls_file.parse('Nurses')
NurseSkillTable = nurse_xls_file.parse('NurseSkills')
NurseVacationTable = nurse_xls_file.parse('NurseVacations')
NurseAssociationTable = nurse_xls_file.parse('NurseAssociations')
NurseIncompatibilityTable = nurse_xls_file.parse('NurseIncompatibilities')
定数定義
以降で使用する定数を定義します。
skills
スキルの種類の配列です。
skills = [SkillTable["name"][i] for i in range(len(SkillTable))]
depts
部門名の配列です。ただし、サンプルコードでは以降で使用されません。
depts = [DeptTable["name"][i] for i in range(len(DeptTable))]
nurses
ナースの名前(name)と時給(pay_rate)のNamedTupleの配列です。
ここで元データのナースのリストから給与(pay_rate)が15以上25未満のナースに絞り込んでいます。
ナースの人数を絞り込んでいることからサンプルコードでは制約条件を厳密に満たすことはできず、解を得るためには制約条件の緩和が必要となります。
nurses =[TNurse(NurseTable["name"][i], NurseTable["pay_rate"][i]) for i in range(len(NurseTable))]
nurses = [nurse for nurse in nurses if nurse.pay_rate <25 and nurse.pay_rate >15]
nurse_skills
ナースの名前をkey, ナースの保有するスキルの配列をvalueとする辞書です。
nurse_skills = {}
for nsk in NurseSkillTable.itertuples(index=False):
nskt= TNurseSkill(*nsk)
nurse_skills.setdefault(nskt.nurse, []).append(nskt.skill)
shifts
NamedTuple: TShiftの配列です。
TShiftには部門名(department), 曜日(day), 開始時間(stat_time), 終了時間(end_time), 最小必要人数(min_requirement), 最大必要人数(max_requirement)が含まれます。
shifts = [TShift(*shift_row) for shift_row in ShiftTable.itertuples(index=False)]
skill_requirements
NamedTuple: TSkillRequirementの配列です。
TSkillRequirementには部門名(department), スキル(skill), 必須(required)が含まれます。
skill_requirements = [TSkillRequirement(*skill_requirement_row) for skill_requirement_row in
SkillRequirementTable.itertuples(index=False)]
vacations
NamedTuple: TVacationの配列です。
TVacationにはナース名(nurse), 曜日(day)が含まれます。
vacations = [TVacation(*vacation_row) for vacation_row in NurseVacationTable.itertuples(index=False)]
nurse_associations
NamedTuple: TNursePairの配列です。
TNurcePairにはナース1(firstNurse), ナース2(secondNurse)が含まれます。
nurse_associations = [TNursePair(*na) for na in NurseAssociationTable.itertuples(index=False)]
nurse_incompatibilities
NamedTuple: TNursePairの配列です。
TNurcePairにはナース1(firstNurse), ナース2(secondNurse)が含まれます。
nurse_incompatibilities = [TNursePair(*na) for na in NurseIncompatibilityTable.itertuples(index=False)]
shift_activities
NamedTuple: TShiftをkey, Class: ShiftActivityをvalueとする辞書です。
Class: ShiftActivityからシフトが同じ時間帯に重複しているか(overlaps)やシフトが何時間勤務であるか(duration)を取得することができます。
shift_activities = {s: ShiftActivity(s.day, s.start_time, s.end_time) for s in shifts}
nurse_by_id
ナース名(name)をkey, IDをvalueとする辞書です。
nurses_by_id = {n.name: n for n in nurses}
work_rules
NamedTuple: TWorkRulesです。
サンプルコードではwork_time_maxを40時間と定義しています。
work_rules = TWorkRules(40)
決定変数
次に決定変数を定義します。
決定変数1: nurse_assignment_vars
各ナースが各シフトにアサインされているかを表す変数です。
binary_var_matrixでバイナリ変数の行列として表現しています。
# One binary variable for each pair (nurse, shift) equal to 1 if nurse n is assigned to shift s
nurse_assignment_vars = mdl.binary_var_matrix(nurses, shifts, 'NurseAssigned')
決定変数2: nurse_work_time_vars
各ナースの勤務時間合計を表す変数です。
couninuous_var_dictで各ナースがkey, 勤務時間合計がvalueの辞書で、勤務時間合計の下限は0として定義しています。
# For each nurse, allocate one variable for worktime
nurse_work_time_vars = mdl.continuous_var_dict(nurses, lb=0, name='NurseWorkTime')
決定変数3: nurse_over_average_time_vars, nurse_under_average_time_vars
ナース全体の平均勤務時間からの各ナースの超過勤務時間(nurse_over_average_time_vars)と過小勤務時間(nurse_under_average_time_vars)を表す変数です。
couninuous_var_dictで各ナースがkey, 勤務時間合計がvalueの辞書で、いずれの変数も下限は0として定義しています。
# And two variables for over_average and under-average work time
nurse_over_average_time_vars = mdl.continuous_var_dict(nurses, lb=0, name='NurseOverAverageWorkTime')
nurse_under_average_time_vars = mdl.continuous_var_dict(nurses, lb=0, name='NurseUnderAverageWorkTime')
決定変数4: average_nurse_work_time
ナース全体の平均勤務時間を表す変数です。
countinuous_varで下限は0として定義しています。
# Finally the global average work time
average_nurse_work_time = mdl.continuous_var(lb=0, name='AverageWorkTime')
制約条件
モデルに対してadd_constraintで制約条件を追加していきます。
各制約条件に対してadd_constraintの引数で名前を設定していて、制約条件の中にはhigh, mediumを含む名前を設定しているものが含まれます。
また、mandatoryとが含まれる名前の制約条件、名前にhigh, medium, lowが含まれない制約条件は緩和されません。
サンプルコードでは問題の解を得る際に、制約条件を満たす解が見つからない場合は制約条件を緩和する処理を実行します。
制約条件はlow > medium > highの優先度の順で緩和されます。
制約条件1: ナースの勤務時間平均の定義
ナースの勤務時間平均を各ナースの勤務時間の合計として定義します。
mdl.add_constraint(len(nurses) * average_nurse_work_time ==
mdl.sum(nurse_work_time_vars[n] for n in nurses), "average")
制約条件2: 各ナースに対する勤務時間の制約
各ナースの勤務時間合計は以下の制約を満たす必要があります。
- 各ナースの勤務時間合計は各ナースがアサインされたシフトの勤務時間の合計であること
- 各ナースの勤務時間合計はナースが超過勤務をしている場合はナース全体の平均勤務時間に各ナースの超過勤務時間を足したものであること
- 各ナースの勤務時間合計はナースが過小勤務をしている場合は時間をナース全体の平均勤務時間から各ナースの過小勤務時間を引いたものであること
- 各ナースの勤務時間合計は最大勤務可能時間を下回ること
for n in nurses:
work_time_var = nurse_work_time_vars[n]
mdl.add_constraint(
work_time_var == mdl.sum(nurse_assignment_vars[n, s] * shift_activities[s].duration for s in shifts),
"work_time_{0!s}".format(n))
mdl.add_constraint(
work_time_var == average_nurse_work_time + nurse_over_average_time_vars[n] - nurse_under_average_time_vars[n],
"average_work_time_{0!s}".format(n))
mdl.add_constraint(work_time_var <= work_rules.work_time_max, "max_time_{0!s}".format(n))
1個目の制約条件は以下の図のように表せます。
図はあるナースのシフト・アサインを表しており、赤い箱はアサインされたシフト、箱の中の数字はシフトの勤務時間です。
この時シフトの勤務時間の合計は6+8+4+6+6+4+8=40となり、このナースの勤務時間合計は40時間になります。
2個目と3個目の制約条件は図のように表せます。
図はナースAとナースBの勤務時間合計を表しており、ナースAの勤務時間合計はナース全体の平均勤務時間にナースAの平均からの超過勤務時間を足したもの、ナースBの勤務時間合計はナース全体の平均勤務時間からナースBの平均からの過小勤務時間を引いたものであることを示しています。
1,2,3個目の制約条件はナースの勤務時間合計を違う形で表現したものであり、2,3個目の制約条件はフェアネスの計算に使用するために定義しています。
制約条件3: 各ナースの休暇にシフトを入れないための制約
各ナースの休暇の曜日にはシフトを入れないようにします。
この制約条件は優先度mediumで緩和されます。
for vac_nurse_id, vac_day in vacations:
vac_n = nurses_by_id.get(vac_nurse_id, -1)
if vac_n != -1:
for shift in (s for s in shifts if s.day == vac_day):
mdl.add_constraint(nurse_assignment_vars[vac_n, shift] == 0,
"medium_vacations_{0!s}_{1!s}_{2!s}".format(vac_n, vac_day, shift))
制約条件4: 各ナースが同じ時間帯に複数のシフトをアサインされないようにする制約
各ナースが同じ時間帯に複数のシフトにアサインされないようにします。
この制約条件は優先度highで緩和されます。
number_of_overlaps = 0
# Post only one constraint per couple(s1, s2)
nb_shifts = len(shifts)
for i1 in range(nb_shifts):
for i2 in range(i1 + 1, nb_shifts):
s1 = shifts[i1]
s2 = shifts[i2]
if shift_activities[s1].overlaps(shift_activities[s2]):
number_of_overlaps += 1
for n in nurses:
mdl.add_constraint(nurse_assignment_vars[n, s1] + nurse_assignment_vars[n, s2] <= 1,
"high_overlapping_{0!s}_{1!s}_{2!s}".format(s1, s2, n))
ここでは異なるシフトs1とs2を抜き出してShiftActivityクラスのoverlapsメソッドにより、s1とs2の時間帯が重複しているかを確認しています。
もし、s1とs2の時間帯が重複している場合（shift_activities[s1].overlaps(shift_activities[s2]) == true）は各ナースnに対してs1とs2の同時にシフトにアサインされないようにnurse_assignment_vars[n, s1] + nurse_assignment_vars[n, s2] <= 1と制約条件を追加しています。
制約条件5: 各シフトに対する最小ナース数および最大ナース数の制約
各シフトの割り当て人数は以下の制約を満たす必要があります。
- 各シフトの割り当て人数は各シフトの最低必要人数以上であること
- 各シフトの割り当て人数は各シフトの最大必要人数以下であること
- 各シフトの割り当て人数は1以上であること
1個目の制約条件は優先度highで緩和され、2個目と3個目の制約条件は優先度mediumで緩和されます。
for s in shifts:
demand_min = s.min_requirement
demand_max = s.max_requirement
total_assigned = mdl.sum(nurse_assignment_vars[n, s] for n in nurses)
mdl.add_constraint(total_assigned >= demand_min,
"high_req_min_{0!s}_{1}".format(s, demand_min))
mdl.add_constraint(total_assigned <= demand_max,
"medium_req_max_{0!s}_{1}".format(s, demand_max))
mdl.add_constraint(total_assigned >= 1,
"mandatory_presence_{0!s}".format(s))
制約条件6: 各シフトに対して必須スキルを持つナースを必要人数アサインするための制約
この制約条件は優先度highで緩和されます。
for (dept, skill, required) in skill_requirements:
if required > 0:
for dsh in (s for s in shifts if dept == s.department):
mdl.add_constraint(mdl.sum(nurse_assignment_vars[skilled_nurse, dsh] for skilled_nurse in
(n for n in nurses if n.name in nurse_skills.keys() and
skill in nurse_skills[n.name])) >= required,
"high_required_{0!s}_{1!s}_{2!s}_{3!s}".format(dept, skill, required, dsh))
ここでは必要スキルのデータに記載されている行でrequired>0となっている部門について、(s for s n shifts if dept == s.department)でその部門についてのシフトを抜き出しています。
そしてナース・スキルに記載されたナース(n.name in nurse_skills.keys())かつ必要とするスキルを持つナース(skill in nurse_skills[n.name])がそのシフトにアサインされた数(nurse_assignment_vars)の和がrequired以上であることを制約条件に追加しています。
制約条件7: 一緒に勤務させたいペアに対する制約
一緒に勤務させたいペアのシフトアサイン変数（nurse_assignment_vars）が一致するようにします。
これにより、一緒に勤務させたいペアは同じシフトにアサインされます。
この制約条件は優先度mediumで緩和されます。
c = 0
for (nurse_id1, nurse_id2) in nurse_associations:
if nurse_id1 in nurses_by_id and nurse_id2 in nurses_by_id:
nurse1 = nurses_by_id[nurse_id1]
nurse2 = nurses_by_id[nurse_id2]
for s in shifts:
c += 1
ctname = 'medium_ct_nurse_assoc_{0!s}_{1!s}_{2:d}'.format(nurse_id1, nurse_id2, c)
mdl.add_constraint(nurse_assignment_vars[nurse1, s] == nurse_assignment_vars[nurse2, s], ctname)
制約条件8: 一緒に勤務させたくないペアに対する制約
一緒に勤務させたくないペアのシフトアサイン変数（nurse_assignment_vars）の合計が1以下になるようにします。
これにより、一緒に勤務させたくないペアは同じシフトにアサインされなくなります。
この制約条件は優先度mediumで緩和されます。
c = 0
for (nurse_id1, nurse_id2) in nurse_incompatibilities:
if nurse_id1 in nurses_by_id and nurse_id2 in nurses_by_id:
nurse1 = nurses_by_id[nurse_id1]
nurse2 = nurses_by_id[nurse_id2]
for s in shifts:
c += 1
ctname = 'medium_ct_nurse_incompat_{0!s}_{1!s}_{2:d}'.format(nurse_id1, nurse_id2, c)
mdl.add_constraint(nurse_assignment_vars[nurse1, s] + nurse_assignment_vars[nurse2, s] <= 1, ctname)
目的関数
目的関数を定義します。
目的関数は"給与コスト + フェアネス + アサイン総数"として定義し、これを最小化することを目指します。
サンプルコードでは目的関数とは別にKPIを定義します。
KPI
KPIとして以下の項目を定義します。
- 給与コスト総額
- アサイン総数
- 全ナースの勤務時間平均
- 全ナースの超過勤務時間合計
- 全ナースの過少勤務時間合計
- フェアネス（超過勤務と過少勤務時間の差分）
KPIは解を求めた後に評価される指標で任意の指標を定義できます。
今回のコードでは目的関数が給与コストとフェアネスとアサイン総数の和になっており、求めた解の値は意味ある値にはなりません。
ここで目的関数を最小化した際に給与コスト総額がどうなったのかなどの意味ある指標の値を知りたい場合にKPIを使用します。
KPIを定義すると解を求めた際にKPIの値がどうなったかが出力されます。
KPIはモデルに対してadd_kpiで追加することできます。
給与コスト総額
給与コスト総額は各ナースのアサインされたシフトの勤務時間(duration)にナースの時給(pay_rate)をかけたものを全ナースで合計したものと定義します。
nurse_costs = [nurse_assignment_vars[n, s] * n.pay_rate * shift_activities[s].duration for n in nurses for s in shifts]
total_salary_cost = mdl.sum(nurse_costs)
mdl.add_kpi(total_salary_cost, "Total salary cost")
アサイン総数
アサイン総数は各ナースのシフトアサイン数の合計として定義します。
total_number_of_assignments = mdl.sum(nurse_assignment_vars[n,s] for n in nurses for s in shifts)
mdl.add_kpi(total_number_of_assignments, "Total number of assignments")
全ナースの勤務時間平均
全ナースの勤務時間平均は決定変数average_nurse_work_timeとして定義します。
average_nurse_work_timeは制約条件より全ナースの勤務時間合計の和をナースの数で割ったものとして設定されています。
mdl.add_kpi(average_nurse_work_time)
全ナースの超過勤務時間合計
全ナースの超過勤務時間合計は各ナースの超過勤務時間の全ナースの和として定義します。
total_over_average_worktime = mdl.sum(nurse_over_average_time_vars[n] for n in nurses)
mdl.add_kpi(total_over_average_worktime, "Total over-average worktime")
全ナースの過小勤務時間合計
全ナースの過小勤務時間合計は各ナースの過小勤務時間の全ナースの和として定義します。
total_under_average_worktime = mdl.sum(nurse_under_average_time_vars[n] for n in nurses)
mdl.add_kpi(total_under_average_worktime, "Total under-average worktime")
フェアネス
フェアネスは全ナースの超過勤務時間合計と全ナースの過小勤務時間合計の和として定義します。
total_fairness = total_over_average_worktime + total_under_average_worktime
mdl.add_kpi(total_fairness, "Total fairness")
目的関数
KPIで指定した給与コスト総額とアサイン総数とフェアネスの和を目的関数として定義し、これを最小化します。
mdl.minimize(total_salary_cost + total_fairness + total_number_of_assignments)
求解結果
上記で定義した制約条件、目的関数に対して解を求めます。
ここで解が整数であること、解の大きさが〜100000のオーダーであることから、mip.tolerances.mipgapを1e-5と設定します。
これにより、小数点以下の正確な解を求めることはせずにsolve処理は停止します。
さらにことでRelaxerを定義して制約条件を緩和しています。
緩和する条件はRelaxerにおいてprioritizerにmatchを指定することで、名前にhigh, medium, lowが含まれる制約条件を緩和の対象とします。
緩和される優先順位はlow, medium, highの順です。
mdl.parameters.mip.tolerances.mipgap = 1e-5
s = mdl.solve(log_output=True)
if not s:
from docplex.mp.relaxer import Relaxer
relaxer = Relaxer(prioritizer='match', verbose=True)
relaxed_sol = relaxer.relax(mdl)
relaxed_ok = relaxed_sol is not None
assert relaxed_ok, "relaxation failed"
relaxer.print_information()
mdl.report()
解を求めると以下のような出力が得られます。
relaxed: と表示されているものは緩和された制約条件です。
ここではシフトの最低必要人数とEmergency部門における必要スキルを持つナースのアサインの制約条件が緩和されています。
出力の最後の方でtotal absolute relaxationという出力があり、ここでは101個の制約条件が緩和されたことを表しています。
また、出力の最後にmodel nurses solved with objectiveという出力があり、これが最適化された目的関数の値を示しています。
その下にKPI: と表示されている値が7個出力されていますが、これが得られた目的関数の解におけるKPIの出力を示しています。
- relaxed: high_req_min_EMER_Sat_20_12, with relaxation: -4.0
- relaxed: high_req_min_EMER_Sun_02_5, with relaxation: -4.0
- relaxed: high_req_min_EMER_Sun_12_7, with relaxation: -5.0
- relaxed: high_req_min_EMER_Sun_20_8, with relaxation: -2.0
- relaxed: high_required_Emergency_Cardiac Care_1_EMER_Tue_02, with relaxation: -1.0
- relaxed: high_required_Emergency_Cardiac Care_1_EMER_Wed_02, with relaxation: -1.0
- relaxed: high_required_Emergency_Cardiac Care_1_EMER_Wed_18, with relaxation: -1.0
- relaxed: high_required_Emergency_Cardiac Care_1_EMER_Thu_12, with relaxation: -1.0
- relaxed: high_required_Emergency_Cardiac Care_1_EMER_Thu_18, with relaxation: -1.0
- relaxed: high_required_Emergency_Cardiac Care_1_EMER_Fri_12, with relaxation: -1.0
- relaxed: high_required_Emergency_Cardiac Care_1_EMER_Fri_18, with relaxation: -1.0
- relaxed: high_required_Emergency_Cardiac Care_1_EMER_Sat_02, with relaxation: -1.0
- relaxed: high_required_Emergency_Cardiac Care_1_EMER_Sat_12, with relaxation: -1.0
- relaxed: high_required_Emergency_Cardiac Care_1_EMER_Sun_02, with relaxation: -1.0
- relaxed: high_required_Emergency_Cardiac Care_1_EMER_Sun_12, with relaxation: -1.0
- relaxed: high_required_Emergency_Cardiac Care_1_EMER_Sun_20, with relaxation: -1.0
* total absolute relaxation: 101.0
* model nurses solved with objective = 14520.000
* KPI: Total salary cost = 14366.000
* KPI: Total number of assignments = 134.000
* KPI: AverageWorkTime = 38.889
* KPI: Total over-average worktime = 10.000
* KPI: Total under-average worktime = 10.000
* KPI: Total fairness = 20.000
結果分析
各部門のシフトアサイン状況
得られた解はシフトの最低必要人数の制約条件が緩和されているため、最低必要人数を満たしていないシフトが存在します。
図で表すと以下のようになります。
塗りつぶされた棒グラフは得られた解の各シフトの人数で、枠線のみの棒グラフはシフトの最低必要人数です。
これを見るとほとんどのシフトで最低必要人数を満たせておらず、シフトの最低必要人数の制約条件を緩和しないと解が存在しないことがわかります。
シフトの最低必要人数の制約条件を満たせない原因として元データからナースを特定の給与レンジのナースに絞り込んでいることがあります。
そこで、試しにナースのデータとして以下のようにナースを絞り込んでいる処理をコメントアウトして解を求めてみます。
nurses =[TNurse(NurseTable["name"][i], NurseTable["pay_rate"][i]) for i in range(len(NurseTable))]
#nurses = [nurse for nurse in nurses if nurse.pay_rate <25 and nurse.pay_rate >15]
解を求めた際の出力は以下のようになります。
ナースの絞り込んでいる場合の出力と比較するとrelaxed: という行の出力はなく、total absolute relaxationも出力されていないため、すべての制約条件が緩和されずに満たされた状態で解が求まったことを示しています。
Root node processing (before b&c):
Real time = 0.28 sec. (248.72 ticks)
Parallel b&c, 8 threads:
Real time = 0.56 sec. (583.07 ticks)
Sync time (average) = 0.07 sec.
Wait time (average) = 0.00 sec.
------------
Total (root+branch&cut) = 0.84 sec. (831.79 ticks)
* model nurses solved with objective = 29756.375
* KPI: Total salary cost = 29508.000
* KPI: Total number of assignments = 223.000
* KPI: AverageWorkTime = 39.562
* KPI: Total over-average worktime = 12.688
* KPI: Total under-average worktime = 12.688
* KPI: Total fairness = 25.375
この時の各部門のシフトアサイン状況は以下の図のようになり、各シフトで最低必要人数を満たすことができ、枠線のみの棒グラフは見えなくなっています。
各ナースのシフトアサイン状況
各ナースのアサイン状況は以下の図のようになります。
ここで緑色のシフトはGeriatrics, 赤色のシフトはEmergency, 青色のシフトはConsultation, 紫色のシフトはCardiac Careへのシフトを表します。
ナースの名前の右のworkedは各ナースの勤務時間合計を表します。
図を見ると同じ時間帯に複数のシフトにアサインされていることはないことがわかります。
よって、各ナースが同じ時間帯に複数のシフトをアサインされないようにする制約は満たされています。
一緒に勤務させたいペアに対する制約はナースが絞り込まれた結果、IsabelleとDeeの1ペアのみになっており、このペアは必ず同じシフトで勤務していることがわかります。
一緒に勤務させたくないペアに対する制約は多数ありますが、例としてJaneとJanelleのペア、DeeとJemmaのペアを見てみます。
どちらのペアも同一時間帯の同じ部門のシフトにアサインされることはなく、この制約条件が満たされていることがわかります。
各ナースの休暇にシフトを入れないための制約についても多数ありますが、例としてCathyとJoanを見てみます。
Cathyは火曜日(Tuesday)と木曜日(Thursday)が休暇となっていて、これらの曜日にシフトが入っていないことがわかります。
Joanは木曜日(Thursday)と土曜日(Saturday)が休暇となっていて、同様にこれらの曜日にシフトが入っていないことがわかります。
給与コスト, フェアネス, アサイン総数
目的関数は"給与コスト + フェアネス + アサイン総数"で定義されていたが、これらの項目についてどのような解が得られたかを確認してみます。
KPIより最適解では給与コストは14366, フェアネスは20, アサイン総数は134となっています。
また、KPIよりナース全体の平均勤務時間は38.889で各ナースのシフトアサイン状況の図から各ナースの勤務時間合計は36〜40時間となっていることがわかります。
ナースの勤務時間の標準偏差を計算してみると1.20となり、ナースの勤務時間合計はナースによって大きなばらつきはありません。
よって、各ナース間の勤務時間のばらつきを抑えるための項目であるフェアネスは有効に働いていたと言えます。
これは目的関数の各項目の値をみると少し意外な結果であると言えます。
目的関数の各項目のオーダーを考慮するとフェアネスを減らすよりも給与コストを減らす方が目的関数の最小化には有効のように見えます。
つまり、フェアネスを無視して時給の低いナースを多くアサインした方が目的関数を最小化できそうですが、結果としてはフェアネスは守られた結果となっています。
これは各ナースに対する勤務時間の制約で各ナースの最大勤務時間が40時間に制限されているためと考えられます。
考察
- サンプルコードではナースの連続勤務時間とシフト間のインターバルが考慮されていません。各ナースのシフトアサイン状況をみると24時間勤務や12時間勤務の後、6時間のインターバル後に10時間勤務が入ったりするなど、現実的には問題のあるシフト・アサインが見られるので連続勤務時間に対する制約や勤務と勤務の間のインターバルに対する制約の追加が必要だと考えます
- サンプルコードでは目的関数のフェアネスが有効に働かないと考えられるので、目的関数の各項のバランスを取るためのパラメータを導入して各項のオーダーを揃えるという改良が考えられます
参考
CPLEXサンプルを読み解く記事一覧