Graphics problem with tsne algorithm

[JKL-STAR]

class TSNEJAMES:

def __init__(self, perplexity, dim, n_iter, learning_rate, random_state, early_exaggeration):
    self.perplexity = perplexity
    self.dim = dim
    self.n_iter = n_iter
    self.learning_rate = learning_rate
    self.random_state = random_state
    self.early_exaggeration = early_exaggeration

@staticmethod
def grid_search(diff_i, i, perplexity):
    result = np.inf
    norm = np.linalg.norm(diff_i, axis=1)
    std_norm = np.std(norm)
    sigma = 0
    for sigma_search in np.linspace(0.01 * std_norm, 5 * std_norm, 200):
        p = np.exp(-norm**2 / (2 * sigma_search**2))
        p[i] = 0
        epsilon = np.nextafter(0, 1)
        p_new = np.maximum(p / np.sum(p), epsilon)
        H = -np.sum(p_new * np.log2(p_new))
        if np.abs(np.log(perplexity) - H * np.log(2)) < np.abs(result):
            result = np.log(perplexity) - H * np.log(2)
            sigma = sigma_search
    return sigma

@staticmethod
def P_ij(X, sigmas):
    n = len(X)
    pij = np.zeros((n, n))
    for i in range(n):
        diff = X[i, :] - X
        pij[i, :] = np.exp(-np.square(np.linalg.norm(diff, axis=1)) / (2 * np.square(sigmas[i])))
        np.fill_diagonal(pij, 0)
        pij[i, :] = pij[i, :] / pij[i, :].sum()
    pij = (pij + pij.T) / (2 * n)
    epsilon = np.nextafter(0, 1)
    pij = np.maximum(pij, epsilon)
    return pij

@staticmethod
def gradient(P, Q, y):
    n, m = len(P), y.shape[1]
    gradient = np.zeros((n, m))
    C = P - Q
    for i in range(n):
        diff = y[i, :] - y
        gradient[i] = 4*C[i,:]@(diff*((1+np.square(np.linalg.norm(diff,axis=1)))**(-1)).reshape((-1,1)))
    return gradient

@staticmethod
def Qij(y):
    n = y.shape[0]
    Qij = np.zeros((n, n))
    for i in range(n):
        diff = y[i, :] - y
        Qij[i, :] = (1 + np.square(np.linalg.norm(diff, axis=1))) ** -1
    np.fill_diagonal(Qij, 0)
    Qij = Qij / Qij.sum()
    epsilon = np.nextafter(0, 1)
    Qij = np.maximum(Qij, epsilon)
    return Qij

def tsne(self, X):
    n = len(X)
    sigmas = np.zeros(n)
    for i in range(n):
        diff = X[i, :] - X
        sigmas[i] = self.grid_search(diff, i, self.perplexity)

    Pij = self.P_ij(X, sigmas)

    def m(t):
        return 0.5 if t < 250 else 0.8

    np.random.seed(self.random_state)
    y = np.random.normal(loc=0.0, scale=1e-4, size=(n, self.dim))
    Y = [y, y]
    for t in range(self.n_iter):
        early_exaggeration = self.early_exaggeration if t < 250 else 1
        Q = self.Qij(Y[-1])
        y = Y[-1] - self.learning_rate * self.gradient(early_exaggeration * Pij, Q, Y[-1]) + \
            m(t) * (Y[-1] - Y[-2])
        Y.append(y)
        if t % 100 == 0 or t == 1:
            cost = np.sum(Pij * np.log(Pij / Q))
            print(f"Iteration {t}: Value of Cost Function is {cost}")
    return y