<jupyter_start><jupyter_text>How to use CGA2M+ package!<jupyter_code>#import CGA2M+
from cga2m_plus.cga2m import Constraint_GA2M
from cga2m_plus.visualize import *

#data processing
import numpy as np
from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_squared_error

#Other
import itertools


np.random.seed(2021)<jupyter_output><empty_output><jupyter_text>Make sample dataset(Example data set)$$\displaystyle y = (1+0.2\epsilon_1)x_0 + \exp(-\frac{1}{3}x_1^2)+ \sin(2\pi x_2) - \frac{1}{3}x_2x_4  +\frac{1}{3}\exp(\frac{1}{8}(x_3 + x_4 - x_5)^2) + \epsilon_2$$$$\epsilon_1,\epsilon_2 \sim N(0,1)$$<jupyter_code>def sample_data(x):
    y =  (1+0.2*np.random.normal(size=1).item())*x[0] +  np.exp(- 1/3 * x[1]**2) +  np.sin(2* np.pi * x[2]) - 1/3 * x[2]*x[4] + 1/3 * np.exp(1/8*(x[3] + x[4] - x[5])**2) 
    y = y + np.random.normal(size=1).item()#noise
    return y
#n_samples： sample size
#n_features: the number of features

num_features = 7 #meaningful features x_0 ~x_5, meaningless feature x_6
n_samples = 30000
data = 4 * np.random.random_sample((n_samples,num_features)) - 2
target = list()
for i in range(data.shape[0]):
    y = sample_data(data[i])
    target.append(y)

target = np.array(target)
print(data.shape)
print(target.shape)


#Split the data into 3 sets(train, validation and test)
X_train, X_test, y_train, y_test = train_test_split(data,target,train_size=0.7,random_state=2021)
X_eval, X_test, y_eval, y_test = train_test_split(X_test,y_test,train_size=0.5,random_state=2021)#random_stateは変更しない

#Check the sample size of 3 sets(train, validation and test)
print('X_train',X_train.shape)
print('y_train',y_train.shape)

print('X_eval',X_eval.shape)
print('y_eval',y_eval.shape)

print('X_test',X_test.shape)
print('y_test',y_test.shape)
print(np.mean(y_train))
print(np.min(y_train))
print(np.max(y_train))<jupyter_output>1.443274768558458
-5.007631785246939
22.337004822595667<jupyter_text>Training CGA2M+ and visualization In this example, the number of backfitting iterations (=20) may be excessive.<jupyter_code># parameters of LightGBM in CGA2M+
lgbm_params = {'objective':'regression',
               'metric':'rmse',
               'verbose': -1}

#Input the data into the model.
#all_interaction_features is List which contains the interactions to be taken into account. Note that the interactions are set to (i<j) for uniqueness.
cga2m_no1 = Constraint_GA2M(X_train,
                            y_train,
                            X_eval,
                            y_eval,
                            lgbm_params,
                            monotone_constraints=None,
                            all_interaction_features=list(itertools.combinations(range(X_test.shape[1]), 2)))
# Train (about 5min)

# max_outer_iteration is the number of iterations in GA2M Framework(==Maximum number of interactions)
## the interactions refer to [(1,2),(1,3),...]
max_outer_iteration = len(list(itertools.combinations(range(X_test.shape[1]), 2)))
print("max_outer_iteration=",max_outer_iteration)

# backfitting_iteration ： The number of iterations of backfitting
# threshold ： If the importance of the variable pair is below this threshold, we finish training.
cga2m_no1.train(max_outer_iteration=max_outer_iteration,backfitting_iteration=20,threshold=0.05)

#Delete ineffective shape functions(e.g. f_i,f_{i,j}) and train again.
## threshold ： If the importance of the variable is below threshold, we consider that it is ineffective.
print('START prune_and_retrain')
cga2m_no1.prune_and_retrain(threshold=0.05,backfitting_iteration=30)

# Training higher-order terms (if you do not need this, comment out!)
cga2m_no1.higher_order_train()
# Pridict
## higher_mode ： True if you want to use the result of the higher order term (False by default).
preds_higher = cga2m_no1.predict(X_test,higher_mode=True)
preds = cga2m_no1.predict(X_test)
# Calculate RMSE for test data
mse_higher = mean_squared_error(y_test, preds_higher)
mse = mean_squared_error(y_test, preds)

