잡담
이 포스트를 작성하면서, WMF 논문을 다시 천천히 읽어봤는데, 이 논문은 정말 잘 쓴 논문이라는 생각이 문득 들었다. (1) 사람들이 잘 인지하지 못했던 기존 방법(explicit feedback methods)의 문제점을 발견하고, (2) 문제를 이해하는 새로운 도구(preference-confidence split)를 제시하고, (3) 그 도구를 이용한 새로운 해결 방법(WMF)을 제안했으며, (4) 이 해결 방법은 너무나도 아름답다(ALS, explaining recommendations). 게다가, (5) 성능도 우수하다(performance, scalability 양쪽 다…).
저자는 얼마나 대단한지, 심지어 내가 잘 쓴 논문이라고 생각한 이유를 논문에 이미 적어놓기까지 했다!
We provide a latent factor algorithm that directly addresses the preference-confidence paradigm. Unlike explicit datasets, here the model should take all user-item preferences as an input, including those which are not related to any input observation (thus hinting to a zero preference). This is crucial, as the given observations are inherently biased towards a positive preference, and thus do not reflect well the user profile. However, taking all user-item values as an input to the model raises serious scalability issues – the number of all those pairs tends to significantly exceed the input size since a typical user would provide feedback only on a small fraction of the available items. We address this by exploiting the algebraic structure of the model, leading to an algorithm that scales linearly with the input size while addressing the full scope of user-item pairs without resorting to any sub-sampling.
Y. Hu, Y. Koren, and C. Volinsky. Collaborative filtering for implicit feedback datasets. In IEEE International Conference on Data Mining (ICDM 2008), pages 263–272, 2008.
처음 봤을 때에는 수식과 알고리즘을 이해하기에 벅찼는데, 지금 다시 보니깐, 내가 이해할 수 있었던 논문 중에서는 가장 아름다운 글 중 하나였다고 생각한다. 일을 하며, 연구실에서 계속 WMF를 사용해왔기 때문에, 조금 이해가 깊어진 걸까, 아니면 그냥 내가 조금 더 성장한 걸까. 잘 모르겠다. 나도 언젠가, 꼭 논문이 아니더라도 아름다운 무엇인가를 남기고 싶다는 생각까지 들게 만드는 페이퍼였다.
note: 이 포스트는 읽는 분께서 Implicit matrix factorization, 혹은 다른 이름으로는 Alternating Least Square(ALS) method가 무엇인지 알고 있다고 가정합니다. 또한, 간단한 선형 대수에 관한 지식이 있다고 가정합니다. 추천 시스템, 혹은 implicit matrix factorization에 대해서는 http://sanghyukchun.github.io/73/ 이 포스트에 자세히 설명이 되어 있는 것 같으니 이를 참조해주세요.
Introduction
Implicit Matrix Factorization의 loss function을 Conjugate gradient method를이용해 Optimize하는 방법에 대해 소개합니다.
기계 학습, 데이터 마이닝 분야에서 어떠한 문제를 풀기 위해 상당히 많은 방법이 제안되었고, 실제로 활용되고 있습니다. 지금은 deep learning이 엄청난 인기를 끌고 있어서, 딥러닝에서 가장 많이 이용되는, Gradient descent 기반의 methods만을 알고, 계신 분도 많을 것 같습니다. deep learning에서 이용하는 loss function의 형태는, 다른 방식의 optimization을 이용하기 어렵거나 불가능한 형태라 gradient descent를 사용합니다(정확히 이게 유일한 이유인가는 잘 모르겠습니다.). 반면, 어떠한 형태의 objective function은 gradient descent를 이용하는 것보다, 더 효율적이게 optimize하는 방법들이 알려 있습니다. 이러한 방법들 중 큰 부분을 차지하는 것은 convex한 objective function을 최적화하는 방법인 convex optimization입니다.
Solving ALS
ALS의 objective function 은 다음과 같다.
\[\sum_{u}\sum_{i}{c_{ui} (p_{ui} - x_{u}^{T}y_{i}^{T} )^2} + \lambda(\sum_u{x_u^Tx_u} +\sum_i{y_i^Ty_i})\]Weighted Matrix Factorization을 solving하는 방법에 대해 알고 계신분은, 다음 부분으로 넘어가 주셔도 괜찮습니다.
\(x_u, y_i\)는 각각 유저 $u$와 아이템 $i$에 관한 latent factor이며, $p_{ui}$는 유저 $u$가 item $i$와 interaction이 있었는지 아닌지를 나타내는 Boolean indicator 그리고, $c_{ui}$는 $p_{ui}$에 관한 confidence를 나타낸다. $p_{ui} =1$일 때의 $c_{ui}$가 $p_{ui} = 0$일 때의 $c_{ui}$보다 큽니다. $\lambda$로 시작하는 term은 regularization term입니다.
Let define \(Y = [y_0|y_1,|...,|y_m]\)이라 정의하고, \(x_u\)와 관계가 없는 식을 모두 정리합니다. (objective function을 \(x_u\)로 미분하면 빠지는 term들). 또한, \(C^u=\text{diag}(c_{ui}), \text{ } P_u = [p_{u0}, p_{u1,}, ... p_{um}]^T\)으로 정의하면, objective function을 다음과 같이 reduce할 수 있습니다. \(\sum_ic_{ui}(p_{ui}-x_u^Ty_i)^2+\lambda x_u^Tx_u\) \(= (P_u -Y^Tx_u)^TC^u(P_u - Y^Tx_u)+ \lambda x_u^Tx_u\)
\(=x_u^T(YC^yY^T+\lambda I)x_u- 2{ YC^uP_u }{ x_u }+P_u^TC^UP_u\) $(YC^yY^T+\lambda I)$을 $A$, $YC^uP_u$을 $b$,$P_u^TC^UP_u$을 $c$로 표현한다면, 이 식은 $x_u$에 대한 Quadratic Equation 입을 알 수 있습니다. \(f(x_u) := x_u^TAx_u + 2bx_u + c\) $A$는 Positive definite이며, Symmetric 이다. 따라서, \(\frac{\partial}{\partial x_u} f(x_u) = 2Ax_u - 2bx_u = 0\) 이 equation을 풀면 \(x_u = A^{-1}b = (YC^uY^T+\lambda I)^{-1}{YC^uP_u}\)이 된다.
같은 방식으로 item $i$에 대한 latent factors는 \(y_i = (XC^iX^T+\lambda I)^{-1}{XC^iP_i}\)으로 표현할 수 있다. (\(X, C^i, P_i\)는 이전과 같은 방식으로 정의된다.)
naive update
\(x_u\)를 구하는 과정에서의 시간 복잡도를 분석해보자. 우선, ${YC^iP_u}$는 vector $P_i$에 행렬곱 1번을 한 것이며, 이는 $O(nf)$이다.
- $YC^uY^T$의 계산
- $C^u$는 Diagonal matrix이므로, $YY^T$를 곱하는 데 드는 복잡도 $O(f^2n)$
- $(YC^iY^T+\lambda I)^-1$의 계산
- $(YC^iY^T+\lambda I)$는 $f \times f$ matrix이므로, 역행렬을 구하는 데 드는 복잡도 $O(f^3)$ -[using Cholesky decomposition] (https://en.wikipedia.org/wiki/Cholesky_decomposition#Matrix_inversion)
즉, $x_u$를 업데이트하기 위해, $O(f^2n + f^3)$의 시간 복잡도가 필요하며, 일반적으로, $f « n$이므$O(f^2n)$의 시간 복잡도가 필요하다.
즉, ALS 알고리즘은 다음과 같다.
// ALS algorithm
for k in [0, num_iteration]
for u in [0,n]:
update x_u ///위에 나온 식대로 x_u를 업데이트
for i in [0,m]
update y_i /// 위에 나온 식대로 y_i를 업데이트
게다가, 위 pseudocode에서 알 수 있듯이, $O(f^2n)$의 계산을 m 번, $O(f^2m)$의 계산을 n번 해주어야 하고, 이게 한 iteration이니까, 이 계산을 k번 더 한다고 하면, 알고리즘의 실제 시간 복잡도는 $O(f^2mnk)$가 된다..! 상당히 느리다.
somewhat faster version 1.
위의 naive solution의 문제점은, $YC^uY^T$를 구하는 데에 시간이 너무 많이 걸린다는 점이다. ALS를 제안한 논문에서, 이 $O(f^2n)$의 시간 복잡도를 줄일 수 있는, 아주 우아한 방법을 제안했는데, 이 아이디어는 다음과 같다. \(YC^uY^T\)는 다음과 같이 나타낼 수 있다.
\(YC^uY^T = YY^T + Y(C_u - I)Y^T\) 이 식에서, 눈치를 채야 할 포인트는 이 두가지이다.
- $YY^T$는 유저 $u$가 변하더라도 변하지 않는다.
- 이 부분은 따라서 한번만 계산하면 된다.
- $(C_u - I)$는 유저 $u$와 interaction이 있는 $n_u$개의 item만 value가 있는 sparse diagonal matrix라는 점이다.
-
$\sum_{i \in u(i)}(c_{ui}-1)y_iy_i^T$, where $u$는 유저가 클릭한 아이템의 집합이고, 이 시간 복잡도는 $O(f^2 u )$이며, 를 모든 유저에 대해 더하면 $O(f^2 \nu )$이다. -
$ \sum u = \nu $이기 때문이다. $u$ . 즉, 역행렬을 구하는 것까지 계산에 포함시킨다면, 실제 시간 복잡도는 $O(f^2 (\nu + \mu + f)k) $가 되며, 엄청나게 트레이닝 속도가 빨라진다..!
-
Using Conjugate Gradient
특정한 조건을 만족시키는 Quadratic Equation \(f(x) = x^TAx - 2bAx^T + c\)가 있을 때, 이를 최소화시키는 켤레기울기법이라는 방법이 존재한다. 이는, matrix \(A\)가, symmetric하며, positive definite일 때 적용 가능하다.
Solving Quadratic Equation using Conjugate Gradient
$A$가 symmetric positive definite이기 때문에, $f’(x) = 0$이 되는 값에서 $f(x)$는 최소이다. 또한, $f’(x) = Ax - b = 0$, thus $Ax = b$와 같은 형태가 되며, 이 식을 만족시키는 $x_$가 존재함을 보장한다. 우선, 일반적으로 $Ax-b=0$을 만족시키는 $x_$를 찾는 방법은 $x = A^{-1}b$이렇게, 역행렬을 구하는 방법이 있다. 이는, 위에서 얘기했던 Cholesky decomposition을 이용해 구하는 것이 일반적이다. 또한, gradient based method를 사용할 수도 있다! 위에서 언급했듯이, WMF 알고리즘에서는, $A$와 $b$가 유저마다, 아이템마다 변하지만 $A$의 변화하지 않는 부분($Y^TY$)을 미리 계산해두어 computation efficiency를 높이는 방법을 이용한다.
이런 특수한 케이스(strictly convex한 quadratic equation, 혹은 ridge regression)의 경우를 효율적으로 푸는 방법 중, Conjugate Gradient Method을 이용해 위 식을 최적화하는 방법과, 그리고 이 방법을 어떻게 Weighted Matrix Factorization에 적용하는지에 대해 소개하는 것이 이 이 글의 목표이다.
Conjugate Gradient Method
Note: Conjugate Gradient에 대해 상당히 길게 글을 적었었는데, 적고 나서 느낀 게, Conjugate Gradient의 유도 및 증명을 저도 잘 이해하고 있지 않다는 생각이 들었습니다. Conjugate Gradient 방법을 자세히 설명한 tutorial의 링크를 남깁니다. 여기서는 알고리즘의 목표와, 간단한 과정만 적도록 하겠습니다. (사실 제가 정확히 이해하지 못 한 것 같습니다.)
Conjugate Gradient Method는 $f(x) = Ax - b =0$이 되게 하는 $x^* \in R^n$ 을 $n$ step만에 construct하는 방법이며, $k « n$인 $k$ step만을 진행할 경우, $\hat{x} \simeq x^*$ 인 $\hat{x}$를 찾는다.
Implementations
제가 전산학에서 가장 좋아하는 부분은, 컴퓨터만 있으면(사실 요즘은 모든 집에 컴퓨터 하나씩은 다 있죠…!) 자신이 공부한 것을 언제든 구현해 볼 수 있다는 점인 것 같아요. 특히, 데이터 사이언스에서 다루는 내용들은 Hardware, OS에 independent하게 테스트해 볼 수 있기 때문에 배운 내용을 구현해보는 것을 추천합니다…! 구현된 코드
WMF model의 skeleton은 다음과 같다.
구현의 편의상 positive interaction의 confidence $c_{ui}$는 한 값으로 fix했다. warning: 약간의 변수명 차이가 있습니다.
- 편의상 자주 쓰는 위에 수식을 전개할 때에는, 선형대수, 수치해석학에서 자주 쓰는 notation을 사용했습니다. 실제 구현시에는, Python을 이용한 numerical computing에서 자주 쓰이는 노테이션을 사용했습니다.
- $P \rightarrow X, X \rightarrow U, Y \rightarrow V$
- 글의 윗부분에서 column-major notation을 사용했지만, 컴퓨터 언어에서는 row-major notation으로 vector 및 matrix를 표현하는 것이 편리하기 때문에 row-major notation을 사용했습니다.
def solve_wmf(X, n_factors, confidence, lamb, n_iters, solver, calc_loss = True):
n_users, n_items = X.shape
XT = X.T
U = np.random.normal(0, 0.01, size = (n_users, n_factors))
V = np.random.normal(0, 0.01, size = (n_items, n_factors))
#for conveinience, I ignored weight(confidence) term
def loss(X, U, V, lamb):
X_pred = np.dot(U,V.T)
error = (X_pred - X)
loss = np.sum( np.multiply(error, error))
reg = lamb * ( np.sum(np.multiply(U, U)) + np.sum(np.multiply(V, V)))
return loss + reg
for iter in range(n_iters):
solver(X, U, V, confidence, lamb)
solver(XT, V, U, confidence, lamb)
if calc_loss == True:
print("at iteration %d\t[loss : %f]" %(iter, loss(X, U, V, lamb) ))
user factor와 item factor을 계산하는 solver 함수의 구현이 변경됩니다.
Naive implementation
def _solve_wmf_naive(X, U, V, confidence, lamb):
n_users, n_factors = U.shape
n_users, n_items = X.shape
VT = V.T
for j in range(n_users):
nonzero_indices = get_nonzero_indices(X,j)
# calculate A = (YCuY^T + lamb * I)
# calculate b = (VC_up(j))
A = lamb * np.eye(n_factors)
b = np.zeros(shape=(n_factors,))
for i in range(n_items):
factor = V[i]
A += confidence * np.outer(factor, factor)
for i in nonzero_indices:
factor = V[i]
b += confidence * factor
Precalculating YYT
def _solve_wmf_precalculate(X, U, V, confidence, lamb):
n_users, n_factors = U.shape
VTV = np.dot(V.T, V)
for j in range(n_users):
# calculate A with precalculated VTV
# calculate B as same way
A = VTV + lamb * np.eye(n_factors)
b = np.zeros(shape = (n_factors,))
nonzero_indices = get_nonzero_indices(X, j)
for i in nonzero_indices:
item_factor = V[i]
A += (confidence - 1) * np.outer(item_factor, item_factor)
#we can bundle two loop into one
#for i in nozero_indices:
b += confidence * item_factor
U[j] = np.linalg.solve(A, b)
Approximation using Conjugate Gradients
def _solve_wmf_cg(X, U, V, confidence, lamb, n_steps=3):
n_users, n_factors = U.shape
VTV = np.dot(V.T, V) + lamb * np.eye(n_factors)
print(VTV.shape)
for j in range(n_users):
nonzero_indices = get_nonzero_indices(X, j)
user_factor = U[j]
## calculate residual r
r = -np.dot(VTV, user_factor)
for i in nonzero_indices:
item_factor = V[i].copy()
r += (confidence - (confidence - 1) * item_factor.dot(user_factor)) * item_factor
p = r.copy()
old_r_square = np.dot(r, r)
for step in range(n_steps):
# calculate Ap without calculating computationaly heavy things
Ap = np.dot(VTV, p)
for i in nonzero_indices:
Ap += (confidence - 1) * np.dot(item_factor, p) * item_factor
# Conjugate Gradient update
alpha = old_r_square / np.dot(p, Ap)
user_factor += alpha * p
r -= alpha * Ap
new_r_square = np.dot(r, r)
p = r + (new_r_square / old_r_square) * p
old_r_square = new_r_square
U[j] = user_factor
Written with StackEdit.