言語処理と医療の狭間で

AtCoderや日記、自然言語処理などについて書きます。

論文まとめ① [ClinicalNLP]

論文まとめ

概要

英語以外の言語でのClinical NLPサーベイ論文

  1. NLPシステムの構築
  2. 英語のNLPシステムの言語拡張
  3. アプリケーション

の3つの観点から研究を俯瞰している.

内容

日本語についての研究,個人的に面白そうな研究を列挙します.

英語のNLPシステムの言語拡張

MetaMapなどの英語NLPシステムの言語拡張を目指した研究

cTakesのドイツ語拡張
Extraction of UMLS® Concepts Using Apache cTAKES™ for German Language. - PubMed - NCBI

アプリケーション

上述のNLPシステムを使ったclinicalなアプリケーションの構築を目指した研究

エピソード記憶の分析,自動分類
Unraveling the linguistic nature of specific autobiographical memories using a computerized classification algorithm. - PubMed - NCBI

議論

英語以外の言語でのClinical NLPでは,主に「英語ですでにあるツールを自国語でも作りたい」というモチベーションで研究されているように感じました.
UMLSの翻訳などがされてない言語では,そもそも医療概念を定義・抽出すること自体が難しく,リソースの自動構築や英語リソースへのリンクなども行われているようです.

AtCoder:ABC140振り返り

AtCoderのABC140に参加したのでその振り返りです.今回はなんとかDまで解けたという感じで,D問題の実装にてこずりました...ただ,一応全部一発で通せたので,そこはよかったです.

A問題

f:id:kuroneko1259:20190908134548p:plain

3桁のうち,それぞれの桁が1からNまでN個の可能性があるので,単純に{N^3}です.

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>
#include <string>

using namespace std;

typedef long long ll;

int main() {
    int N; cin >> N;
    cout << pow(N, 3) << endl;
    return 0;
}

B問題

f:id:kuroneko1259:20190908134859p:plain

i番目の料理{A_i}を食べた時の満足度{B_{A_i}}{A_{i+1} = A_i + 1}なら追加の満足度{C_{A_i}}を足して行きます.これも複雑なアルゴリズムは必要なく,言われた通りに実装すれば良い問題でした.ここ1ヶ月くらいABCに参加した感じでは,基本B問題までは素直に実装すれば良い感じですね.

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>
#include <string>

using namespace std;

typedef long long ll;

int main() {
    int N; cin >> N;
    vector<int> a(N), b(N), c(N-1);
    for (int i = 0; i < N; i++) {
        cin >> a[i];
        a[i]--;
    }
    for (int i = 0; i < N; i++) cin >> b[i];
    for (int i = 0; i < N-1; i++) cin >> c[i];

    int res = 0;
    for (int i = 0; i < N; i++) {
        res += b[a[i]];
        if (i < N-1 && a[i+1] == a[i] + 1) res += c[a[i]];
    }

    cout << res << endl;
    return 0;
}

C問題

f:id:kuroneko1259:20190908135333p:plain
数列{A_i}に関する2項間の関係{B_i}がわかっていて,{\sum A_i}の最大値を求める問題でした.{B_i \geq \max (A_i, A_{i+1})}を読み替えると{A_i \leq min(B_{i-1}, B_i)}になります.{A_i}をそれぞれ最大にするには,{A_i = min(B_{i-1}, B_i)}にすればよいことになります.ただし,{A_1}{A_n}{B_0, B_n}がないので,{A_1 = B_1}{A_n = B_{n-1}}になります.

こういうところで{B_0, B_n = inf}(infは十分大きい数)とおければ場合分けなしでかっこよくコードが書けるんですが,なかなかとっさに思いつかないですね...

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>
#include <string>

using namespace std;

typedef long long ll;

