diff --git a/test/test_adaptative_permutation_threshold_SVR.py b/test/test_adaptative_permutation_threshold_SVR.py index b175b51..144e1f1 100644 --- a/test/test_adaptative_permutation_threshold_SVR.py +++ b/test/test_adaptative_permutation_threshold_SVR.py @@ -4,9 +4,11 @@ import numpy as np from numpy.testing import assert_almost_equal +from sklearn.svm import SVR from hidimstat.adaptative_permutation_threshold_SVR import ada_svr, ada_svr_pvalue from hidimstat.scenario import multivariate_1D_simulation +from hidimstat.permutation_test import permutation_test def test_ada_svr(): @@ -42,3 +44,54 @@ def test_ada_svr(): assert_almost_equal(pval[:support_size], expected[:support_size], decimal=1) assert_almost_equal(pval_corr[support_size:], expected[support_size:], decimal=1) + + +def test_ada_svr_rcond(): + """ + Testing the effect of rcond + """ + # create dataset + X, y, beta, _ = multivariate_1D_simulation( + n_samples=20, + n_features=50, + support_size=3, + sigma=0.1, + shuffle=False, + seed=42, + ) + X[:10] *= 1e-5 + beta_hat, scale = ada_svr(X, y) + beta_hat_2, scale_2 = ada_svr(X, y, rcond=1e-15) + assert np.max(np.abs(beta_hat - beta_hat_2)) > 1 + assert np.max(np.abs(scale - scale_2)) > 1 + + +def test_ada_svr_vs_permutation(): + """ + Validate the adaptive permutation threshold procedure against a permutation + test. The adaptive permutation threshold procedure should good approciation + of the proba of the permutation test. + """ + # create dataset + X, y, beta, _ = multivariate_1D_simulation( + n_samples=10, + n_features=100, + support_size=1, + sigma=0.1, + shuffle=False, + seed=42, + ) + beta_hat, scale = ada_svr(X, y) + # fit a SVR to get the coefficients + estimator = SVR(kernel="linear", epsilon=0.0, gamma="scale", C=1.0) + estimator.fit(X, y) + beta_hat_svr = estimator.coef_ + + # compare that the coefficiants are the same that the one of SVR + assert np.max(np.abs(beta_hat - beta_hat_svr.T[:, 0])) < 2e-4 + + proba = permutation_test( + X, y, estimator=estimator, n_permutations=10000, n_jobs=8, seed=42, proba=True + ) + assert np.max(np.abs(np.mean(proba, axis=0))) < 1e-3 + assert np.max(np.abs(scale - np.std(proba, axis=0)) / scale) < 1e-1