グラフィカルラッソやってみた - 星を求める蛾の願い

けっこう前の話になるけど、
ふと思い立ったのでグラフィカルラッソのアルゴリズムを実装してみた。

入力として、なにかの分散共分散行列 $S$ を貰った時に、
$\min_\Lambda \mathrm{tr}\left(\Lambda S \right) -\log \Lambda S + || \Lambda ||_1$
をそこそこ効率的にもとめるアルゴリズムのこと。

数理情報工学特論第一【機械学習とデータマイニング】
http://www.geocities.co.jp/Technopolis/5893/course_machine_learning_2010.html

にあるスライドあたりを参考にしてます。

#coding: utf-8

from numpy import *
from numpy.linalg import *

def GraphicalLasso(S,l):
    (dim,m) = S.shape

    #対角成分をまずセットする
    L = diag(diag(S)) + 2*l

    sdiag = diag(S)
    for k in xrange(dim):
        if k==0:
            continue
        print k

        beta=zeros(k)
        sk = S[k,:k]
        sdiagk = sdiag[:k]

        #列ごとの最適化。
        epsilon=1e-8
        oldbeta=beta
        import itertools
        for time in itertools.count():
            beta = (sk - dot(S[:k,:k],beta) + beta*sdiagk) / sdiagk

            #shrinkage
            gamma = 2.0 * l / sdiagk
            big = where(abs(beta)>gamma)
            small = where(abs(beta)<=gamma)
            print big
            beta[big] -= (gamma*sign(beta))[big]
            beta[small] = 0

            #収束条件は適当に。
            if norm(beta-oldbeta) < epsilon:
                break
            oldbeta[:]=beta

        Lnew = -L[k,k]*beta
        L[k,:k] = Lnew
        L[:k,k] = Lnew

    return L

これでそこそこ動いているようだ。