αダイバージェンスやβダイバージェンスでのNMFを一般化するために,αβダイバージェンスというものを定義します.入力行列と出力低ランク行列間のαβダイバージェンスが小さくなるようなNMFアルゴリズムをJuliaで実装しました.αβダイバージェンスやアルゴリズムは以下の論文の通りに実装しました.
Entropy | Free Full-Text | Generalized Alpha-Beta Divergences and Their Application to Robust Nonnegative Matrix Factorization
あんまりちゃんとテストしてないので,使うときは気を付けてね.
とりあえず,αβダイバージェンスD_abを定義に従って計算できるようにD_ab関数を用意しておく.α=β=1で通常のフロベニウスノルムになる.
using LinearAlgebra
function D_KL(P, Q)
n, m = size(P)
d_kl = 0.0
for i = 1:n
for j = 1:m
d_kl += P[i,j] * log( P[i,j] / Q[i,j] ) - P[i,j] + Q[i,j]
end
end
return d_kl
end
function D_IS(P, Q)
n, m = size(P)
d_is = 0.0
for i = 1:n
for j = 1:m
d_is += log( Q[i,j] / P[i,j] ) + P[i,j] / Q[i,j] + - 1.0
end
end
return d_is
end
function D_logeuc(P, Q)
n, m = size(P)
log_P = log.(P)
log_Q = log.(Q)
return sum( (log_P - log_Q).^2 )
end
function D_ab(P, Q, alpha, beta)
if alpha != 0 && beta == 0
return 1.0 / alpha^2 * D_KL(P.^alpha, Q.^beta)
end
if alpha + beta == 0 && alpha !=0 && beta !=0
return 1.0 / alpha^2 * D_IS(P.^alpha, Q.^alpha)
end
if alpha == 0 && beta != 0
return 1.0 / beta^2 * D_KL(Q.^beta, P.^beta)
end
if alpha == 0 && beta == 0
return 1.0 / 2.0 * D_logeuc(P, Q)
end
n, m = size(P)
d_ab = 0.0
for i = 1:n
for t = 1:m
d_ab += P[i,t]^alpha * Q[i,t]^beta
d_ab += -alpha / (alpha + beta) * P[i,t]^(alpha+beta)
d_ab += -beta / (alpha + beta) * Q[i,t]^(alpha+beta)
end
end
return -1.0 * d_ab / (alpha*beta)
end
あとは論文の更新ルールをそのまま実装する.
function nmf(P, r ; alpha=1, beta=1, max_iter=200, tol = 1.0E-3, verbose = false)
n, m = size(P)
A = rand(n, r)
X = rand(r, m)
one_mn = ones(n, m)
error_at_init = D_ab(X, A*X, alpha, beta)
previous_error = error_at_init
for iter in 1:max_iter
Q = A*X
Z = P.^alpha .* Q.^( beta-1 )
X .= X .* ( (A'*Z) ./ (A'*Q.^(alpha+beta-1))).^(1.0/alpha)
Q = A*X
Z = P.^alpha .* Q.^( beta-1 )
A .= A .* ( (Z*X') ./ (Q.^(alpha+beta-1) * X')).^(1.0/alpha)
if tol > 0 && iter % 10 == 0
error = D_ab(P, A*X, alpha, beta)
if verbose
println("iter: $iter cost: $error")
end
if (previous_error - error) / error_at_init < tol
break
end
previous_error = error
end
end
return A, X
end
A*X
が入力行列P
を近似するr
ランク行列となる.
追記
私の記憶が正しければKL-NMFは行和列和を保持する.しかし,このアルゴリズムでα=1,β=0にしても行和列和が保持されない.何か間違えていると思う.そもそも論文の式(71)なんかおかしい気がするのは気のせいだろうか.
追記2
以前アップロードしたNMFの記事のコードと,上のコードとで同じ行列を入れても別の答えが返ってくる.やはりこの実装あやしい...(さすがに論文が間違っている可能性は低いと思うので,私の実装がおかしいのだろう...)
genkaiphd.hatenablog.com
追記3
α=1, β=0で更新式が以前の記事で引用した文献にある式にならない.Normalize Factor分ずれている気がする.めんどくさいが,今度,自分で更新式を導いて,どこが怪しいか調べる.