int main() {
    int N; cin >> N;
    vector<ll> b(N-1);
    for (int i = 0; i < N-1; i++) cin >> b[i];

    ll res = b[0];
    for (int i = 0; i < N - 2; i++) {
        if (b[i] <= b[i+1]) res += b[i];
        else res += b[i+1];
    }
    res += b[N-2];
    cout << res << endl;
}

D問題

f:id:kuroneko1259:20190908140424p:plain
一応ACになりましたが,無理やり通した感...前にもLとRの文字列の問題がありましたが,規則はわかりますが実装が難しいです...

この問題では

  • LRLLを反転するとRRLR
  • LRLRを反転するとLRLR

みたいな感じで,両端以外は反転しても幸福な人は増えません.しかも両端がLRRLみたいに揃っていないと,ひっくり返しても幸福な人は増えないことになります(LLRRをひっくり返してもLLRR).

なので,基本的な方針としてLLRRLLのRのように別の文字に囲まれている文字を一気にひっくり返すことにしました.この操作1回ずつに幸福な人は2人(両端の人)増えることになります.また,両端を別の文字で囲まれて,同じのみで構成される部分列の個数はLとRで1個違いになります.

ここまで考えて,両端の文字が揃っていない場合にうまくいかない(最後の一回は1人しか幸福な人が増えない)ことに気づいて,慌てて場合分けを増やそうとして,あってるのかわからないけどサンプルは全部通ったし出しちゃえ!でACでした.

今思えば,両端が揃っていなくても最後の1回しか処理は変わらない +各操作で2人ずつ増えていくので,最初の幸福な人の人数をfとすると{f + 2 \times K}が幸福な人の最大値N-1以上になるかどうかで場合分けすれば良いだけの話でした.つまり,K回の操作すべてで2人増えることにしても差が生じるのは最後の操作だけであり,その場合は幸福な人の最大値を超えてしまうので最大値を出力すればよかったです.うまく言えませんが...

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>
#include <string>

using namespace std;

typedef long long ll;

int main() {
    int N, K; cin >> N >> K;
    string s; cin >> s;

    int res = 0;
    for (int i = 1; i < int(s.size()); i++) {
        if (s[i-1] == s[i]) res++;
    }

    vector<int> R, L;
    int idx = 0, r=0, l=0;
    while (idx != N && s[idx] != 'R') idx++;
    while (idx != N) {
        while (idx != N && s[idx] != 'L') idx++;
        while (idx != N && s[idx] != 'R') idx++;
        if (idx != N) l++;
    }
    idx = 0;
    while (idx != N && s[idx] != 'R') idx++;
    while (idx != N) {
        while (idx != N && s[idx] != 'R') idx++;
        while (idx != N && s[idx] != 'L') idx++;
        if (idx != N) r++;
    }

    if (r==0 || l==0) {
        cout << N - 1 << endl;
        return 0;
    }
    int sr=0, sl=0;
    if (s[0] == 'L') sl++;
    else sr++;
    if (s[int(s.size())-1] == 'L') sl++;
    else sr++;

    if (K >= min(r+sr, l+sl)) {
        cout << N-1 << endl;
        return 0;
    } else if (K >= min(r, l)){
        cout << res + min(r, l)*2 << endl;
    } else {
        cout << res + K*2 << endl;
    }
    return 0;
}

まとめ

今回はE問題にまわす時間はなかったです.見た感じE問題は時間があっても解けなかったと思いますが...D問題のような問題はなんの情報を保持しておけば良いのか(LとR,RとLの変わり目を分けるのかどうか,とか),それをどう処理すれば良いのか難しいです.

AtCoder:ABC139の振り返り

ABC139に参加したのでその振り返りです.今回はE問題まで解くことができました!ゴリゴリの競プロ的な問題ではなかったのがよかったのかな??と思います.まだまだ典型的なアルゴリズムや基本的な発想が身についていないので,今回のD問題やE問題のような数学的な問題の方が解きやすく感じます.

A問題

f:id:kuroneko1259:20190903091254p:plain
単純に3文字を比較するだけでした.

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>
#include <string>

using namespace std;

