import numpy as np
import time
def pca_powit(D):
print("Starting power iteration...")
t_start = time.time()
# Reshape and center
X = D.reshape(len(D), -1).T
X -= X.mean(axis=1, keepdims=True)
# Initial random vector
rng = np.random.mtrand.RandomState(57)
w = rng.normal(0, 1, [X.shape[0]])
w /= np.sqrt((w**2).sum())
# Covariance matrix
S = np.dot(X, X.T) / X.shape[1]
# Power iteration loop
J_old = float('nan')
for it in range(25):
Sw = S @ w
J = w @ Sw
print(f"iteration {it:2d} J(w)={J:10.3f}")
if not np.isnan(J_old) and abs(J - J_old) < 0.01:
print("Stopping criterion satisfied")
break
w = Sw / np.sqrt((Sw**2).sum())
J_old = J
t_end = time.time()
print(f"Time: {t_end - t_start:.3f} seconds")
# Comment render while debugging
# utils.render(w, 1, 1)
return w