function [cv_err, cv_nu, cv_err_nu, cv_spe] = tlasso_CVTest(X, y, nu, percent_train, niter)

n = size(X,1);

% Augment 'nu's with a large value (i.e., Gaussian)
nu = [nu, 1e4];
n_nu = length(nu);
cv_nu = zeros(niter, 1);
cv_err = zeros(niter, 2);
cv_err_nu = zeros(niter, n_nu);
cv_spe = zeros(niter, 2);
spe_nu = zeros(niter, n_nu);

b = cell(1, n_nu);
b0 = cell(1, n_nu);
retval = cell(1, n_nu);

%% Peform the cross-validationt tests
for i = 1:niter
    fprintf('(Iter = %3d/%3d)\n', i, niter);
    
    % Split the data into training and testing sets
    cv = cvpartition(n,'HoldOut',1-percent_train);
    
    % Try the different nu's on the training dataset
    for j = 1:n_nu
        [b{j}, b0{j}, retval{j}] = tlasso(X(cv.training,:), y(cv.training), nu(j), 'CV', 10 );
    end
    
    % Compute CV errors on testing dataset
    e = zeros(1, n_nu);
    s = zeros(1, n_nu);
    best_err = zeros(1, n_nu);
    for j = 1:n_nu
        I = retval{j}.IndexMinCV;
        best_err(j) = retval{j}.CVErr(I);
        e(j) = tlasso_negll(X(cv.test,:), y(cv.test), b0{j}(I), b{j}(:,I), retval{j}.sigma2(I), nu(j));
        s(j) = mean((y(cv.test) - X(cv.test,:)*b{j}(:,I) + b0{j}(I)).^2);
    end
    
    % Store the results
    [~,cv_nu(i)] = min(best_err);
    cv_err(i,1) = e(cv_nu(i));
    cv_err(i,2) = e(end);
    cv_spe(i,1) = s(cv_nu(i));
    cv_spe(i,2) = s(end);
    
    cv_err_nu(i,:) = e;
    spe_nu(i,:) = s;
end

%% Done -- display results
fprintf('\n');
fprintf(' Results\n');
fprintf(' ------------------------------------------\n');
fprintf('%16s  %12s %12s\n', 'nu', 'CV error ', '% selected');
%nu = [1, 2, 10, 1e4, inf];
nu = [nu,inf];
for j = 1:length(nu)
    if (j == n_nu)
        s = 'Gaussian';
    else
        s = sprintf('%d', nu(j));
    end
    
    if (j <= (n_nu))
        fprintf('%16s %12.2f %12d\n', s, mean(cv_err_nu(:,j)), sum(cv_nu==j));
    else
        fprintf('%16s %12.2f %12s\n', 'Best nu (by CV)', mean(cv_err(:,1)), '-');
    end
end
fprintf(' ------------------------------------------\n\n');

return