typedef long long ll;

int main() {
    string s, t; cin >> s >> t;
    int res = 0;
    if (s[0] == t[0]) res++;
    if (s[1] == t[1]) res++;
    if (s[2] == t[2]) res++;
    cout << res << endl;
}

B問題

f:id:kuroneko1259:20190903091416p:plain

最初,1口のときを考えていなくてWAが出てしまいました.もったいない!
場合分けしない方が良いと思いますが,へんに書き直すより場合分けで追加した方が良いかなと思ったのでこんな感じのコードになりました.

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>

using namespace std;

typedef long long ll;

int main() {
    int a, b; cin >> a >> b;
    if (b == 1) {
        cout << 0 << endl;
        return 0;
    }
    int tap = a, res = 1;

    while (tap < b) {
        tap += a - 1;
        res++;
    }
    cout << res << endl;
    return 0;
}

C問題

f:id:kuroneko1259:20190903091843p:plain

一つ右の階段での最大移動回数がわかれば,今の最大移動回数がわかるので(高さが右以上であれば+1,右より下であれば0)後ろからDPみたいな感じで回していけばO(N),その中で最大のものを探索してもO(N)だと思ってその方針で行きました.比較的すんなり解けて,この問題終わった時点で15分くらいだったので,次以降の問題につながったと思います.

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>

using namespace std;

typedef long long ll;

int main() {
    int N; cin >> N;
    vector<int> v(N);
    for (int i = 0; i < N; i++) cin >> v[i];

    int dp[N];
    for (int i = 0; i < N; i++) dp[i] = 0;

    for (int i = N-2; i >= 0; i--) {
        if (v[i] >= v[i+1]) dp[i] = dp[i+1] + 1;
    }
    
    int res = 0;
    for (int i = 0; i < N; i++) {
        if (dp[i] > res) res = dp[i];
    }
    
    cout << res << endl;
}

D問題

f:id:kuroneko1259:20190903092335p:plain

もはや数学...?競プロはこんな数学の問題もでるんですね...コードが最終的な計算式かくだけで謎の罪悪感に包まれました.

{\sum_{i}^n i \ mod\ P_i}を求めたいのですが,{i \ mod\ P_i}はi以下であり{P_i}より小さいので,{n\ mod\ P_n}はnよりだいぶ小さい値になるなーというのを出発点にして考えました.そうなるとできるだけ{i = n-1}以下の組み合わせについて最大化したくて,ひとつずつ考えていくと例えば{n-1\ mod\ P_{n-1}}を最大化するのは{P_{n-1}=n}{n-1\ mod\ P_{n-1} = n-1}になります.このように後ろから見ていくとiより1大きい数を{P_i}として選べば,{\sum_i^{n-1} i = \frac{n(n-1)}{2}}で最大化されるのかな??と思って実装したらACでした.

正直厳密な証明に裏付けされた答えじゃなくて,直感に基づく実装だったので実力かと言われると微妙ですが...

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>

using namespace std;

typedef long long ll;

int main() {
    ll N; cin >> N;
    cout << (N-1) * N / ll(2)<< endl;
}

E問題

f:id:kuroneko1259:20190903093430p:plain
Dまでで30分くらいだったので,E問題も解いて見たら意外と解けました.

  1. 各選手の試合順が決まっているので,各試合における最短日数は一意にきまりそう(本当にそうなのかは不明)
  2. しかも各試合の最短日数は,出場選手が前の試合を行った日にのみ依存しそう(その前にしなければならない試合があれば,その試合の最短日数が決まるまで待てば良い)
  3. あれ??でもそれだとforで各選手を回して,その外側にwhile回さないといけない??時間足りる??
  4. あ,でも最短日の更新はたかだか試合数分しか起こらないから大丈夫かな?更新がなくなってwhile抜ければ問題なさそう

みたいな思考回路でDPで解きました.実装にもわりと時間がかかって,入力は列が試合順のインデックスにしてるのに,dpは選手のインデックスでこんがらがったりしてました.

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>

