The One with ...

授業・研究用のメモなどを公開しています

DAアルゴリズムにかんするメモ

RによるDAアルゴリズムの再現

DAアルゴリズム(GSアルゴリズムとも呼ばれます)のコードをRで書きなおして簡単な実験をやってみました. オリジナルのコードは奥村晴彦先生の『C言語による最新アルゴリズム事典』技術評論社を参照しました.

Cではprefix increment,つまりループの中に ++ v を入れることでコードが簡略化されています. 同じ記法はMathematicaでも可能だったので,Mathematicaへの翻訳は簡単でした. しかしRには仕様上prefix incrementがないらしく,その部分で多少の変更を余儀なくされました. まずコードの概略は以下のようなものです.DAアルゴリズムそのものの説明は,ここでは省略します (近々刊行予定の本で説明していますので、よろしければそちらを参照してください(宣伝)).

『その問題、数理モデルが解決します』 https://www.amazon.co.jp/dp/4860645685

###################################################
# GS Algorithm 概略
###################################################
#関数gsは引数として,男性の選好(proposer)と女性の選好(receiver)を
#入力すると,GS Algorithmに基づく男女のマッチングを出力します

#二つの引数proposerとreceiverは男女である必要はありませんが,
#先行研究にならって,以下男女の喩えを使います.

# はじめにproposer, receiverの選好を,集団別にlist型で入力します.
# 例えば男性が求婚する時,proposerとして
# list(c(2, 3, 1), c(1, 2, 3), c(3, 1, 2)) を入力する
# これは
#男性1の選好は  女2>女3>女1
#男性2の選好は  女1>女2>女3
#男性3の選好は  女3>女1>女2
#を意味します.
#同様にreceiverは,女性が持つ男性に対する選好のリストです.

関数内のローカル変数の意味は以下の通りです.カウンタがたくさんあるため,構造は結構複雑です.

########################################
# gs内のローカル変数の意味
########################################
#prefb:男性が持つ,女性に対する選好. prefg:女性が持つ,男性に対する選好.
#nextgirl:プロポーズする順番を記録.各男性が独立に持つ.振られる度に1増加
#tempboy:仮マッチした男性のリスト,最初はn+1(番兵)が入る.
#gid:男性がプロポーズする相手女性のid, 
#b:for用のカウンタ.男性1~nに対応
#s:現在プロポーズ中の男性id
#t:振られた男id.この値をsに代入する事で,
#「今回振られた男t」が自動的に次のsになる.

結果以下のようなコードになりました *1

##################################
# 関数gsの定義
##################################
gs130<- function(proposer,receiver)
{ n<-length(proposer)
prefb = proposer #男性の好みは引数proposerをそのまま入力
z=receiver  #一旦,女性の選好リストを別のzに代入する
for(j in 1:n){
    for(i in 1:n){  
        z[[j]][i]=which(receiver[[j]]==i)}#for j
}#for i
prefg =z #女性の選好を順位表記にする
print("男性の選好") 
print(paste0(matrix(prefb))) 
print(paste0("女性の選好") )
print(paste0(matrix(receiver)) )
for(i in 1:n){
    prefg[[i]] = append(prefg[[i]], n + 1)
}#女性の選好に番兵n+1を追加
nextgirl = rep(1, n)  #個人別のプロポーズ用カウンタ.1位からスタートして,2位,3位と順番に進む
tempboy = rep(n + 1, n)  # while条件用 番人sentinel  
for(b in 1:n){
    s = b
    while(s != n + 1 ){ #sが1,2,3である間ループ*)
        k <- nextgirl[s]
        gid =prefb[[s]][k]#k=1からスタート
        if(prefg[[gid]][s] < prefg[[gid]][tempboy[gid]])#現在の相手順位<仮の相手順位
        {   t = tempboy[gid]#仮の相手を捨てる
        print(c("振られた男=",t))#初回はt=n+1,入れ替わるとt=以前に仮マッチした男性ID 
        tempboy[gid] = s # 仮の相手を入れ替え
        s = t  # 次のs(proposer)を指定する.
        } #if 条件が真の場合の実行文終わり
        else{#もしsが拒否されたらnextgirl[s]に1を足す
        nextgirl[s] <- nextgirl[s]+1
        }
    } # whileここで終わり
}# for ここで終わり

for(i in 1:n)#結果の表示
{print(paste0("男性",tempboy[i], " と 女性",i))}
}#gs ここで終わり

さて,ランダムな選好を関数

lg <- function(n){replicate(n,list(sample(1:n)))}

で発生させると,関数gsは次のように実行できます.

gs130(lg(10),lg(10))

これで10対10の,マッチング結果が出てきます.n=3の場合の計算結果を確認してみましょう.

gs130(lg(3),lg(3))
[1] "男性の選好"
[1] "c(3, 1, 2)" "c(3, 1, 2)" "c(1, 3, 2)"
[1] "女性の選好"
[1] "c(1, 3, 2)" "c(1, 3, 2)" "c(3, 1, 2)"
[1] "振られた男" "4"         
[1] "振られた男" "4"         
[1] "振られた男" "2"         
[1] "振られた男" "4"         
[1] "男性3 と 女性1"
[1] "男性2 と 女性2"
[1] "男性1 と 女性3"

"振られた男" "4" は番兵(センチネル)です.whileループを抜け出す条件として使います.

このアルゴリズムは相手集団に対する完全な選好が必要なのですが,現実的には全員に順位をつけるのは難しいでしょう. 一方でアルゴリズムの原理から,選好によっては10人中5人くらいまで順位をつけておけば事足りる場合もあることは想像できます.

そこで10対10マッチングを実行した場合に,選好リストの下位何番までの情報が必要になるのかを調べてみました.

f:id:hamada7418:20181107180016p:plain

この図は,最も望ましくない相手とペアになったマッチングがどのくらいの確率で生じるのかを表しています. 500回ランダムな選好を発生させてGSアルゴリズムでマッチングを作り,それを100回計測した分布です. 例えば一番左は,proposerから見て一番順位の低い人とのペアが実現したマッチングの割合です.全体の3%程度このようなマッチングが存在します. 逆に言えば,97%は最後まで選好を表明しなくても安定マッチングができるという意味です.

図の横軸-2は1~2番目に低い順位とのマッチング,-3は1~3番目に順位の低い人とのマッチングが生じた割合を示しています.

DAアルゴリズムはとても単純なアルゴリズムなのに,結果は強力です(提案側最適な安定マッチングを実現する).

コードを書くよい練習になるので,みなさんも一度再現に挑戦してみてください.

*1:ちょっと細かい話ですが,関数の内部でreceiverの選好をidではなく順位に変換しています. 例えば,ある女性の選好がc(2,3,1) 男2>男3>男1 のとき,この選好を (男1は)3番目, (男2は)1番目, (男3は)2番目  という順位のベクトルに変換します. つまり c(2, 3, 1)] というreceiverの選好の入力をいったん  c(3, 1, 2) と変換します. これは receiver側で相手を留保するかどうかを順位を参照して決定するためです