[PYTHON] Inversely analyze a machine learning model

What is inverse analysis?

--In a broad sense, estimating the input from the output or finding the solution of some equation. ――In a narrow sense, in fields such as material design and chemistry, it is an analysis that first determines the desired physical properties and finds the conditions of the material that realizes them.

In general, finding the physical properties of a synthesized substance from the conditions of the synthesis source is called a ** forward problem **, and in the opposite direction, it is sometimes expressed as ** solving an inverse problem **.

Purpose of this article

** Perform inverse analysis of the machine learning model. ** **

If machine learning makes it possible to predict physical properties, we should be able to search for input values that give the output of the model a predetermined value.

However, as the number of input dimensions (the number of explanatory variables) increases, the space to be searched becomes enormous, and there should be cases where inverse analysis cannot be performed depending on the search time and computer performance.

1. Therefore, first create a prediction model for simple data and confirm that inverse analysis is possible.
2. After that, while increasing the number of explanatory variables, we will investigate how the accuracy of the inverse analysis decreases. (Scheduled to be added later)

basic design

Without doing anything complicated, we try to ** search for the input value that minimizes the output for the regression model **.

Random forest is used as a regression model for the time being.

SMBO (Sequential Model-based Global Optimization) is used as the search algorithm, and hyperopt is used as the library. (There are other various methods).

1. Inverse analysis for toy model

environment

• Python 3.6.10
• scikit-learn 0.22.0
• hyperopt 0.2.4

Setting

As a simple model

y= x_1 {}^2 + x_2 {}^2, \qquad (x_1, x_2) \in \mathbb{R} ^2

I think about the correspondence. Obviously the minimum value is $ 0 $, and the point that gives this is $ (x_1, x_2) = (0,0) $.

import matplotlib.pyplot as plt
import numpy as np
from mpl_toolkits.mplot3d import Axes3D

def true_function(x, y):
    """True function"""
    return x ** 2 + y ** 2

X, Y = np.mgrid[-100:100, -100:100]
Z = true_function(X, Y)

plt.rcParams["font.size"] = 10  #Increase the font size
fig = plt.figure(figsize = (12, 9))
ax = fig.add_subplot(111, projection="3d", facecolor="w")
ax.plot_surface(X, Y, Z, cmap="rainbow", rstride=3, cstride=3)
ax.set_xlabel('x1', fontsize=15)
ax.set_ylabel('x2', fontsize=15)  
ax.set_zlabel('y', fontsize=15) 
plt.show()

1.1 Generation of training test data

Based on the above correspondence, input / output sample groups are generated.

Draw the generated training data.

from sklearn.model_selection import train_test_split

def true_model(X):
    return true_function(X[:,[0]], X[:,[1]])

X = np.random.uniform(low=-100,high=100,size=(3000,2))
Y = true_model(X)

test_size = 0.3  #Split ratio

x_train, x_test, y_train, y_test = train_test_split(X, Y, test_size=test_size, random_state=0)
  
fig = plt.figure(figsize = (12, 9)) 
ax = plt.axes(projection ="3d")
sctt = ax.scatter3D(x_train[:,0], x_train[:,1], y_train[:,0], c=y_train[:,0], s=8, alpha = 0.6,
                    cmap = plt.get_cmap('rainbow'), marker ='^') 

plt.title("x_train, y_train") 
ax.set_xlabel('x1', fontsize=15)
ax.set_ylabel('x2', fontsize=15)  
ax.set_zlabel('y', fontsize=15) 
plt.show() 

1.2 Learning / Reasoning

Random forest is trained with the training data mentioned earlier, and the result of inference for the test data is drawn.

It seems that the correct value can be estimated.

1.3 Searching for the minimum value for the regression model

Try to search for the minimum value by hyperopt.

After defining the function to be minimized, search for the point that gives the minimum value, and draw the obtained point on top of the previous figure.

from hyperopt import hp
from hyperopt import fmin
from hyperopt import tpe

def objective_hyperopt_by_reg(args):
    """Objective function for hyperopt"""
    global reg
    x, y = args
    return float(reg.predict([[x,y]]))

def hyperopt_exe():
    """Performing optimization with hyperopt"""
    #Search space settings
    space = [
        hp.uniform('x', -100, 100),
        hp.uniform('y', -100, 100)
    ]

    #Start exploration
    best = fmin(objective_hyperopt_by_reg, space, algo=tpe.suggest, max_evals=500)
    return best

best = hyperopt_exe()
print(f"best: {best}")

fig = plt.figure(figsize = (12, 9)) 
ax = plt.axes(projection ="3d")
sctt = ax.scatter3D(x_test[:,0], x_test[:,1], y_test[:,0], c=y_test[:,0], s=6, alpha = 0.5,
                    cmap = plt.get_cmap('rainbow'), marker ='^')
ax.scatter3D([best["x"]], [best["y"]], [objective_hyperopt_by_reg((best["x"], best["y"]))], 
                    c="red", s=250, marker="*", label="minimum") 

plt.title("x_test, y_pred", fontsize=18) 
ax.set_xlabel('x1', fontsize=15)
ax.set_ylabel('x2', fontsize=15)  
ax.set_zlabel('y', fontsize=15) 
plt.legend(fontsize=15)
plt.show() 

output

100%|██████████████████████████████████████████████| 500/500 [00:09<00:00, 52.54trial/s, best loss: 27.169204190118908]
best: {'x': -0.6924078319870626, 'y': -1.1731945130395605}

A point close to the minimum point was obtained.

Summary and challenges

In this article, we performed the procedure of regression model learning ⇒ inverse analysis on simple data.

Since the input dimension was small this time, we were able to search for the minimum value by chance, but it is expected that various problems will occur when applying it to actual data.

-** There is too little data to learn . - Because you have selected an inappropriate regression model **, you will derive a solution that is impossible in reality. --The dimension of the explanatory variable is too high to search well, or it is trapped in the polar solution. --The distribution of data is biased, and learning in sparse areas is insufficient. ――Excessive noise is learned.

Also, ** inverse problem is an extremely nonsense analysis if not handled properly. It is expected that there is a risk of doing **.

--There is no solution in the first place --The solution is not the only one --The true distribution is discontinuous and the solution is not stable

I definitely want to avoid wasting effort because I have set the inverse problem in spite of such a problem.

reference

-New inductive approach in materials science -Inverse problem analysis using optimization -What is hyperopt? -Try function optimization using Hyperopt