def METHOD_NAME(self, x: DNDarray, y: DNDarray) -> None:
    """
    Fit lasso model with coordinate descent
    Parameters
    ----------
    x : DNDarray
        Input data, Shape = (n_samples, n_features)
    y : DNDarray
        Labels, Shape = (n_samples,)
    """
    # Get number of model parameters
    _, n = x.shape
    if y.ndim > 2:
        raise ValueError(f"y.ndim must <= 2, currently: {y.ndim}")
    if x.ndim != 2:
        raise ValueError(f"X.ndim must == 2, currently: {x.ndim}")
    if len(y.shape) == 1:
        y = ht.expand_dims(y, axis=1)
    # Initialize model parameters
    theta = ht.zeros((n, 1), dtype=float, device=x.device)
    # Looping until max number of iterations or convergence
    for i in range(self.max_iter):
        theta_old = theta.copy()
        # Looping through each coordinate
        for j in range(n):
            X_j = ht.array(x.larray[:, j : j + 1], is_split=0)
            y_est = x @ theta
            theta_j = theta.larray[j].item()
            rho = (X_j * (y - y_est + theta_j * X_j)).mean()
            # Intercept parameter theta[0] not be regularized
            if j == 0:
                theta[j] = rho
            else:
                theta[j] = self.soft_threshold(rho)
        diff = self.rmse(theta, theta_old)
        if self.tol is not None and diff < self.tol:
            self.n_iter = i + 1
            self.__theta = theta
            break
    self.n_iter = i + 1
    self.__theta = theta