In [2]:
import numpy as np
import pandas as pd
from matplotlib import pylab as plt
import time
%matplotlib inline
In [3]:
## Y = A*B
M = 3
N = 5
K = 2
A_true = np.mat(np.random.randn(M, K))
B_true = np.mat(np.random.randn(K, N))
print(A_true)
print(B_true)
In [4]:
Y = A_true * B_true
Y
Out[4]:
In [5]:
B = A_true.I * Y
B
Out[5]:
In [6]:
In [7]:
In [8]:
U * M
Out[8]:
In [9]:
U
Out[9]:
In [10]:
M
Out[10]:
In [11]:
U = A * M.I
In [12]:
M = U.I * A
In [13]:
U
Out[13]:
In [14]:
M
Out[14]:
In [15]:
U * M
Out[15]:
In [234]:
A = np.mat([
# colums: words, rows: documents
# the, a, sparse, love, war, gun
[21, 22, 0, 70, 10, 6],
[22, 21, 0, 70, 10, 6],
[20, 23, 0, 0, 10, 6],
[22, 24, 65, 0, 10, 6],
[21, 22, 65, 0, 10, 6],
])
In [17]:
np.linalg.svd(A)[1]
Out[17]:
In [18]:
np.linalg.svd(A)[1]
Out[18]:
In [19]:
np.linalg.svd(A)[1]
Out[19]:
In [20]:
np.linalg.svd(A)[0][0] np.linalg.svd(A)[1][0]
In [21]:
np.linalg.svd(A)[0][:,0] * np.linalg.svd(A)[1][0] * np.linalg.svd(A)[2][0]
Out[21]:
In [22]:
def reconstruct(A, k=1):
return np.linalg.svd(A)[0][:,:k] * np.diag(np.linalg.svd(A)[1][:k]) * np.linalg.svd(A)[2][:k]
In [23]:
reconstruct(A, 1)
Out[23]:
In [24]:
reconstruct(A, 2)
Out[24]:
In [25]:
reconstruct(A, 3)
Out[25]:
In [235]:
A = np.mat([
[5, 0, 5, 3, 0, 0, 1, 0, 0, 1, 0],
[4, 5, 0, 4, 3, 0, 2, 1, 0, 0, 2],
[5, 4, 4, 4, 0, 0, 1, 0, 1, 0, 1],
[0, 1, 0, 0, 0, 3, 0, 0, 0, 3, 0],
[1, 0, 0, 1, 1, 5, 4, 0, 3, 5, 0],
[0, 1, 1, 0, 1, 4, 5, 5, 4, 0, 4]
])
shape = A.shape
print("rows, cols: ", shape)
In [27]:
U, S, V = np.linalg.svd(A)
(U, S, V)
Out[27]:
In [28]:
def reconstruct(U, S, V, rank):
return U[:,:rank] * np.diag(S[:rank]) * V[:rank]
In [29]:
reconstruct(U, S, V, 1)
Out[29]:
In [30]:
def calcError(A, rA):
return np.sum(np.power(A - rA, 2))
In [31]:
S_sum = np
for i in range(1, len(S)):
rA = reconstruct(U, S, V, i)
percentage = S[:i].sum() / S.sum()
error = calcError(A, rA)
print("with rank {}, coverrage: {:0.4f}, error= {:0.4f}".format(i, percentage, error))
In [32]:
reconstruct(U, S, V, 4)
Out[32]:
In [33]:
A
Out[33]:
In [34]:
A.shape
Out[34]:
In [35]:
mask = np.ones(A.shape)
mask[A==0] = 0
mask
Out[35]:
In [36]:
A_true = np.mat(np.random.randn(6, 5))
B_true = np.mat(np.random.randn(5, 11))
In [37]:
B_true = A.I * A
A_true = A * B_true.I
calcError(A, A_true * B_true)
Out[37]:
In [38]:
B_true = A.I * A
A_true = A * B_true.I
calcError(A, A_true * B_true)
Out[38]:
In [39]:
A
Out[39]:
In [224]:
In [269]:
A = np.random.random((10, 20))
In [270]:
rank = 4
U1 = np.mat(np.random.randn(A.shape[0], rank))
M1 = np.mat(np.random.randn(rank, A.shape[1]))
U1.shape, M1.shape
Out[270]:
In [424]:
U = np.copy(U1)
M = np.copy(M1)
errors = []
start = time.time()
for i in range(200):
U = np.linalg.solve(M.dot(M.T), M.dot(A.T)).T
M = np.linalg.solve(U.T.dot(U), U.T.dot(A))
if i % 50 == 0:
e = np.sqrt(np.power(A - U.dot(M), 2).sum()) / (A.shape[0] + A.shape[1])
errors.append(e)
elapsed = time.time() - start
e = np.sqrt(np.power(A - U.dot(M), 2).sum()) / (A.shape[0] + A.shape[1])
errors.append(e)
print(elapsed, errors[-1])
plt.plot(errors)
plt.show()
In [451]:
U = np.copy(U1)
M = np.copy(M1)
eta = 0.2
MAX_ITERATION = 20000
errors = []
start = time.time()
for i in range(MAX_ITERATION):
if i % 200 == 0:
eta = eta * 0.9
for j in range(A.shape[0]):
# u_ind = np.random.randint(0, A.shape[0])
u_ind = j
m_ind = np.random.randint(0, A.shape[1])
a_i = A[u_ind, m_ind]
u_i = U[u_ind,:]
m_i = M[:,m_ind]
err = a_i - u_i.dot(m_i)
dU = - err * m_i.T
dM = - u_i.T.dot(err)
U[u_ind,:] = U[u_ind,:] - eta * dU
M[:,m_ind] = M[:,m_ind] - eta * dM
if i % 400 == 0:
e = np.sqrt(np.power(A - U.dot(M), 2).sum()) / (A.shape[0] + A.shape[1])
errors.append(e)
elapsed = time.time() - start
e = np.sqrt(np.power(A - U.dot(M), 2).sum()) / (A.shape[0] + A.shape[1])
errors.append(e)
print(elapsed, errors[-1])
plt.plot(errors)
plt.show()
In [452]:
U = np.copy(U1)
M = np.copy(M1)
eta = 0.14
MAX_ITERATION = 40000
block_size = 4
errors = []
start = time.time()
for i in range(MAX_ITERATION):
# for j in range(0, A.shape[0] - block_size):
# u_ind = j
u_ind = np.random.randint(0, A.shape[0]-block_size)
m_ind = np.random.randint(0, A.shape[1]-block_size)
a_i = A[u_ind:u_ind+block_size, m_ind:m_ind+block_size]
u_i = U[u_ind:u_ind+block_size,:]
m_i = M[:,m_ind:m_ind+block_size]
err = a_i - u_i.dot(m_i)
dU = - err.dot(m_i.T)
dM = - u_i.T.dot(err)
U[u_ind:u_ind+block_size,:] = U[u_ind:u_ind+block_size,:] - eta * dU
M[:,m_ind:m_ind+block_size] = M[:,m_ind:m_ind+block_size] - eta * dM
if i % 100 == 0:
e = np.sqrt(np.power(A - U.dot(M), 2).sum()) / (A.shape[0] + A.shape[1])
errors.append(e)
elapsed = time.time() - start
print('k')
e = np.sqrt(np.power(A - U.dot(M), 2).sum()) / (A.shape[0] + A.shape[1])
errors.append(e)
print(elapsed, errors[-1])
plt.plot(errors)
plt.show()
In [458]:
U = np.copy(U1)
M = np.copy(M1)
eta = 0.01
MAX_ITERATION = 12000
block_size = 1
errors = []
start = time.time()
for i in range(MAX_ITERATION):
for j in range(0, A.shape[0] - block_size):
u_ind = j
m_ind = np.random.randint(0, A.shape[1]-block_size)
a_i = A[u_ind:u_ind+block_size, m_ind:m_ind+block_size]
u_i = U[u_ind:u_ind+block_size,:]
m_i = M[:,m_ind:m_ind+block_size]
# u_i = np.linalg.solve(m_i.dot(m_i.T), m_i.dot(a_i.T)).T
err = a_i - u_i.dot(m_i)
dU = - err.dot(m_i.T)
U[u_ind:u_ind+block_size,:] = U[u_ind:u_ind+block_size,:] - eta * dU
# U[u_ind:u_ind+block_size,:] = u_i
for j in range(0, A.shape[1] - block_size):
u_ind = np.random.randint(0, A.shape[0]-block_size)
m_ind = j
a_i = A[u_ind:u_ind+block_size, m_ind:m_ind+block_size]
u_i = U[u_ind:u_ind+block_size,:]
m_i = M[:,m_ind:m_ind+block_size]
# m_i = np.linalg.solve(u_i.T.dot(u_i), u_i.T.dot(a_i))
err = a_i - u_i.dot(m_i)
dM = - u_i.T.dot(err)
M[:,m_ind:m_ind+block_size] = M[:,m_ind:m_ind+block_size] - eta * dM
if i % 2400 == 0:
e = np.sqrt(np.power(A - U.dot(M), 2).sum()) / (A.shape[0] + A.shape[1])
errors.append(e)
elapsed = time.time() - start
print('k')
e = np.sqrt(np.power(A - U.dot(M), 2).sum()) / (A.shape[0] + A.shape[1])
errors.append(e)
print(elapsed, errors[-1])
plt.plot(errors)
plt.show()
In [462]:
U = np.copy(U1)
M = np.copy(M1)
eta = 0.01
MAX_ITERATION = 12000
block_size = 1
errors = []
start = time.time()
for i in range(MAX_ITERATION):
for j in range(0, A.shape[0] - block_size):
u_ind = j
m_ind = np.random.randint(0, A.shape[1]-block_size)
a_i = A[u_ind:u_ind+block_size, :]
u_i = U[u_ind:u_ind+block_size,:]
err = a_i - u_i.dot(M)
dU = - err.dot(M.T)
U[u_ind:u_ind+block_size,:] = U[u_ind:u_ind+block_size,:] - eta * dU
# U[u_ind:u_ind+block_size,:] = u_i
for j in range(0, A.shape[1] - block_size):
u_ind = np.random.randint(0, A.shape[0]-block_size)
m_ind = j
a_i = A[:, m_ind:m_ind+block_size]
m_i = M[:,m_ind:m_ind+block_size]
err = a_i - U.dot(U)
dM = - U.T.dot(err)
M[:,m_ind:m_ind+block_size] = M[:,m_ind:m_ind+block_size] - eta * dM
if i % 2400 == 0:
e = np.sqrt(np.power(A - U.dot(M), 2).sum()) / (A.shape[0] + A.shape[1])
errors.append(e)
elapsed = time.time() - start
print('k')
e = np.sqrt(np.power(A - U.dot(M), 2).sum()) / (A.shape[0] + A.shape[1])
errors.append(e)
print(elapsed, errors[-1])
plt.plot(errors)
plt.show()
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [431]:
U.dot(M)[3,:5]
Out[431]:
In [457]:
a_i
Out[457]:
In [ ]:
In [408]:
U.dot(M)[3,:5]
Out[408]:
In [441]:
M
Out[441]:
In [ ]: