pca_svds.py

def pca(X, npc):
    n_samples, n_features = X.shape
    Xmean = np.mean(X, axis=0)
    U, s, Vt = scipy.sparse.linalg.svds(X - Xmean, k=npc)
    order = np.argsort(-s) # sort s in descending order

    # svds returns U, s, Vt sorder in ascending order. We want descending
    s = s[order]
    W = Vt[order,:]

    explained_variance = (s**2) / float(n_samples)

    return Xmean, W, explained_variance