using namespace std;

typedef long long ll;

int dp[1005][1005];

int main() {
    int N; cin >> N;
    int A[N][N-1];
    for (int i = 0; i < N; i++) {
        for (int j = 0; j < N-1; j++) {
            cin >> A[i][j];
            A[i][j]--;
        }
    }
    int next[N];
    for (int i = 0; i < N; i++) next[i] = 0;

    bool flag = false;

    while (!flag) {
        for (int i = 0; i < N; i++) {
            int idx = A[i][next[i]];
            if (next[i] >= N-1 || next[idx] >= N-1) continue;
            int prev_i, prev_idx;
            if (next[i] > 0) prev_i = A[i][next[i]-1]+1;
            else prev_i = 0;

            if (next[idx] > 0) prev_idx = A[idx][next[idx]-1]+1;
            else prev_idx = 0;

            if (A[idx][next[idx]] == i) {
                dp[i+1][idx+1] = dp[idx+1][i+1] = max(dp[i+1][prev_i], dp[idx+1][prev_idx]) + 1;
                next[i]++;
                next[idx]++;
                flag = true;
            }
        }
        flag = !flag;
    }

    int res = 0;

    for (int i = 0; i < N; i++) {
        for (int j = 0; j < N; j++) {
            if (i != j && dp[i+1][j+1] == 0) {
                cout << -1 << endl;
                return 0;
            } else if (i != j) {
                res = max(res, dp[i+1][j+1]);
            }
        }
    }

    cout << res << endl;
    return 0;
}

まとめ

今回は初めてE問題まで解けてかなりうれしかったです.レートもかなり上がって茶色の上の方までいきました.最初の10回までは補正がかかるらしいので(今6回目),とりあえずレートを気にせず10回参加することを目標に頑張りたいと思います.

ABC137振り返り

ABC137の振り返りをしていきます.今回はC問題も解けず,,,まだまだ使える典型の種類が足りてないなとおもいます...

A問題

f:id:kuroneko1259:20190811140641p:plain

これはさすがにそのまま書くだけでした.特に何もないです.

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>

using namespace std;

typedef long long ll;

int main() {
    int A, B; cin >> A >> B;
    int p = A + B, s = A - B, m = A * B;

    cout << max(p, max(s, m)) << endl;
    return 0;
}

B問題

f:id:kuroneko1259:20190811140811p:plain

制約から範囲がはみ出す可能性はないので,これも素直に書きました.

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>

using namespace std;

typedef long long ll;

int main() {
    int K, X; cin >> K >> X;

    for (int i = X - K + 1; i < X + K; i++) {
        cout << i;
        if (i != X + K - 1) cout << " ";
    }
    cout << endl;
    return 0;
}

C問題

f:id:kuroneko1259:20190811141023p:plain

解けなかったです...

とりあえず,アナグラムを判定する部分と組み合わせを考える部分で分けて考えました.アナグラムの判定は,あらかじめ文字列をソートしておくと定数時間で実行できるので,問題なかったです.

アナグラムになっている文字列をグループ分けしてあげて,同じグループに{n \geq 2}個の文字列があるときに,{{}_n C_2}の総和を出せばいいことに終了10分前くらいに気づきました.

グループ分けは,文字列の配列をソートしてO(NlogN),前から順にO(N)でみていけばよかったので,あとは{{}_n C_2}の総和をとって提出!!したのですが,最後のテストケースだけWAでした.後から考えてみると,組み合わせを{1e9 + 7}で割った余りを出してしまっていたのでそれが原因だと思います.

アナグラムである文字列の個数の数え方も,単純に足して行って,違う文字列がでてきたらリセットでよかったんですが,配列やらフラグやらをこねくり回してわかりにくくなってます.

