Merge pull request #32 from MaxBalmus/master

GPErks/perks/history_matching.augment_nimp: optimised the code and cr…

GPErks/perks/history_matching.py

@@ -128,6 +128,7 @@ def augment_nimp(
         self,
         n_total_points,
         scaling=0.1,
+        n_max=2000
     ):
         X = np.copy(self.NIMP)
         lbounds = self.Itrain[:, 0]
@@ -158,18 +159,19 @@ def augment_nimp(
             )
 
             temp = np.random.normal(loc=X, scale=scale)
+            count2 = 0
             while True:
+                count2 +=1
                 l = []
-                for i in range(temp.shape[0]):
-                    d1 = temp[i, :] - lbounds
-                    d2 = ubounds - temp[i, :]
-                    if (
-                        np.sum(np.sign(d1)) != temp.shape[1]
-                        or np.sum(np.sign(d2)) != temp.shape[1]
-                    ):
-                        l.append(i)
-                if l:
-                    temp[l, :] = np.random.normal(loc=X[l, :], scale=scale)
+                d1 = temp - lbounds.reshape((1,-1))
+                d2 = ubounds.reshape((1,-1)) - temp
+                flag = np.logical_or(np.sum(np.sign(d1), axis=1) != temp.shape[1], np.sum(np.sign(d2), axis=1) != temp.shape[1])
+                if count2 > n_max:
+                    temp = temp[~flag,:]
+                    break
+                if np.sum(flag) > 0:
+                    temp[flag, :] = np.random.normal(loc=X[flag, :], scale=scale)
+                    continue
                 else:
                     break
 
@@ -229,7 +231,7 @@ def copy(self):
         W.imp_idx = self.imp_idx.copy()
         return W
 
-    def plot_wave(self, xlabels=None, display="impl", filepath=None):
+    def plot_wave(self, xlabels=None, display="impl", filepath=None, figsize=None):
         X = self.reconstruct_tests()
         input_dim = X.shape[1]
 
@@ -268,7 +270,10 @@ def plot_wave(self, xlabels=None, display="impl", filepath=None):
 
         height = 9.36111
         width = 5.91667
