Calculating the PCE coefficients with Gaussian quadrature using existing data

hello everyone
My uncertainty quantification problem doesn’t have a definite function relationship between input and output because the output is calculated by CFD based on corresponding input values. Which means Unable to create a “model” as defined in uqlab. So, I can’t define a model and have to use the manual input data method(existing data). I want to input existing data and still use Gaussian quadrature to calculate the PCE coefficients. I found that there is a significant error in the computation, but theoretically, the accuracy of a fourth-order PCE should be sufficient to predict the uncertainty quantification (UQ) problem in our field. How can I get the correct results?

  1. List item

Research problem Description: Uncertainty quantification for a one-dimensional normal distribution N(0,0.16672) input variable.( Just preliminary research, there are high-dimensional problems to be addressed later.")

My input variables: The input variable is one-dimensional and follows a normal distribution N(0,0.16672).
image
image
I selected a PCE degree of 4, so there are 5 sample points according to Gaussian quadrature .
I noticed in the manual:


For projection-based PCE, the experimental design is determined by the choice of quadrature scheme and degree. Since I chose a PCE degree of 4, there are 5 sample points, all of which are specific points.
For my input variable distribution, using Gaussian quadrature, the experimental design is determined by the choice of quadrature scheme and degree.
image
My input data is in Excel, where X is obtained through the quadrature scheme and degree based on the existing distribution, and Y is obtained through CFD for each data point of X.
image
After calculating like this, the error is very large and cannot achieve good results.
So I chose a one-dimensional model function for testing and found the problem

image

by defining the model.

clearvars
rng(100,'twister')
uqlab

%% 2 - Model Computation
% Define MATLAB expression string
modelopts.mString = '10/sqrt(2*pi) * (exp(-(X+2.7).^2/2) + exp(-(X-2.7).^2/8))'; % One-dimensional function

modelopts.isVectorized = true;
myModel = uq_createModel(modelopts);

%% 3 - Define Input Random Variable Model
InputOpts.Marginals(1).Type = 'Gaussian';
InputOpts.Marginals(1).Parameters = [0 0.1667]; 
% Create an input object based on specified marginal distributions
myInput = uq_createInput(InputOpts);

uq_print(myInput)
uq_display(myInput);

%% 4 - Polynomial Chaos Expansion (PCE) Metamodel
% Choose PCE as the metamodel tool in UQLab:
MetaOpts.Type = 'Metamodel';
MetaOpts.MetaType = 'PCE';
MetaOpts.FullModel = myModel;
MetaOpts.Method = 'Quadrature';
MetaOpts.Degree = 4;
%%
myPCE = uq_createModel(MetaOpts);
% Print the computed PCE report:
uq_print(myPCE)
% Visualize the spectrum of nonzero coefficients obtained
uq_display(myPCE)

%% 5-2-1 NIPC Response vs. True Value Plot Validation
figure;

uq_plot(Yval, YQuadrature, 'ro', 'MarkerSize', 10, 'LineWidth', 2, 'MarkerFaceColor', 'none'); % Plot PCE model response values, using red
hold on
uq_plot([min(Yval) max(Yval)], [min(Yval) max(Yval)],'k')
% Set figure properties
title('True Value vs. PCE Model Response');
xlabel('True Value');
ylabel('PCE Model Response');
grid;
hold off;
% Compute validation errors
fprintf('NIPC Method Accuracy Validation:\n');
fprintf('%10s | Relative Error | Experimental Design Samples\n', 'Method');
fprintf('---------------------------------\n');
fprintf('%10s | %10.8e | %7d\n', MetaOpts.Method, mean((Yval - YQuadrature).^2)/var(Yval), myPCE.ExpDesign.NSamples);

%% 5 - NIPC Method Accuracy Validation
% Create a validation sample of size $10^4$ from the input model:
Xval = uq_getSample(20);
% Evaluate the true values at the sample points:
Yval = uq_evalModel(myModel,Xval);
% Generate the response values of the created PCE model at the sample points:
YQuadrature = uq_evalModel(myPCE,Xval);

%% 5-5-3 NIPC Input vs. Predicted Values Plot (Applicable for One Dimension)
% Plot scatter plot of true values vs. predicted values:
 uq_figure
% 
 uq_plot(Xval, YQuadrature, 'b+')
 hold on
 uq_plot(Xval, Yval, 'r+')
 legend('YPCE', 'Yreal');

From the results, it can be seen that only five samples were used, and the constructed error is very small.



However, using the same computational method and providing only five data points without defining a model have big error.

clearvars
uqlab

%% 2 - Obtaining Dataset
% Read data from Excel
filename = 'E:\UQ\UQLab_Rel2.0.0\uq_test37\uqdata.xlsx';
% Read data from Excel file
data = xlsread(filename);
% Assign data to X and Y
X = data(:, 1); % Save the first column of data to X
Y = data(:, end); % Save the last column of data to Y

% Save data to MAT file
save('angle.mat', 'X', 'Y');
%%
% Read experimental design dataset file and store the content in a matrix:
load(fullfile('angle.mat'), 'X', 'Y');
%%
%% 3 - Input Model
%
% Since PCE requires the selection of polynomial bases,
% a marginal distribution of the probability input model needs to be defined.
% Specify the marginal distribution of the probability input model:

InputOpts.Marginals(1).Type = 'Gaussian';
InputOpts.Marginals(1).Parameters = [0 0.1667]; 
%%
% Create an INPUT object based on the specified marginal distribution:
myInput = uq_createInput(InputOpts);

%% 4 -1 Polynomial Chaos Expansion (PCE) Metamodel
%
% Choose PCE as the metamodel tool:

MetaOpts.Type = 'Metamodel';
MetaOpts.MetaType = 'PCE';
%%
% Use the experimental design loaded from the data file:

MetaOpts.ExpDesign.X = X;
MetaOpts.ExpDesign.Y = Y;
MetaOpts.Method = 'quadrature';
MetaOpts.Degree = 4;
myPCE = uq_createModel(MetaOpts);

%%
% Provide validation dataset for validation error calculation:
Xval = uq_getSample(200);
% 
%%
% Print summary information of the generated PCE metamodel:
uq_print(myPCE)
%% 4 -2 PCE Validation
%
% Evaluate the PCE metamodel at validation set points:
YPCE = uq_evalModel(myPCE,Xval);

uq_figure
% 
 uq_plot(Xval, YPCE, 'b+')
 hold on
 uq_plot(X, Y, 'r+')
 legend('YPCE', 'CFD');
hold off



You can see that the error is very large.
Data display in E:\UQ\UQLab_Rel2.0.0\uq_test37\uqdata.xlsx
image
my questions are:

  1. Can I Calculating the PCE coefficients with Gaussian quadrature using existing data?
  2. If yes, how should I proceed?
    Thanks in advance

Hi

It is not possible to compute PCE coefficients with quadrature from existing data. The whole point of quadrature is to define an optimal set of points (on which the computational model is evaluated) and weights, according to the underlying distribution of the input parameters.
If you have given pre-computed data, this data does usually NOT correspond to quadrature points, so this does make sense to try quadrature. You can use OLS or any of the sparse solvers (LARS, SP, etc.)
Best regards
Bruno

Dear professor, thank you for your reply and guidance.
I have another ideas and questions.

  1. Yes, I already noticed that there are not only points but also weights involved in experimental design using the ‘Quadrature’ method. When using existing data and customizing the experimental design, I might wonder if I can input the weights required by quadrature using MetaOpts.ExpDesign.W = W; . However, I attempted it and it seems to have failed.
    image

  2. The sample points I input Do correspond to quadrature points. I understand that for projection-based Polynomial Chaos Expansion (PCE), the experimental design is determined by the choice of quadrature scheme and degree. By inputting a test model and variables distributions similar to those in my research, I obtained corresponding quadrature points. Then, I used these points to compute my Computational Fluid Dynamics (CFD) results and finally established a PCE model using existing data. Therefore, the sample points(existing data) I input in do correspond to quadrature points. you can see my tests in my origin posts.

  3. I have already tested simultaneously inputting corresponding quadrature points and weights manually using MetaOpts.ExpDesign.X = X; MetaOpts.ExpDesign.W = W; the weights are found using the same method as quadrature points as mentioned in 2. However, I still obtained the same results as in test 2.

  4. So, I want to confirm whether even though my existing data has found corresponding quadrature points and weights based on my input distribution, I still cannot compute like this?

Thanks in advance
Best regards !