In [1]:
import numpy as np
import escape as esc
esc.require("0.9.7")
from escape.utils.widgets import show

Loading material database from /home/dkor/Data/Development/workspace_escape/escape/python/src/escape/scattering/../data/mdb/materials.db

Multi-layer Perceptron regressor¶

In this notebook we demonstrate how to use Multilayer Perceptor Regressor for batch processing of experimental data. If you are not familiar with neural networks and particularly with multi-layer regression, please find a brief overview with some details here

https://scikit-learn.org/stable/modules/neural_networks_supervised.html

There are several very popular frameworks for deep learning, like PyTorch, Tensorflow, Keras, etc... it will be difficult and unnecessary to compete with them in terms of flexibility and optimizations speed. That's why instead of a full implementation of MLP regressor, we provide a wrapper. This wrapper is responsible for preparing of training and testing data in terms of ESCAPE framework and its parameterized objects. Prepared training and testing data is provided for further analysis, optimization and prediction to a wrapped regressor instance, which is a python object. Below we give an example how to use this wrapper with sklearn.MLPRegressor class.

As an interface basis for the wrapped regressor object we took sklearn MLPRegressor, which means, that if you like to use your own regressor, it should have the following methods:

  1. fit(X, y) for optimization
  2. score(X, y) for score calculation
  3. predict(X) for prediction

Additionally there are methods which are responsible for noisifying and normalizing of the training data, noisify(X) and normalize(X). These methods can be a part of regressor instance or can be provided by user when creating a regressor wrapper. See example below.

Tip: Your generated training and experimental data which you are using for prediction of parameters are supposed to have the same coordinates. If this is not the case, you, probably, should use interpolation.

In the example below we use the previous model for Ni/Si/Ni trilayer on a Si substrate. Regressor as an optimizer requires a model or modelstack object, thus, we need a preliminary experimental data, which is generated. The coordinates array of the preliminary experimental data will be used by regressor wrapper to generate a training data and test data.

After optimization, the regressor wrapper can be applied to real experimental data for prediction of parameters and further optimization if required.

As usual, we start with the parameters and model description.

In [2]:
# Mass densities of Ni and Si

thknNi1 = esc.par("Ni1 Thkn", 20, units="nm", userlim=[15, 25])
roughNi1 = esc.par("Ni1 Roughness", 1.5, units="nm", userlim=[1.0, 2.0])
thknNi2 = esc.par("Ni2 Thkn", 20, units="nm", userlim=[15, 25])
roughNi2 = esc.par("Ni2 Roughness", 1.5, units="nm", userlim=[1.0, 2.0])
thknSi = esc.par("Si Thkn", 8, units="nm", userlim=[5, 10])
roughSi = esc.par("Si Roughness", 0.2, units="nm", userlim=[0, 0.5])

#we also add a background parameter
B = esc.par("Background", 0, units="", scale=1e-6, userlim=[0, 10])
In [3]:
LayNi1 = esc.layer("Layer: Ni1", material="Ni", thkn="20+-5nm", rough="1.5+-0.5nm", bydensity=True)
LayNi2 = esc.layer("Layer: Ni2", material="Ni", thkn="20+-5nm", rough="1.5+-0.5nm", bydensity=True)
LaySi = esc.layer("Layer: Si", material="Si", thkn="8+-3nm", rough="0.2+-0.2nm", bydensity=True)
Sub = esc.substrate("Substrate: Si", material="Si", rough=0.0, bydensity=True) # roughness is zero

sample = esc.multilayer("Sample", formula="LayNi1/LaySi/LayNi2//Sub", globals=globals())

Now we create calculation kernel for specular reflectivity

In [4]:
#we generate specrefl functor and add background to it
Qz=esc.var("qz")
src =esc.xrays(wavelength=0.154, units="nm")
Rf = esc.specrefl("Specrefl", Qz, sample, "matrix", source=src)+B
R=esc.kernel("Kernel", Rf, True)

Now we can generate experimental data. We also add some poisson noise.

In [5]:
qz=np.linspace(0.001, 2.5, 1000).copy()
#scaling coefficient for intensities. Smaller values will produce more noisy data 
I0 = 1e8

y=np.empty(qz.shape)
R(qz, y)

y = np.random.poisson(y*I0)/I0
err = np.sqrt(y)/np.sqrt(I0)

dobj = esc.data("Ni/Si/Ni", qz, y, err, copy=True)
mobj = esc.model("Model: Ni/Si/Ni", R, dobj, residuals_scale="none", weight_type="data")

show(mobj, xlabel="Q[nm⁻¹]", ylabel="I/I0", ylog=True)

Below we create a regressor wrapper. As an actual regressor instance we are going to use MLPRegressor from sklearn library. Our generated training data will be noisified and logarithmically normalized. The simulated specular reflectivity curve is usually normalized, thus, scale parameter is required to get a realsitic poisson noise.

We are ready to create regressor instance and its wrapper.

In [7]:
from sklearn.neural_network import MLPRegressor 
rgi=MLPRegressor(max_iter=1000, activation="relu", alpha=1e-5, verbose=False)


def noisify(X):
    X = np.asarray(X)
    X = np.random.poisson(X*I0)
    X = X/np.max(I0)
    return X

def normalize(X):
    X = np.asarray(X)
    X = X*qz**4
    X = X/np.max(X)
    return X

# 25000 training and 150 testing sets

rg = esc.regressor("", [mobj], 30000, 150, impl=rgi, norm_method=normalize,
                   noisify_method=noisify
                  )
#the widget is used for training and makes a prediction of parameters  only for the current experimental data
# see below how to use the resulted neural network for other experimental data
#rg.train(True)
show(rg, ylog=True)

We have trained our regressor. Let's check now how prediction of parameters works for some user provided data. You can run the next cell several times to see how good prediction works for the data out of training set.

In [33]:
#Since we do not have a real experimental data, we are going to generate
#the data object for the randomly changed model parameters
mobj.shake()
costb = mobj(True)
#create data object with simulated array as intensities
y = np.asarray(mobj.simulation)
y = np.random.poisson(y*I0)/I0
err = np.sqrt(y)/np.sqrt(I0)
dobj = esc.data("Ni/Si/Ni Generated", qz, y, err, copy=True)
#run prediction
rg([dobj])
costa = mobj(True)
show(mobj, ylog=True)
In [ ]:
In [ ]: