2015-05-10

ABC023復習

　昨日AtCoder Beginner Contest 023 に参加しました。結果はボロボロ..orzせめてC問題はちゃんと解きたかったなあ。

<a href="http://abc023.contest.atcoder.jp/">Welcome to AtCoder Beginner Contest 023 - AtCoder Beginner Contest 023 | AtCoder</a>abc023.contest.atcoder.jp

珍しくちゃんと全部復習したので上げることにしました。ブログを書くことで復習へのモチベが上がると嬉しいですね。

A：加算王

　まあワンライナーでいきますよね

print sum(map(int, list(raw_input())))

B：手芸王

　解説では素直にシミュレーションしても十分間に合うってことだったんですが、pythonでstringの連結はあんまりやりたくなかったのでこのような形にしました。真ん中から見て行って、不正なブロックをくっつけている場合はその時点で-1を返します。文字列作ってないだけでやってることは一緒ですね。

N = int(raw_input())
S = list(raw_input())
block = ["a", "b", "c"]
if N % 2 == 0:
    print -1
else:
    center = N / 2
    for i in range(center+1):
        if S[center-i] == block[(1-i) % 3] and S[center+i] == block[(1+i) % 3]:
            continue
        print -1
        break
    else:
        print center

C：収集王

　全マスをチェックするのは計算量的に絶対無理って分かってたんですが、良い手が思いつかず部分点だけ取りに行きました。最後の二重カウントを省くのがN回回すだけで済むとなぜ気付かなかったのか..

from collections import defaultdict as ddict
R, C, K = map(int,raw_input().split(" "))
N = int(raw_input())
rc = [0 for i in range(R)]
cc = [0 for i in range(C)]
candy = [None for i in range(N)]
for i in range(N):
    r, c = map(int,raw_input().split(" "))
    candy[i] = (r-1,c-1)
    rc[r-1] += 1
    cc[c-1] += 1
rn = ddict(int)
cn = ddict(int)
for i in range(R):
    rn[rc[i]] += 1
for i in range(C):
    cn[cc[i]] += 1
count = 0
for key in rn:
    count += rn[key] * cn[K-key]
for r, c in candy:
    if rc[r] + cc[c] == K:
        count -= 1
    elif rc[r] + cc[c] == K + 1:
        count += 1
print count

D：撃墜王

　これもオーダーがきつかったですね。最初から部分点狙いでbit-DPしてみたけどそれでも間に合わないっていう(^^;;) chokudaiさんも仰ってましたけど、二分探索をちゃんと使いこなせるようになりたいです。

def check(X, H, S, N):
    time_limit = [0 for i in range(N)]
    for i in range(N):
        tmp = (X - H[i]) / S[i]
        if tmp < 0:
            return False
        if tmp < N:
            time_limit[tmp] += 1
    count = 0
    for i in range(N):
        count += time_limit[i]
        if count > i + 1:
            return False
    else:
        return True

N = int(raw_input())
H = [0 for i in range(N)]
S = [0 for i in range(N)]
for i in range(N):
    H[i], S[i] = map(int,raw_input().split(" "))
left = max(H) - 1
right = max(H[i] + S[i]*(N-1) for i in range(N))
while right - left > 1:
    mid = (left + right) / 2
    if check(mid, H, S, N):
        right = mid
    else:
        left = mid
print right

いつになったらビギナーを脱せるんですかねえ。まあ頑張ります。

2015-05-09

言語処理の為の機械学習入門をつまみ食い(その1)：パラメータの最尤推定

　最近MLaPP(Machine Learning: a Probabilistic Perspective)の輪講が我がラボで始まったのですが、話が大分抽象的+英語で早くもヤバい空気が流れ始めました(^^;)

そこで、ちょっと機械学習で使う確率系の話を復習する為に、

言語処理のための機械学習入門 (自然言語処理シリーズ)

を読み直そうと思います。この本は二年前から所持しており分からなくなったときにちょくちょく読んでいるのですが、そろそろまとめ的なのを作ろうかと思いここに書くことにしました。

あくまで入門なのですが、簡潔に書いてある＋日本語（超重要！）であるところが素晴らしい。もちろんさらに深層に挑むなら、MLaPPもそうだけどPRML読めって話ですけどね。

そんなこんなでちょびちょびやって行こうと思います。

具体例：コイントスはフェアなのか？

　コイントスに使うコインは表と裏が出る確率が同じだからこそ意味がありますが、ふとイカサマを疑うときがあります。そこで、実際に10回コイントスして調べてみましょう。

コインA : 表, 裏, 表, 裏, 表, 表, 表, 裏, 裏, 裏
コインB : 表, 裏, 裏, 裏, 裏, 表, 裏, 裏, 裏, 裏

すごく大雑把ですが、コインAは表と裏が出る確率はそれぞれ0.5だと考えられ、フェアと言えそうです。逆にコインBは表と裏が出る確率はそれぞれ0.2と0.8、コインへの細工を疑いますね笑

このように、実際に試行したデータを用いることで我々は直感的にコインの表の出る確率（と裏の出る確率）を調べる事が出来ます。今回のテーマである最尤推定法は、すごくシンプルに言えばこのコインの表の出る確率を求める手法です。

コイントスの確率分布の記述(読み飛ばし可)

　コイントスは取り得る値(確率変数X)が表か裏かの二つしかなく、確率θで表が、確率(1-θ)で裏が出ると書けます*1(θは確率分布のパラメータ。 0 < θ < 1)。すなわち確率変数Xが値x(表か裏)をとる確率 $p(X=x;\theta)$ は以下のように書けます。

$p(X=x;\theta) = \delta(x,表)\theta + \delta(x,裏)(1-\theta) \\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ = \delta(x,表)\theta + (1-\delta(x,表))(1-\theta) \\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ = \theta^{\delta(x,表)}(1-\theta)^{1-\delta(x,表)}$

ここで $\delta(x,a)$ はデルタ関数と呼ばれ、x=aのときに1を、そうでないときに0を返す関数。

※要するに $p(X=表;\theta) = \theta$ 、 $p(X=裏;\theta) = 1-\theta$ になるってことです。これさえ分かってれば問題ないです。

尤度とは

　さて先ほど $p(X=x;\theta)$ を記述しましたが、これにより実データが生成される確率を記述することが出来ます。ここで、コイントスは一回ごとに結果が独立であり、かつ全てが同一の確率分布に従っている(それぞれの試行の確率を $p(X=x;\theta)$ で書ける)*2ので、実データ $D = [x_1, x_2, ... , x_n$ ]の生成確率 $P(D)$ は

$P(D)=\prod_{i}^{n}p(x_i)$

と書けます。それぞれのコイントスの確率の積になってますね。この実データDの生成確率を尤度と呼びます。通常は積の形だと扱いにくいので、対数を取った値、

$\log P(D) = \log \prod_{i}^{n}p(x_i) \\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ = \sum_{i}^{n} \log p(x_i)$

を用いることが多く、これを対数尤度と呼びます。

最尤推定

　やっと本題に入れます。最尤推定法とは、尤度が最も大きくなるようにパラメータ(今回はθ)を決定する方法です。対数尤度の方が扱いやすいため、通常は対数尤度を最大化することが多いです。尤度（実際のデータが生成される（起こる）確率）が最も大きくなるようにパラメータを決定することで、パラメータを現実のデータにフィットさせている訳です。

さてコイントスの話に戻ると、コインAのデータ $D_A = [表, 裏, 表, 裏, 表, 表, 表, 裏, 裏, 裏$ ] の対数尤度は

$\log P(D_A) = \sum_{i=1}^{10} \log p(x_i) = 5\log p(表) + 5\log p(裏) = 5(\log\theta + \log(1-\theta))$

これを微分すると

$\displaystyle \frac{d log P(D_A)}{d \theta} = 5( \frac{1}{\theta} - \frac{1}{1-\theta} ) = \frac{5(1-2\theta)}{\theta(1-\theta)}$

となります。直感の通り、このデータにとって最も尤もらしい(※駄洒落ではない)表の出る確率(θ)は0.5となります。

同様にコインBのデータ $D_B =[表, 裏, 裏, 裏, 裏, 表, 裏, 裏, 裏, 裏$ ] についても

$\log P(D_B) = \sum_{i=1}^{10} \log p(x_i) = 2\log p(表) + 8\log p(裏) = 2(\log\theta + 4\log(1-\theta))$

$\displaystyle \frac{d log P(D_B)}{d \theta} = 2( \frac{1}{\theta} - \frac{4}{1-\theta} ) = \frac{2(1-5\theta)}{\theta(1-\theta)}$

となり、これもまた直感の通りθ=0.2となっていることがわかります。非常にシンプルな例ではありますが、最尤推定法でパラメータを求めることが出来ました！

終わりに

　繰り返しになりますが、最尤推定法は実際の観測データにフィットするようにパラメータを決定する方法です。今回の例は非常にシンプルですが、もっと複雑な問題を扱う場合も根本の考えは一緒だと思います。

　ところで、たとえイカサマをしていないコインであっても、10回の試行中偶然にも裏を9回出した場合は最尤推定法でθ=0.1が求まります。実データからはこう求まっても、コインは平等だと思っている我々に取っては少し納得がいきません。θが大体0.5位になりそうだと分かっている場合に、それをθに反映させたいと思うときがあるわけです。こんな場合に用いるのが最大事後確率推定(MAP推定)なのですが、これについてはまたの機会に。

*1:このような確率分布をベルヌーイ分布(Bernoulli distribution)と言う

*2:独立に同一の確率分布に従う(independently, identically distributed)。 i.i.d.としばしば略記される

2015-05-08

ブログ始めました

　技術ブログっぽいものをかいて備忘録にしよう！と仲間内で言い続けて早半年あまり。ついに始めることとなりました。

ゆるゆるとプロコンとか勉強関連の話を書いて行けたらなと思っております。よろしければお付き合い下さい。

俵言

しがない社会人が書く、勉強とかのこと。最近は機械学習や kaggle 関連がメイン。