diff --git a/aifeynman/RPN_to_eq.py b/aifeynman/RPN_to_eq.py index 018bc54..c4810fa 100644 --- a/aifeynman/RPN_to_eq.py +++ b/aifeynman/RPN_to_eq.py @@ -17,7 +17,7 @@ def RPN_to_eq(expr): elif i == "0": stack = np.append(stack,"0") elif i == "1": - stack = np.append(stack,"1") + stack = np.append(stack,"1") else: stack = np.append(stack,"x" + str(ord(i)-97)) elif i in operations_2: diff --git a/aifeynman/S_NN_eval.py b/aifeynman/S_NN_eval.py index 93a4d81..30ff66c 100644 --- a/aifeynman/S_NN_eval.py +++ b/aifeynman/S_NN_eval.py @@ -1,4 +1,5 @@ from __future__ import print_function +from typing import Any, Callable, Optional import torch import torch.nn as nn import torch.nn.functional as F @@ -12,6 +13,8 @@ from matplotlib import pyplot as plt import time +from aifeynman.model import DefaultSimpleNet + is_cuda = torch.cuda.is_available() bs = 2048 @@ -40,7 +43,7 @@ def rmse_loss(pred, targ): return torch.sqrt(F.mse_loss(pred, targ))/denom -def NN_eval(pathdir,filename): +def NN_eval(pathdir,filename, torch_model_class: Optional[Callable[[Any], nn.Module]]=None): try: n_variables = np.loadtxt(pathdir+filename, dtype='str').shape[1]-1 variables = np.loadtxt(pathdir+filename, usecols=(0,)) @@ -76,35 +79,20 @@ def NN_eval(pathdir,filename): else: factors_val = factors_val factors_val = factors_val.float() - product_val = torch.from_numpy(f_dependent[int(5*len(variables)/6):int(len(variables))]) + product_val = torch.from_numpy(f_dependent[int(5*len(variables)/6):int(len(variables))]) if is_cuda: product_val = product_val.cuda() else: product_val = product_val product_val = product_val.float() - class SimpleNet(nn.Module): - def __init__(self, ni): - super().__init__() - self.linear1 = nn.Linear(ni, 128) - self.linear2 = nn.Linear(128, 128) - self.linear3 = nn.Linear(128, 64) - self.linear4 = nn.Linear(64,64) - self.linear5 = nn.Linear(64,1) - - def forward(self, x): - x = F.tanh(self.linear1(x)) - x = F.tanh(self.linear2(x)) - x = F.tanh(self.linear3(x)) - x = F.tanh(self.linear4(x)) - x = self.linear5(x) - return x + Net = torch_model_class or DefaultSimpleNet if is_cuda: - model = SimpleNet(n_variables).cuda() + model = Net(n_variables).cuda() else: - model = SimpleNet(n_variables) - + model = Net(n_variables) + model.load_state_dict(torch.load("results/NN_trained_models/models/"+filename+".h5")) model.eval() return(rmse_loss(model(factors_val),product_val),model) diff --git a/aifeynman/S_NN_train.py b/aifeynman/S_NN_train.py index a0d1b7b..40ed982 100644 --- a/aifeynman/S_NN_train.py +++ b/aifeynman/S_NN_train.py @@ -1,18 +1,17 @@ from __future__ import print_function +from typing import Any, Callable, Optional import torch import torch.nn as nn import torch.nn.functional as F import torch.optim as optim -import pandas as pd import numpy as np import torch from torch.utils import data -import pickle -from matplotlib import pyplot as plt import torch.utils.data as utils -import time import os +from aifeynman.model import DefaultSimpleNet + bs = 2048 wd = 1e-2 @@ -40,7 +39,7 @@ def rmse_loss(pred, targ): denom = torch.sqrt(denom.sum()/len(denom)) return torch.sqrt(F.mse_loss(pred, targ))/denom -def NN_train(pathdir, filename, epochs=1000, lrs=1e-2, N_red_lr=4, pretrained_path=""): +def NN_train(pathdir, filename, epochs=1000, lrs=1e-2, N_red_lr=4, pretrained_path="", torch_model_class: Optional[Callable[[Any], nn.Module]]=None): try: os.mkdir("results/NN_trained_models/") except: @@ -54,9 +53,8 @@ def NN_train(pathdir, filename, epochs=1000, lrs=1e-2, N_red_lr=4, pretrained_pa n_variables = np.loadtxt(pathdir+"%s" %filename, dtype='str').shape[1]-1 variables = np.loadtxt(pathdir+"%s" %filename, usecols=(0,)) - # epochs = 200*n_variables + epochs = 200*n_variables if len(variables)<5000: - print('WARNING: tripling epochs since len(variables)<5000...') epochs = epochs*3 if n_variables==0 or n_variables==1: @@ -84,30 +82,15 @@ def NN_train(pathdir, filename, epochs=1000, lrs=1e-2, N_red_lr=4, pretrained_pa product = product product = product.float() - class SimpleNet(nn.Module): - def __init__(self, ni): - super().__init__() - self.linear1 = nn.Linear(ni, 128) - self.linear2 = nn.Linear(128, 128) - self.linear3 = nn.Linear(128, 64) - self.linear4 = nn.Linear(64,64) - self.linear5 = nn.Linear(64,1) - - def forward(self, x): - x = F.tanh(self.linear1(x)) - x = F.tanh(self.linear2(x)) - x = F.tanh(self.linear3(x)) - x = F.tanh(self.linear4(x)) - x = self.linear5(x) - return x + Net = torch_model_class or DefaultSimpleNet my_dataset = utils.TensorDataset(factors,product) # create your datset my_dataloader = utils.DataLoader(my_dataset, batch_size=bs, shuffle=True) # create your dataloader if is_cuda: - model_feynman = SimpleNet(n_variables).cuda() + model_feynman = Net(n_variables).cuda() else: - model_feynman = SimpleNet(n_variables) + model_feynman = Net(n_variables) if pretrained_path!="": model_feynman.load_state_dict(torch.load(pretrained_path)) @@ -120,18 +103,18 @@ def forward(self, x): model_feynman.train() for i, data in enumerate(my_dataloader): optimizer_feynman.zero_grad() - + if is_cuda: fct = data[0].float().cuda() prd = data[1].float().cuda() else: fct = data[0].float() prd = data[1].float() - + loss = rmse_loss(model_feynman(fct),prd) loss.backward() optimizer_feynman.step() - + ''' # Early stopping if epoch%20==0 and epoch>0: @@ -145,7 +128,7 @@ def forward(self, x): torch.save(model_feynman.state_dict(), "results/NN_trained_models/models/" + filename + ".h5") check_es_loss = loss ''' - torch.save(model_feynman.state_dict(), "results/NN_trained_models/models/" + filename + ".h5") + torch.save(model_feynman.state_dict(), "results/NN_trained_models/models/" + filename + ".h5") lrs = lrs/10 return model_feynman diff --git a/aifeynman/S_compositionality.py b/aifeynman/S_compositionality.py index 3ffa182..9c21932 100644 --- a/aifeynman/S_compositionality.py +++ b/aifeynman/S_compositionality.py @@ -62,7 +62,7 @@ def check_compositionality(pathdir,filename,model,express,mu,sigma,nu=10): else: i = i + 1 - + if i==len(data[0:1000]) and np.mean(list_z)1): for j in range(len(possible_vars)-len(temp_list),temp_list[-1]-len(temp_list)+1,-1): exp1 = exp1.replace(possible_vars[j],possible_vars[j+1]) temp_list = np.delete(temp_list,-1) - + # replace variables in bf_eq arr_idx = np.flip(np.arange(0,len(gen_sym_idx),1), axis=0) actual_idx = np.flip(gen_sym_idx, axis=0) diff --git a/aifeynman/S_get_number_DL.py b/aifeynman/S_get_number_DL.py index a69f813..e712397 100644 --- a/aifeynman/S_get_number_DL.py +++ b/aifeynman/S_get_number_DL.py @@ -1,6 +1,6 @@ # Calculates the complexity of a number to be used for the Pareto frontier -import numpy as np +import numpy as np def get_number_DL(n): epsilon = 1e-10 diff --git a/aifeynman/S_gradient_decomposition.py b/aifeynman/S_gradient_decomposition.py index 92c8842..1ba5c8f 100644 --- a/aifeynman/S_gradient_decomposition.py +++ b/aifeynman/S_gradient_decomposition.py @@ -31,7 +31,7 @@ def powerset_atleast_2(iterable, max_subset_size): return r def evaluate_derivatives(model, s, pts): - + pts = pts.clone().detach() try: device = 'cuda' if model.is_cuda else 'cpu' @@ -59,7 +59,7 @@ def evaluate_derivatives_andrew(model, s, pts): pts = pts.cuda() model = model.cuda() grad_weights = grad_weights.cuda() - + pts.requires_grad_(True) outs = model(pts) grad = torch.autograd.grad(outs, pts, grad_outputs=grad_weights, create_graph=True)[0] @@ -74,7 +74,7 @@ def forward(self, X): def draw_samples(X, y, model, s, NUM_SAMPLES, point = None): ''' - Draw samples by sampling each dimension independently, + Draw samples by sampling each dimension independently, keeping the positions at s fixed to given point if exists, sampled point if not. ''' @@ -205,7 +205,7 @@ def filter_decompositions_relative_scoring(X, y, model, max_subset_size=None, vi bench_scores.append(score) snr = signal_to_noise(hypot_scores, bench_scores) # penalizes larger decompositions - snr -= np.log10(2)*len(s) + snr -= np.log10(2)*len(s) results.append((snr, s)) print((snr, s)) if visualize: @@ -254,7 +254,7 @@ def identify_decompositions(pathdir,filename, model, max_subset_size=2, visualiz data = np.loadtxt(pathdir+filename) X = torch.Tensor(data[:, :-1]) y = torch.Tensor(data[:, [-1]]) - # Return best decomposition + # Return best decomposition all_scores = filter_decompositions_relative_scoring(X, y, model, visualize=visualize) assert(all_scores) best_decomposition = all_scores[0][1] diff --git a/aifeynman/S_polyfit_utils.py b/aifeynman/S_polyfit_utils.py index 1f009d0..f8d09f6 100644 --- a/aifeynman/S_polyfit_utils.py +++ b/aifeynman/S_polyfit_utils.py @@ -17,7 +17,7 @@ def as_tall(x): def multipolyfit(xs, y, deg): - + y = asarray(y).squeeze() rows = y.shape[0] xs = asarray(xs) @@ -28,9 +28,9 @@ def multipolyfit(xs, y, deg): xs = np.reshape(xs,(len(xs),1)) xs = hstack((ones((xs.shape[0], 1), dtype=xs.dtype) , xs)) - + generators = [basis_vector(num_covariates+1, i) for i in range(num_covariates+1)] - + # All combinations of degrees powers = map(sum, itertools.combinations_with_replacement(generators, deg)) @@ -38,7 +38,7 @@ def multipolyfit(xs, y, deg): A = hstack(asarray([as_tall((xs**p).prod(1)) for p in powers])) params = lsqr(A, y)[0] # get the best params of the fit rms = lsqr(A, y)[4] # get the rms params of the fit - + return (params, rms) diff --git a/aifeynman/S_remove_input_neuron.py b/aifeynman/S_remove_input_neuron.py index 7ac383a..d0516b4 100644 --- a/aifeynman/S_remove_input_neuron.py +++ b/aifeynman/S_remove_input_neuron.py @@ -19,7 +19,7 @@ def remove_input_neuron(net,n_inp,idx_neuron,ct_median,save_filename): removed_weights = net.linear1.weight[:,idx_neuron] - # Remove the weights associated with the removed input neuron + # Remove the weights associated with the removed input neuron t = torch.transpose(net.linear1.weight,0,1) preserved_ids = torch.LongTensor(np.array(list(set(range(n_inp)) - set([idx_neuron])))) t = nn.Parameter(t[preserved_ids, :]) diff --git a/aifeynman/S_run_aifeynman.py b/aifeynman/S_run_aifeynman.py index 54a80dc..396cc65 100644 --- a/aifeynman/S_run_aifeynman.py +++ b/aifeynman/S_run_aifeynman.py @@ -1,26 +1,21 @@ +from typing import Any, Callable, Optional import numpy as np import matplotlib.pyplot as plt import os from os import path from .get_pareto import Point, ParetoSet -from .RPN_to_pytorch import RPN_to_pytorch -from .RPN_to_eq import RPN_to_eq +from torch.nn import Module from .S_NN_train import NN_train from .S_NN_eval import NN_eval from .S_symmetry import * from .S_separability import * from .S_change_output import * -from .S_brute_force import brute_force from .S_combine_pareto import combine_pareto -from sympy.parsing.sympy_parser import parse_expr -from sympy import preorder_traversal, count_ops -from .S_polyfit import polyfit from .S_get_symbolic_expr_error import get_symbolic_expr_error from .S_add_snap_expr_on_pareto import add_snap_expr_on_pareto from .S_add_sym_on_pareto import add_sym_on_pareto from .S_run_bf_polyfit import run_bf_polyfit from .S_final_gd import final_gd -from .S_add_bf_on_numbers_on_pareto import add_bf_on_numbers_on_pareto from .dimensionalAnalysis import dimensionalAnalysis from .S_NN_get_gradients import evaluate_derivatives from .S_brute_force_comp import brute_force_comp @@ -30,7 +25,8 @@ from .S_gradient_decomposition import identify_decompositions PA = ParetoSet() -def run_AI_all(pathdir,filename,BF_try_time=60,BF_ops_file_type="14ops", polyfit_deg=4, NN_epochs=4000, PA=PA): +def run_AI_all(pathdir,filename,BF_try_time=60,BF_ops_file_type="14ops", polyfit_deg=4, NN_epochs=4000, PA=PA, + torch_model_class: Optional[Callable[[Any], Module]]=None): try: os.mkdir("results/") except: @@ -41,7 +37,7 @@ def run_AI_all(pathdir,filename,BF_try_time=60,BF_ops_file_type="14ops", polyfit # Run bf and polyfit PA = run_bf_polyfit(pathdir,pathdir,filename,BF_try_time,BF_ops_file_type, PA, polyfit_deg) - PA = get_squared(pathdir,"results/mystery_world_squared/",filename,BF_try_time,BF_ops_file_type, PA, polyfit_deg) + PA = get_squared(pathdir,"results/mystery_world_squared/",filename,BF_try_time,BF_ops_file_type, PA, polyfit_deg) # Run bf and polyfit on modified output @@ -64,31 +60,30 @@ def run_AI_all(pathdir,filename,BF_try_time=60,BF_ops_file_type="14ops", polyfit pass elif path.exists("results/NN_trained_models/models/" + filename + ".h5"):# or len(data[0])<3: print("NN already trained \n") - print("NN loss: ", NN_eval(pathdir,filename)[0], "\n") - model_feynman = NN_eval(pathdir,filename)[1] + print("NN loss: ", NN_eval(pathdir,filename, torch_model_class=torch_model_class)[0], "\n") + model_feynman = NN_eval(pathdir,filename, torch_model_class=torch_model_class)[1] elif path.exists("results/NN_trained_models/models/" + filename + "_pretrained.h5"): print("Found pretrained NN \n") - model_feynman = NN_train(pathdir,filename,NN_epochs/2,lrs=1e-3,N_red_lr=3,pretrained_path="results/NN_trained_models/models/" + filename + "_pretrained.h5") - print("NN loss after training: ", NN_eval(pathdir,filename), "\n") + model_feynman = NN_train(pathdir,filename,NN_epochs/2,lrs=1e-3,N_red_lr=3,pretrained_path="results/NN_trained_models/models/" + filename + "_pretrained.h5", torch_model_class=torch_model_class) + print("NN loss after training: ", NN_eval(pathdir,filename, torch_model_class=torch_model_class), "\n") else: print("Training a NN on the data... \n") - model_feynman = NN_train(pathdir,filename,NN_epochs) - print("NN loss: ", NN_eval(pathdir,filename), "\n") + model_feynman = NN_train(pathdir,filename,NN_epochs, torch_model_class=torch_model_class) + print("NN loss: ", NN_eval(pathdir,filename, torch_model_class=torch_model_class), "\n") - # Check which symmetry/separability is the best # Symmetries print("Checking for symmetries...") - symmetry_minus_result = check_translational_symmetry_minus(pathdir,filename) - symmetry_divide_result = check_translational_symmetry_divide(pathdir,filename) - symmetry_multiply_result = check_translational_symmetry_multiply(pathdir,filename) - symmetry_plus_result = check_translational_symmetry_plus(pathdir,filename) + symmetry_minus_result = check_translational_symmetry_minus(pathdir,filename, torch_model_class=torch_model_class) + symmetry_divide_result = check_translational_symmetry_divide(pathdir,filename, torch_model_class=torch_model_class) + symmetry_multiply_result = check_translational_symmetry_multiply(pathdir,filename, torch_model_class=torch_model_class) + symmetry_plus_result = check_translational_symmetry_plus(pathdir,filename, torch_model_class=torch_model_class) print("") print("Checking for separabilities...") # Separabilities - separability_plus_result = check_separability_plus(pathdir,filename) - separability_multiply_result = check_separability_multiply(pathdir,filename) + separability_plus_result = check_separability_plus(pathdir,filename, torch_model_class=torch_model_class) + separability_multiply_result = check_separability_multiply(pathdir,filename, torch_model_class=torch_model_class) if symmetry_plus_result[0]==-1: idx_min = -1 @@ -113,7 +108,7 @@ def run_AI_all(pathdir,filename,BF_try_time=60,BF_ops_file_type="14ops", polyfit print("Checking for compositionality...") # Save the gradients for compositionality try: - succ_grad = evaluate_derivatives(pathdir,filename,model_feynman) + succ_grad = evaluate_derivatives(pathdir, filename, model_feynman) except: succ_grad = 0 @@ -140,7 +135,7 @@ def run_AI_all(pathdir,filename,BF_try_time=60,BF_ops_file_type="14ops", polyfit else: idx_comp = 0 print("") - + if idx_comp==1: idx_min = 6 @@ -154,7 +149,7 @@ def run_AI_all(pathdir,filename,BF_try_time=60,BF_ops_file_type="14ops", polyfit decomp_idx = identify_decompositions(pathdir,filename, model_feynman) brute_force_gen_sym("results/","gradients_gen_sym_%s" %filename,600,"14ops.txt") bf_all_output = np.loadtxt("results_gen_sym.dat", dtype="str") - + for bf_i in range(len(bf_all_output)): idx_gen_sym_temp = 0 try: @@ -175,54 +170,54 @@ def run_AI_all(pathdir,filename,BF_try_time=60,BF_ops_file_type="14ops", polyfit # Apply the best symmetry/separability and rerun the main function on this new file if idx_min == 0: print("Translational symmetry found for variables:", symmetry_plus_result[1],symmetry_plus_result[2]) - new_pathdir, new_filename = do_translational_symmetry_plus(pathdir,filename,symmetry_plus_result[1],symmetry_plus_result[2]) + new_pathdir, new_filename = do_translational_symmetry_plus(pathdir,filename,symmetry_plus_result[1],symmetry_plus_result[2], torch_model_class=torch_model_class) PA1_ = ParetoSet() - PA1 = run_AI_all(new_pathdir,new_filename,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_) + PA1 = run_AI_all(new_pathdir,new_filename,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_, torch_model_class=torch_model_class) PA = add_sym_on_pareto(pathdir,filename,PA1,symmetry_plus_result[1],symmetry_plus_result[2],PA,"+") return PA elif idx_min == 1: print("Translational symmetry found for variables:", symmetry_minus_result[1],symmetry_minus_result[2]) - new_pathdir, new_filename = do_translational_symmetry_minus(pathdir,filename,symmetry_minus_result[1],symmetry_minus_result[2]) + new_pathdir, new_filename = do_translational_symmetry_minus(pathdir,filename,symmetry_minus_result[1],symmetry_minus_result[2], torch_model_class=torch_model_class) PA1_ = ParetoSet() - PA1 = run_AI_all(new_pathdir,new_filename,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_) + PA1 = run_AI_all(new_pathdir,new_filename,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_, torch_model_class=torch_model_class) PA = add_sym_on_pareto(pathdir,filename,PA1,symmetry_minus_result[1],symmetry_minus_result[2],PA,"-") return PA elif idx_min == 2: print("Translational symmetry found for variables:", symmetry_multiply_result[1],symmetry_multiply_result[2]) - new_pathdir, new_filename = do_translational_symmetry_multiply(pathdir,filename,symmetry_multiply_result[1],symmetry_multiply_result[2]) + new_pathdir, new_filename = do_translational_symmetry_multiply(pathdir,filename,symmetry_multiply_result[1],symmetry_multiply_result[2], torch_model_class=torch_model_class) PA1_ = ParetoSet() - PA1 = run_AI_all(new_pathdir,new_filename,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_) + PA1 = run_AI_all(new_pathdir,new_filename,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_, torch_model_class=torch_model_class) PA = add_sym_on_pareto(pathdir,filename,PA1,symmetry_multiply_result[1],symmetry_multiply_result[2],PA,"*") return PA elif idx_min == 3: print("Translational symmetry found for variables:", symmetry_divide_result[1],symmetry_divide_result[2]) - new_pathdir, new_filename = do_translational_symmetry_divide(pathdir,filename,symmetry_divide_result[1],symmetry_divide_result[2]) + new_pathdir, new_filename = do_translational_symmetry_divide(pathdir,filename,symmetry_divide_result[1],symmetry_divide_result[2], torch_model_class=torch_model_class) PA1_ = ParetoSet() - PA1 = run_AI_all(new_pathdir,new_filename,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_) + PA1 = run_AI_all(new_pathdir,new_filename,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_, torch_model_class=torch_model_class) PA = add_sym_on_pareto(pathdir,filename,PA1,symmetry_divide_result[1],symmetry_divide_result[2],PA,"/") return PA elif idx_min == 4: print("Additive separability found for variables:", separability_plus_result[1],separability_plus_result[2]) - new_pathdir1, new_filename1, new_pathdir2, new_filename2, = do_separability_plus(pathdir,filename,separability_plus_result[1],separability_plus_result[2]) + new_pathdir1, new_filename1, new_pathdir2, new_filename2, = do_separability_plus(pathdir,filename,separability_plus_result[1],separability_plus_result[2], torch_model_class=torch_model_class) PA1_ = ParetoSet() - PA1 = run_AI_all(new_pathdir1,new_filename1,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_) + PA1 = run_AI_all(new_pathdir1,new_filename1,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_, torch_model_class=torch_model_class) PA2_ = ParetoSet() - PA2 = run_AI_all(new_pathdir2,new_filename2,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA2_) + PA2 = run_AI_all(new_pathdir2,new_filename2,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA2_, torch_model_class=torch_model_class) combine_pareto_data = np.loadtxt(pathdir+filename) PA = combine_pareto(combine_pareto_data,PA1,PA2,separability_plus_result[1],separability_plus_result[2],PA,"+") return PA elif idx_min == 5: print("Multiplicative separability found for variables:", separability_multiply_result[1],separability_multiply_result[2]) - new_pathdir1, new_filename1, new_pathdir2, new_filename2, = do_separability_multiply(pathdir,filename,separability_multiply_result[1],separability_multiply_result[2]) + new_pathdir1, new_filename1, new_pathdir2, new_filename2, = do_separability_multiply(pathdir,filename,separability_multiply_result[1],separability_multiply_result[2], torch_model_class=torch_model_class) PA1_ = ParetoSet() - PA1 = run_AI_all(new_pathdir1,new_filename1,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_) + PA1 = run_AI_all(new_pathdir1,new_filename1,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_, torch_model_class=torch_model_class) PA2_ = ParetoSet() - PA2 = run_AI_all(new_pathdir2,new_filename2,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA2_) + PA2 = run_AI_all(new_pathdir2,new_filename2,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA2_, torch_model_class=torch_model_class) combine_pareto_data = np.loadtxt(pathdir+filename) PA = combine_pareto(combine_pareto_data,PA1,PA2,separability_multiply_result[1],separability_multiply_result[2],PA,"*") return PA @@ -231,7 +226,7 @@ def run_AI_all(pathdir,filename,BF_try_time=60,BF_ops_file_type="14ops", polyfit print("Compositionality found") new_pathdir, new_filename = do_compositionality(pathdir,filename,math_eq_comp) PA1_ = ParetoSet() - PA1 = run_AI_all(new_pathdir,new_filename,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_) + PA1 = run_AI_all(new_pathdir,new_filename,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_, torch_model_class=torch_model_class) PA = add_comp_on_pareto(PA1,PA,math_eq_comp) return PA @@ -239,13 +234,13 @@ def run_AI_all(pathdir,filename,BF_try_time=60,BF_ops_file_type="14ops", polyfit print("Generalized symmetry found") new_pathdir, new_filename = do_gen_sym(pathdir,filename,decomp_idx,math_eq_gen_sym) PA1_ = ParetoSet() - PA1 = run_AI_all(new_pathdir,new_filename,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_) + PA1 = run_AI_all(new_pathdir,new_filename,BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA1_, torch_model_class=torch_model_class) PA = add_gen_sym_on_pareto(PA1,PA, decomp_idx, math_eq_gen_sym) return PA else: return PA # this runs snap on the output of aifeynman -def run_aifeynman(pathdir,filename,BF_try_time,BF_ops_file_type, polyfit_deg=4, NN_epochs=4000, vars_name=[],test_percentage=20): +def run_aifeynman(pathdir,filename,BF_try_time,BF_ops_file_type, polyfit_deg=4, NN_epochs=4000, vars_name=[],test_percentage=20, torch_model_class: Optional[Callable[[Any], Module]]=None): # If the variable names are passed, do the dimensional analysis first filename_orig = filename try: @@ -271,7 +266,7 @@ def run_aifeynman(pathdir,filename,BF_try_time,BF_ops_file_type, polyfit_deg=4, PA = ParetoSet() # Run the code on the train data - PA = run_AI_all(pathdir,filename+"_train",BF_try_time,BF_ops_file_type, polyfit_deg, NN_epochs, PA=PA) + PA = run_AI_all(pathdir, filename+"_train", BF_try_time, BF_ops_file_type, polyfit_deg, NN_epochs, PA=PA, torch_model_class=torch_model_class) PA_list = PA.get_pareto_points() ''' @@ -297,8 +292,7 @@ def run_aifeynman(pathdir,filename,BF_try_time,BF_ops_file_type, polyfit_deg=4, # Run gradient descent on the data one more time for i in range(len(PA_list)): try: - dt = np.loadtxt(pathdir+filename) - gd_update = final_gd(dt,PA_list[i][-1]) + gd_update = final_gd(pathdir,filename,PA_list[i][-1]) PA.add(Point(x=gd_update[1],y=gd_update[0],data=gd_update[2])) except: continue diff --git a/aifeynman/S_separability.py b/aifeynman/S_separability.py index 47a5016..8f6fc96 100644 --- a/aifeynman/S_separability.py +++ b/aifeynman/S_separability.py @@ -1,44 +1,24 @@ from __future__ import print_function +from typing import Any, Callable, Optional import torch import os import torch.nn as nn import torch.nn.functional as F -import torch.optim as optim -import pandas as pd import numpy as np import torch -from torch.utils import data -import pickle -from torch.optim.lr_scheduler import CosineAnnealingLR -from matplotlib import pyplot as plt from itertools import combinations -import time + +from aifeynman.model import DefaultSimpleNet is_cuda = torch.cuda.is_available() -class SimpleNet(nn.Module): - def __init__(self, ni): - super().__init__() - self.linear1 = nn.Linear(ni, 128) - self.linear2 = nn.Linear(128, 128) - self.linear3 = nn.Linear(128, 64) - self.linear4 = nn.Linear(64,64) - self.linear5 = nn.Linear(64,1) - - def forward(self, x): - x = F.tanh(self.linear1(x)) - x = F.tanh(self.linear2(x)) - x = F.tanh(self.linear3(x)) - x = F.tanh(self.linear4(x)) - x = self.linear5(x) - return x def rmse_loss(pred, targ): denom = targ**2 denom = torch.sqrt(denom.sum()/len(denom)) return torch.sqrt(F.mse_loss(pred, targ))/denom -def check_separability_plus(pathdir, filename): +def check_separability_plus(pathdir, filename, torch_model_class: Optional[Callable[[Any], nn.Module]]=None): try: pathdir_weights = "results/NN_trained_models/models/" @@ -54,12 +34,12 @@ def check_separability_plus(pathdir, filename): for j in range(1,n_variables): v = np.loadtxt(pathdir+filename, usecols=(j,)) variables = np.column_stack((variables,v)) - + f_dependent = np.loadtxt(pathdir+filename, usecols=(n_variables,)) f_dependent = np.reshape(f_dependent,(len(f_dependent),1)) - factors = torch.from_numpy(variables) + factors = torch.from_numpy(variables) if is_cuda: factors = factors.cuda() else: @@ -73,11 +53,12 @@ def check_separability_plus(pathdir, filename): product = product product = product.float() + Net = torch_model_class or DefaultSimpleNet # load the trained model and put it in evaluation mode if is_cuda: - model = SimpleNet(n_variables).cuda() + model = Net(n_variables).cuda() else: - model = SimpleNet(n_variables) + model = Net(n_variables) model.load_state_dict(torch.load(pathdir_weights+filename+".h5")) model.eval() @@ -122,13 +103,13 @@ def check_separability_plus(pathdir, filename): best_mu = mu best_sigma = sigma return min_error, best_i, best_j, best_mu, best_sigma - + except Exception as e: print(e) - return (-1,-1,-1,-1,-1) - - -def do_separability_plus(pathdir, filename, list_i,list_j): + return (-1,-1,-1,-1,-1) + + +def do_separability_plus(pathdir, filename, list_i,list_j, torch_model_class: Optional[Callable[[Any], nn.Module]]=None): try: pathdir_weights = "results/NN_trained_models/models/" @@ -144,12 +125,12 @@ def do_separability_plus(pathdir, filename, list_i,list_j): for j in range(1,n_variables): v = np.loadtxt(pathdir+filename, usecols=(j,)) variables = np.column_stack((variables,v)) - + f_dependent = np.loadtxt(pathdir+filename, usecols=(n_variables,)) f_dependent = np.reshape(f_dependent,(len(f_dependent),1)) - factors = torch.from_numpy(variables) + factors = torch.from_numpy(variables) if is_cuda: factors = factors.cuda() else: @@ -163,18 +144,20 @@ def do_separability_plus(pathdir, filename, list_i,list_j): product = product product = product.float() + Net = torch_model_class or DefaultSimpleNet + # load the trained model and put it in evaluation mode if is_cuda: - model = SimpleNet(n_variables).cuda() + model = Net(n_variables).cuda() else: - model = SimpleNet(n_variables) + model = Net(n_variables) model.load_state_dict(torch.load(pathdir_weights+filename+".h5")) model.eval() # make some variables at the time equal to the median of factors models_one = [] models_rest = [] - + fact_vary = factors.clone() for k in range(len(factors[0])): fact_vary[:,k] = torch.full((len(factors),),torch.median(factors[:,k])) @@ -192,7 +175,7 @@ def do_separability_plus(pathdir, filename, list_i,list_j): data_sep_1 = variables data_sep_1 = np.delete(data_sep_1,list_j,axis=1) data_sep_1 = np.column_stack((data_sep_1,model(fact_vary_one).cpu())) - # save the second half + # save the second half data_sep_2 = variables data_sep_2 = np.delete(data_sep_2,list_i,axis=1) data_sep_2 = np.column_stack((data_sep_2,model(fact_vary_rest).cpu()-model(fact_vary).cpu())) @@ -209,8 +192,8 @@ def do_separability_plus(pathdir, filename, list_i,list_j): print(e) return (-1,-1) - -def check_separability_multiply(pathdir, filename): + +def check_separability_multiply(pathdir, filename, torch_model_class: Optional[Callable[[Any], nn.Module]]=None): try: pathdir_weights = "results/NN_trained_models/models/" @@ -226,11 +209,11 @@ def check_separability_multiply(pathdir, filename): for j in range(1,n_variables): v = np.loadtxt(pathdir+filename, usecols=(j,)) variables = np.column_stack((variables,v)) - + f_dependent = np.loadtxt(pathdir+filename, usecols=(n_variables,)) - # Pick only data which is close enough to the maximum value (5 times less or higher) + # Pick only data which is close enough to the maximum value (5 times less or higher) max_output = np.max(abs(f_dependent)) use_idx = np.where(abs(f_dependent)>=max_output/5) f_dependent = f_dependent[use_idx] @@ -251,11 +234,13 @@ def check_separability_multiply(pathdir, filename): product = product product = product.float() + Net = torch_model_class or DefaultSimpleNet + # load the trained model and put it in evaluation mode if is_cuda: - model = SimpleNet(n_variables).cuda() + model = Net(n_variables).cuda() else: - model = SimpleNet(n_variables) + model = Net(n_variables) model.load_state_dict(torch.load(pathdir_weights+filename+".h5")) model.eval() @@ -299,14 +284,12 @@ def check_separability_multiply(pathdir, filename): best_mu = mu best_sigma = sigma return min_error, best_i, best_j, best_mu, best_sigma - + except Exception as e: print(e) - return (-1,-1,-1,-1,-1) + return (-1,-1,-1,-1,-1) - - -def do_separability_multiply(pathdir, filename, list_i,list_j): +def do_separability_multiply(pathdir, filename, list_i,list_j, torch_model_class: Optional[Callable[[Any], nn.Module]]=None): try: pathdir_weights = "results/NN_trained_models/models/" @@ -322,7 +305,7 @@ def do_separability_multiply(pathdir, filename, list_i,list_j): for j in range(1,n_variables): v = np.loadtxt(pathdir+filename, usecols=(j,)) variables = np.column_stack((variables,v)) - + f_dependent = np.loadtxt(pathdir+filename, usecols=(n_variables,)) f_dependent = np.reshape(f_dependent,(len(f_dependent),1)) @@ -341,17 +324,19 @@ def do_separability_multiply(pathdir, filename, list_i,list_j): product = product product = product.float() + Net = torch_model_class or DefaultSimpleNet + # load the trained model and put it in evaluation mode if is_cuda: - model = SimpleNet(n_variables).cuda() + model = Net(n_variables).cuda() else: - model = SimpleNet(n_variables) + model = Net(n_variables) model.load_state_dict(torch.load(pathdir_weights+filename+".h5")) model.eval() # make some variables at the time equal to the median of factors models_one = [] - models_rest = [] + models_rest = [] fact_vary = factors.clone() for k in range(len(factors[0])): @@ -361,8 +346,8 @@ def do_separability_multiply(pathdir, filename, list_i,list_j): for t1 in list_j: fact_vary_one[:,t1] = torch.full((len(factors),),torch.median(factors[:,t1])) for t2 in list_i: - fact_vary_rest[:,t2] = torch.full((len(factors),),torch.median(factors[:,t2])) - + fact_vary_rest[:,t2] = torch.full((len(factors),),torch.median(factors[:,t2])) + with torch.no_grad(): str1 = filename+"-mult_a" str2 = filename+"-mult_b" @@ -370,7 +355,7 @@ def do_separability_multiply(pathdir, filename, list_i,list_j): data_sep_1 = variables data_sep_1 = np.delete(data_sep_1,list_j,axis=1) data_sep_1 = np.column_stack((data_sep_1,model(fact_vary_one).cpu())) - # save the second half + # save the second half data_sep_2 = variables data_sep_2 = np.delete(data_sep_2,list_i,axis=1) data_sep_2 = np.column_stack((data_sep_2,model(fact_vary_rest).cpu()/model(fact_vary).cpu())) @@ -387,4 +372,4 @@ def do_separability_multiply(pathdir, filename, list_i,list_j): print(e) return (-1,-1) - + diff --git a/aifeynman/S_snap.py b/aifeynman/S_snap.py index 5272d93..cd0fe80 100644 --- a/aifeynman/S_snap.py +++ b/aifeynman/S_snap.py @@ -19,7 +19,7 @@ def float2contfrac(x,nmax): y = y - i k = k + 1 return c - + def contfrac2frac(seq): ''' Convert the simple continued fraction in `seq` into a fraction, num / den @@ -28,27 +28,27 @@ def contfrac2frac(seq): for u in reversed(seq): num, den = den + num*u, num return num, den - + def contFracRationalApproximations(c): return np.array(list(contfrac2frac(c[:i+1]) for i in range(len(c)))) - + def contFracApproximations(c): q = contFracRationalApproximations(c) return q[:,0] / float(q[:,1]) - + def truncateContFrac(q,imax): k = 0 while k < len(q) and np.maximum(np.abs(q[k,0]), q[k,1]) <= imax: k = k + 1 return q[:k] - + def pval(p): p = p.astype(float) return 1 - np.exp(-p ** 0.87 / 0.36) - + xsign = np.sign(x) q = truncateContFrac(contFracRationalApproximations(float2contfrac(abs(x),20)),imax) - + if len(q) > 0: p = np.abs(q[:,0] / q[:,1] - abs(x)).astype(float) * (1 + np.abs(q[:,0])) * q[:,1] p = pval(p) @@ -75,10 +75,10 @@ def rationalSnap(p, top=1): """Snap to nearest rational number using continued fraction.""" p = np.array(p) snaps = np.array(list(bestApproximation(x,10) for x in p)) - chosen = np.argsort(snaps[:, 3])[:top] + chosen = np.argsort(snaps[:, 3])[:top] d = dict(list(zip(chosen, snaps[chosen, 1:3]))) d = {k: f"{val[0]}/{val[1]}" for k,val in d.items()} - + return d diff --git a/aifeynman/S_symmetry.py b/aifeynman/S_symmetry.py index 06fd210..87a185e 100644 --- a/aifeynman/S_symmetry.py +++ b/aifeynman/S_symmetry.py @@ -1,39 +1,19 @@ # checks for symmetries in the data from __future__ import print_function +from typing import Any, Callable, Optional import torch import os import torch.nn as nn import torch.nn.functional as F -import torch.optim as optim -import pandas as pd import numpy as np import torch -from torch.utils import data -import pickle -from torch.optim.lr_scheduler import CosineAnnealingLR -from matplotlib import pyplot as plt + +from aifeynman.model import DefaultSimpleNet from .S_remove_input_neuron import remove_input_neuron -import time is_cuda = torch.cuda.is_available() -class SimpleNet(nn.Module): - def __init__(self, ni): - super().__init__() - self.linear1 = nn.Linear(ni, 128) - self.linear2 = nn.Linear(128, 128) - self.linear3 = nn.Linear(128, 64) - self.linear4 = nn.Linear(64,64) - self.linear5 = nn.Linear(64,1) - - def forward(self, x): - x = F.tanh(self.linear1(x)) - x = F.tanh(self.linear2(x)) - x = F.tanh(self.linear3(x)) - x = F.tanh(self.linear4(x)) - x = self.linear5(x) - return x def rmse_loss(pred, targ): denom = targ**2 @@ -41,7 +21,7 @@ def rmse_loss(pred, targ): return torch.sqrt(F.mse_loss(pred, targ))/denom # checks if f(x,y)=f(x-y) -def check_translational_symmetry_minus(pathdir, filename): +def check_translational_symmetry_minus(pathdir, filename, torch_model_class: Optional[Callable[[Any], nn.Module]]=None): try: pathdir_weights = "results/NN_trained_models/models/" @@ -57,12 +37,12 @@ def check_translational_symmetry_minus(pathdir, filename): for j in range(1,n_variables): v = np.loadtxt(pathdir+"/%s" %filename, usecols=(j,)) variables = np.column_stack((variables,v)) - + f_dependent = np.loadtxt(pathdir+"/%s" %filename, usecols=(n_variables,)) f_dependent = np.reshape(f_dependent,(len(f_dependent),1)) - factors = torch.from_numpy(variables) + factors = torch.from_numpy(variables) if is_cuda: factors = factors.cuda() else: @@ -76,18 +56,18 @@ def check_translational_symmetry_minus(pathdir, filename): product = product product = product.float() + Net = torch_model_class or DefaultSimpleNet + # load the trained model and put it in evaluation mode if is_cuda: - model = SimpleNet(n_variables).cuda() + model = Net(n_variables).cuda() else: - model = SimpleNet(n_variables) + model = Net(n_variables) model.load_state_dict(torch.load(pathdir_weights+filename+".h5")) model.eval() - models_one = [] - models_rest = [] - with torch.no_grad(): + with torch.no_grad(): # make the shift x->x+a for 2 variables at a time (different variables) min_error = 1000 best_i = -1 @@ -116,8 +96,9 @@ def check_translational_symmetry_minus(pathdir, filename): except Exception as e: print(e) return (-1,-1,-1,-1,-1) - -def do_translational_symmetry_minus(pathdir, filename, i,j): + + +def do_translational_symmetry_minus(pathdir, filename, i,j, torch_model_class: Optional[Callable[[Any], nn.Module]]=None): try: pathdir_weights = "results/NN_trained_models/models/" @@ -128,11 +109,11 @@ def do_translational_symmetry_minus(pathdir, filename, i,j): for k in range(1,n_variables): v = np.loadtxt(pathdir+"/%s" %filename, usecols=(k,)) variables = np.column_stack((variables,v)) - + f_dependent = np.loadtxt(pathdir+"/%s" %filename, usecols=(n_variables,)) f_dependent = np.reshape(f_dependent,(len(f_dependent),1)) - factors = torch.from_numpy(variables) + factors = torch.from_numpy(variables) if is_cuda: factors = factors.cuda() else: @@ -146,18 +127,20 @@ def do_translational_symmetry_minus(pathdir, filename, i,j): product = product product = product.float() + Net = torch_model_class or DefaultSimpleNet + # load the trained model and put it in evaluation mode if is_cuda: - model = SimpleNet(n_variables).cuda() + model = Net(n_variables).cuda() else: - model = SimpleNet(n_variables) + model = Net(n_variables) model.load_state_dict(torch.load(pathdir_weights+filename+".h5")) model.eval() models_one = [] models_rest = [] - with torch.no_grad(): + with torch.no_grad(): file_name = filename + "-translated_minus" ct_median = torch.median(torch.from_numpy(variables[:,j])) data_translated = variables @@ -175,10 +158,10 @@ def do_translational_symmetry_minus(pathdir, filename, i,j): except Exception as e: print(e) return (-1,-1) - + # checks if f(x,y)=f(x/y) -def check_translational_symmetry_divide(pathdir, filename): +def check_translational_symmetry_divide(pathdir, filename, torch_model_class: Optional[Callable[[Any], nn.Module]]=None): try: pathdir_weights = "results/NN_trained_models/models/" @@ -194,12 +177,12 @@ def check_translational_symmetry_divide(pathdir, filename): for j in range(1,n_variables): v = np.loadtxt(pathdir+"/%s" %filename, usecols=(j,)) variables = np.column_stack((variables,v)) - + f_dependent = np.loadtxt(pathdir+"/%s" %filename, usecols=(n_variables,)) f_dependent = np.reshape(f_dependent,(len(f_dependent),1)) - factors = torch.from_numpy(variables) + factors = torch.from_numpy(variables) if is_cuda: factors = factors.cuda() else: @@ -213,16 +196,16 @@ def check_translational_symmetry_divide(pathdir, filename): product = product product = product.float() + Net = torch_model_class or DefaultSimpleNet + # load the trained model and put it in evaluation mode if is_cuda: - model = SimpleNet(n_variables).cuda() + model = Net(n_variables).cuda() else: - model = SimpleNet(n_variables) + model = Net(n_variables) model.load_state_dict(torch.load(pathdir_weights+filename+".h5")) model.eval() - models_one = [] - models_rest = [] with torch.no_grad(): a = 1.2 @@ -253,9 +236,9 @@ def check_translational_symmetry_divide(pathdir, filename): except Exception as e: print(e) return (-1,-1,-1,-1,-1) - - -def do_translational_symmetry_divide(pathdir, filename, i,j): + + +def do_translational_symmetry_divide(pathdir, filename, i,j, torch_model_class: Optional[Callable[[Any], nn.Module]]=None): try: pathdir_weights = "results/NN_trained_models/models/" @@ -266,11 +249,11 @@ def do_translational_symmetry_divide(pathdir, filename, i,j): for k in range(1,n_variables): v = np.loadtxt(pathdir+"/%s" %filename, usecols=(k,)) variables = np.column_stack((variables,v)) - + f_dependent = np.loadtxt(pathdir+"/%s" %filename, usecols=(n_variables,)) f_dependent = np.reshape(f_dependent,(len(f_dependent),1)) - factors = torch.from_numpy(variables) + factors = torch.from_numpy(variables) if is_cuda: factors = factors.cuda() else: @@ -284,18 +267,18 @@ def do_translational_symmetry_divide(pathdir, filename, i,j): product = product product = product.float() + Net = torch_model_class or DefaultSimpleNet + # load the trained model and put it in evaluation mode if is_cuda: - model = SimpleNet(n_variables).cuda() + model = Net(n_variables).cuda() else: - model = SimpleNet(n_variables) + model = Net(n_variables) model.load_state_dict(torch.load(pathdir_weights+filename+".h5")) model.eval() - models_one = [] - models_rest = [] - with torch.no_grad(): + with torch.no_grad(): file_name = filename + "-translated_divide" data_translated = variables ct_median =torch.median(torch.from_numpy(variables[:,j])) @@ -315,7 +298,7 @@ def do_translational_symmetry_divide(pathdir, filename, i,j): return (-1,1) # checks if f(x,y)=f(x*y) -def check_translational_symmetry_multiply(pathdir, filename): +def check_translational_symmetry_multiply(pathdir, filename, torch_model_class: Optional[Callable[[Any], nn.Module]]=None): try: pathdir_weights = "results/NN_trained_models/models/" @@ -331,12 +314,12 @@ def check_translational_symmetry_multiply(pathdir, filename): for j in range(1,n_variables): v = np.loadtxt(pathdir+"/%s" %filename, usecols=(j,)) variables = np.column_stack((variables,v)) - + f_dependent = np.loadtxt(pathdir+"/%s" %filename, usecols=(n_variables,)) f_dependent = np.reshape(f_dependent,(len(f_dependent),1)) - factors = torch.from_numpy(variables) + factors = torch.from_numpy(variables) if is_cuda: factors = factors.cuda() else: @@ -349,17 +332,16 @@ def check_translational_symmetry_multiply(pathdir, filename): else: product = product product = product.float() + Net = torch_model_class or DefaultSimpleNet # load the trained model and put it in evaluation mode if is_cuda: - model = SimpleNet(n_variables).cuda() + model = Net(n_variables).cuda() else: - model = SimpleNet(n_variables) + model = Net(n_variables) model.load_state_dict(torch.load(pathdir_weights+filename+".h5")) model.eval() - models_one = [] - models_rest = [] with torch.no_grad(): a = 1.2 @@ -391,7 +373,7 @@ def check_translational_symmetry_multiply(pathdir, filename): print(e) return (-1,-1,-1,-1,-1) -def do_translational_symmetry_multiply(pathdir, filename, i,j): +def do_translational_symmetry_multiply(pathdir, filename, i,j, torch_model_class: Optional[Callable[[Any], nn.Module]]=None): try: pathdir_weights = "results/NN_trained_models/models/" @@ -402,11 +384,11 @@ def do_translational_symmetry_multiply(pathdir, filename, i,j): for k in range(1,n_variables): v = np.loadtxt(pathdir+"/%s" %filename, usecols=(k,)) variables = np.column_stack((variables,v)) - + f_dependent = np.loadtxt(pathdir+"/%s" %filename, usecols=(n_variables,)) f_dependent = np.reshape(f_dependent,(len(f_dependent),1)) - factors = torch.from_numpy(variables) + factors = torch.from_numpy(variables) if is_cuda: factors = factors.cuda() else: @@ -420,18 +402,18 @@ def do_translational_symmetry_multiply(pathdir, filename, i,j): product = product product = product.float() + Net = torch_model_class or DefaultSimpleNet + # load the trained model and put it in evaluation mode if is_cuda: - model = SimpleNet(n_variables).cuda() + model = Net(n_variables).cuda() else: - model = SimpleNet(n_variables) + model = Net(n_variables) model.load_state_dict(torch.load(pathdir_weights+filename+".h5")) model.eval() - models_one = [] - models_rest = [] - with torch.no_grad(): + with torch.no_grad(): file_name = filename + "-translated_multiply" data_translated = variables ct_median =torch.median(torch.from_numpy(variables[:,j])) @@ -451,7 +433,7 @@ def do_translational_symmetry_multiply(pathdir, filename, i,j): return (-1,1) # checks if f(x,y)=f(x+y) -def check_translational_symmetry_plus(pathdir, filename): +def check_translational_symmetry_plus(pathdir, filename, torch_model_class: Optional[Callable[[Any], nn.Module]]=None): try: pathdir_weights = "results/NN_trained_models/models/" @@ -467,12 +449,12 @@ def check_translational_symmetry_plus(pathdir, filename): for j in range(1,n_variables): v = np.loadtxt(pathdir+"/%s" %filename, usecols=(j,)) variables = np.column_stack((variables,v)) - + f_dependent = np.loadtxt(pathdir+"/%s" %filename, usecols=(n_variables,)) f_dependent = np.reshape(f_dependent,(len(f_dependent),1)) - factors = torch.from_numpy(variables) + factors = torch.from_numpy(variables) if is_cuda: factors = factors.cuda() else: @@ -486,16 +468,16 @@ def check_translational_symmetry_plus(pathdir, filename): product = product product = product.float() + Net = torch_model_class or DefaultSimpleNet + # load the trained model and put it in evaluation mode if is_cuda: - model = SimpleNet(n_variables).cuda() + model = Net(n_variables).cuda() else: - model = SimpleNet(n_variables) + model = Net(n_variables) model.load_state_dict(torch.load(pathdir_weights+filename+".h5")) model.eval() - models_one = [] - models_rest = [] with torch.no_grad(): min_error = 1000 @@ -525,8 +507,8 @@ def check_translational_symmetry_plus(pathdir, filename): except Exception as e: print(e) return (-1,-1,-1,-1,-1) - -def do_translational_symmetry_plus(pathdir, filename, i,j): + +def do_translational_symmetry_plus(pathdir, filename, i,j, torch_model_class: Optional[Callable[[Any], nn.Module]]=None): try: pathdir_weights = "results/NN_trained_models/models/" @@ -537,11 +519,11 @@ def do_translational_symmetry_plus(pathdir, filename, i,j): for k in range(1,n_variables): v = np.loadtxt(pathdir+"/%s" %filename, usecols=(k,)) variables = np.column_stack((variables,v)) - + f_dependent = np.loadtxt(pathdir+"/%s" %filename, usecols=(n_variables,)) f_dependent = np.reshape(f_dependent,(len(f_dependent),1)) - factors = torch.from_numpy(variables) + factors = torch.from_numpy(variables) if is_cuda: factors = factors.cuda() else: @@ -555,18 +537,18 @@ def do_translational_symmetry_plus(pathdir, filename, i,j): product = product product = product.float() + Net = torch_model_class or DefaultSimpleNet + # load the trained model and put it in evaluation mode if is_cuda: - model = SimpleNet(n_variables).cuda() + model = Net(n_variables).cuda() else: - model = SimpleNet(n_variables) + model = Net(n_variables) model.load_state_dict(torch.load(pathdir_weights+filename+".h5")) model.eval() - models_one = [] - models_rest = [] - with torch.no_grad(): + with torch.no_grad(): file_name = filename + "-translated_plus" data_translated = variables ct_median =torch.median(torch.from_numpy(variables[:,j])) @@ -580,7 +562,7 @@ def do_translational_symmetry_plus(pathdir, filename, i,j): np.savetxt("results/translated_data_plus/"+file_name , data_translated) remove_input_neuron(model,n_variables,j,ct_median,"results/NN_trained_models/models/"+filename + "-translated_plus_pretrained.h5") return ("results/translated_data_plus/", file_name) - + except Exception as e: print(e) return (-1,-1) diff --git a/aifeynman/get_pareto.py b/aifeynman/get_pareto.py index 145a2de..0eceb22 100644 --- a/aifeynman/get_pareto.py +++ b/aifeynman/get_pareto.py @@ -70,15 +70,15 @@ def _input_check(self, p): return Point(x=p[0], y=p[1], data=None) else: raise TypeError("Must be instance of Point or 2-tuple.") - - + + def get_id_list(self): id_list = [] for point in self: id_list.append(point.id) return id_list - + def add(self, p): """Insert Point into set if minimal in first two indices. @@ -220,14 +220,14 @@ def to_array(self): def get_pareto_points(self): """Returns the x, y and data for each point in the pareto frontier - + """ pareto_points = [] for i, p in enumerate(self): pareto_points = pareto_points + [[p.x, p.y, p.data]] - + return pareto_points - + def from_list(self, A): """Convert iterable of Points into ParetoSet. @@ -241,8 +241,8 @@ def from_list(self, A): """ for a in A: self.add(a) - - + + def plot(self): """Plotting the Pareto frontier.""" array = self.to_array() @@ -254,14 +254,14 @@ def plot(self): if __name__ == "__main__": PA = ParetoSet() A = np.zeros((40, 2)) - + for i in range(40): x = np.random.rand() y = np.random.rand() - + A[i, 0] = x A[i, 1] = y - + PA.add(Point(x=x, y=y, data=None)) paretoA = PA.to_array() diff --git a/aifeynman/model.py b/aifeynman/model.py new file mode 100644 index 0000000..f17a432 --- /dev/null +++ b/aifeynman/model.py @@ -0,0 +1,20 @@ +import torch.nn as nn +import torch.nn.functional as F + + +class DefaultSimpleNet(nn.Module): + def __init__(self, ni): + super().__init__() + self.linear1 = nn.Linear(ni, 128) + self.linear2 = nn.Linear(128, 128) + self.linear3 = nn.Linear(128, 64) + self.linear4 = nn.Linear(64, 64) + self.linear5 = nn.Linear(64, 1) + + def forward(self, x): + x = F.tanh(self.linear1(x)) + x = F.tanh(self.linear2(x)) + x = F.tanh(self.linear3(x)) + x = F.tanh(self.linear4(x)) + x = self.linear5(x) + return x diff --git a/aifeynman/test_points.py b/aifeynman/test_points.py index 76853c0..25e1de9 100644 --- a/aifeynman/test_points.py +++ b/aifeynman/test_points.py @@ -34,7 +34,7 @@ def project_plane(normal, point, guess): class TestPoints: def __init__(self, eq, symbols, X, y, bitmask, mode='general'): ''' - mode is in {'general', 'add', 'minus', 'times', 'divide'} + mode is in {'general', 'add', 'minus', 'times', 'divide'} eq, symbols are None if mode is not general ''' self.mode = mode @@ -63,7 +63,7 @@ def general_init(self, eq, symbols, X, y, bitmask): self.lambda_grads = [sympy.lambdify(self.symbols, x) for x in self.symp_grads] self.low_range = np.percentile(self.X[:, self.bitmask], 0, axis=0) self.high_range = np.percentile(self.X[:, self.bitmask], 100, axis=0) - + self.init_median_projection() self.init_scatter() @@ -84,7 +84,7 @@ def find_initial_guess_median_projection(self, pt): num_grads = self.evaluate_gradients(pt) num_grads /= np.linalg.norm(num_grads) return [project_plane(num_grads, pt, self.median_point)] - + def find_initial_guess_scatter(self, pt, low=2, high=3): guess = self.find_initial_guess_median_projection(pt) @@ -93,7 +93,7 @@ def find_initial_guess_scatter(self, pt, low=2, high=3): candidates = list(range(max(0, target_index-low), min(self.X.shape[0], target_index+high))) results = [self.X[self.hindices[guess], self.bitmask] for guess in candidates] return results - + def find_initial_guess_random(self, pt, num=2): return [self.X[t, self.bitmask] for t in np.random.randint(0, self.X.shape[0], num)] @@ -110,8 +110,8 @@ def optimize_fmin(self, guess, target): def optimize_bfgs(self, guess, target): MAXITER = 20 res = scipy.optimize.minimize(lambda x: .5 * (self.lambda_eq(*x) - target)**2, - guess, method='BFGS', - jac=lambda x: self.evaluate_gradients(x) * (self.lambda_eq(*x) - target), + guess, method='BFGS', + jac=lambda x: self.evaluate_gradients(x) * (self.lambda_eq(*x) - target), options={'maxiter':MAXITER, 'disp': False}) if res.success: return res.x @@ -146,8 +146,8 @@ def find_reference_point(self, pt): if result is not None and self.in_domain(result): results.append(result) return max(results, key=lambda x: np.linalg.norm(pt - x), default=None) - - + + def analyze_reference_point(self, pt, opt, disp): reference_point_rel_error = relative_error(self.lambda_eq(*opt), self.lambda_eq(*pt)) reference_point_distance = np.linalg.norm(opt - pt) diff --git a/requirements.txt b/requirements.txt index 2c2955f..dbb9395 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,8 +1,6 @@ -python >= 3.6 -torch >= 1.4.0 -matplotlib -sympy = 1.4 +torch>=1.4.0 +sympy==1.4 +sortedcontainers pandas scipy -sortedcontainers gfortran diff --git a/setup.cfg b/setup.cfg index b88034e..a5eb9f6 100644 --- a/setup.cfg +++ b/setup.cfg @@ -1,2 +1,4 @@ [metadata] description-file = README.md +license=MIT +license_files=("LICENSE",) diff --git a/setup.py b/setup.py index 85835ab..2997094 100644 --- a/setup.py +++ b/setup.py @@ -41,14 +41,11 @@ include_package_data=True, ext_modules=[sr1, sr2, sr3, sr_mdl_mult, sr_mdl_plus, sr_mdl4, sr_mdl5], python_requires='>3.6', - install_requires=['matplotlib', - 'numpy', - 'seaborn', - 'sklearn', + install_requires=['numpy', 'sortedcontainers', + 'scikit-learn==0.24.1', 'sympy >= 1.4', 'torch >= 1.4.0', - 'torchvision', ], entry_points={ 'console_scripts': [