Merge branch 'main' of github.com:NFFT/pyANOVAapprox

Author: passscoed

simpleTest/exampleCheb.py

@@ -11,12 +11,7 @@
 
 
 def TestFunction(x):
-    return (
-        x[0] * x[4]
-        + 2
-        - np.exp(x[3])
-        + np.sqrt(x[5] + 3 + x[1])
-    )
+    return x[0] * x[4] + 2 - np.exp(x[3]) + np.sqrt(x[5] + 3 + x[1])
 
 
 rng = np.random.default_rng(1234)
@@ -39,12 +34,15 @@ def TestFunction(x):
 b = M / (
     math.log10(M) * num
 )  # number for the number of frequencies if we use logarithmic oversampling and distribute it evenly to all subsets
-bw = [math.floor(b / 2) * 2, math.floor(math.sqrt(b) / 2) * 2]  # bandwidths (use even numbers)
+bw = [
+    math.floor(b / 2) * 2,
+    math.floor(math.sqrt(b) / 2) * 2,
+]  # bandwidths (use even numbers)
 # Use all subsets up to ds and use bw[1] many frequences in the the subsets with one element, b[2]^2 many for subsets with two elements and so on
 #
 ########### Variant 2:
 # used subsets:
-# U = [(), (0,), (1,), (2,), (3,), (4,), (5,), 
+# U = [(), (0,), (1,), (2,), (3,), (4,), (5,),
 #      (0, 1), (0, 2), (0, 3), (0, 4), (0, 5), (1, 2), (1, 3), (1, 4), (1, 5), (2, 3), (2, 4), (2, 5), (3, 4), (3, 5), (4, 5)]
 # Bandwidths for these subsets:
 # N = [0 ,  100,  100,  100,  100,  100,  100,
@@ -53,7 +51,7 @@ def TestFunction(x):
 #
 ########### Variant 3:
 # used subsets:
-# U = [(),   (0,),   (1,),   (2,),   (3,),   (4,),   (5,), 
+# U = [(),   (0,),   (1,),   (2,),   (3,),   (4,),   (5,),
 #        (0, 1),  (0, 2),  (0, 3),  (0, 4),  (0, 5),  (1, 2),  (1, 3),  (1, 4),  (1, 5),  (2, 3),  (2, 4),  (2, 5),  (3, 4),  (3, 5),  (4, 5)]
 # Bandwidths for these subsets:
 # N = [(), (100,), (100,), (100,), (100,), (100,), (100,),
@@ -156,7 +154,9 @@ def TestFunction(x):
 Umask = np.append(np.array([True]), gsis > 1e-2)
 U = [ads.U[i] for i in np.arange(0, len(Umask))[Umask]]  # get important subsets
 bws = M / (math.log10(M) * (len(U) - 1))  # calculate frequencies per subset
-N = [math.floor(bws ** (1 / max(1, len(u))) / 2) * 2 for u in U]  # distribute the frequencies evenly and make them even
+N = [
+    math.floor(bws ** (1 / max(1, len(u))) / 2) * 2 for u in U
+]  # distribute the frequencies evenly and make them even
 N[0] = 0
 
 a = ANOVAapprox.approx(

simpleTest/exampleClassification.py

@@ -9,10 +9,14 @@
 
 import pyANOVAapprox as ANOVAapprox
 
+
 def TestFunction(x):
-    e = (abs(x[1]+1.0j*x[2])+(np.angle(x[1]+1.0j*x[2])/(math.pi*8))) % 0.25>0.125
-    return e*2-1
-    
+    e = (
+        abs(x[1] + 1.0j * x[2]) + (np.angle(x[1] + 1.0j * x[2]) / (math.pi * 8))
+    ) % 0.25 > 0.125
+    return e * 2 - 1
+
+
 rng = np.random.default_rng(1234)
 
 ##################################
@@ -85,8 +89,8 @@ def TestFunction(x):
 ## get classification accuracy ##
 #################################
 
-y_approx = ads.evaluate(X=X_test, lam=0.0) # evaluate the classification
-acc = sum(np.sign(y_approx) == y_test)/M_test # calculate the accuracity
+y_approx = ads.evaluate(X=X_test, lam=0.0)  # evaluate the classification
+acc = sum(np.sign(y_approx) == y_test) / M_test  # calculate the accuracity
 print("accuracity = " + str(acc))
 
 ###############################################
@@ -156,30 +160,31 @@ def TestFunction(x):
     X, y, U, N, "cos", classification=True
 )  # generate the data structure for the classification
 a.approximate(
-    lam=lambdas,
-    max_iter=max_iter
+    lam=lambdas, max_iter=max_iter
 )  # do the approximation for all specified regularisation parameters
 