解説によれば,同じ文字列の個数を数える際は,ハッシュを使うと高速にできるみたいです(検索も追加も平均O(1) ??全てのデータが衝突したとしてO(N)).C++ではunordered_mapで使えます.

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>
#include <string>

using namespace std;

typedef long long ll;

int main() {
    int N; cin >> N;
    vector<string> s(N);

    for (int i = 0; i < N; i++) {
        string st; cin >> st;
        sort(st.begin(), st.end());
        s[i] = st;
    }

    sort(s.begin(), s.end());

    bool prev = false;
    vector<int> res;
    for (int i = 0; i < N-1; i++) {
        bool flag = true;
        for (int j = 0; j < 10; j++) {
            if (s[i][j] != s[i+1][j]) flag = false;
        }
        if (flag) {
            if (prev) {
                res[int(res.size())-1] += 1;
            } else {
                res.push_back(2);
            }

            prev = true;
        } else {
            prev = false;
        }
    }
    ll ans = 0;
    for (int i = 0; i < res.size(); i++) {
        ans += ll(res[i]) * (res[i] - 1) / 2;
    }

    cout << ans << endl;
    return 0;
}

D問題

f:id:kuroneko1259:20190811143619p:plain

C問題が解けそうで解けなかったのでこっちも目を通しましたが,こっちも解けなかったです...

時間がたつにつれてできる仕事の選択肢は減っていく(制約が厳しくなっていく)ので,こういう場合は制約が厳しい方から考えるのが定石みたいです.例えば時刻tのときを考えると,選べる仕事はM - t以内に報酬が振り込まれるものになるので,その中で報酬が最大のバイトを選べば良いみたいです.

これに気づいても,愚直に最大値を調べるとO(N^2)かかると思うのですが,毎回選ぶのは最大値なのでpriority_queueを使うと,要素の追加がO(N),最大値を取り出すのがO(M),それぞれの操作はO(logN)かかるのでO((M+N)logN)で計算できます.

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>
#include <queue>

using namespace std;

typedef long long ll;

int main() {
    int N, M; cin >> N >> M;
    vector<pair<int, int>> v(N);

    for (int i = 0; i < N; i++) {
        int a, b; cin >> a >> b;
        v[i].first = a;
        v[i].second = b;
    }
    sort(v.begin(), v.end());

    priority_queue<int> q;

    int idx = 0;
    int ans = 0;
    for (int t = 0; t <= M; t++) {
        while (idx != N && v[idx].first <= t) {
            q.push(v[idx].second);
            idx++;
        }

        if (!q.empty()) {
            ans += q.top();
            q.pop();
        }
    }

    cout << ans << endl;
    return 0;
}

まとめ

まずは最低でもCまでは確実に解けるようにしていきたいです

AtCoder ABC136の振り返り

ABC136に参加したので,復習がてら振り返ります.

A問題

f:id:kuroneko1259:20190806192219p:plain

普通に{C - (A - B)}をしました.負になる場合だけ気をつけました.こういう簡単な時にさらっと三項演算子使えるようになりたいです.

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>

using namespace std;

typedef long long ll;

int main() {
    int A, B, C; cin >> A >> B >> C;

    if (A - B > C) cout << 0 << endl;
    else cout << C - A + B << endl;

    return 0;
}

B問題

f:id:kuroneko1259:20190806192555p:plain

一度文字列になおせばsize()の判定でいけるじゃん!と思いましたが,C++で数字から文字列に直す方法を知らなかったので(調べればよかったですが),数字のまま処理しました.1から10倍していって,奇数回目の分だけ数えました.

文字列に直す場合,to_string()なる便利なものがあるみたいです.常識かもしれませんが...

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>
#include <sstream>

using namespace std;

typedef long long ll;

int main() {
    int N; cin >> N;

    int cnt = 0, i = 1, j = 1;
    bool flag = true;

    while (1) {
        i *= 10;
        if (i <= N) {
            if (flag) cnt += i - j;
        } else {
            if (flag) {
                cnt += N - j + 1;
            }
            break;
        }
        j = i;
        flag = !flag;
    }

    cout << cnt << endl;
}