-        fig = plt.figure(figsize=(2 * width, 1.2 * 2 * height / 3))
+        if figsize is None:
+            fig = plt.figure(figsize=(2 * width, 1.2 * 2 * height / 3))
+        else:
+            fig = plt.figure(figsize=figsize)
         gs = grsp.GridSpec(
             input_dim - 1,
             input_dim,

examples/example_1.py

@@ -7,9 +7,10 @@ def main():
     import numpy as np
     import torch
 
+    from GPErks.constants import DEFAULT_RANDOM_SEED
+
     # enforce reproducibility
     from GPErks.utils.random import set_seed
-    from GPErks.constants import DEFAULT_RANDOM_SEED
     seed = DEFAULT_RANDOM_SEED
     set_seed(seed)  # reproducible sampling
 

examples/example_10.py

@@ -4,13 +4,14 @@
 #
 def main():
     # import main libraries
-    import numpy as np
     import matplotlib.pyplot as plt
+    import numpy as np
     import torch
-
+
+    from GPErks.constants import DEFAULT_RANDOM_SEED
+
     # enforce reproducibility
     from GPErks.utils.random import set_seed
-    from GPErks.constants import DEFAULT_RANDOM_SEED
     seed = DEFAULT_RANDOM_SEED
     set_seed(seed)
 
@@ -64,6 +65,7 @@ def f(x):  # elliptic paraboloid: y = x1^2/a^2 + x2^2/b^2, with a=b=1
     # choose likelihood, mean function and covariance function
     from gpytorch.kernels import MaternKernel, ScaleKernel
     from gpytorch.likelihoods import GaussianLikelihood
+
     from GPErks.gp.mean import LinearMean
     likelihood = GaussianLikelihood()
     mean = LinearMean(degree=2, input_size=dataset.input_size, bias=True)

examples/example_11.py

@@ -2,29 +2,27 @@
 #
 # 11. Bayesian History Matching technique (advanced use)
 #
+import json
 import os
 from pathlib import Path
 
 import matplotlib.pyplot as plt
 import numpy as np
-import json
-
-from GPErks.constants import DEFAULT_TMP_OUTFILE_DIR
-from GPErks.perks.history_matching import Wave
-from GPErks.serialization.path import posix_path
-from GPErks.utils.array import get_minmax
-from GPErks.utils.plotting import interp_col, get_col
-from GPErks.utils.sampling import Sampler
 from gpytorch.kernels import MaternKernel, ScaleKernel
 from gpytorch.likelihoods import GaussianLikelihood
 from torchmetrics import MeanSquaredError, R2Score
 
-from GPErks.constants import DEFAULT_RANDOM_SEED
+from GPErks.constants import DEFAULT_RANDOM_SEED, DEFAULT_TMP_OUTFILE_DIR
 from GPErks.gp.data.dataset import Dataset
 from GPErks.gp.experiment import GPExperiment
 from GPErks.gp.mean import LinearMean
+from GPErks.perks.history_matching import Wave
+from GPErks.serialization.path import posix_path
 from GPErks.train.emulator import GPEmulator
+from GPErks.utils.array import get_minmax
+from GPErks.utils.plotting import get_col, interp_col
 from GPErks.utils.random import set_seed
+from GPErks.utils.sampling import Sampler
 
 
 def main():

examples/example_2.py

@@ -6,9 +6,10 @@ def main():
     # import main libraries
     import torch
 
+    from GPErks.constants import DEFAULT_RANDOM_SEED
+
     # enforce reproducibility
     from GPErks.utils.random import set_seed
-    from GPErks.constants import DEFAULT_RANDOM_SEED
     seed = DEFAULT_RANDOM_SEED
     set_seed(seed)  # reproducible sampling
 
@@ -44,7 +45,10 @@ def main():
 
     # choose metrics
     from torchmetrics import R2Score
-    from GPErks.utils.metrics import IndependentStandardError  # not available from torchmetrics -> we implemented it
+
+    from GPErks.utils.metrics import (
+        IndependentStandardError,  # not available from torchmetrics -> we implemented it
+    )
     metrics = [IndependentStandardError(), R2Score()]
 
     # define experiment

examples/example_3.py

@@ -5,12 +5,14 @@
 def main():
     # import main libraries
     import os
-    import torch
     from pathlib import Path
 
+    import torch
+
+    from GPErks.constants import DEFAULT_RANDOM_SEED
+
     # enforce reproducibility
     from GPErks.utils.random import set_seed
-    from GPErks.constants import DEFAULT_RANDOM_SEED
     seed = DEFAULT_RANDOM_SEED
     set_seed(seed)  # reproducible sampling
 

examples/example_4.py

@@ -5,13 +5,15 @@
 def main():
     # import main libraries
     import os
+    from pathlib import Path
+
     import numpy as np
     import torch
-    from pathlib import Path
+
+    from GPErks.constants import DEFAULT_RANDOM_SEED
 
     # enforce reproducibility
     from GPErks.utils.random import set_seed
-    from GPErks.constants import DEFAULT_RANDOM_SEED
     seed = DEFAULT_RANDOM_SEED
     set_seed(seed)
 
@@ -60,9 +62,10 @@ def main():
     dataset.plot_pairwise()
 
     # define experiment
-    from gpytorch.likelihoods import GaussianLikelihood
     from gpytorch.kernels import MaternKernel, ScaleKernel
+    from gpytorch.likelihoods import GaussianLikelihood
     from torchmetrics import MeanSquaredError, R2Score
+
     from GPErks.gp.experiment import GPExperiment
     from GPErks.gp.mean import LinearMean
 

examples/example_5.py

@@ -5,13 +5,15 @@
 def main():
     # import main libraries
     import os
+    from pathlib import Path
+
     import numpy as np
     import torch
-    from pathlib import Path
+
+    from GPErks.constants import DEFAULT_RANDOM_SEED
 
     # enforce reproducibility
     from GPErks.utils.random import set_seed
-    from GPErks.constants import DEFAULT_RANDOM_SEED
     seed = DEFAULT_RANDOM_SEED
     set_seed(seed)
 
@@ -47,9 +49,10 @@ def main():
     dataset.summary()
 
     # define experiment
-    from gpytorch.likelihoods import GaussianLikelihood
     from gpytorch.kernels import MaternKernel, ScaleKernel
+    from gpytorch.likelihoods import GaussianLikelihood
     from torchmetrics import MeanSquaredError, R2Score
+
     from GPErks.gp.experiment import GPExperiment
     from GPErks.gp.mean import LinearMean
 
@@ -122,13 +125,15 @@ def main():
 def main():
     # import main libraries
     import os
+    from pathlib import Path
+
     import numpy as np
     import torch
-    from pathlib import Path
+
+    from GPErks.constants import DEFAULT_RANDOM_SEED
 
     # enforce reproducibility
     from GPErks.utils.random import set_seed
-    from GPErks.constants import DEFAULT_RANDOM_SEED
     seed = DEFAULT_RANDOM_SEED
     set_seed(seed)
 
@@ -159,9 +164,10 @@ def main():
     )
 
     # define experiment
-    from gpytorch.likelihoods import GaussianLikelihood
     from gpytorch.kernels import MaternKernel, ScaleKernel
+    from gpytorch.likelihoods import GaussianLikelihood
     from torchmetrics import MeanSquaredError, R2Score