rmse_higher = np.sqrt(mse_higher)
rmse = np.sqrt(mse)

print('RMSE(test) with the higher-order term',rmse_higher)
print('RMSE(test) without the higher-order term',rmse)


preds_higher = cga2m_no1.predict(X_train,higher_mode=True)
preds = cga2m_no1.predict(X_train)
# Calculate RMSE for training data
mse_higher = mean_squared_error(y_train, preds_higher)
mse = mean_squared_error(y_train, preds)

rmse_higher = np.sqrt(mse_higher)
rmse = np.sqrt(mse)

print('RMSE(train) with the higher-order term',rmse_higher)
print('RMSE(train) without the higher-order term',rmse)
# Visualize the effect of features on the target variables
plot_main(cga2m_no1,X_train)
#Visualize (3d) the effect of pairs of features on the target variables
plot_interaction(cga2m_no1,X_train,mode = "3d")
#Visualize (2d) the effect of pairs of features on the target variables
plot_interaction(cga2m_no1,X_train,mode = "2d")
# Feature Importance
## after_prune： If True, display the feature importance after pruning.
## higher_mode: If True, display the feature importance which contain the effect of the higher-order term.
show_importance(cga2m_no1,after_prune=True,higher_mode=True)
# Selected features
print("Selected features")
print(cga2m_no1.use_main_features,"\n")

# Pairs of selected features (interaction)
print("Pairs of selected features (interaction)")
print(cga2m_no1.use_interaction_features,"\n")

# Feature Importance after pruning
print("Feature Importance after pruning")
print(cga2m_no1.after_feature_importance_,"\n")

# Feature Importance before pruning
print("Feature Importance before pruning")
print(cga2m_no1.before_feature_importance_)<jupyter_output>Selected features
[0, 1, 2, 3, 4, 5, 6] 

Pairs of selected features (interaction)
[(3, 5), (3, 4), (4, 5), (2, 4), (2, 3), (1, 3)] 

Feature Importance after pruning
{0: 0.633392043866198, 1: 0.1375802338586143, 2: 0.40978674472942883, 3: 0.12813562977042145, 4: 0.1307073499774606, 5: 0.12063237512840794, 6: 0.056437720552546415, (3, 5): 0.19483815505085672, (3, 4): 0.20685491267673345, (4, 5): 0.19268456099874567, (2, 4): 0.22713819687295386, (2, 3): 0.054471658396229986, (1, 3): 0.04786235992822891, 'higher': 0.06875713697261573} 

Feature Importance before pruning
{0: 0.6295763014463029, 1: 0.1374524920845846, 2: 0.410635617796357, 3: 0.13027530481618385, 4: 0.13592643124116868, 5: 0.12444194050845682, 6: 0.05621094219527024, (3, 5): 0.19441289105214876, (3, 4): 0.2093836690781018, (4, 5): 0.1916523847390596, (2, 4): 0.22907313106379773, (2, 3): 0.051719230408024223, (1, 3): 0.05301057504551453, (0, 3): 0.04755690895513956, (0, 4): 0.048329654535105245}<jupyter_text>Apply monotonic constraints to some shape functions of the modelYou have a lot of knowledge about the data. (For example, in real estate prices, you know that the larger the room, the higher the rent.)In this case, you may impose monotonic constrains on some features.In this example, we apply a monotonic constraint to x_0.   How to set monotonic constraintsThe monotonic constraints are set as `monotone_constraints=[1,0,-1]`.Monotonically increasing if `1`, monotonically decreasing if `-1`, unconstrained(No setting) if `0`.Assign 1, -1, or 0 to all indices for each feature in training data.  Ex)Four features. The second feature is monotonically increasing. The fourth feature is monotonically decreasing.  `monotone_constraints=[0,1,0,-1]`<jupyter_code># Training model with constraint

cga2m_no2 = Constraint_GA2M(X_train,
                            y_train,
                            X_eval,
                            y_eval,
                            lgbm_params,
                            monotone_constraints=[1,0,0,0,0,0,0],# monotonic constraints
                            all_interaction_features=list(itertools.combinations(range(X_test.shape[1]), 2)))

