<jupyter_start><jupyter_text>Polynomial Chaos Expansion Example 2 Author: <NAME> \Date: December 8, 2020 In this example, PCE is used to generate a surrogate model for a given set of 2D data.  Sphere function $$ f(\textbf{x}) = x_1^2 + x_2^2$$**Description:**  Dimensions: 2**Input Domain:**  This function is evaluated on the hypercube $x_i \in [-5.12, 5.12]$ for all i = 1,2.**Global minimum:** $f(x^*)=0,$ at $x^* = (0,0)$.**Reference:**  <NAME>., & <NAME>. (1978). The global optimization problem: an introduction. Towards global optimization, 2, 1-15. Import necessary libraries.<jupyter_code>import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm
from matplotlib.ticker import LinearLocator, FormatStrFormatter
from UQpy.Surrogates import *
from UQpy.Distributions import Uniform, JointInd<jupyter_output><empty_output><jupyter_text>Define the function.<jupyter_code>def function(x,y):
    return x**2 + y**2<jupyter_output><empty_output><jupyter_text>Create a distribution object, generate samples and evaluate the function at the samples.<jupyter_code>np.random.seed(1)

dist_1 = Uniform(loc=-5.12, scale=10.24)   
dist_2 = Uniform(loc=-5.12, scale=10.24)    

marg = [dist_1, dist_2]
joint = JointInd(marginals=marg)

n_samples = 100
x = joint.rvs(n_samples)
y = function(x[:,0], x[:,1])<jupyter_output><empty_output><jupyter_text>Visualize the 2D function.<jupyter_code>xmin, xmax = -6,6
ymin, ymax = -6,6
X1 = np.linspace(xmin, xmax, 50)
X2 = np.linspace(ymin, ymax, 50)
X1_, X2_ = np.meshgrid(X1, X2) # grid of points
f = function(X1_, X2_)

fig = plt.figure(figsize=(10,6))
ax = fig.gca(projection='3d')
surf = ax.plot_surface(X1_, X2_, f, rstride=1, cstride=1, cmap='gnuplot2', linewidth=0, antialiased=False)
ax.set_title('True function')
ax.set_xlabel('$x_1$', fontsize=15)
ax.set_ylabel('$x_2$', fontsize=15)
ax.zaxis.set_major_locator(LinearLocator(10))
ax.zaxis.set_major_formatter(FormatStrFormatter('%.02f'))
ax.view_init(20, 140)
fig.colorbar(surf, shrink=0.5, aspect=7)

plt.show()<jupyter_output><empty_output><jupyter_text>Visualize training data.<jupyter_code>fig = plt.figure(figsize=(10,6))
ax = fig.gca(projection='3d')
ax.scatter(x[:,0], x[:,1], y, s=20, c='r')

ax.set_title('Training data')
ax.zaxis.set_major_locator(LinearLocator(10))
ax.zaxis.set_major_formatter(FormatStrFormatter('%.02f'))
ax.view_init(20,140)
ax.set_xlabel('$x_1$', fontsize=15)
ax.set_ylabel('$x_2$', fontsize=15)
#ax.set_xlim(-10,10)
#ax.set_ylim(-6,6)
#ax.set_zlim(-1,1.5)
plt.show()<jupyter_output><empty_output><jupyter_text>Create an object from the PCE class.<jupyter_code>max_degree = 3
polys = Polynomials(dist_object=joint, degree=max_degree)<jupyter_output><empty_output><jupyter_text>Compute PCE coefficients using least squares regression.<jupyter_code>lstsq = PolyChaosLstsq(poly_object=polys)
pce = PCE(method=lstsq) 
pce.fit(x,y)<jupyter_output><empty_output><jupyter_text>Compute PCE coefficients using LASSO.<jupyter_code>lasso = PolyChaosLasso(poly_object=polys, learning_rate=0.1, iterations=1000, penalty=0.5)
pce2 = PCE(method=lasso) 
pce2.fit(x,y)<jupyter_output><empty_output><jupyter_text>Compute PCE coefficients with Ridge regression.<jupyter_code>ridge = PolyChaosRidge(poly_object=polys, learning_rate=0.01, iterations=1000, penalty=0.5)
pce3 = PCE(method=ridge) 
pce3.fit(x,y)<jupyter_output><empty_output><jupyter_text>PCE surrogate is used to predict the behavior of the function at new samples.<jupyter_code>n_test_samples = 10000
x_test = joint.rvs(n_test_samples)
y_test = pce.predict(x_test)<jupyter_output><empty_output><jupyter_text>Plot PCE prediction.<jupyter_code>fig = plt.figure(figsize=(10,6))
ax = fig.gca(projection='3d')
ax.scatter(x_test[:,0], x_test[:,1], y_test, s=1)

ax.set_title('PCE predictor')
ax.zaxis.set_major_locator(LinearLocator(10))
ax.zaxis.set_major_formatter(FormatStrFormatter('%.02f'))
ax.view_init(20,140)
ax.set_xlim(-6,6)
ax.set_ylim(-6,6)
ax.set_xlabel('$x_1$', fontsize=15)
ax.set_ylabel('$x_2$', fontsize=15)
#ax.set_zlim(0,136)
plt.show()<jupyter_output><empty_output><jupyter_text>Error Estimation Validation error.<jupyter_code>n_samples = 150
x_val = joint.rvs(n_samples)
y_val = function(x_val[:,0], x_val[:,1])

error = ErrorEstimation(surr_object=pce)
error2 = ErrorEstimation(surr_object=pce2)
error3 = ErrorEstimation(surr_object=pce3)

print('Error from least squares regression is: ', error.validation(x_val, y_val))
print('Error from LASSO regression is: ', error2.validation(x_val, y_val))
print('Error from Ridge regression is: ', error3.validation(x_val, y_val))<jupyter_output>Error from least squares regression is:  0.0
Error from LASSO regression is:  1e-07
Error from Ridge regression is:  0.0001055<jupyter_text>Moment Estimation Returns mean and variance of the PCE surrogate.<jupyter_code>n_mc = 1000000
x_mc = joint.rvs(n_mc)
y_mc = function(x_mc[:,0], x_mc[:,1])    
mu = np.mean(y_mc)

print('Moments from least squares regression :', MomentEstimation(surr_object=pce).get())
print('Moments from LASSO regression :', MomentEstimation(surr_object=pce2).get())
print('Moments from Ridge regression :', MomentEstimation(surr_object=pce3).get())
print('Moments from Monte Carlo integration: ', (round((1/n_mc)*np.sum(y_mc),6), round((1/n_mc)*np.sum((y_mc-mu)**2),6)))<jupyter_output>Moments from least squares regression : (17.4763, 122.168)
Moments from LASSO regression : (17.4756, 122.0909)
Moments from Ridge regression : (17.4582, 120.4658)
Moments from Monte Carlo integration:  (17.479458, 122.282943)
{'fraction_non_alphanumeric': 0.14793536804308796, 'fraction_numerical': 0.05170556552962298, 'mean_word_length': 7.6105100463678514, 'pattern_counts': {'":': 0, '<': 54, '<?xml version=': 0, '>': 54, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '13', 'dir': 'jupyter-structured-clean-dedup', 'id': '591574', 'license': 'MIT', 'n_tokens_mistral': 2358, 'n_tokens_neox': 2185, 'n_words': 417, 'path': 'example/Surrogates/PCE/PCE_Example2.ipynb', 'repo_name': 'shellshocked2003/UQpy', 'size': '5594'}