C問題

f:id:kuroneko1259:20190806193204p:plain

最初は,{H_i - H_{i+1} = 1 かつ H_{i-1} \geq H_i - 1}のときに{H_i}から1をひく方針でいきましたが,それだと1332とかが無理になるので,少し悩みました.

最初からみていく場合,1引ける場合は引いておいて損はない(引く操作が最適??)ので,その方針でいけました.回答は後ろからみていってました.確かに後ろからだと1332みたいな状況でもいけますね.

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>

using namespace std;

typedef long long ll;

int main() {
    int N; cin >> N;
    vector<int> a(N);
    for (int i = 0; i < N; i++) cin >> a[i];

    a[0] -= 1;
    for (int i = 1; i < N - 1; i++) {
        if (a[i] - a[i+1] > 1 || (a[i] - a[i+1] == 1 && a[i] - a[i-1] <= 0)) {
            cout << "No" << endl;
            return 0;
        } else if (a[i] - a[i+1] <= 1 && a[i] - a[i-1] > 0){
            a[i] -= 1;
        }
    }

    cout << "Yes" << endl;
    return 0;
}

D問題

f:id:kuroneko1259:20190806194304p:plain

今回は本番中に初めてD問題が解けました!といっても時間ギリギリですが...

子供の人数に比べて移動回数が十分に大きいので,RLのようにRとLが切り替わるところに左右の子供がたまっていくことはすぐにわかりました.切り替わりの位置との相対距離の偶奇によってRLのどっちにいるのか判断できることも割とさくっと気づきました.が,実装力が足りず,実装にとても時間がかかりました.

RからLに切り替わる位置を覚えておいて,前から順番に子供のいる位置を決めていきましたが,もっといい実装の仕方があると思います.解説のように(R, 3), (L, 5)...みたいに保持するやり方を試しに実装したいですが,時間がある時&気が向いたらやりたいです.

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>

using namespace std;

typedef long long ll;

int main() {
    string s; cin >> s;

    vector<int> ans(int(s.size()));

    vector<int> ch;
    for (int i = 1; i < s.size(); i++) {
        if (s[i] == 'L' && s[i-1] == 'R') ch.push_back(i);
    }

    int tmp = 0, i = 0;
    while (i != ch.size()) {
        while (tmp != s.size()) {
            if (s[tmp] == 'R' && tmp < ch[i]) {
                if ((ch[i] - tmp) % 2 == 0) ans[ch[i]]++;
                else ans[ch[i]-1]++;
            } else if (s[tmp] == 'L' && tmp >= ch[i]){
                if ((tmp - ch[i]) % 2 == 0) ans[ch[i]]++;
                else ans[ch[i]-1]++;
            } else {
                i++;
                break;
            }
            tmp++;
        }
        if (tmp == int(s.size())) break;
    }


    for (int i = 0; i < s.size(); i++) {
        cout << ans[i];
        if (i != int(s.size())-1) cout << " ";
    }
    cout << endl;
    return 0;
}

AtCoder Educational DP Contest:O問題(bit DP)

DPコンテストのO問題に引っかかったので,まとめます.

f:id:kuroneko1259:20190801183729p:plain

方針

男女N人ずつのマッチングですが、例えば男性iを昇順に固定した場合、女性N人の順列を試せば全通り考えられます。しかし,O(N!)の時間がかかり間に合わないです。一般に,順列などO(N!)をO(2^N)に落とすテクニックとしてbit DPなるものがあるようなのでそれを使っていきます。

bit DP

N個の要素の順列で条件に合うものを数え上げたりするときに、順列全ての情報を保持する必要がなく,どの要素をつかっているかを保持したので十分だという場合に使えるみたいです。どの要素を使っているかは,N個の要素の部分集合として表せ,bit列を用いることで割と簡単に表現できます。

考え方