print(cga2m_no2.use_main_features)
print(cga2m_no2.remaining_interaction_features)

cga2m_no2.train(max_outer_iteration=max_outer_iteration,backfitting_iteration=20,threshold=0.05)

print('START prune_and_retrain')
cga2m_no2.prune_and_retrain(threshold=0.05,backfitting_iteration=30)

cga2m_no2.higher_order_train()
#Pridict
## higher_mode ： True if you want to use the result of the higher order term (False by default).

preds_higher = cga2m_no2.predict(X_test,higher_mode=True)
preds = cga2m_no2.predict(X_test)
# Calculate RMSE for test data
mse_higher = mean_squared_error(y_test, preds_higher)
mse = mean_squared_error(y_test, preds)

rmse_higher = np.sqrt(mse_higher)
rmse = np.sqrt(mse)

print('RMSE(test) with the higher-order term',rmse_higher)
print('RMSE(test) without the higher-order term',rmse)


preds_higher = cga2m_no2.predict(X_train,higher_mode=True)
preds = cga2m_no2.predict(X_train)
# Calculate RMSE for training data
mse_higher = mean_squared_error(y_train, preds_higher)
mse = mean_squared_error(y_train, preds)

rmse_higher = np.sqrt(mse_higher)
rmse = np.sqrt(mse)

print('RMSE(train) with the higher-order term',rmse_higher)
print('RMSE(train) without the higher-order term',rmse)
# Visualize the effect of features on the target variables
plot_main(cga2m_no2,X_train)
#Visualize the effect of pairs of features on the target variables
plot_interaction(cga2m_no2,X_train)<jupyter_output><empty_output><jupyter_text>Delete some features or pairs of featuresWe will show you how to remove a feature (or a pair of features) from the model you have just trained, and how to retrain the model.There are two ways to do this.1. From the features and pairs of features used in the previously trained model, delete the features and pairs which are meaningless. Then, you train the model with the features and pairs of features determined.2. Decide which features and feature pairs are not used, and train the model from the beginning.  1. From the features and pairs of features used in the previously trained model, delete the features and pairs which are meaningless. Then, you train the model with the features and pairs of features determined.You have examined your data and your model(```cga2m_no1```) and determined that (x_1,x_3) is not needed among the features and feature pairs.  In this example, you remove (x_1,x_3) from your model(```cga2m_no1```) and retrain the model.<jupyter_code># Remove (x_1,x_3) from the model(cga2m_no1).
cga2m_no1.use_interaction_features.remove((1,3))

# if you remove x_5 from the model, type "cga2m.use_main_features.remove(5)".

print("use_main_features",cga2m_no1.use_main_features)
print("use_interaction_features",cga2m_no1.use_interaction_features)
# Training again.
cga2m_no1.retrain(backfitting_iteration=30)
# Train higher-order term
cga2m_no1.higher_order_train()<jupyter_output><empty_output><jupyter_text>2. Decide which features and feature pairs are not used, and train the model from the beginning.<jupyter_code>use_interaction_features = list(itertools.combinations(range(X_test.shape[1]), 2))
use_interaction_features.remove((1,3))# remove (x_1,x_3)

cga2m_no1 = Constraint_GA2M(X_train,
                            y_train,
                            X_eval,
                            y_eval,
                            lgbm_params,
                            monotone_constraints=None,
                            all_interaction_features=use_interaction_features)

cga2m_no1.use_main_features.remove(5)# remove x_5
max_outer_iteration = len(use_interaction_features)
print("max_outer_iteration=",max_outer_iteration)
cga2m_no1.train(max_outer_iteration=max_outer_iteration,backfitting_iteration=20,threshold=0.05)
cga2m_no1.prune_and_retrain(threshold=0.05,backfitting_iteration=30)
cga2m_no1.higher_order_train()<jupyter_output>max_outer_iteration= 20
START 1ST ITERATION
START 2ST ITERATION
START 3ST ITERATION
START 4ST ITERATION
START 5ST ITERATION
START 6ST ITERATION
START 7ST ITERATION
START 8ST ITERATION
START 9ST ITERATION
START 10ST ITERATION
End of training