-y_approx = a.evaluate(X=X_test, lam=0.0) # evaluate the classification
-acc = sum(np.sign(y_approx) == y_test)/M_test # calculate the accuracity
+y_approx = a.evaluate(X=X_test, lam=0.0)  # evaluate the classification
+acc = sum(np.sign(y_approx) == y_test) / M_test  # calculate the accuracity
 print("accuracity = " + str(acc))
 
 ########################
 ## Evaluate the model ##
 ########################
 
-#y_approx = a.evaluate(X=X_test) # evaluate the classification at the training points for the regularisation λ_min
-y_approx = a.evaluate(X=X_test, lam=0.0) # evaluate the classification at the points X_test for the regularisation λ_min
+# y_approx = a.evaluate(X=X_test) # evaluate the classification at the training points for the regularisation λ_min
+y_approx = a.evaluate(
+    X=X_test, lam=0.0
+)  # evaluate the classification at the points X_test for the regularisation λ_min
 
 plt.figure()
 scatter = plt.scatter(
-    X_test[:,1],
-    X_test[:,2],
-    c=np.sign(y_approx), 
-    cmap='winter',      
-    s=100,               
-    alpha=0.8,          
-    edgecolors='black',  
-    linewidths=0.5     
+    X_test[:, 1],
+    X_test[:, 2],
+    c=np.sign(y_approx),
+    cmap="winter",
+    s=100,
+    alpha=0.8,
+    edgecolors="black",
+    linewidths=0.5,
 )
 plt.show()

simpleTest/exampleNonPeriodic.py

@@ -11,12 +11,7 @@
 
 
 def TestFunction(x):
-    return (
-        x[0] * x[4]
-        + 2
-        - np.exp(x[3])
-        + np.sqrt(x[5] + 3 + x[1])
-    )
+    return x[0] * x[4] + 2 - np.exp(x[3]) + np.sqrt(x[5] + 3 + x[1])
 
 
 rng = np.random.default_rng(1234)
@@ -39,12 +34,15 @@ def TestFunction(x):
 b = M / (
     math.log10(M) * num
 )  # number for the number of frequencies if we use logarithmic oversampling and distribute it evenly to all subsets
-bw = [math.floor(b / 2) * 2, math.floor(math.sqrt(b) / 2) * 2]  # bandwidths (use even numbers)
+bw = [
+    math.floor(b / 2) * 2,
+    math.floor(math.sqrt(b) / 2) * 2,
+]  # bandwidths (use even numbers)
 # Use all subsets up to ds and use bw[1] many frequences in the the subsets with one element, b[2]^2 many for subsets with two elements and so on
 #
 ########### Variant 2:
 # used subsets:
-# U = [(), (0,), (1,), (2,), (3,), (4,), (5,), 
+# U = [(), (0,), (1,), (2,), (3,), (4,), (5,),
 #      (0, 1), (0, 2), (0, 3), (0, 4), (0, 5), (1, 2), (1, 3), (1, 4), (1, 5), (2, 3), (2, 4), (2, 5), (3, 4), (3, 5), (4, 5)]
 # Bandwidths for these subsets:
 # N = [0 ,  100,  100,  100,  100,  100,  100,
@@ -53,7 +51,7 @@ def TestFunction(x):
 #
 ########### Variant 3:
 # used subsets:
-# U = [(),   (0,),   (1,),   (2,),   (3,),   (4,),   (5,), 
+# U = [(),   (0,),   (1,),   (2,),   (3,),   (4,),   (5,),
 #        (0, 1),  (0, 2),  (0, 3),  (0, 4),  (0, 5),  (1, 2),  (1, 3),  (1, 4),  (1, 5),  (2, 3),  (2, 4),  (2, 5),  (3, 4),  (3, 5),  (4, 5)]
 # Bandwidths for these subsets:
 # N = [(), (100,), (100,), (100,), (100,), (100,), (100,),
@@ -154,7 +152,9 @@ def TestFunction(x):
 Umask = np.append(np.array([True]), gsis > 1e-2)
 U = [ads.U[i] for i in np.arange(0, len(Umask))[Umask]]  # get important subsets
 bws = M / (math.log10(M) * (len(U) - 1))  # calculate frequencies per subset
