def fit(self, X, y=None):
        N = X.shape[0]
        if y is None:
            y = np.zeros(N)
        self.classes = list(set(y))
        self.classes.sort()
        self.n_classes = len(self.classes)
        if not self.sigma:
            self.sigma = median_kneighbour_distance(X)
            self.gamma = self.sigma**-2
        if not self.gamma:
            self.gamma = self.sigma**-2
        if not self.rho:
            self.rho = 0.1
        if self.kernel_pos is None:
            B = min(self.n_kernels_max, N)
            kernel_idx = np.random.permutation(N)
            self.kernel_pos = X[kernel_idx[:B]]
        else:
            B = self.kernel_pos.shape[0]
        Phi = metrics.pairwise.rbf_kernel(X, self.kernel_pos, self.gamma)
        theta = {}
        Phi_PhiT = np.dot(Phi.T, Phi)
        inverse_term = np.linalg.inv(Phi_PhiT + self.rho*np.eye(B))
        for c in self.classes:
            m = (y == c).astype(int)
            theta[c] = np.dot(inverse_term, np.dot(Phi.T, m))
        self.theta = theta