Version bump to 0.1.4;

Made the __INCLUDE_SE_KERNEL act less harshly (i.e. do not remove SE at the root, but just in the production rules) so that tests would pass untouched;
Noticed and fixed WN not acting idempotently

GPy_ABCD/KernelExpansion/grammar.py

@@ -1,6 +1,6 @@
 from GPy_ABCD.KernelExpansion.kernelExpressionOperations import *
 from GPy_ABCD.Util.genericUtil import *
-from GPy_ABCD.Kernels.baseKernels import __USE_LIN_KERNEL_HORIZONTAL_OFFSET
+from GPy_ABCD.Kernels.baseKernels import __INCLUDE_SE_KERNEL, __USE_LIN_KERNEL_HORIZONTAL_OFFSET
 
 
 ## Expansion functions
@@ -26,15 +26,17 @@ def deep_apply(operator, S, *args): # Deepcopy the tree and connect the new node
 
 ## Production Rules
 
-def plus_base(S): return [deep_apply(add, S, B) for B in base_kerns]
-def times_base(S): return [deep_apply(multiply, S, B) for B in base_kerns - {'C'}]
-def replace_base(S): return [deep_apply(swap_base, S, B) for B in base_kerns]
-def change_new_base(S): return [deep_apply(both_changes, S, B) for B in base_kerns - {'C'}] # Not in original ABCD
+base_kerns_for_prod = base_kerns if __INCLUDE_SE_KERNEL else base_kerns - {'SE'}
+
+def plus_base(S): return [deep_apply(add, S, B) for B in base_kerns_for_prod]
+def times_base(S): return [deep_apply(multiply, S, B) for B in base_kerns_for_prod - {'C'}]
+def replace_base(S): return [deep_apply(swap_base, S, B) for B in base_kerns_for_prod]
+def change_new_base(S): return [deep_apply(both_changes, S, B) for B in base_kerns_for_prod - {'C'}] # Not in original ABCD
 def change_same(S): return [deep_apply(both_changes, S)]
 def change_window_constant(S): return [deep_apply(one_change, S, 'CW', 'C')]
 def change_point_linear(S): return [deep_apply(one_change, S, 'CP', 'LIN')] # Not in original ABCD
-def times_shifted_base(S): return [deep_apply(multiply, S, SumKE([B, 'C'])) for B in base_kerns - {'C'}]
-def replace_with_singleton(S): return [deep_apply(replace_node, S, SumKE([B])) for B in base_kerns]
+def times_shifted_base(S): return [deep_apply(multiply, S, SumKE([B, 'C'])) for B in base_kerns_for_prod - {'C'}]
+def replace_with_singleton(S): return [deep_apply(replace_node, S, SumKE([B])) for B in base_kerns_for_prod]
 def remove_some_term(S): return [deep_apply(remove_a_term, S)]
 def try_higher_curves(S): return [deep_apply(higher_curves, S)] # Not in original ABCD
 
@@ -70,7 +72,6 @@ def pseudo_to_real_kex(pkex):
 def make_simple_kexs(pseudo_kexs): return [pseudo_to_real_kex(pkex)._initialise() for pkex in pseudo_kexs]
 
 
-# TODO: SumKE(['LIN', 'C']) instead of 'LIN'? Only for non-horizontal-offset-including version?
 standard_start_kernels = make_simple_kexs(list(base_kerns - {'SE'}) + # Base Kernels without SE
                                           [ProductKE(['LIN', 'LIN']), ProductKE(['LIN', 'LIN', 'LIN']), SumKE(['PER', 'C'])] + # More generic LIN and PER
                                           both_changes('LIN')) # To catch a possible changepoint or changewindow with simple enough shapes
