'Fit the model
Parameters
X : numpy array of shape [n_samples,n_features]
Training data
y : numpy array of shape [n_samples]
Target values. Will be cast to X\'s dtype if necessary
Returns
self : returns an instance of self.'

def fit(self, X, y):

    (X, y) = check_X_y(X, y, dtype=np.float64, y_numeric=True)
    (X, y, X_offset_, y_offset_, X_scale_) = self._preprocess_data(X, y, self.fit_intercept, self.normalize, self.copy_X)
    self.X_offset_ = X_offset_
    self.X_scale_ = X_scale_
    (n_samples, n_features) = X.shape
    alpha_ = (1.0 / np.var(y))
    lambda_ = 1.0
    verbose = self.verbose
    lambda_1 = self.lambda_1
    lambda_2 = self.lambda_2
    alpha_1 = self.alpha_1
    alpha_2 = self.alpha_2
    self.scores_ = list()
    coef_old_ = None
    XT_y = np.dot(X.T, y)
    (U, S, Vh) = linalg.svd(X, full_matrices=False)
    eigen_vals_ = (S ** 2)
    for iter_ in range(self.n_iter):
        if (n_samples > n_features):
            coef_ = np.dot(Vh.T, (Vh / (eigen_vals_ + (lambda_ / alpha_))[:, np.newaxis]))
            coef_ = np.dot(coef_, XT_y)
            if self.compute_score:
                logdet_sigma_ = (- np.sum(np.log((lambda_ + (alpha_ * eigen_vals_)))))
        else:
            coef_ = np.dot(X.T, np.dot((U / (eigen_vals_ + (lambda_ / alpha_))[None, :]), U.T))
            coef_ = np.dot(coef_, y)
            if self.compute_score:
                logdet_sigma_ = (lambda_ * np.ones(n_features))
                logdet_sigma_[:n_samples] += (alpha_ * eigen_vals_)
                logdet_sigma_ = (- np.sum(np.log(logdet_sigma_)))
        self.alpha_ = alpha_
        self.lambda_ = lambda_
        rmse_ = np.sum(((y - np.dot(X, coef_)) ** 2))
        gamma_ = np.sum(((alpha_ * eigen_vals_) / (lambda_ + (alpha_ * eigen_vals_))))
        lambda_ = ((gamma_ + (2 * lambda_1)) / (np.sum((coef_ ** 2)) + (2 * lambda_2)))
        alpha_ = (((n_samples - gamma_) + (2 * alpha_1)) / (rmse_ + (2 * alpha_2)))
        if self.compute_score:
            s = ((lambda_1 * log(lambda_)) - (lambda_2 * lambda_))
            s += ((alpha_1 * log(alpha_)) - (alpha_2 * alpha_))
            s += (0.5 * ((((((n_features * log(lambda_)) + (n_samples * log(alpha_))) - (alpha_ * rmse_)) - (lambda_ * np.sum((coef_ ** 2)))) - logdet_sigma_) - (n_samples * log((2 * np.pi)))))
            self.scores_.append(s)
        if ((iter_ != 0) and (np.sum(np.abs((coef_old_ - coef_))) < self.tol)):
            if verbose:
                print('Convergence    after    ', str(iter_), '    iterations')
            break
        coef_old_ = np.copy(coef_)
    self.coef_ = coef_
    sigma_ = np.dot(Vh.T, (Vh / (eigen_vals_ + (lambda_ / alpha_))[:, np.newaxis]))
    self.sigma_ = ((1.0 / alpha_) * sigma_)
    self._set_intercept(X_offset_, y_offset_, X_scale_)
    return self