-N = [math.floor(bws ** (1 / max(1, len(u))) / 2) * 2 for u in U]  # distribute the frequencies evenly and make them even
+N = [
+    math.floor(bws ** (1 / max(1, len(u))) / 2) * 2 for u in U
+]  # distribute the frequencies evenly and make them even
 N[0] = 0
 
 a = ANOVAapprox.approx(

simpleTest/exampleWavelet.py

@@ -13,12 +13,9 @@
 # for 'cos' the samples have to be in [0,1]^d
 # for a periodic function use 'per' or wavelets  'chui2', 'chui3', 'chui4' (samples have to be in [-0.5,0.5]^d here, 'chuim' are the Chui-Wang wavelets of order m)
 
+
 def TestFunction(x):  # this function is of the form f_0 + f_1 + f_2 + f_3 + f_4 + f_2,3
-    return (
-        2 * abs(x[0])
-        + abs(math.sin(math.pi * x[1] * x[2]))
-        + np.cos(3 + x[3])
-    )
+    return 2 * abs(x[0]) + abs(math.sin(math.pi * x[1] * x[2])) + np.cos(3 + x[3])
 
 
 rng = np.random.default_rng(1234)
@@ -28,10 +25,10 @@ def TestFunction(x):  # this function is of the form f_0 + f_1 + f_2 + f_3 + f_4
 ##################################
 
 d = 8  # dimension
-q = 2      # superposition dimension
+q = 2  # superposition dimension
 M = 10000  # number of used evaluation points to train the model
 M_test = 10000  # number of used evaluation points to test the accuracity the model
-N = [5,2]  # number of parameters, should be vector of length q: 
+N = [5, 2]  # number of parameters, should be vector of length q:
 # for wavelets the total number of parameters scales exponentially, i.e.:
 # for q = 1 and N = [N1] the total number of parameters scales like ~O(d*2^N1)
 # for q = 2 and N = [N1 , N2] the total number of parameters scales like ~O(d*2^N1) + O(d^2 * N2*2^N2)
@@ -42,18 +39,28 @@ def TestFunction(x):  # this function is of the form f_0 + f_1 + f_2 + f_3 + f_4
 ## Generation of the data ##
 ############################
 
-if basis =="chui2" or basis =="chui3" or basis =="chui4" or basis =="per":
-    X = rng.random((M, d)) -0.5        # for perioidic approximation samples have to be in [-0.5,0.5]^d
+if basis == "chui2" or basis == "chui3" or basis == "chui4" or basis == "per":
+    X = (
+        rng.random((M, d)) - 0.5
+    )  # for perioidic approximation samples have to be in [-0.5,0.5]^d
 elif basis == "cos":
     X = rng.random((M, d))
 y = np.array(
     [TestFunction(X[i, :].T) for i in range(M)]
 )  # evaluate the function at these points
 
-if basis == "chui1" or basis =="chui2" or basis =="chui3" or basis =="chui4" or basis =="per":
-    X_test = rng.random((M_test, d)) - 0.5      # for perioidic approximation samples have to be in [-0.5,0.5]^d
+if (
+    basis == "chui1"
+    or basis == "chui2"
+    or basis == "chui3"
+    or basis == "chui4"
+    or basis == "per"
+):
+    X_test = (
+        rng.random((M_test, d)) - 0.5
+    )  # for perioidic approximation samples have to be in [-0.5,0.5]^d
 elif basis == "cos":
-    X_test = rng.random((M_test, d))     
+    X_test = rng.random((M_test, d))
 y_test = np.array(
     [TestFunction(X_test[i, :].T) for i in range(M_test)]
 )  # the same for the test points
@@ -74,8 +81,10 @@ def TestFunction(x):  # this function is of the form f_0 + f_1 + f_2 + f_3 + f_4
 #######################
 
 ### Do sensitivity analysis ####
-gsis = ANOVAapprox.get_GSI(anova_model,0.0)                     #calculates indices for importance of terms (gsis is vector, with indices belonging to terms in anova_model.U)
-#gsis_as_dict = ANOVAapprox.get_GSI(anova_model,0.0,dict=true) 
+gsis = ANOVAapprox.get_GSI(
+    anova_model, 0.0
+)  # calculates indices for importance of terms (gsis is vector, with indices belonging to terms in anova_model.U)
+# gsis_as_dict = ANOVAapprox.get_GSI(anova_model,0.0,dict=true)
 
 y_min_calc = 10 ** (np.min(np.log10(gsis)) - 0.5)
 label = list(anova_model.U[1:])
@@ -104,22 +113,22 @@ def TestFunction(x):  # this function is of the form f_0 + f_1 + f_2 + f_3 + f_4
 ################################
 
 ### error analysis ###
-mse_train = ANOVAapprox.get_mse(anova_model,lam=0.0)
-mse_test = ANOVAapprox.get_mse(anova_model,X_test,y_test, lam=0.0)
+mse_train = ANOVAapprox.get_mse(anova_model, lam=0.0)
+mse_test = ANOVAapprox.get_mse(anova_model, X_test, y_test, lam=0.0)
 
-print("MSE on test points: " + str(mse_test)) 
+print("MSE on test points: " + str(mse_test))
 
 ################################################
 ## Approximation with better suited index set ##
 ################################################
 
 U = ANOVAapprox.get_ActiveSet(anova_model, [0.01, 0.01], lam=0.0)
-print("Found index-set U: " + str(U) )
-anova_model = ANOVAapprox.approx(X, y, U=U, N=[i+2 for i in N] , basis=basis)               # increase number of paramers in N for the important terms
+print("Found index-set U: " + str(U))
+anova_model = ANOVAapprox.approx(
+    X, y, U=U, N=[i + 2 for i in N], basis=basis
+)  # increase number of paramers in N for the important terms
 anova_model.approximate(lam=lambdas)
 print("Total number of used parameters = " + str(len(anova_model.fc[lambdas[0]].vec())))
-mse_train = ANOVAapprox.get_mse(anova_model,lam=0.0)
-mse_test = ANOVAapprox.get_mse(anova_model,X_test,y_test, lam=0.0)
-print("MSE on test points after ANOVA truncation: " + str(mse_test)) 
-
-
+mse_train = ANOVAapprox.get_mse(anova_model, lam=0.0)
+mse_test = ANOVAapprox.get_mse(anova_model, X_test, y_test, lam=0.0)
+print("MSE on test points after ANOVA truncation: " + str(mse_test))

src/pyANOVAapprox/approx.py

@@ -1,4 +1,4 @@
-#from pyGroupedTransforms.GroupedTransform import *  # TODO: Kann wahrscheinlich weg sobald in pyGroupedTransform GreoupedTransform exportiert wird
+# from pyGroupedTransforms.GroupedTransform import *  # TODO: Kann wahrscheinlich weg sobald in pyGroupedTransform GreoupedTransform exportiert wird
 
 from pyANOVAapprox import *
 from pyANOVAapprox.fista import *

src/pyANOVAapprox/fista.py

@@ -1,4 +1,4 @@
-#from pyGroupedTransforms.GroupedCoefficients import *  # TODO: Kann wahrscheinlich weg sobald in pyGroupedTransform GreoupedTransform exportiert wird
+# from pyGroupedTransforms.GroupedCoefficients import *  # TODO: Kann wahrscheinlich weg sobald in pyGroupedTransform GreoupedTransform exportiert wird
 
 from pyANOVAapprox import *
 

tests/wav_lsqr.py

@@ -5,8 +5,8 @@
 sys.path.insert(0, src_aa)
 
 import numpy as np
-from TestFunctionPeriodic import *
 from pyGroupedTransforms import *
+from TestFunctionPeriodic import *
 
 import pyANOVAapprox as ANOVAapprox
 
@@ -26,8 +26,8 @@
 ads = ANOVAapprox.approx(X, y, ds=ds, N=bw, basis="chui2")
 ads.approximate(lam=lambdas, solver="lsqr")
 
-print( "AR: "+ str(sum(ANOVAapprox.get_AttributeRanking(ads, 0.0))) )
-assert abs( sum(ANOVAapprox.get_AttributeRanking(ads, 0.0)) - 1 ) < 0.0001
+print("AR: " + str(sum(ANOVAapprox.get_AttributeRanking(ads, 0.0))))
+assert abs(sum(ANOVAapprox.get_AttributeRanking(ads, 0.0)) - 1) < 0.0001
 
 bw = ANOVAapprox.get_orderDependentBW(AS, [4, 4])
 aU = ANOVAapprox.approx(X, y, U=AS, N=bw, basis="chui2")
@@ -45,6 +45,6 @@
 print("l2 rand U: ", err_l2_rand_U)
 
 assert err_l2_ds < 0.01
-assert err_l2_U < 0.01 # maybe restrict to 0.005
+assert err_l2_U < 0.01  # maybe restrict to 0.005
 assert err_l2_rand_ds < 0.01
 assert err_l2_rand_U < 0.01