cpals.py

import numpy as np

class CPALS:
    def __init__(self,k,lamb,max_iter):
        self.k=k
        self.lamb=lamb
        self.max_iter=max_iter

    def _calc_data_shape(self,X):
        max_i = -1
        max_j = -1
        max_k = -1
        for i,j,k in X:
            if max_i<i: max_i=i
            if max_j<j: max_j=j
            if max_k<k: max_k=k
        return (max_i+1,max_j+1,max_k+1)

    def _init_latent_vectors(self):
        A = {}
        for m in range(3): # 3rd-order tensor
            A[m] = np.random.normal(loc=0, scale=0.1, size=(self.data_shape[m],self.k))
        return A

    def _loss(self,X):
        L = 0
        count = 0
        for indices in X:
            L += 0.5 * (X[indices]-self.predict(indices))**2
            count += 1
        return L/count

    def _update(self,X,mode):
        X_bar = np.zeros_like(self.A[mode])
        for indices in X: # indices = (i,j,k)
            X_bar[indices[mode]] += X[indices] * np.prod([self.A[m][indices[m]] for m in range(3) if m!=mode],axis=0)
        K = np.prod([self.A[m].T.dot(self.A[m]) for m in range(3) if m!=mode],axis=0) + self.lamb*np.identity(self.k)
        return np.dot(X_bar, np.linalg.inv(K))

    def fit(self,X):
        self.data_shape = self._calc_data_shape(X)
        self.A = self._init_latent_vectors()
        self._losses = []
        remained_iter = self.max_iter
        while remained_iter>0:
            for m in range(3):
                self.A[m] = self._update(X,m)
            remained_iter-=1
            self._losses.append(self._loss(X))
        return self

    def predict(self,indices):
        return np.prod([self.A[m][indices[m]] for m in range(3)],axis=0).sum()