+
     from GPErks.gp.experiment import GPExperiment
     from GPErks.gp.mean import LinearMean
 
@@ -218,8 +224,8 @@ def main():
     print(best_epochs)
 
     # train emulator on all the available points
-    from GPErks.train.emulator import GPEmulator
     from GPErks.train.early_stop import NoEarlyStoppingCriterion
+    from GPErks.train.emulator import GPEmulator
 
     # build a new dataset now including a testing set
     del dataset

examples/example_6.py

@@ -5,13 +5,15 @@
 def main():
     # import main libraries
     import os
+    from pathlib import Path
+
     import numpy as np
     import torch
-    from pathlib import Path
+
+    from GPErks.constants import DEFAULT_RANDOM_SEED
 
     # enforce reproducibility
     from GPErks.utils.random import set_seed
-    from GPErks.constants import DEFAULT_RANDOM_SEED
     seed = DEFAULT_RANDOM_SEED
     set_seed(seed)
 
@@ -42,9 +44,10 @@ def main():
     )
 
     # define experiment
-    from gpytorch.likelihoods import GaussianLikelihood
     from gpytorch.kernels import MaternKernel, ScaleKernel
+    from gpytorch.likelihoods import GaussianLikelihood
     from torchmetrics import MeanSquaredError, R2Score
+
     from GPErks.gp.experiment import GPExperiment
     from GPErks.gp.mean import LinearMean
 
@@ -101,8 +104,8 @@ def main():
     print(best_epochs)
 
     # train emulator on all the available points
-    from GPErks.train.emulator import GPEmulator
     from GPErks.train.early_stop import NoEarlyStoppingCriterion
+    from GPErks.train.emulator import GPEmulator
 
     # build a new dataset now including a testing set
     del dataset

examples/example_7.py

@@ -4,12 +4,14 @@
 #
 def main():
     # import main libraries
-    import numpy as np
     from functools import partial
 
+    import numpy as np
+
+    from GPErks.constants import DEFAULT_RANDOM_SEED
+
     # enforce reproducibility
     from GPErks.utils.random import set_seed
-    from GPErks.constants import DEFAULT_RANDOM_SEED
     seed = DEFAULT_RANDOM_SEED
     set_seed(seed)
 
@@ -83,9 +85,10 @@ def main():
     # let's estimate Sobol' indices using an emulator of the SobolGstar function
 
     # define experiment
-    from gpytorch.likelihoods import GaussianLikelihood
     from gpytorch.kernels import MaternKernel, ScaleKernel
+    from gpytorch.likelihoods import GaussianLikelihood
     from torchmetrics import MeanSquaredError, R2Score
+
     from GPErks.gp.experiment import GPExperiment
     from GPErks.gp.mean import LinearMean
     likelihood = GaussianLikelihood()
@@ -104,8 +107,9 @@ def main():
 
     # train emulator
     import torch
-    from GPErks.train.emulator import GPEmulator
+
     from GPErks.train.early_stop import GLEarlyStoppingCriterion
+    from GPErks.train.emulator import GPEmulator
     device = torch.device("cuda" if torch.cuda.is_available() else "cpu")
     emulator = GPEmulator(experiment, device)
     # emulator.train_auto()

examples/example_8.py

@@ -3,20 +3,20 @@
 # 8. GPE auto-training + GSA using external dataset (from publication) loaded from json file
 #
 import os
-import torch
 
+import torch
+from gpytorch.kernels import MaternKernel, ScaleKernel
 from gpytorch.likelihoods import GaussianLikelihood
 from gpytorch.means import LinearMean
-from gpytorch.kernels import MaternKernel, ScaleKernel
 
 from GPErks.constants import DEFAULT_RANDOM_SEED
 from GPErks.gp.data.dataset import Dataset
 from GPErks.gp.experiment import GPExperiment
 from GPErks.gp.mean import LinearMean
+from GPErks.perks.gsa import SobolGSA
 from GPErks.serialization.path import posix_path
 from GPErks.train.emulator import GPEmulator
 from GPErks.utils.random import set_seed
-from GPErks.perks.gsa import SobolGSA
 
 
 def main():

examples/example_9.py

@@ -7,9 +7,10 @@ def main():
     # import main libraries
     import torch
 
+    from GPErks.constants import DEFAULT_RANDOM_SEED
+
     # enforce reproducibility
     from GPErks.utils.random import set_seed
-    from GPErks.constants import DEFAULT_RANDOM_SEED
     seed = DEFAULT_RANDOM_SEED
     set_seed(seed)
 
@@ -46,6 +47,7 @@ def main():
 
     # choose metrics
     from torchmetrics import MeanSquaredError, R2Score
+
     from GPErks.utils.metrics import IndependentStandardError
     metrics = [MeanSquaredError(), R2Score(), IndependentStandardError()]