@@ -79,7 +80,7 @@ def make_simple_kexs(pseudo_kexs): return [pseudo_to_real_kex(pkex)._initialise(
                                           [ProductKE(['LIN', 'LIN']), ProductKE(['LIN', 'LIN', 'LIN']), SumKE(['PER', 'C'])] + # More generic LIN and PER
                                           [SumKE(['C'], [ck]) for ck in both_changes('LIN')]) # To catch a possible changepoint or changewindow with simple enough shapes
 
-test_start_kernels = make_simple_kexs(list(base_kerns) + # Base Kernels
+test_start_kernels = make_simple_kexs(list(base_kerns - {'SE'}) + # Base Kernels without SE
                                           [ProductKE(['LIN', 'LIN']), ProductKE(['LIN', 'LIN', 'LIN']), SumKE(['PER', 'C'])] + # More generic LIN and PER
                                           both_changes('LIN')) # To catch a possible changepoint or changewindow with simple enough shapes
 

GPy_ABCD/KernelExpansion/kernelExpression.py

@@ -322,6 +322,7 @@ def bracket_if_needed(kex):
 
     def simplify_base_terms(self):
         if self.base_terms['WN'] > 0: # WN acts as a multiplicative zero for all stationary kernels, i.e. all but LIN and sigmoids
+            self.base_terms['WN'] = 1 # It is also idempotent
             for bt in list(self.base_terms.keys()):
                 if bt not in ['WN', 'LIN'] + list(base_sigmoids): del self.base_terms[bt]
         else:

GPy_ABCD/Kernels/baseKernels.py

@@ -5,7 +5,7 @@
 
 
 #### CORE CONFIGURATION OF BASE KERNELS ####
-__INCLUDE_SE_KERNEL = True # The most generic kernel; always a bargain in terms of parameters
+__INCLUDE_SE_KERNEL = False # The most generic kernel; always a bargain in terms of parameters
 __USE_LIN_KERNEL_HORIZONTAL_OFFSET = True # Identifies the polynomial roots; more accurate but one extra parameter per degree
 __USE_NON_PURELY_PERIODIC_PER_KERNEL = False # Full standard periodic kernel [MacKay (1998)] instead of only its purely periodic part
 __FIX_SIGMOIDAL_KERNELS_SLOPE = True # Hence one parameter fewer for each sigmoidal and related kernel
@@ -15,7 +15,7 @@
 
 ## Useful Kernel Sets
 
-base_kerns = frozenset(['WN', 'C', 'LIN', 'SE', 'PER']) if __INCLUDE_SE_KERNEL else frozenset(['WN', 'C', 'LIN', 'PER'])
+base_kerns = frozenset(['WN', 'C', 'LIN', 'SE', 'PER']) #if __INCLUDE_SE_KERNEL else frozenset(['WN', 'C', 'LIN', 'PER'])
 
 stationary_kerns = frozenset(['WN', 'C', 'SE', 'PER'])
 # addition_idempotent_kerns = frozenset(['WN', 'C'])
@@ -31,63 +31,43 @@
 
 ## Base Kernels
 
-# WN = _Gk.White(1)
 def WN(): return _Gk.White(1)
 
-# C = _Gk.Bias(1)
 def C(): return _Gk.Bias(1)
 
-# LIN = _Gk.Linear(1)
 # def LIN(): return _Gk.Linear(1) # Not the same as ABCD's; missing horizontal offset
-# LIN = _Lk.Linear(1)
 def LIN(): return _Lk.Linear(1) # Not the same as ABCD's; missing horizontal offset
 if __USE_LIN_KERNEL_HORIZONTAL_OFFSET: # The version in ABCD; not sure if a good idea; the horizontal offset is the same as a vertical one, which is just kC
-    # LIN = _Lk.LinearWithOffset(1)
     def LIN(): return _LOk.LinearWithOffset(1)
 
-# SE = _Gk.RBF(1)
 def SE(): return _Gk.RBF(1)
 
-# PER = _Pk.PureStdPeriodicKernel(1)
 def PER(): return _Pk.PureStdPeriodicKernel(1)
 if __USE_NON_PURELY_PERIODIC_PER_KERNEL: # Not the same as ABCD's
-    # PER = _Gk.StdPeriodic(1)
     def PER(): return _Gk.StdPeriodic(1)
 
 
 ## Sigmoidal Kernels
 
-# S = _Sk.SigmoidalKernel(1, False)
 def S(): return _Sk.SigmoidalKernel(1, False, fixed_slope = __FIX_SIGMOIDAL_KERNELS_SLOPE)
-# Sr = _Sk.SigmoidalKernel(1, True)
 def Sr(): return _Sk.SigmoidalKernel(1, True, fixed_slope = __FIX_SIGMOIDAL_KERNELS_SLOPE)
 
 
-# SIW = _Sk.SigmoidalIndicatorKernelWithWidth(1, False)
 def SI(): return _Sk.SigmoidalIndicatorKernelWithWidth(1, False, fixed_slope = __FIX_SIGMOIDAL_KERNELS_SLOPE)
-# SIWr = _Sk.SigmoidalIndicatorKernelWithWidth(1, True)
 def SIr(): return _Sk.SigmoidalIndicatorKernelWithWidth(1, True, fixed_slope = __FIX_SIGMOIDAL_KERNELS_SLOPE)
 
-# SI = _Sk.SigmoidalIndicatorKernel(1, False)
 def SIT(): return _Sk.SigmoidalIndicatorKernel(1, False, fixed_slope = __FIX_SIGMOIDAL_KERNELS_SLOPE)
-# SIr = _Sk.SigmoidalIndicatorKernel(1, True)
 def SITr(): return _Sk.SigmoidalIndicatorKernel(1, True, fixed_slope = __FIX_SIGMOIDAL_KERNELS_SLOPE)
 
-# SIO = _Sk.SigmoidalIndicatorKernelOneLocation(1, False)
 def SIO(): return _Sk.SigmoidalIndicatorKernelOneLocation(1, False, fixed_slope = __FIX_SIGMOIDAL_KERNELS_SLOPE)
-# SIOr = _Sk.SigmoidalIndicatorKernelOneLocation(1, True)
 def SIOr(): return _Sk.SigmoidalIndicatorKernelOneLocation(1, True, fixed_slope = __FIX_SIGMOIDAL_KERNELS_SLOPE)
 
 
 # Change-Operator Kernels
 
-CP = _Cs.ChangePointKernel
 def CP(left, right): return _Cs.ChangePointKernel(left, right, fixed_slope = __FIX_SIGMOIDAL_KERNELS_SLOPE)
-# # CW = _Cs.ChangeWindowKernel
 # def CW(left, right): return _Cs.ChangeWindowKernel(left, right, fixed_slope = __FIX_SIGMOIDAL_KERNELS_SLOPE)
-# # CW = _Cs.ChangeWindowKernelOneLocation
 # def CW(left, right): return _Cs.ChangeWindowKernelOneLocation(left, right, fixed_slope = __FIX_SIGMOIDAL_KERNELS_SLOPE)
-# CW = _Cs.ChangeWindowKernelWithWidth
 def CW(left, right): return _Cs.ChangeWindowKernelWithWidth(left, right, fixed_slope = __FIX_SIGMOIDAL_KERNELS_SLOPE)
 
 # CP = _CFs.kCP

Tests/TempTests/testData.py

@@ -1,7 +1,8 @@
-from modelSearch import *
 import numpy as np
 import scipy.io as sio
 
+from GPy_ABCD.Models.modelSearch import *
+
 
 # TODO:
 #   Solar data: the changewindow is not appropriately captured; need to investigate how to best address this;
@@ -21,12 +22,12 @@
 
 
     best_mods, all_mods, all_exprs = find_best_model(X, Y, start_kernels = standard_start_kernels, p_rules = production_rules_all,
-                                                     restarts = 2, utility_function = 'BIC', rounds = 1, buffer = 3, verbose= True)
+                                                     restarts = 2, utility_function = 'BIC', rounds = 2, buffer = 2, dynamic_buffer = True, verbose = True, parallel = True)
 
-    for mod_depth in all_mods: print(', '.join([str(mod.kernel_expression) for mod in mod_depth]))
+    for mod_depth in all_mods: print(', '.join([str(mod.kernel_expression) for mod in mod_depth]) + f'\n{len(mod_depth)}')
 
     from matplotlib import pyplot as plt
-    for bm in best_mods:
+    for bm in best_mods[:3]:
         print(bm.kernel_expression)
         print(bm.model.kern)
         print(bm.model.log_likelihood())

Tests/checkModelSearch.py

@@ -9,15 +9,19 @@
     # np.seterr(all='raise') # Raise exceptions instead of RuntimeWarnings. The exceptions can then be caught by the debugger
 
     # X, Y = generate_data(lambda x: x * np.cos( (x - 5) / 2 )**2, np.linspace(-15, 15, 101), 2, 1)
-    X, Y = generate_data(lambda x: (x + 3) * (x - 5) * (x - 7), np.linspace(-15, 15, 101), 1, 10)
-    # X, Y = generate_changepoint_data(np.linspace(-15, 15, 101), lambda x: 0.1 * x, lambda x: 2 + 3 * np.sin(x*3), 0, 1, 0.3)
-    # X, Y = generate_changewindow_data(np.linspace(-20, 20, 101), lambda x: 0.1 * x, lambda x: 4 * np.sin(x*2), -5, 5, 1, 0.3, True)
+    # kernel = LIN() * (PER() + C())
+    # X, Y = generate_data(lambda x: (x + 30) * (x - 5) * (x - 7), np.linspace(-15, 15, 101), 1, 30)
+    # kernel = LIN() * LIN() * LIN()
+    X, Y = generate_changepoint_data(np.linspace(-20, 30, 101), lambda x: 0.1 * x, lambda x: 2 * np.sin(x), 3, 1, 0.3, True)
+    kernel = CP(LIN(), PER() + C())
+    # X, Y = generate_changewindow_data(np.linspace(-30, 30, 101), lambda x: 0.1 * x, lambda x: 4 * np.sin(x), -10, 10, 1, 0.3, True)
+    # kernel = CW(LIN(), PER() + C())
 
     # print(gg_plot(X, Y))
 
 
     # from GPy_ABCD.Util.kernelUtil import doGPR
-    # mod = doGPR(X, Y, LIN() * (PER() + C()), 10)
+    # mod = doGPR(X, Y, kernel, 10)
     # predict_X = np.linspace(10, 15, 50)[:, None]
     # (p_mean, p_var) = mod.predict(predict_X)
     # print(predicted)
@@ -32,9 +36,9 @@
 
 
     # best_mods, all_mods, all_exprs = find_best_model(X, Y, start_kernels = ['WN'], p_rules = production_rules_all,
-    best_mods, all_mods, all_exprs = find_best_model(X, Y, start_kernels = test_start_kernels, p_rules = production_rules_all,
+    # best_mods, all_mods, all_exprs = find_best_model(X, Y, start_kernels = test_start_kernels, p_rules = production_rules_all,
     # best_mods, all_mods, all_exprs = find_best_model(X, Y, start_kernels = extended_start_kernels, p_rules = production_rules_all,
-    # best_mods, all_mods, all_exprs = find_best_model(X, Y, start_kernels = standard_start_kernels, p_rules = production_rules_all,
+    best_mods, all_mods, all_exprs = find_best_model(X, Y, start_kernels = standard_start_kernels, p_rules = production_rules_all,
                                                      restarts = 2, utility_function = 'BIC', rounds = 2, buffer = 2, dynamic_buffer = True, verbose = True, parallel = True)
 
     for mod_depth in all_mods: print(', '.join([str(mod.kernel_expression) for mod in mod_depth]) + f'\n{len(mod_depth)}')

Tests/regressChangeOperators.py

@@ -24,4 +24,4 @@
 # # kernel = CW(PER + C, LIN)
 
 
-m = doGPR(X, Y, kernel, 5)
+m = doGPR(X, Y, kernel, 10)

Tests/regressLinearKernels.py

@@ -15,7 +15,7 @@
 # kernel = LinearWithOffset(1)
 
 
-X, Y = generate_data(lambda x: (x - 2) * (x  + 3), np.linspace(-20, 20, 201), 1, 1)
+X, Y = generate_data(lambda x: - (x - 2) * (x + 3), np.linspace(-20, 20, 201), 1, 20)
 # kernel = Linear(1) * Linear(1)
 # kernel = Linear(1) * Linear(1) + C()
 # kernel = ProductKE(['LIN', 'LIN']).to_kernel()

Tests/test_kernelExpressions.py

@@ -13,7 +13,7 @@
     # Base terms
     (SumKE(['WN', 'WN', 'C', 'C', 'LIN', 'LIN', 'SE', 'SE', 'PER', 'PER']), SumKE(['WN', 'C', 'LIN', 'LIN', 'SE', 'SE', 'PER', 'PER'])),
     (ProductKE(['C', 'LIN', 'LIN', 'SE', 'SE', 'PER', 'PER']), ProductKE(['LIN', 'LIN', 'SE', 'PER', 'PER'])),
-    (ProductKE(['C', 'LIN', 'LIN', 'SE', 'SE', 'PER', 'PER', 'WN']), ProductKE(['LIN', 'LIN', 'WN'])),
+    (ProductKE(['C', 'LIN', 'LIN', 'SE', 'SE', 'PER', 'PER', 'WN', 'WN']), ProductKE(['LIN', 'LIN', 'WN'])),
     # Nested Singleton Extractions
         # Base term singletons
     (ChangeKE('CP', ProductKE(['PER', 'C'], [SumKE(['WN', 'C', 'C'])]), SumKE([], [ProductKE(['WN', 'C'])]))._initialise(), ChangeKE('CP', ProductKE(['PER'], [SumKE(['WN', 'C'])]), 'WN')._initialise()),

setup.py

@@ -14,7 +14,7 @@ def read(filename):
 
 setup(
     name="GPy-ABCD",
-    version="0.1.3",
+    version="0.1.4",
     url="https://github.com/T-Flet/GPy-ABCD",
     license='BSD 3-Clause',