例えば,3人の男女のマッチングで3組のペアが何通りできるかを考えます。このとき,どの女性がマッチングに参加しているかをbit列で101のように表現することにします。男は昇順に固定されていることを考えると,この例でいえば女1と女3と,男1,男2のマッチングを考えることになります。また,bit列で表される女M人がマッチングに参加した場合に,M組みのペアがdp[bit]通り考えられるとします.

3人全員がマッチングに参加するとき,bit列は111になります.漸化式を立てたいので、ペアを1組固定することで2人でのマッチングとの関係を考えます。ペアを1つ固定する上で,男は昇順に固定されているので当然男3になりますが、女は3人のうち1人を選ぶことになります。

仮に女1が選ぶとしましょう。このとき女1以外の参加者のbit列は011になります。男3と女1が実際にペアになれるとき、女1を除いた2人のマッチングにこのペアを追加するだけで3人のマッチングが完成することから,男3と女1がペアになる3人のマッチングはdp[011]通りあることになります(女2と女3での2組のマッチングが全くない場合は男3と女1がペアになっても0通りですが,dp[011]=0が成り立っているので問題ないです)。

女1が選ばれたのに男3と女1がペアになれない場合、女2と女3がペアになれたところで3人のマッチングは不可能なので0通りです。

これらをマッチングに参加している女性すべてでやっていくと,以下の漸化式が立てられます。このとき男iと女jがペアになれるとき1,なれないとき0になる関数として{a[i][j]}を定義します。

{dp[111] = dp[011] * a[3][1] + dp[101] * a[3][2] + dp[110] * a[3][3]}

これを一般化すると,bit列において1の数をi個としたとき(参加している女・男の人数はbit列の1の数に対応),参加しているすべての女と男iについて上記のことをやっていけばよいわけです。

この漸化式を使って動的計画法で実装します。bitは10進数で{2^N(1 << N) - 1}まで1ずつ足していけばすべてのbit列を網羅できます。

このとき、

  • 1 << i:右からi+1番目の桁だけ1をたてる
  • bit xor (1 << i):bitのi+1番目の桁を反転させる

ということに注意すると

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>

using namespace std;

typedef long long ll;

#define MOD 1000000007LL

ll dp[1<<25];
int a[21][21];

int main() {
  int N; cin >> N;
  for (int i = 0; i < N; i++) {
    for (int j = 0; j < N; j++) {
      cin >> a[i][j];
    }
  }
  dp[0] = 1;
  for (int bit = 0; bit < (1<<N); bit++) {
    int pos = __builtin_popcount(bit) - 1;
    
    for (int i = 0; i < N; i++) {
      if ((bit & (1<<i)) && a[pos][i]) {
        dp[bit] += dp[bit ^ (1<<i)];
        dp[bit] %= MOD;
      }
    }
  }
  cout << dp[(1<<N) - 1] << endl;
  return 0;
}

みたいになります。

まとめ

見返してて、自己満な解説になってしまった感がありますが、とりあえずこれでいいや...同じようなレベルの人の少しでも参考になれば嬉しいです。

巡回セールスマン問題とかいろいろな場面でbit DPが使えるみたいです。

AtCoder Educational DP Contest:J問題

AtCoderをやっていて、Educational DP ContestのJ問題に詰まったので復習がてら自分なりにまとめて見ます。数学的に厳密なものではないのでご了承ください。

f:id:kuroneko1259:20190731115328p:plain

文字の定義

  • {(i, j, k)}:1個の寿司が乗っている皿がi個、2個の寿司が乗っている皿がj個、3個の寿司が乗っている皿がk個ある状態
  • {dp[i][j][k]}:(i, j, k)の状態からすべて食べる場合の試行回数の期待値
  • N:皿の数
  • {sm = i + j + k}:寿司が乗っている皿の数

めんどうなので,ある状態からすべて寿司を食べるまでに必要な試行回数の期待値のことを、ある状態の期待値と表現します。

考え方

