From 1ba7867de523fb1f0a25aaf9efc762613a77714b Mon Sep 17 00:00:00 2001
From: faridm69 <faridmohammadi69@gmail.com>
Date: Sat, 14 Nov 2020 18:09:14 +0100
Subject: [PATCH] [MCMC] added KDEMove as the default.
---
BayesValidRox/BayesInference/BayesInference.py | 15 +++++++++++----
BayesValidRox/BayesInference/MCMC.py | 14 ++++++--------
2 files changed, 17 insertions(+), 12 deletions(-)
diff --git a/BayesValidRox/BayesInference/BayesInference.py b/BayesValidRox/BayesInference/BayesInference.py
index 801cd3523..1b291e03f 100644
--- a/BayesValidRox/BayesInference/BayesInference.py
+++ b/BayesValidRox/BayesInference/BayesInference.py
@@ -59,7 +59,8 @@ class BayesInference:
self.perturbedData = []
self.MeasuredData = None
self.MeasurementError = []
-
+ self.selectedIndices = None
+
self.ModelPriorPred = []
self.ModelPosteriorPred = []
self.PCEAllPred_df = []
@@ -183,12 +184,18 @@ class BayesInference:
np.fill_diagonal(covMatrix_PCE, varPCE)
covMatrix += covMatrix_PCE
+ # Select the data points to compare
+ if self.selectedIndices is not None:
+ TotalOutputs = TotalOutputs[:,self.selectedIndices]
+ Data = Data[self.selectedIndices]
+ covMatrix = covMatrix[self.selectedIndices,self.selectedIndices]
+
# Compute loglikelihood
try:
logL = multivariate_normal.logpdf(TotalOutputs, mean=Data, cov=covMatrix)
except:
logL = -np.inf
- # print("logL: {0:.3f}".format(logL))
+
return logL
#--------------------------------------------------------------------------------------------------------
@@ -719,9 +726,9 @@ class BayesInference:
figPosterior.set_size_inches((24,16))
if self.emulator:
- plotname = '/Poisterior_Dist_'+Model.Name+'_emulator'
+ plotname = '/Posterior_Dist_'+Model.Name+'_emulator'
else:
- plotname = '/Poisterior_Dist_'+Model.Name
+ plotname = '/Posterior_Dist_'+Model.Name
figPosterior.savefig('./'+OutputDir+ plotname + '.pdf',
bbox_inches='tight')
diff --git a/BayesValidRox/BayesInference/MCMC.py b/BayesValidRox/BayesInference/MCMC.py
index 11b7fef81..96f4afee6 100755
--- a/BayesValidRox/BayesInference/MCMC.py
+++ b/BayesValidRox/BayesInference/MCMC.py
@@ -64,10 +64,8 @@ class MCMC():
initsamples = priorDist.sample(self.nwalkers).T
except:
# when aPCE selected - gaussian kernel distribution
- # initsamples = priorDist.resample(self.nwalkers).T
- randIdx = np.random.randint(0, len(PCEModel.ExpDesign.Input_distributions[0]),
- self.nwalkers)
- initsamples = PCEModel.ExpDesign.Input_distributions.T[randIdx]
+ theta = PCEModel.ExpDesign.Input_distributions.mean(axis=1)
+ initsamples = [theta + 1e-1*np.multiply(np.random.randn(ndim), theta) for i in range(self.nwalkers)]
else:
# Pick samples based on a uniform dist between min and max of each dim
@@ -82,7 +80,7 @@ class MCMC():
PCEModel.ExpDesign.BoundTuples = BoundTuples
- # TODO: Check if bias term needs to be inferred
+ # Check if bias term needs to be inferred
if Discrepancy.optSigma != 'B':
sigma2Samples = Discrepancy.get_Sample(self.nwalkers)
@@ -135,7 +133,7 @@ class MCMC():
# Burn-in
print("\n Burn-in period is starting:")
pos = sampler.run_mcmc(initsamples, self.nburn, progress=True)
-
+
# Reset sampler
sampler.reset()
@@ -328,7 +326,7 @@ class MCMC():
# Evaluate Model/PCEModel at theta
meanPCEPriorPred, BayesOpts.stdPCEPriorPred = self.eval_model(theta)
-
+
# Likelihood
return BayesOpts.normpdf(meanPCEPriorPred, self.Observation, self.TotalSigma2, Sigma2)
@@ -361,7 +359,7 @@ class MCMC():
else:
# Compute log prior
log_prior = self.log_prior(theta)
- # Compute loLikelihood
+ # Compute log Likelihood
log_likelihood = self.log_likelihood(theta)
return log_prior + log_likelihood
else:
--
GitLab