ビンゴの当選者を決めよう

## 問題

会社の忘年会でビンゴをします。あなたは、ビンゴをするためのシステム開発を頼まれました。ただし、以下のような要件があります（フィクションです）。

- 「賞品は7個」（→当選者数は7）
- 「なるべく男女で当選者の比率を一定にしてね」
- 「なるべく部署で当選者の比率を一定にしてね」
- 「若者、ベテラン、役員の当選確率（重み）を 1, 2, 3にしてね」

（参考：[よく似た話](https://el.jibun.atmarkit.co.jp/pressenter/2019/12/2.html)）

## 考え方

組合せ最適化を使うと、色々な条件を考慮した組み合わせを求められます。
当選確率に従って当選者を順次選び並べます。その並びに従って$i$番目のスコアを$2^{1-i}$とします。そのスコアの和を最大化することで、なるべく当選確率に従った当選者を選択します。
男女の比率や部署の比率は、比率から外れる量をペナルティーにして目的関数から引きます。

## Pythonで解く

参加者の情報をpandas.DataFrameに読み込みます（コードが実行できない場合、URLからCSVをダウンロードして読み込んでください）。

```py
import numpy as np, pandas as pd, urllib
from pulp import lpDot, lpSum
from ortoolpy import model_max, addbinvars, addvars, addvals

url = 'https://www.dropbox.com/s/refo0vfj5wbmv2h/bingo.csv?dl=1'
with urllib.request.urlopen(url) as fp:
    df = pd.read_csv(fp)
```

### スコアと変数作成

DataFrame（`df`）の各行が参加者になります。
当選確率（`Rate`列）を比率にして `numpy.random.choice` で順番に当選させます。その順番をindexにしたスコアを計算し、`Score`列にします。
最適化では、スコアの総和を最大化します。

当選かどうかを表す変数の列`Var`を追加します。
また、`pena1`を男女の比率から外れる量のペナルティーとします。`pena2`を所属の比率から外れる量のペナルティーとします。

```py
nn = len(df)  # 人数
rnd = np.random.choice(df.index, nn, False, df.Rate / df.Rate.sum())
df['Score'] = pd.Series(2.0**-np.arange(nn), index=rnd)  # スコア
addbinvars(df)  # Var列として変数（当選かどうか）追加
pena1, pena2 = addvars(2)  # 比率から外れる分のペナルティー2つ
print(df)
```

|     | Name   | IsMale | Div | Rate | Score    | Var     |
| --- | ------ | ------ | --- | ---- | -------- | ------- |
| 0  | 朝倉 英仁  | True   | 総務部 | 2    | 0.000488 | v000001 |
| 1  | 宮川 亮晴  | True   | 経理部 | 1    | 0.001953 | v000002 |
| 2  | 古橋 謙次  | True   | 総務部 | 2    | 0.000977 | v000003 |
| 3  | 本村 文   | False  | 総務部 | 3    | 0.000015 | v000004 |
| 4  | 大坪 慶二  | True   | 経理部 | 1    | 0.000002 | v000005 |
| 5  | 畠山 一隆  | True   | 総務部 | 2    | 0.000031 | v000006 |
| 6  | 竹内 和佳奈 | False  | 経理部 | 2    | 0.000008 | v000007 |
| 7  | 依田 香名恵 | False  | 人事部 | 2    | 0.015625 | v000008 |
| 8  | 古田 景輔  | True   | 人事部 | 3    | 0.125000 | v000009 |
| 9  | 石川 慶汰  | True   | 人事部 | 1    | 0.000061 | v000010 |
| 10 | 冨永 英莉子 | False  | 人事部 | 2    | 1.000000 | v000011 |
| 11 | 飛田 沙都  | False  | 総務部 | 1    | 0.003906 | v000012 |
| 12 | 神山 誠斗  | True   | 人事部 | 2    | 0.250000 | v000013 |
| 13 | 田端 龍一朗 | True   | 総務部 | 3    | 0.500000 | v000014 |
| 14 | 亀井 利恵  | False  | 経理部 | 1    | 0.000004 | v000015 |
| 15 | 曽根 温香  | False  | 総務部 | 1    | 0.000244 | v000016 |
| 16 | 春日 清志朗 | True   | 経理部 | 2    | 0.007812 | v000017 |
| 17 | 中井 明紗  | False  | 経理部 | 3    | 0.062500 | v000018 |
| 18 | 金城 幸帆  | False  | 人事部 | 2    | 0.031250 | v000019 |
| 19 | 柿沼 義克  | True   | 総務部 | 1    | 0.000122 | v000020 |

男女別部署別の人数を確認します。

```py
print(df.groupby('IsMale').Div.value_counts())
```

```
IsMale  Div
False   人事部    3
        経理部    3
        総務部    3
True    総務部    5
        人事部    3
        経理部    3
Name: Div, dtype: int64
```

### 数理モデルを作成して求解

- 目的関数：スコアの総和 - 人数 ×（男女比率外れ量 ＋ 部署比率外れ量）
- 制約条件
  - 当選者数は7
  - 男女別の比を一定に
  - 部署別の比を一定に

※ ペナルティーの係数は適当です。

```py
ns = 7  # 当選者数

m = model_max()  # 数理モデル
m += lpDot(df.Score, df.Var) - nn * (pena1 + pena2)  # 目的関数
m += lpSum(df.Var) == ns  # 当選者数
for _, gr in df.groupby('IsMale'):  # 男女別の比を一定に
    m += lpSum(gr.Var) <= ns * len(gr) / nn + pena1
for _, gr in df.groupby('Div'):  # 部署別の比を一定に
    m += lpSum(gr.Var) <= ns * len(gr) / nn + pena2
m.solve()  # 求解
addvals(df)  # Val列として結果追加
print(df[df.Val > 0])  # 当選者
```

|     | Name   | IsMale | Div | Rate | Score    | Var     | Val |
| --- | ------ | ------ | --- | ---- | -------- | ------- | --- |
| 10 | 冨永 英莉子 | False  | 人事部 | 2    | 1.000000 | v000011 | 1.0 |
| 11 | 飛田 沙都  | False  | 総務部 | 1    | 0.003906 | v000012 | 1.0 |
| 12 | 神山 誠斗  | True   | 人事部 | 2    | 0.250000 | v000013 | 1.0 |
| 13 | 田端 龍一朗 | True   | 総務部 | 3    | 0.500000 | v000014 | 1.0 |
| 16 | 春日 清志朗 | True   | 経理部 | 2    | 0.007812 | v000017 | 1.0 |
| 17 | 中井 明紗  | False  | 経理部 | 3    | 0.062500 | v000018 | 1.0 |
| 19 | 柿沼 義克  | True   | 総務部 | 1    | 0.000122 | v000020 | 1.0 |

当選者が決まりました。
男女別部署別の人数を見てみましょう。

```py
print(df[df.Val > 0].groupby('IsMale').Div.value_counts())
```

```
IsMale  Div
False   人事部    1
        経理部    1
        総務部    1
True    総務部    2
        人事部    1
        経理部    1
Name: Div, dtype: int64
```

良さそうです。

## まとめ

- 色々な条件を満たす組み合わせは、組合せ最適化という手法で求められます。
- 当選確率に従って当選する順番を決め、その順番ごとに`[1, 0.5, 0.25, ...]`をスコアとします。このスコアの和を最大化することで、予め決めた順番に当選者を選択しようとします。

（個人的には、この記事の需要がないことを祈ります）

以上