寿司を1個食べるための試行回数の期待値を使って帰納的にときます。
状態(i, j, k)を考えると、

{dp[i][j][k] = ((i, j, k) から1個寿司を食べた後の状態の期待値) + ((i, j, k)から1個寿司を食べるための試行回数の期待値)}

という漸化式が成り立ちます。これを元に、右辺をそれぞれ求めれば動的計画法で解けそうです。

右辺第1項の計算

(i, j, k)からの遷移先には(i-1, j, k), (i+1, j-1, k), (i, j+1, k-1)の3つがあります。食べる皿は等確率で選ぶので、(i, j, k)の遷移先のうち、それぞれの状態である確率は{\frac{i}{sm}, \frac{j}{sm}, \frac{k}{sm}}となります。よって(i, j, k)から1個寿司を食べていける状態の期待値は

{dp[i-1][j][k] \times \frac{i}{sm} + dp[i+1][j-1][k] \times \frac{j}{sm} + dp[i][j+1][k-1] \times \frac{k}{sm}}

となります。

右辺第2項の計算

{1 \times \frac{sm}{N} + 2 \times (1 - \frac{sm}{N}) \frac{sm}{N} + 3 \times (1 - \frac{sm}{N})^2 \frac{sm}{N}  \cdots}

となります。(1個の寿司を食べるためにn回かかるということは、n-1回はすべて寿司が乗っていない皿を選び続けなければならない)

一般にqが1未満の時、{1 \times p + 2 \times q p + 3 \times q^2 p \cdots = \frac{p}{(1 - q)^2}}が成り立つので、(i, j, k)から1個寿司を食べるための試行回数の期待値は{\frac{N}{sm}}となります。これが漸化式の右辺第二項になります。

漸化式

以上を踏まえると漸化式は

{dp[i][j][k] = dp[i-1][j][k] \times \frac{i}{sm} + dp[i+1][j-1][k] \times \frac{j}{sm} + dp[i][j+1][k-1] \times \frac{k}{sm} + \frac{N}{sm}}

となります。あとはこれを使って動的計画法で解けると思います。私はメモ化再帰で書きました。直感的にメモ化再帰がわかりやすいものと、ループ回す方がわかりやすいものがありませんか...?この問題はループ回すやり方の直感的なイメージが湧かなかったのでメモ化再帰をつかいました。

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>

using namespace std;

typedef long long ll;

#define INF 1000000009LL

double dp[305][305][305];//dp[i][j][k]:1個の寿司がi個,2個の寿司がj個,3個の寿司がk個あるときすべて食べ終わるまでの試行回数の期待値
int N;

double dfs(int i, int j, int k) {
    if (dp[i][j][k] > -1) return dp[i][j][k];

    double res = 1.0 * N / (i + j + k);
    if (i > 0) res += 1.0 * i / (i + j + k) * dfs(i-1, j, k);
    if (j > 0) res += 1.0 * j / (i + j + k) * dfs(i+1, j-1, k);
    if (k > 0) res += 1.0 * k / (i + j + k) * dfs(i, j+1, k-1);

    return dp[i][j][k] = res;
}

int main() {
    cin >> N;
    int one = 0, two = 0, three = 0;
    for (int i = 0; i < N; i++) {
        int a; cin >> a;
        if (a == 1) one++;
        else if (a == 2) two++;
        else three++;
    }

    for (int i = 0; i <= N; i++) {
        for (int j = 0; j <= N; j++) {
            for (int k = 0; k <= N; k++) {
                dp[i][j][k] = -1;
            }
        }
    }

    dp[0][0][0] = 0;
    cout << fixed;
    cout << setprecision(10) << dfs(one, two, three) << endl;

    return 0;
}

まとめ

Educational DP Contestをアルファベット順に解いていますが、J問題で大きくつまずきました。中級者以上の方はこの問題はさくっと解けてしまうんでしょうか...地道に典型的な問題や考え方を身につけていきたいです。