After creating myBayesianAnalysis, we checked the maximum and minimum values of myBayesianAnalysis.Results.Sample which gave ~ [0.9,0]. The problem is that these values are outside the input range, is this normal?
We also noticed that in myBayesianAnalysis.PriorDist or myBayesianAnalysis.Internal.FullPrior the bounds do not appear anymore. It is not in fact included in the uq_remove_constants_from_input which writes only ‘Name’, ‘Parameters’, ‘Type’.
Consider the simply-supported beam example in the documentation of the Bayesian module :
% Initialisation
clc % Nettoye la console
clear all % Supprime les variables actuelles
uqlab % Demarre l'environnement UQLAB
%% Model definition
addpath('C:\Program Files\UQLab\UQLab_Rel2.0.0\Examples\SimpleTestFunctions')
ModelOpts.Name = 'Forward model';
ModelOpts.mFile = 'uq_SimplySupportedBeam';
myModel = uq_createModel(ModelOpts);
clear ModelOpts
%% Input definition
PriorOpts.Name = 'Model parameters prior';
PriorOpts.Marginals(1).Name = 'b';
PriorOpts.Marginals(1).Type = 'Constant';
PriorOpts.Marginals(1).Parameters = 0.15; % (m)
PriorOpts.Marginals(2).Name = 'h';
PriorOpts.Marginals(2).Type = 'Constant';
PriorOpts.Marginals(2).Parameters = 0.3; % (m)
PriorOpts.Marginals(3).Name = 'L';
PriorOpts.Marginals(3).Type = 'Constant';
PriorOpts.Marginals(3).Parameters = 5; % (m)
PriorOpts.Marginals(4).Name = 'E';
PriorOpts.Marginals(4).Type = 'Lognormal';
PriorOpts.Marginals(4).Moments = [30 4.5]*1e9 ; % (N/mˆ2)
PriorOpts.Marginals(4).Bounds = [21,30]*1e9;
PriorOpts.Marginals(5).Name = 'p';
PriorOpts.Marginals(5).Type = 'Constant';
PriorOpts.Marginals(5).Parameters = 12000; % (N/m)
myPriorInput = uq_createInput(PriorOpts);
clear PriorOpts
%% Measurments
V_mid = [12.84; 13.12; 12.13; 12.19; 12.67]/1000; % (m)
myData.Name = 'Beam mid-span deflection';
myData.y = V_mid;
%% Inversion
BayesOpts.Type = 'Inversion';
BayesOpts.Data = myData;
clear myData
BayesOpts.Prior = myPriorInput;
BayesOpts.Solver.MCMC.Steps = 1000;
BayesOpts.Solver.MCMC.NChains = 100;
BayesOpts.Solver.MCMC.Sampler = 'AM';
BayesOpts.Discrepancy.Type = 'Gaussian';
BayesOpts.Discrepancy.Parameters = 1e-7;
myBayesian = uq_createAnalysis(BayesOpts);
uq_postProcessInversion(myBayesian,'gelmanRubin',true)
clear BayesOpts
%% Check the bounds of the sampling
% Samples
Sampling_E_Sample=myBayesian.Results.Sample;
Sampling_E_Sample=reshape(Sampling_E_Sample,1,[]);
Sampling_Bayes_Sample=[max(Sampling_E_Sample)/10^9,min(Sampling_E_Sample)/10^9]; % take the max and min value from Bayes' Sample
% Prior samples
Sampling_E_Prior=myBayesian.Results.PostProc.PriorSample;
Sampling_E_Prior=reshape(Sampling_E_Prior,1,[]);
Sampling_Bayes_Prior=[max(Sampling_E_Prior)/10^9,min(Sampling_E_Prior)/10^9]; % take the max and min value from Bayes' sampling
% Input distribution samples
Sampling_test=uq_getSample(myPriorInput,100000);
Sampling_test=[max(Sampling_test(:,4))/10^9,min(Sampling_test(:,4))/10^9]; % take the max and min value from Input sampling
As you will notice, the maximum and minimum value of the two sets “myBayesian.Results.PostProc.PriorSample” and “myBayesian.Results.Sample” are not contained within the bounds defined in the input, whereas those from uq_getSample() do. On my machine I got :
Sampling_test = [29.9999 ; 20.002] %OK with the Input
Sampling_Bayes_Prior = [49.5657 ; 17.9774] %Outside the bounds
Sampling_Bayes_Sample= [40.9093 ; 20.65] %Outside the bounds
I must admit I have not been in depth with the theory behind the MCMC algorithm but the fact that values are outside the bounds in the prior’s sampling seems odd to me.
Sorry for the late reply - I checked your input script and it looks fine to me. Unfortunately, I cannot debug it at the moment, but I have a suspect .
First of all, is there a small typo in the bounds you report for Sampling_test? The lower bound is also outside the set lower bound of 21.
Secondly, I suspect the problem is inside the uq_remove_constants_from_input function that unintentionally removes the bounds of the original prior. I couldn’t test the following, so please backup the file before proceeding:
Add the following lines after line 81 and after line 102 of the existing file:
if isfield(Marginals(rv), 'Bounds')
newOpts.Marginals(nrMargs+kk).Bounds= Marginals(rv).Bounds;
end
After doing this, add a debug point in line 370 of uq_initialize_uq_inversion and verify that the bounds you specified are still present in the generated ModelPrior_noConstants object.
I tested your solution in uq_remove_constants_from_input, it works well. However, be careful because this seems to have an impact on other modules like the Polynomial Choas Kriging module. Unfortunately I couldn’t simulate this effect on the minimal example I gave you, so I’ll give you the original code:
If everything is OK, you should get an LOO of ~6.0e-06, which is the right value.
Now let’s change these two input variables by adding boundaries (adding these boundaries to a Lognormal distribution should not change anything right ?) :
I don’t think that these issues are related. In fact, (I think that) PCE in UQLab numerically constructs the polynomial basis functions if bounds are explicitly given by the user - independent of whether they match the theoretical bounds of the marginals or not. The observed decrease in accuracy can probably be attributed to this.
.