From 80c3b0333143025d4a5adcbdd6e080d79d2b3e07 Mon Sep 17 00:00:00 2001 From: Dimitris Tsapetis Date: Mon, 21 Nov 2022 11:54:22 -0500 Subject: [PATCH 1/6] Removes old kriging files --- src/UQpy/reliability/taylor_series/FORM.py | 2 +- src/UQpy/surrogates/__init__.py | 1 - src/UQpy/surrogates/kriging/Kriging.py | 354 ------------------ src/UQpy/surrogates/kriging/__init__.py | 4 - .../correlation_models/CubicCorrelation.py | 28 -- .../ExponentialCorrelation.py | 20 - .../correlation_models/GaussianCorrelation.py | 21 -- .../correlation_models/LinearCorrelation.py | 35 -- .../SphericalCorrelation.py | 28 -- .../correlation_models/SplineCorrelation.py | 58 --- .../kriging/correlation_models/__init__.py | 7 - .../baseclass/Correlation.py | 47 --- .../correlation_models/baseclass/__init__.py | 1 - .../regression_models/ConstantRegression.py | 10 - .../regression_models/LinearRegression.py | 12 - .../regression_models/QuadraticRegression.py | 37 -- .../kriging/regression_models/__init__.py | 4 - .../regression_models/baseclass/Regression.py | 14 - .../regression_models/baseclass/__init__.py | 1 - .../sampling/test_adaptive_kriging.py | 1 - .../sampling/test_refined_stratified.py | 1 - tests/unit_tests/surrogates/test_kriging.py | 198 ---------- 22 files changed, 1 insertion(+), 883 deletions(-) delete mode 100755 src/UQpy/surrogates/kriging/Kriging.py delete mode 100644 src/UQpy/surrogates/kriging/__init__.py delete mode 100644 src/UQpy/surrogates/kriging/correlation_models/CubicCorrelation.py delete mode 100644 src/UQpy/surrogates/kriging/correlation_models/ExponentialCorrelation.py delete mode 100644 src/UQpy/surrogates/kriging/correlation_models/GaussianCorrelation.py delete mode 100644 src/UQpy/surrogates/kriging/correlation_models/LinearCorrelation.py delete mode 100644 src/UQpy/surrogates/kriging/correlation_models/SphericalCorrelation.py delete mode 100644 src/UQpy/surrogates/kriging/correlation_models/SplineCorrelation.py delete mode 100644 src/UQpy/surrogates/kriging/correlation_models/__init__.py delete mode 100644 src/UQpy/surrogates/kriging/correlation_models/baseclass/Correlation.py delete mode 100644 src/UQpy/surrogates/kriging/correlation_models/baseclass/__init__.py delete mode 100644 src/UQpy/surrogates/kriging/regression_models/ConstantRegression.py delete mode 100644 src/UQpy/surrogates/kriging/regression_models/LinearRegression.py delete mode 100644 src/UQpy/surrogates/kriging/regression_models/QuadraticRegression.py delete mode 100644 src/UQpy/surrogates/kriging/regression_models/__init__.py delete mode 100644 src/UQpy/surrogates/kriging/regression_models/baseclass/Regression.py delete mode 100644 src/UQpy/surrogates/kriging/regression_models/baseclass/__init__.py delete mode 100644 tests/unit_tests/surrogates/test_kriging.py diff --git a/src/UQpy/reliability/taylor_series/FORM.py b/src/UQpy/reliability/taylor_series/FORM.py index a3cd250fe..dee48a4fc 100644 --- a/src/UQpy/reliability/taylor_series/FORM.py +++ b/src/UQpy/reliability/taylor_series/FORM.py @@ -312,7 +312,7 @@ def run(self, seed_x: Union[list, np.ndarray] = None, else: k = k + 1 - self.logger.error("Error: %s", error_record[-1]) + self.logger.info("Error: %s", error_record[-1]) if converged is True or k > self.n_iterations: break diff --git a/src/UQpy/surrogates/__init__.py b/src/UQpy/surrogates/__init__.py index 168ebdbfa..fe76e50e2 100644 --- a/src/UQpy/surrogates/__init__.py +++ b/src/UQpy/surrogates/__init__.py @@ -1,6 +1,5 @@ from UQpy.surrogates.polynomial_chaos import * from UQpy.surrogates.stochastic_reduced_order_models import * -from UQpy.surrogates.kriging import * from UQpy.surrogates.gaussian_process import * from UQpy.surrogates.baseclass import * diff --git a/src/UQpy/surrogates/kriging/Kriging.py b/src/UQpy/surrogates/kriging/Kriging.py deleted file mode 100755 index 19cbfaf1d..000000000 --- a/src/UQpy/surrogates/kriging/Kriging.py +++ /dev/null @@ -1,354 +0,0 @@ -import logging -from typing import Callable - -import numpy as np -from scipy.linalg import cholesky -import scipy.stats as stats -from beartype import beartype - -from UQpy.utilities import MinimizeOptimizer -from UQpy.utilities.Utilities import process_random_state -from UQpy.surrogates.baseclass.Surrogate import Surrogate -from UQpy.utilities.ValidationTypes import RandomStateType -from UQpy.surrogates.kriging.correlation_models.baseclass.Correlation import Correlation -from UQpy.surrogates.kriging.regression_models.baseclass.Regression import Regression - - -class Kriging(Surrogate): - @beartype - def __init__( - self, - regression_model: Regression, - correlation_model: Correlation, - correlation_model_parameters: list, - optimizer, - bounds: list = None, - optimize: bool = True, - optimizations_number: int = 1, - normalize: bool = True, - random_state: RandomStateType = None, - ): - """ - Κriging generates an Gaussian process regression-based surrogate model to predict the model output at new sample - points. - - :param regression_model: `regression_model` specifies and evaluates the basis functions and their coefficients, - which defines the trend of the model. Built-in options: :class:`Constant`, :class:`Linear`, :class:`Quadratic` - :param correlation_model: `correlation_model` specifies and evaluates the correlation function. - Built-in options: :class:`Exponential`, :class:`Gaussian`, :class:`Linear`, :class:`Spherical`, - :class:`Cubic`, :class:`Spline` - :param correlation_model_parameters: List or array of initial values for the correlation model - hyperparameters/scale parameters. - :param bounds: Bounds on the hyperparameters used to solve optimization problem to estimate maximum likelihood - estimator. This should be a closed bound. - Default: :math:`[0.001, 10^7]` for each hyperparameter. - :param optimize: Indicator to solve MLE problem or not. If :any:'True' corr_model_params will be used as initial - solution for optimization problem. Otherwise, correlation_model_parameters will be directly use as the - hyperparamters. - Default: :any:`True`. - :param optimizations_number: Number of times MLE optimization problem is to be solved with a random starting - point. Default: :math:`1`. - :param normalize: Boolean flag used in case data normalization is required. - :param optimizer: Object of the :class:`Optimizer` optimizer used during the Kriging surrogate. - Default: :class:`.MinimizeOptimizer`. - :param random_state: Random seed used to initialize the pseudo-random number generator. If an :any:`int` is - provided, this sets the seed for an object of :class:`numpy.random.RandomState`. Otherwise, the - object itself can be passed directly. - """ - self.regression_model = regression_model - self.correlation_model = correlation_model - self.correlation_model_parameters = np.array(correlation_model_parameters) - self.bounds = bounds - self.optimizer = optimizer - self.optimizations_number = optimizations_number - self.optimize = optimize - self.normalize = normalize - self.logger = logging.getLogger(__name__) - self.random_state = random_state - - # Variables are used outside the __init__ - self.samples = None - self.values = None - self.sample_mean, self.sample_std = None, None - self.value_mean, self.value_std = None, None - self.rmodel, self.cmodel = None, None - self.beta: list = None - """Regression coefficients.""" - self.gamma = None - self.err_var: float = None - """Variance of the Gaussian random process.""" - self.F_dash = None - self.C_inv = None - self.G = None - self.F, self.R = None, None - - if isinstance(self.optimizer, str): - raise ValueError("The optimization function provided a input parameter cannot be None.") - - if optimizer._bounds is None: - optimizer.update_bounds([[0.001, 10 ** 7]] * self.correlation_model_parameters.shape[0]) - - self.jac = optimizer.supports_jacobian() - self.random_state = process_random_state(random_state) - - def fit( - self, - samples, - values, - optimizations_number: int = None, - correlation_model_parameters: list = None, - ): - """ - Fit the surrogate model using the training samples and the corresponding model values. - - The user can run this method multiple time after initiating the :class:`.Kriging` class object. - - This method updates the samples and parameters of the :class:`.Kriging` object. This method uses - `correlation_model_parameters` from previous run as the starting point for MLE problem unless user provides a - new starting point. - - :param samples: :class:`numpy.ndarray` containing the training points. - :param values: :class:`numpy.ndarray` containing the model evaluations at the training points. - :param optimizations_number: number of optimization iterations - :param correlation_model_parameters: List or array of initial values for the correlation model - hyperparameters/scale parameters. - - The :meth:`fit` method has no returns, although it creates the :py:attr:`beta`, :py:attr:`err_var` and - :py:attr:`C_inv` attributes of the :class:`.Kriging` class. - """ - self.logger.info("UQpy: Running kriging.fit") - - if optimizations_number is not None: - self.optimizations_number = optimizations_number - if correlation_model_parameters is not None: - self.correlation_model_parameters = np.array(correlation_model_parameters) - self.samples = np.array(samples) - - # Number of samples and dimensions of samples and values - nsamples, input_dim = self.samples.shape - output_dim = int(np.size(values) / nsamples) - - self.values = np.array(values).reshape(nsamples, output_dim) - - # Normalizing the data - if self.normalize: - self.sample_mean, self.sample_std = np.mean(self.samples, 0), np.std(self.samples, 0) - self.value_mean, self.value_std = np.mean(self.values, 0), np.std(self.values, 0) - s_ = (self.samples - self.sample_mean) / self.sample_std - y_ = (self.values - self.value_mean) / self.value_std - else: - s_ = self.samples - y_ = self.values - - self.F, jf_ = self.regression_model.r(s_) - - # Maximum Likelihood Estimation : Solving optimization problem to calculate hyperparameters - if self.optimize: - starting_point = self.correlation_model_parameters - - minimizer, fun_value = np.zeros([self.optimizations_number, input_dim]),\ - np.zeros([self.optimizations_number, 1]) - for i__ in range(self.optimizations_number): - p_ = self.optimizer.optimize(function=Kriging.log_likelihood, - initial_guess=starting_point, - args=(self.correlation_model, s_, self.F, y_, self.jac), - jac=self.jac) - print(p_.success) - # print(self.kwargs_optimizer) - minimizer[i__, :] = p_.x - fun_value[i__, 0] = p_.fun - # Generating new starting points using log-uniform distribution - if i__ != self.optimizations_number - 1: - starting_point = stats.reciprocal.rvs([j[0] for j in self.optimizer._bounds], - [j[1] for j in self.optimizer._bounds], 1, - random_state=self.random_state) - print(starting_point) - - if min(fun_value) == np.inf: - raise NotImplementedError("Maximum likelihood estimator failed: Choose different starting point or " - "increase nopt") - t = np.argmin(fun_value) - self.correlation_model_parameters = minimizer[t, :] - - # Updated Correlation matrix corresponding to MLE estimates of hyperparameters - self.R = self.correlation_model.c(x=s_, s=s_, params=self.correlation_model_parameters) - - self.beta, self.gamma, tmp = self._compute_additional_parameters(self.R) - self.C_inv, self.F_dash, self.G, self.err_var = tmp[1], tmp[3], tmp[2], tmp[5] - - self.logger.info("UQpy: kriging fit complete.") - - def _compute_additional_parameters(self, correlation_matrix): - if self.normalize: - y_ = (self.values - self.value_mean) / self.value_std - else: - y_ = self.values - # Compute the regression coefficient (solving this linear equation: F * beta = Y) - # Eq: 3.8, DACE - c = cholesky(correlation_matrix + (10 + self.samples.shape[0]) * 2 ** (-52) * np.eye(self.samples.shape[0]), - lower=True, check_finite=False) - c_inv = np.linalg.inv(c) - f_dash = np.linalg.solve(c, self.F) - y_dash = np.linalg.solve(c, y_) - q_, g_ = np.linalg.qr(f_dash) # Eq: 3.11, DACE - # Check if F is a full rank matrix - if np.linalg.matrix_rank(g_) != min(np.size(self.F, 0), np.size(self.F, 1)): - raise NotImplementedError("Chosen regression functions are not sufficiently linearly independent") - # Design parameters (beta: regression coefficient) - beta = np.linalg.solve(g_, np.matmul(np.transpose(q_), y_dash)) - - # Design parameter (R * gamma = Y - F * beta = residual) - gamma = np.linalg.solve(c.T, (y_dash - np.matmul(f_dash, beta))) - - # Computing the process variance (Eq: 3.13, DACE) - err_var = np.zeros(self.values.shape[1]) - for i in range(self.values.shape[1]): - err_var[i] = (1 / self.samples.shape[0]) * (np.linalg.norm(y_dash[:, i] - - np.matmul(f_dash, beta[:, i])) ** 2) - - return beta, gamma, (c, c_inv, g_, f_dash, y_dash, err_var) - - def predict(self, points: np.ndarray, return_std: bool = False, correlation_model_parameters: list = None): - """ - Predict the model response at new points. - - This method evaluates the regression and correlation model at new sample points. Then, it predicts the function - value and standard deviation. - - :param points: Points at which to predict the model response. - :param return_std: Indicator to estimate standard deviation. - :param correlation_model_parameters: Hyperparameters for correlation model. - :return: Predicted values at the new points, Standard deviation of predicted values at the new points - """ - x_ = np.atleast_2d(points) - if self.normalize: - x_ = (x_ - self.sample_mean) / self.sample_std - s_ = (self.samples - self.sample_mean) / self.sample_std - else: - s_ = self.samples - fx, jf = self.regression_model.r(x_) - if correlation_model_parameters is None: - correlation_model_parameters = self.correlation_model_parameters - rx = self.correlation_model.c( - x=x_, s=s_, params=correlation_model_parameters - ) - if correlation_model_parameters is None: - beta, gamma = self.beta, self.gamma - c_inv, f_dash, g_, err_var = self.C_inv, self.F_dash, self.G, self.err_var - else: - beta, gamma, tmp = self._compute_additional_parameters( - self.correlation_model.c(x=s_, s=s_, params=correlation_model_parameters)) - c_inv, f_dash, g_, err_var = tmp[1], tmp[3], tmp[2], tmp[5] - y = np.einsum("ij,jk->ik", fx, beta) + np.einsum( - "ij,jk->ik", rx, gamma - ) - if self.normalize: - y = self.value_mean + y * self.value_std - if x_.shape[1] == 1: - y = y.flatten() - if return_std: - r_dash = np.matmul(c_inv, rx.T) - u = np.matmul(f_dash.T, r_dash) - fx.T - norm1 = np.linalg.norm(r_dash, 2, 0) - norm2 = np.linalg.norm(np.linalg.solve(g_, u), 2, 0) - mse = np.sqrt(err_var * np.atleast_2d(1 + norm2 - norm1).T) - if self.normalize: - mse = self.value_std * mse - if x_.shape[1] == 1: - mse = mse.flatten() - return y, mse - else: - return y - - def jacobian(self, points: np.ndarray): - """ - Predict the gradient of the model at new points. - - This method evaluates the regression and correlation model at new sample point. Then, it predicts the gradient - using the regression coefficients and the training second_order_tensor. - - :param points: Points at which to evaluate the gradient. - :return: Gradient of the surrogate model evaluated at the new points. - """ - x_ = np.atleast_2d(points) - if self.normalize: - x_ = (x_ - self.sample_mean) / self.sample_std - s_ = (self.samples - self.sample_mean) / self.sample_std - else: - s_ = self.samples - - fx, jf = self.regression_model.r(x_) - rx, drdx = self.correlation_model.c( - x=x_, s=s_, params=self.correlation_model_parameters, dx=True - ) - y_grad = np.einsum("ikj,jm->ik", jf, self.beta) + np.einsum( - "ijk,jm->ki", drdx.T, self.gamma - ) - if self.normalize: - y_grad = y_grad * self.value_std / self.sample_std - if x_.shape[1] == 1: - y_grad = y_grad.flatten() - return y_grad - - @staticmethod - def log_likelihood(p0, cm, s, f, y, return_grad): - # Return the log-likelihood function and it's gradient. Gradient is calculate using Central Difference - m = s.shape[0] - n = s.shape[1] - r__, dr_ = cm.c(x=s, s=s, params=p0, dt=True) - try: - cc = cholesky(r__ + 2 ** (-52) * np.eye(m), lower=True) - except np.linalg.LinAlgError: - return np.inf, np.zeros(n) - - # Product of diagonal terms is negligible sometimes, even when cc exists. - if np.prod(np.diagonal(cc)) == 0: - return np.inf, np.zeros(n) - - cc_inv = np.linalg.inv(cc) - r_inv = np.matmul(cc_inv.T, cc_inv) - f__ = cc_inv.dot(f) - y__ = cc_inv.dot(y) - - q__, g__ = np.linalg.qr(f__) # Eq: 3.11, DACE - - # Check if F is a full rank matrix - if np.linalg.matrix_rank(g__) != min(np.size(f__, 0), np.size(f__, 1)): - raise NotImplementedError( - "Chosen regression functions are not sufficiently linearly independent" - ) - - # Design parameters - beta_ = np.linalg.solve(g__, np.matmul(np.transpose(q__), y__)) - - # Computing the process variance (Eq: 3.13, DACE) - sigma_ = np.zeros(y.shape[1]) - - ll = 0 - for out_dim in range(y.shape[1]): - sigma_[out_dim] = (1 / m) * ( - np.linalg.norm(y__[:, out_dim] - np.matmul(f__, beta_[:, out_dim])) ** 2) - # Objective function:= log(det(sigma**2 * R)) + constant - ll = (ll + ( np.log(np.linalg.det(sigma_[out_dim] * r__)) + m * (np.log(2 * np.pi) + 1)) / 2) - - # Gradient of loglikelihood - # Reference: C. E. Rasmussen & C. K. I. Williams, Gaussian Processes for Machine Learning, the MIT Press, - # 2006, ISBN 026218253X. (Page 114, Eq.(5.9)) - residual = y - np.matmul(f, beta_) - gamma = np.matmul(r_inv, residual) - grad_mle = np.zeros(n) - for in_dim in range(n): - r_inv_derivative = np.matmul(r_inv, np.matmul(dr_[:, :, in_dim], r_inv)) - tmp = np.matmul(residual.T, np.matmul(r_inv_derivative, residual)) - for out_dim in range(y.shape[1]): - alpha = gamma / sigma_[out_dim] - tmp1 = np.matmul(alpha, alpha.T) - r_inv / sigma_[out_dim] - cov_der = sigma_[out_dim] * dr_[:, :, in_dim] + tmp * r__ / m - grad_mle[in_dim] = grad_mle[in_dim] - 0.5 * np.trace( - np.matmul(tmp1, cov_der) - ) - - if return_grad: - return ll, grad_mle - else: - return ll diff --git a/src/UQpy/surrogates/kriging/__init__.py b/src/UQpy/surrogates/kriging/__init__.py deleted file mode 100644 index 55a50199d..000000000 --- a/src/UQpy/surrogates/kriging/__init__.py +++ /dev/null @@ -1,4 +0,0 @@ -from UQpy.surrogates.kriging.Kriging import Kriging - -from UQpy.surrogates.kriging.regression_models import * -from UQpy.surrogates.kriging.correlation_models import * diff --git a/src/UQpy/surrogates/kriging/correlation_models/CubicCorrelation.py b/src/UQpy/surrogates/kriging/correlation_models/CubicCorrelation.py deleted file mode 100644 index 3e909507f..000000000 --- a/src/UQpy/surrogates/kriging/correlation_models/CubicCorrelation.py +++ /dev/null @@ -1,28 +0,0 @@ -from UQpy.surrogates.kriging.correlation_models.baseclass.Correlation import * - - -class CubicCorrelation(Correlation): - def c(self, x, s, params, dt=False, dx=False): - zeta_matrix, dtheta_derivs, dx_derivs = Correlation.derivatives( - x_=x, s_=s, params=params - ) - # Initial matrices containing derivates for all values in array. Note since - # dtheta_s and dx_s already accounted for where derivative should be zero, all - # that must be done is multiplying the |dij| or thetaj matrix on top of a - # matrix of derivates w.r.t zeta (in this case, dzeta = -6zeta+6zeta**2) - drdt = (-6 * zeta_matrix + 6 * zeta_matrix ** 2) * dtheta_derivs - drdx = (-6 * zeta_matrix + 6 * zeta_matrix ** 2) * dx_derivs - # Also, create matrix for values of equation, 1 - 3zeta**2 + 2zeta**3, for loop - zeta_function_cubic = 1 - 3 * zeta_matrix ** 2 + 2 * zeta_matrix ** 3 - rx = np.prod(zeta_function_cubic, 2) - if dt: - # Same as previous example, loop over zeta matrix by shifting index - for i in range(len(params) - 1): - drdt = drdt * np.roll(zeta_function_cubic, i + 1, axis=2) - return rx, drdt - if dx: - # Same as previous example, loop over zeta matrix by shifting index - for i in range(len(params) - 1): - drdx = drdx * np.roll(zeta_function_cubic, i + 1, axis=2) - return rx, drdx - return rx diff --git a/src/UQpy/surrogates/kriging/correlation_models/ExponentialCorrelation.py b/src/UQpy/surrogates/kriging/correlation_models/ExponentialCorrelation.py deleted file mode 100644 index 94702760c..000000000 --- a/src/UQpy/surrogates/kriging/correlation_models/ExponentialCorrelation.py +++ /dev/null @@ -1,20 +0,0 @@ -from UQpy.surrogates.kriging.correlation_models.baseclass.Correlation import * - - -class ExponentialCorrelation(Correlation): - def c(self, x, s, params, dt=False, dx=False): - stack = Correlation.check_samples_and_return_stack(x, s) - rx = np.exp(np.sum(-params * abs(stack), axis=2)) - if dt: - drdt = -abs(stack) * np.transpose( - np.tile(rx, (np.size(x, 1), 1, 1)), (1, 2, 0) - ) - return rx, drdt - if dx: - drdx = ( - -params - * np.sign(stack) - * np.transpose(np.tile(rx, (np.size(x, 1), 1, 1)), (1, 2, 0)) - ) - return rx, drdx - return rx diff --git a/src/UQpy/surrogates/kriging/correlation_models/GaussianCorrelation.py b/src/UQpy/surrogates/kriging/correlation_models/GaussianCorrelation.py deleted file mode 100644 index 05ce09830..000000000 --- a/src/UQpy/surrogates/kriging/correlation_models/GaussianCorrelation.py +++ /dev/null @@ -1,21 +0,0 @@ -from UQpy.surrogates.kriging.correlation_models.baseclass.Correlation import * - - -class GaussianCorrelation(Correlation): - def c(self, x, s, params, dt=False, dx=False): - stack = Correlation.check_samples_and_return_stack(x, s) - rx = np.exp(np.sum(-params * (stack ** 2), axis=2)) - if dt: - drdt = -(stack ** 2) * np.transpose( - np.tile(rx, (np.size(x, 1), 1, 1)), (1, 2, 0) - ) - return rx, drdt - if dx: - drdx = ( - -2 - * params - * stack - * np.transpose(np.tile(rx, (np.size(x, 1), 1, 1)), (1, 2, 0)) - ) - return rx, drdx - return rx diff --git a/src/UQpy/surrogates/kriging/correlation_models/LinearCorrelation.py b/src/UQpy/surrogates/kriging/correlation_models/LinearCorrelation.py deleted file mode 100644 index 69d7f1506..000000000 --- a/src/UQpy/surrogates/kriging/correlation_models/LinearCorrelation.py +++ /dev/null @@ -1,35 +0,0 @@ -from UQpy.surrogates.kriging.correlation_models.baseclass.Correlation import * - - -class LinearCorrelation(Correlation): - def c(self, x, s, params, dt=False, dx=False): - stack = Correlation.check_samples_and_return_stack(x, s) - # Taking stack and turning each d value into 1-theta*dij - after_parameters = 1 - params * abs(stack) - # Define matrix of zeros to compare against (not necessary to be defined separately, - # but the line is bulky if this isn't defined first, and it is used more than once) - comp_zero = np.zeros((np.size(x, 0), np.size(s, 0), np.size(s, 1))) - # Compute matrix of max{0,1-theta*d} - max_matrix = np.maximum(after_parameters, comp_zero) - rx = np.prod(max_matrix, 2) - # Create matrix that has 1s where max_matrix is nonzero - # -Essentially, this acts as a way to store the indices of where the values are nonzero - ones_and_zeros = max_matrix.astype(bool).astype(int) - # Set initial derivatives as if all were positive - first_dtheta = -abs(stack) - first_dx = np.negative(params) * np.sign(stack) - # Multiply derivs by ones_and_zeros...this will set the values where the - # derivative should be zero to zero, and keep all other values the same - drdt = np.multiply(first_dtheta, ones_and_zeros) - drdx = np.multiply(first_dx, ones_and_zeros) - if dt: - # Loop over parameters, shifting max_matrix and multiplying over derivative matrix with each iter - for i in range(len(params) - 1): - drdt = drdt * np.roll(max_matrix, i + 1, axis=2) - return rx, drdt - if dx: - # Loop over parameters, shifting max_matrix and multiplying over derivative matrix with each iter - for i in range(len(params) - 1): - drdx = drdx * np.roll(max_matrix, i + 1, axis=2) - return rx, drdx - return rx diff --git a/src/UQpy/surrogates/kriging/correlation_models/SphericalCorrelation.py b/src/UQpy/surrogates/kriging/correlation_models/SphericalCorrelation.py deleted file mode 100644 index 1f6b8173d..000000000 --- a/src/UQpy/surrogates/kriging/correlation_models/SphericalCorrelation.py +++ /dev/null @@ -1,28 +0,0 @@ -from UQpy.surrogates.kriging.correlation_models.baseclass.Correlation import * - - -class SphericalCorrelation(Correlation): - def c(self, x, s, params, dt=False, dx=False): - zeta_matrix, dtheta_derivs, dx_derivs = Correlation.derivatives( - x_=x, s_=s, params=params - ) - # Initial matrices containing derivates for all values in array. Note since - # dtheta_s and dx_s already accounted for where derivative should be zero, all - # that must be done is multiplying the |dij| or thetaj matrix on top of a - # matrix of derivates w.r.t zeta (in this case, dzeta = -1.5+1.5zeta**2) - drdt = (-1.5 + 1.5 * zeta_matrix ** 2) * dtheta_derivs - drdx = (-1.5 + 1.5 * zeta_matrix ** 2) * dx_derivs - # Also, create matrix for values of equation, 1 - 1.5zeta + 0.5zeta**3, for loop - zeta_function = 1 - 1.5 * zeta_matrix + 0.5 * zeta_matrix ** 3 - rx = np.prod(zeta_function, 2) - if dt: - # Same as previous example, loop over zeta matrix by shifting index - for i in range(len(params) - 1): - drdt = drdt * np.roll(zeta_function, i + 1, axis=2) - return rx, drdt - if dx: - # Same as previous example, loop over zeta matrix by shifting index - for i in range(len(params) - 1): - drdx = drdx * np.roll(zeta_function, i + 1, axis=2) - return rx, drdx - return rx diff --git a/src/UQpy/surrogates/kriging/correlation_models/SplineCorrelation.py b/src/UQpy/surrogates/kriging/correlation_models/SplineCorrelation.py deleted file mode 100644 index 0aa6282d1..000000000 --- a/src/UQpy/surrogates/kriging/correlation_models/SplineCorrelation.py +++ /dev/null @@ -1,58 +0,0 @@ -from UQpy.surrogates.kriging.correlation_models.baseclass.Correlation import * - - -class SplineCorrelation(Correlation): - def c(self, x, s, params, dt=False, dx=False): - # x_, s_ = np.atleast_2d(x_), np.atleast_2d(s_) - # # Create stack matrix, where each block is x_i with all s - # stack = np.tile(np.swapaxes(np.atleast_3d(x_), 1, 2), (1, np.size(s_, 0), 1)) - np.tile(s_, ( - # np.size(x_, 0), - # 1, 1)) - stack = Correlation.check_samples_and_return_stack(x, s) - # In this case, the zeta value is just abs(stack)*parameters, no comparison - zeta_matrix = abs(stack) * params - # So, dtheta and dx are just |dj| and theta*sgn(dj), respectively - dtheta_derivs = abs(stack) - # dx_derivs = np.ones((np.size(x,0),np.size(s,0),np.size(s,1)))*parameters - dx_derivs = np.sign(stack) * params - - # Initialize empty sigma and dsigma matrices - sigma = np.ones( - (zeta_matrix.shape[0], zeta_matrix.shape[1], zeta_matrix.shape[2]) - ) - dsigma = np.ones( - (zeta_matrix.shape[0], zeta_matrix.shape[1], zeta_matrix.shape[2]) - ) - - # Loop over cases to create zeta_matrix and subsequent dR matrices - for i in range(zeta_matrix.shape[0]): - for j in range(zeta_matrix.shape[1]): - for k in range(zeta_matrix.shape[2]): - y = zeta_matrix[i, j, k] - if 0 <= y <= 0.2: - sigma[i, j, k] = 1 - 15 * y ** 2 + 30 * y ** 3 - dsigma[i, j, k] = -30 * y + 90 * y ** 2 - elif 0.2 < y < 1.0: - sigma[i, j, k] = 1.25 * (1 - y) ** 3 - dsigma[i, j, k] = 3.75 * (1 - y) ** 2 * -1 - elif y >= 1: - sigma[i, j, k] = 0 - dsigma[i, j, k] = 0 - - rx = np.prod(sigma, 2) - - if dt: - # Initialize derivative matrices incorporating chain rule - drdt = dsigma * dtheta_derivs - # Loop over to create proper matrices - for i in range(len(params) - 1): - drdt = drdt * np.roll(sigma, i + 1, axis=2) - return rx, drdt - if dx: - # Initialize derivative matrices incorporating chain rule - drdx = dsigma * dx_derivs - # Loop over to create proper matrices - for i in range(len(params) - 1): - drdx = drdx * np.roll(sigma, i + 1, axis=2) - return rx, drdx - return rx diff --git a/src/UQpy/surrogates/kriging/correlation_models/__init__.py b/src/UQpy/surrogates/kriging/correlation_models/__init__.py deleted file mode 100644 index 10f39dafc..000000000 --- a/src/UQpy/surrogates/kriging/correlation_models/__init__.py +++ /dev/null @@ -1,7 +0,0 @@ -from UQpy.surrogates.kriging.correlation_models.baseclass import * -from UQpy.surrogates.kriging.correlation_models.CubicCorrelation import CubicCorrelation -from UQpy.surrogates.kriging.correlation_models.ExponentialCorrelation import ExponentialCorrelation -from UQpy.surrogates.kriging.correlation_models.GaussianCorrelation import GaussianCorrelation -from UQpy.surrogates.kriging.correlation_models.LinearCorrelation import LinearCorrelation -from UQpy.surrogates.kriging.correlation_models.SphericalCorrelation import SphericalCorrelation -from UQpy.surrogates.kriging.correlation_models.SplineCorrelation import SplineCorrelation diff --git a/src/UQpy/surrogates/kriging/correlation_models/baseclass/Correlation.py b/src/UQpy/surrogates/kriging/correlation_models/baseclass/Correlation.py deleted file mode 100644 index 703461b5f..000000000 --- a/src/UQpy/surrogates/kriging/correlation_models/baseclass/Correlation.py +++ /dev/null @@ -1,47 +0,0 @@ -from abc import ABC, abstractmethod -import numpy as np - - -class Correlation(ABC): - """ - Abstract base class of all Correlations. Serves as a template for creating new Kriging correlation - functions. - """ - - @abstractmethod - def c(self, x, s, params, dt=False, dx=False): - """ - Abstract method that needs to be implemented by the user when creating a new Correlation function. - """ - pass - - @staticmethod - def check_samples_and_return_stack(x, s): - x_, s_ = np.atleast_2d(x), np.atleast_2d(s) - # Create stack matrix, where each block is x_i with all s - stack = np.tile( - np.swapaxes(np.atleast_3d(x_), 1, 2), (1, np.size(s_, 0), 1) - ) - np.tile(s_, (np.size(x_, 0), 1, 1)) - return stack - - @staticmethod - def derivatives(x_, s_, params): - stack = Correlation.check_samples_and_return_stack(x_, s_) - # Taking stack and creating array of all thetaj*dij - after_parameters = params * abs(stack) - # Create matrix of all ones to compare - comp_ones = np.ones((np.size(x_, 0), np.size(s_, 0), np.size(s_, 1))) - # zeta_matrix has all values min{1,theta*dij} - zeta_matrix_ = np.minimum(after_parameters, comp_ones) - # Copy zeta_matrix to another matrix that will used to find where derivative should be zero - indices = zeta_matrix_.copy() - # If value of min{1,theta*dij} is 1, the derivative should be 0. - # So, replace all values of 1 with 0, then perform the .astype(bool).astype(int) - # operation like in the linear example, so you end up with an array of 1's where - # the derivative should be caluclated and 0 where it should be zero - indices[indices == 1] = 0 - # Create matrix of all |dij| (where non zero) to be used in calculation of dR/dtheta - dtheta_derivs_ = indices.astype(bool).astype(int) * abs(stack) - # Same as above, but for matrix of all thetaj where non-zero - dx_derivs_ = indices.astype(bool).astype(int) * params * np.sign(stack) - return zeta_matrix_, dtheta_derivs_, dx_derivs_ diff --git a/src/UQpy/surrogates/kriging/correlation_models/baseclass/__init__.py b/src/UQpy/surrogates/kriging/correlation_models/baseclass/__init__.py deleted file mode 100644 index e8cf1815d..000000000 --- a/src/UQpy/surrogates/kriging/correlation_models/baseclass/__init__.py +++ /dev/null @@ -1 +0,0 @@ -from UQpy.surrogates.kriging.correlation_models.baseclass.Correlation import Correlation diff --git a/src/UQpy/surrogates/kriging/regression_models/ConstantRegression.py b/src/UQpy/surrogates/kriging/regression_models/ConstantRegression.py deleted file mode 100644 index 0e4f9e984..000000000 --- a/src/UQpy/surrogates/kriging/regression_models/ConstantRegression.py +++ /dev/null @@ -1,10 +0,0 @@ -import numpy as np -from UQpy.surrogates.kriging.regression_models.baseclass.Regression import Regression - - -class ConstantRegression(Regression): - def r(self, s): - s = np.atleast_2d(s) - fx = np.ones([np.size(s, 0), 1]) - jf = np.zeros([np.size(s, 0), np.size(s, 1), 1]) - return fx, jf diff --git a/src/UQpy/surrogates/kriging/regression_models/LinearRegression.py b/src/UQpy/surrogates/kriging/regression_models/LinearRegression.py deleted file mode 100644 index 118d8d73c..000000000 --- a/src/UQpy/surrogates/kriging/regression_models/LinearRegression.py +++ /dev/null @@ -1,12 +0,0 @@ -import numpy as np -from UQpy.surrogates.kriging.regression_models.baseclass.Regression import Regression - - -class LinearRegression(Regression): - def r(self, s): - s = np.atleast_2d(s) - fx = np.concatenate((np.ones([np.size(s, 0), 1]), s), 1) - jf_b = np.zeros([np.size(s, 0), np.size(s, 1), np.size(s, 1)]) - np.einsum("jii->ji", jf_b)[:] = 1 - jf = np.concatenate((np.zeros([np.size(s, 0), np.size(s, 1), 1]), jf_b), 2) - return fx, jf diff --git a/src/UQpy/surrogates/kriging/regression_models/QuadraticRegression.py b/src/UQpy/surrogates/kriging/regression_models/QuadraticRegression.py deleted file mode 100644 index fdddefbb5..000000000 --- a/src/UQpy/surrogates/kriging/regression_models/QuadraticRegression.py +++ /dev/null @@ -1,37 +0,0 @@ -import numpy as np -from UQpy.surrogates.kriging.regression_models.baseclass.Regression import Regression - - -class QuadraticRegression(Regression): - def r(self, s): - s = np.atleast_2d(s) - fx = np.zeros( - [np.size(s, 0), int((np.size(s, 1) + 1) * (np.size(s, 1) + 2) / 2)] - ) - jf = np.zeros( - [ - np.size(s, 0), - np.size(s, 1), - int((np.size(s, 1) + 1) * (np.size(s, 1) + 2) / 2), - ] - ) - for i in range(np.size(s, 0)): - temp = np.hstack((1, s[i, :])) - for j in range(np.size(s, 1)): - temp = np.hstack((temp, s[i, j] * s[i, j::])) - fx[i, :] = temp - # definie H matrix - h_ = 0 - for j in range(np.size(s, 1)): - tmp_ = s[i, j] * np.eye(np.size(s, 1)) - t1 = np.zeros([np.size(s, 1), np.size(s, 1)]) - t1[j, :] = s[i, :] - tmp = tmp_ + t1 - if j == 0: - h_ = tmp[:, j::] - else: - h_ = np.hstack((h_, tmp[:, j::])) - jf[i, :, :] = np.hstack( - (np.zeros([np.size(s, 1), 1]), np.eye(np.size(s, 1)), h_) - ) - return fx, jf diff --git a/src/UQpy/surrogates/kriging/regression_models/__init__.py b/src/UQpy/surrogates/kriging/regression_models/__init__.py deleted file mode 100644 index e6da265b3..000000000 --- a/src/UQpy/surrogates/kriging/regression_models/__init__.py +++ /dev/null @@ -1,4 +0,0 @@ -from UQpy.surrogates.kriging.regression_models.baseclass import * -from UQpy.surrogates.kriging.regression_models.ConstantRegression import ConstantRegression -from UQpy.surrogates.kriging.regression_models.LinearRegression import LinearRegression -from UQpy.surrogates.kriging.regression_models.QuadraticRegression import QuadraticRegression diff --git a/src/UQpy/surrogates/kriging/regression_models/baseclass/Regression.py b/src/UQpy/surrogates/kriging/regression_models/baseclass/Regression.py deleted file mode 100644 index 3d435e51d..000000000 --- a/src/UQpy/surrogates/kriging/regression_models/baseclass/Regression.py +++ /dev/null @@ -1,14 +0,0 @@ -from abc import ABC, abstractmethod - - -class Regression(ABC): - """ - Abstract base class of all Regressions. Serves as a template for creating new Kriging regression - functions. - """ - @abstractmethod - def r(self, s): - """ - Abstract method that needs to be implemented by the user when creating a new Regression function. - """ - pass diff --git a/src/UQpy/surrogates/kriging/regression_models/baseclass/__init__.py b/src/UQpy/surrogates/kriging/regression_models/baseclass/__init__.py deleted file mode 100644 index 004e95e23..000000000 --- a/src/UQpy/surrogates/kriging/regression_models/baseclass/__init__.py +++ /dev/null @@ -1 +0,0 @@ -from UQpy.surrogates.kriging.regression_models.baseclass.Regression import Regression diff --git a/tests/unit_tests/sampling/test_adaptive_kriging.py b/tests/unit_tests/sampling/test_adaptive_kriging.py index ef932fb7f..ab0c28714 100644 --- a/tests/unit_tests/sampling/test_adaptive_kriging.py +++ b/tests/unit_tests/sampling/test_adaptive_kriging.py @@ -3,7 +3,6 @@ from UQpy.run_model.model_execution.PythonModel import PythonModel from UQpy.utilities.MinimizeOptimizer import MinimizeOptimizer -from UQpy.surrogates.kriging.Kriging import Kriging from UQpy.sampling import MonteCarloSampling, AdaptiveKriging from UQpy.run_model.RunModel import RunModel from UQpy.distributions.collection import Normal diff --git a/tests/unit_tests/sampling/test_refined_stratified.py b/tests/unit_tests/sampling/test_refined_stratified.py index 15a3ccd12..6915add90 100644 --- a/tests/unit_tests/sampling/test_refined_stratified.py +++ b/tests/unit_tests/sampling/test_refined_stratified.py @@ -9,7 +9,6 @@ from UQpy.sampling.stratified_sampling.refinement.RandomRefinement import * from UQpy.sampling.stratified_sampling.strata.VoronoiStrata import * from UQpy.run_model.RunModel import * -from UQpy.surrogates.kriging.Kriging import Kriging def test_rss_simple_rectangular(): diff --git a/tests/unit_tests/surrogates/test_kriging.py b/tests/unit_tests/surrogates/test_kriging.py deleted file mode 100644 index d25084b02..000000000 --- a/tests/unit_tests/surrogates/test_kriging.py +++ /dev/null @@ -1,198 +0,0 @@ -import pytest -from UQpy.utilities.MinimizeOptimizer import MinimizeOptimizer -from beartype.roar import BeartypeCallHintPepParamException - -from UQpy.surrogates.kriging.Kriging import Kriging -import numpy as np -from UQpy.surrogates.kriging.regression_models import LinearRegression, ConstantRegression -from UQpy.surrogates.kriging.correlation_models import GaussianCorrelation - - -samples = np.linspace(0, 5, 20).reshape(-1, 1) -values = np.cos(samples) -optimizer = MinimizeOptimizer(method="L-BFGS-B") -krig = Kriging(regression_model=LinearRegression(), correlation_model=GaussianCorrelation(), optimizer=optimizer, - correlation_model_parameters=[0.14], optimize=False, random_state=1) -krig.fit(samples=samples, values=values, correlation_model_parameters=[0.3]) - -optimizer = MinimizeOptimizer(method="L-BFGS-B") -krig2 = Kriging(regression_model=ConstantRegression(), correlation_model=GaussianCorrelation(), optimizer=optimizer, - correlation_model_parameters=[0.3], bounds=[[0.01, 5]], - optimize=False, optimizations_number=100, normalize=False, - random_state=2) -krig2.fit(samples=samples, values=values) - - -# Using the in-built linear regression model as a function -linear_regression_model = Kriging(regression_model=LinearRegression(), correlation_model=GaussianCorrelation(), optimizer=optimizer, - correlation_model_parameters=[1]).regression_model -optimizer = MinimizeOptimizer(method="L-BFGS-B") -gaussian_corrleation_model = Kriging(regression_model=LinearRegression(), correlation_model=GaussianCorrelation(), optimizer=optimizer, - correlation_model_parameters=[1]).correlation_model - -optimizer = MinimizeOptimizer(method="L-BFGS-B") -krig3 = Kriging(regression_model=linear_regression_model, correlation_model=gaussian_corrleation_model, - optimizer=optimizer, - correlation_model_parameters=[1], optimize=False, normalize=False, random_state=0) -krig3.fit(samples=samples, values=values) - - -def test_predict(): - prediction = np.round(krig.predict([[1], [np.pi/2], [np.pi]], True), 3) - expected_prediction = np.array([[0.54, 0., -1.], [0., 0., 0.]]) - assert (expected_prediction == prediction).all() - - -def test_predict1(): - prediction = np.round(krig2.predict([[1], [2*np.pi], [np.pi]], True), 3) - expected_prediction = np.array([[0.54, 1.009, -1.], [0., 0.031, 0.]]) - assert (expected_prediction == prediction).all() - - -def test_predict2(): - prediction = np.round(krig3.predict([[1], [np.pi/2], [np.pi]]), 3) - expected_prediction = np.array([[0.54, -0., -1.]]) - assert (expected_prediction == prediction).all() - - -def test_jacobian(): - jacobian = np.round(krig.jacobian([[np.pi], [np.pi/2]]), 3) - expected_jacobian = np.array([-0., -1.]) - assert (expected_jacobian == jacobian).all() - - -def test_jacobian1(): - jacobian = np.round(krig3.jacobian([[np.pi], [np.pi/2]]), 3) - expected_jacobian = np.array([0., -1.]) - assert (expected_jacobian == jacobian).all() - - -def test_regression_models(): - from UQpy.surrogates.kriging.regression_models import ConstantRegression, LinearRegression, QuadraticRegression - krig.regression_model = ConstantRegression() - tmp = krig.regression_model.r([[0], [1]]) - tmp_test1 = (tmp[0] == np.array([[1.], [1.]])).all() and (tmp[1] == np.array([[[0.]], [[0.]]])).all() - - krig.regression_model = LinearRegression() - tmp = krig.regression_model.r([[0], [1]]) - tmp_test2 = (tmp[0] == (np.array([[1., 0.], [1., 1.]]))).all() and \ - (tmp[1] == np.array([[[0., 1.]], [[0., 1.]]])).all() - - krig.regression_model = QuadraticRegression() - tmp = krig.regression_model.r([[-1, 1], [2, -0.5]]) - tmp_test3 = (tmp[0] == np.array([[1., -1., 1., 1., -1., 1.], [1., 2., -0.5, 4., -1., 0.25]])).all() and \ - (tmp[1] == np.array([[[0., 1., 0., -2., 1., 0.], - [0., 0., 1., 0., -1., 2.]], - [[0., 1., 0., 4., -0.5, 0.], - [0., 0., 1., 0., 2., -1.]]])).all() - - assert tmp_test1 and tmp_test2 and tmp_test3 - - -def test_correlation_models(): - from UQpy.surrogates.kriging.correlation_models import ExponentialCorrelation, LinearCorrelation, SphericalCorrelation, CubicCorrelation, SplineCorrelation - krig.correlation_model = ExponentialCorrelation() - rx_exponential = (np.round(krig.correlation_model.c([[0], [1], [2]], [[2]], np.array([1])), 3) == - np.array([[0.135], [0.368], [1.]])).all() - drdt_exponential = (np.round(krig.correlation_model.c([[0], [1], [2]], [[2]], np.array([1]), dt=True)[1], 3) == - np.array([[[-0.271]], [[-0.368]], [[0.]]])).all() - drdx_exponential = (np.round(krig.correlation_model.c([[0], [1], [2]], [[2]], np.array([1]), dx=True)[1], 3) == - np.array([[[0.135]], [[0.368]], [[0.]]])).all() - expon = rx_exponential and drdt_exponential and drdx_exponential - - krig.correlation_model = LinearCorrelation() - rx_linear = (np.round(krig.correlation_model.c([[0], [1], [2]], [[2]], np.array([1])), 3) == - np.array([[0.], [0.], [1.]])).all() - drdt_linear = (np.round(krig.correlation_model.c([[0.4], [0.5], [0.6]], [[0.5]], np.array([1]), dt=True)[1], 3) == - np.array([[[-0.1]], [[-0.]], [[-0.1]]])).all() - drdx_linear = (np.round(krig.correlation_model.c([[0.4], [0.5], [0.6]], [[0.5]], np.array([1]), dx=True)[1], 3) == - np.array([[[1.]], [[-0.]], [[-1.]]])).all() - linear = rx_linear and drdt_linear and drdx_linear - - krig.correlation_model = SphericalCorrelation() - rx_spherical = (np.round(krig.correlation_model.c([[0], [1], [2]], [[2]], np.array([1])), 3) == - np.array([[0.], [0.], [1.]])).all() - drdt_spherical = (np.round(krig.correlation_model.c([[0.4], [0.5], [0.6]], [[0.5]], np.array([1]), dt=True)[1], 3) - == np.array([[[-0.148]], [[-0.]], [[-0.148]]])).all() - drdx_spherical = (np.round(krig.correlation_model.c([[0.4], [0.5], [0.6]], [[0.5]], np.array([1]), dx=True)[1], 3) - == np.array([[[1.485]], [[-0.]], [[-1.485]]])).all() - spherical = rx_spherical and drdt_spherical and drdx_spherical - - krig.correlation_model = CubicCorrelation() - rx_cubic = (np.round(krig.correlation_model.c([[0.2], [0.5], [1]], [[0.5]], np.array([1])), 3) == - np.array([[0.784], [1.], [0.5]])).all() - drdt_cubic = (np.round(krig.correlation_model.c([[0.4], [0.5], [0.6]], [[0.5]], np.array([1]), dt=True)[1], 3) == - np.array([[[-0.054]], [[0.]], [[-0.054]]])).all() - drdx_cubic = (np.round(krig.correlation_model.c([[0.4], [0.5], [0.6]], [[0.5]], np.array([1]), dx=True)[1], 3) == - np.array([[[0.54]], [[0.]], [[-0.54]]])).all() - cubic = rx_cubic and drdt_cubic and drdx_cubic - - krig.correlation_model = SplineCorrelation() - rx_spline = (np.round(krig.correlation_model.c([[0.2], [0.5], [1]], [[0.5]], np.array([1])), 3) == - np.array([[0.429], [1.], [0.156]])).all() - drdt_spline = (np.round(krig.correlation_model.c([[0.4], [0.5], [0.6]], [[0.5]], np.array([1]), dt=True)[1], 3) == - np.array([[[-0.21]], [[0.]], [[-0.21]]])).all() - drdx_spline = (np.round(krig.correlation_model.c([[0.4], [0.5], [0.6]], [[0.5]], np.array([1]), dx=True)[1], 3) == - np.array([[[2.1]], [[0.]], [[-2.1]]])).all() - spline = rx_spline and drdt_spline and drdx_spline - - assert expon and linear and spherical and cubic and spline - - -def test_wrong_regression_model(): - """ - Raises an error if reg_model is not callable or a string of an in-built model. - """ - with pytest.raises(BeartypeCallHintPepParamException): - Kriging(regression_model='A', correlation_model=GaussianCorrelation(), correlation_model_parameters=[1]) - - -def test_wrong_correlation_model(): - """ - Raises an error if corr_model is not callable or a string of an in-built model. - """ - with pytest.raises(BeartypeCallHintPepParamException): - Kriging(regression_model=LinearRegression(), correlation_model='A', correlation_model_parameters=[1]) - - -def test_missing_correlation_model_parameters(): - """ - Raises an error if corr_model_params is not defined. - """ - with pytest.raises(TypeError): - Kriging(regression_model=LinearRegression(), correlation_model=GaussianCorrelation(), bounds=[[0.01, 5]], - optimizations_number=100, random_state=1) - - -def test_optimizer(): - """ - Raises an error if corr_model_params is not defined. - """ - with pytest.raises(ValueError): - Kriging(regression_model=LinearRegression(), correlation_model=GaussianCorrelation(), - correlation_model_parameters=[1], optimizer='A') - - -def test_random_state(): - """ - Raises an error if type of random_state is not correct. - """ - with pytest.raises(BeartypeCallHintPepParamException): - Kriging(regression_model=LinearRegression(), correlation_model=GaussianCorrelation(), - correlation_model_parameters=[1], random_state='A') - - -def test_log_likelihood(): - prediction = np.round(krig3.log_likelihood(np.array([1.5]), - krig3.correlation_model, np.array([[1], [2]]), - np.array([[1], [1]]), np.array([[1], [2]]), return_grad=False), 3) - expected_prediction = 1.679 - assert (expected_prediction == prediction).all() - - -def test_log_likelihood_derivative(): - prediction = np.round(krig3.log_likelihood(np.array([1.5]), krig3.correlation_model, np.array([[1], [2]]), - np.array([[1], [1]]), np.array([[1], [2]]), return_grad=True)[1], 3) - expected_prediction = np.array([-0.235]) - assert (expected_prediction == prediction).all() - From 08b5f33f821a270361aeb8666be4dc7b95572dd8 Mon Sep 17 00:00:00 2001 From: Dimitris Tsapetis Date: Mon, 21 Nov 2022 12:25:17 -0500 Subject: [PATCH 2/6] Fixes AKMCS tests --- src/UQpy/utilities/MinimizeOptimizer.py | 4 +- .../sampling/test_adaptive_kriging.py | 152 ++++++++---------- 2 files changed, 72 insertions(+), 84 deletions(-) diff --git a/src/UQpy/utilities/MinimizeOptimizer.py b/src/UQpy/utilities/MinimizeOptimizer.py index aa76a99ec..d15904b01 100644 --- a/src/UQpy/utilities/MinimizeOptimizer.py +++ b/src/UQpy/utilities/MinimizeOptimizer.py @@ -24,11 +24,11 @@ def optimize(self, function, initial_guess, args=(), jac=False): return minimize(function, initial_guess, args=args, method=self.method, bounds=self._bounds, constraints=self.constraints, jac=jac, - options={'disp': True, 'maxiter': 10000, 'catol': 0.002}) + options={'disp': False, 'maxiter': 10000, 'catol': 0.002}) else: return minimize(function, initial_guess, args=args, method=self.method, bounds=self._bounds, jac=jac, - options={'disp': True, 'maxiter': 10000, 'catol': 0.002}) + options={'disp': False, 'maxiter': 10000, 'catol': 0.002}) def apply_constraints(self, constraints): if self.method.lower() in ['cobyla', 'slsqp', 'trust-constr']: diff --git a/tests/unit_tests/sampling/test_adaptive_kriging.py b/tests/unit_tests/sampling/test_adaptive_kriging.py index ab0c28714..73070dba5 100644 --- a/tests/unit_tests/sampling/test_adaptive_kriging.py +++ b/tests/unit_tests/sampling/test_adaptive_kriging.py @@ -1,5 +1,6 @@ import pytest +from UQpy import GaussianProcessRegression, RBF, LinearRegression from UQpy.run_model.model_execution.PythonModel import PythonModel from UQpy.utilities.MinimizeOptimizer import MinimizeOptimizer @@ -7,172 +8,159 @@ from UQpy.run_model.RunModel import RunModel from UQpy.distributions.collection import Normal from UQpy.sampling.adaptive_kriging_functions import * -import shutil def test_akmcs_weighted_u(): - from UQpy.surrogates.kriging.regression_models.LinearRegression import LinearRegression - from UQpy.surrogates.kriging.correlation_models.ExponentialCorrelation import ExponentialCorrelation - marginals = [Normal(loc=0., scale=4.), Normal(loc=0., scale=4.)] x = MonteCarloSampling(distributions=marginals, nsamples=20, random_state=0) model = PythonModel(model_script='series.py', model_object_name="series") rmodel = RunModel(model=model) - regression_model = LinearRegression() - correlation_model = ExponentialCorrelation() - K = Kriging(regression_model=regression_model, correlation_model=correlation_model, - optimizer=MinimizeOptimizer('l-bfgs-b'), - optimizations_number=10, correlation_model_parameters=[1, 1], random_state=1) + + kernel1 = RBF() + bounds_1 = [[10 ** (-4), 10 ** 3], [10 ** (-3), 10 ** 2], [10 ** (-3), 10 ** 2]] + optimizer1 = MinimizeOptimizer(method='L-BFGS-B', bounds=bounds_1) + gpr = GaussianProcessRegression(kernel=kernel1, hyperparameters=[1, 10 ** (-3), 10 ** (-2)], optimizer=optimizer1, + optimizations_number=10, noise=False, regression_model=LinearRegression(), + random_state=1) # OPTIONS: 'U', 'EFF', 'Weighted-U' learning_function = WeightedUFunction(weighted_u_stop=2) - a = AdaptiveKriging(distributions=marginals, runmodel_object=rmodel, surrogate=K, - learning_nsamples=10**3, n_add=1, learning_function=learning_function, + a = AdaptiveKriging(distributions=marginals, runmodel_object=rmodel, surrogate=gpr, + learning_nsamples=10 ** 3, n_add=1, learning_function=learning_function, random_state=2) a.run(nsamples=25, samples=x.samples) - assert a.samples[23, 0] == 1.083176685073489 - assert a.samples[20, 1] == 0.20293978126855253 - + assert a.samples[23, 0] == -0.48297825309989356 + assert a.samples[20, 1] == 0.39006110248010434 def test_akmcs_u(): - from UQpy.surrogates.kriging.regression_models.LinearRegression import LinearRegression - from UQpy.surrogates.kriging.correlation_models.ExponentialCorrelation import ExponentialCorrelation marginals = [Normal(loc=0., scale=4.), Normal(loc=0., scale=4.)] x = MonteCarloSampling(distributions=marginals, nsamples=20, random_state=1) model = PythonModel(model_script='series.py', model_object_name="series") rmodel = RunModel(model=model) - regression_model = LinearRegression() - correlation_model = ExponentialCorrelation() - K = Kriging(regression_model=regression_model, correlation_model=correlation_model, - optimizer=MinimizeOptimizer('l-bfgs-b'), - optimizations_number=10, correlation_model_parameters=[1, 1], random_state=0) + kernel1 = RBF() + bounds_1 = [[10 ** (-4), 10 ** 3], [10 ** (-3), 10 ** 2], [10 ** (-3), 10 ** 2]] + optimizer1 = MinimizeOptimizer(method='L-BFGS-B', bounds=bounds_1) + gpr = GaussianProcessRegression(kernel=kernel1, hyperparameters=[1, 10 ** (-3), 10 ** (-2)], optimizer=optimizer1, + optimizations_number=10, noise=False, regression_model=LinearRegression(), + random_state=0) # OPTIONS: 'U', 'EFF', 'Weighted-U' learning_function = UFunction(u_stop=2) - a = AdaptiveKriging(distributions=marginals, runmodel_object=rmodel, surrogate=K, - learning_nsamples=10**3, n_add=1, learning_function=learning_function, + a = AdaptiveKriging(distributions=marginals, runmodel_object=rmodel, surrogate=gpr, + learning_nsamples=10 ** 3, n_add=1, learning_function=learning_function, random_state=2) a.run(nsamples=25, samples=x.samples) - assert a.samples[23, 0] == -4.141979058326188 - assert a.samples[20, 1] == -1.6476534435429009 - + assert a.samples[23, 0] == -3.781937137406927 + assert a.samples[20, 1] == 0.17610325620498946 def test_akmcs_expected_feasibility(): - from UQpy.surrogates.kriging.regression_models.LinearRegression import LinearRegression - from UQpy.surrogates.kriging.correlation_models.ExponentialCorrelation import ExponentialCorrelation marginals = [Normal(loc=0., scale=4.), Normal(loc=0., scale=4.)] x = MonteCarloSampling(distributions=marginals, nsamples=20, random_state=1) model = PythonModel(model_script='series.py', model_object_name="series") rmodel = RunModel(model=model) - regression_model = LinearRegression() - correlation_model = ExponentialCorrelation() - K = Kriging(regression_model=regression_model, correlation_model=correlation_model, - optimizations_number=10, correlation_model_parameters=[1, 1], - optimizer=MinimizeOptimizer('l-bfgs-b'),) + kernel1 = RBF() + bounds_1 = [[10 ** (-4), 10 ** 3], [10 ** (-3), 10 ** 2], [10 ** (-3), 10 ** 2]] + optimizer1 = MinimizeOptimizer(method='L-BFGS-B', bounds=bounds_1) + gpr = GaussianProcessRegression(kernel=kernel1, hyperparameters=[1, 10 ** (-3), 10 ** (-2)], optimizer=optimizer1, + optimizations_number=20, noise=False, regression_model=LinearRegression(), + random_state=0) # OPTIONS: 'U', 'EFF', 'Weighted-U' learning_function = ExpectedFeasibility(eff_a=0, eff_epsilon=2, eff_stop=0.001) - a = AdaptiveKriging(distributions=marginals, runmodel_object=rmodel, surrogate=K, - learning_nsamples=10**3, n_add=1, learning_function=learning_function, + a = AdaptiveKriging(distributions=marginals, runmodel_object=rmodel, surrogate=gpr, + learning_nsamples=10 ** 3, n_add=1, learning_function=learning_function, random_state=2) a.run(nsamples=25, samples=x.samples) - assert a.samples[23, 0] == 1.366058523912817 - assert a.samples[20, 1] == -12.914668932772358 - + assert a.samples[23, 0] == 5.423754197908594 + assert a.samples[20, 1] == 2.0355505295053384 def test_akmcs_expected_improvement(): - from UQpy.surrogates.kriging.regression_models.LinearRegression import LinearRegression - from UQpy.surrogates.kriging.correlation_models.ExponentialCorrelation import ExponentialCorrelation marginals = [Normal(loc=0., scale=4.), Normal(loc=0., scale=4.)] x = MonteCarloSampling(distributions=marginals, nsamples=20, random_state=1) model = PythonModel(model_script='series.py', model_object_name="series") rmodel = RunModel(model=model) - regression_model = LinearRegression() - correlation_model = ExponentialCorrelation() - K = Kriging(regression_model=regression_model, correlation_model=correlation_model, - optimizations_number=10, correlation_model_parameters=[1, 1], - optimizer=MinimizeOptimizer('l-bfgs-b'),) + kernel1 = RBF() + bounds_1 = [[10 ** (-4), 10 ** 3], [10 ** (-3), 10 ** 2], [10 ** (-3), 10 ** 2]] + optimizer1 = MinimizeOptimizer(method='L-BFGS-B', bounds=bounds_1) + gpr = GaussianProcessRegression(kernel=kernel1, hyperparameters=[1, 10 ** (-3), 10 ** (-2)], optimizer=optimizer1, + optimizations_number=50, noise=False, regression_model=LinearRegression(), + random_state=0) # OPTIONS: 'U', 'EFF', 'Weighted-U' learning_function = ExpectedImprovement() - a = AdaptiveKriging(distributions=marginals, runmodel_object=rmodel, surrogate=K, - learning_nsamples=10**3, n_add=1, learning_function=learning_function, + a = AdaptiveKriging(distributions=marginals, runmodel_object=rmodel, surrogate=gpr, + learning_nsamples=10 ** 3, n_add=1, learning_function=learning_function, random_state=2) a.run(nsamples=25, samples=x.samples) - assert a.samples[23, 0] == 4.553078100499578 - assert a.samples[20, 1] == -3.508949564718469 - + assert a.samples[21, 0] == 6.878734574049913 + assert a.samples[20, 1] == -6.3410533857909215 def test_akmcs_expected_improvement_global_fit(): - from UQpy.surrogates.kriging.regression_models.LinearRegression import LinearRegression - from UQpy.surrogates.kriging.correlation_models.ExponentialCorrelation import ExponentialCorrelation marginals = [Normal(loc=0., scale=4.), Normal(loc=0., scale=4.)] x = MonteCarloSampling(distributions=marginals, nsamples=20, random_state=1) model = PythonModel(model_script='series.py', model_object_name="series") rmodel = RunModel(model=model) - regression_model = LinearRegression() - correlation_model = ExponentialCorrelation() - K = Kriging(regression_model=regression_model, correlation_model=correlation_model, - optimizations_number=10, correlation_model_parameters=[1, 1], - optimizer=MinimizeOptimizer('l-bfgs-b'),) + kernel1 = RBF() + bounds_1 = [[10 ** (-4), 10 ** 3], [10 ** (-3), 10 ** 2], [10 ** (-3), 10 ** 2]] + optimizer1 = MinimizeOptimizer(method='L-BFGS-B', bounds=bounds_1) + gpr = GaussianProcessRegression(kernel=kernel1, hyperparameters=[1, 10 ** (-3), 10 ** (-2)], optimizer=optimizer1, + optimizations_number=50, noise=False, regression_model=LinearRegression(), + random_state=0) # OPTIONS: 'U', 'EFF', 'Weighted-U' learning_function = ExpectedImprovementGlobalFit() - a = AdaptiveKriging(distributions=marginals, runmodel_object=rmodel, surrogate=K, - learning_nsamples=10**3, n_add=1, learning_function=learning_function, + a = AdaptiveKriging(distributions=marginals, runmodel_object=rmodel, surrogate=gpr, + learning_nsamples=10 ** 3, n_add=1, learning_function=learning_function, random_state=2) a.run(nsamples=25, samples=x.samples) - assert a.samples[23, 0] == 11.939859785098493 - assert a.samples[20, 1] == -8.429899469300118 + assert a.samples[23, 0] == -10.24267076486663 + assert a.samples[20, 1] == -11.419510366469687 def test_akmcs_samples_error(): - from UQpy.surrogates.kriging.regression_models.LinearRegression import LinearRegression - from UQpy.surrogates.kriging.correlation_models.ExponentialCorrelation import ExponentialCorrelation marginals = [Normal(loc=0., scale=4.), Normal(loc=0., scale=4.)] x = MonteCarloSampling(distributions=marginals, nsamples=20, random_state=0) model = PythonModel(model_script='series.py', model_object_name="series") rmodel = RunModel(model=model) - regression_model = LinearRegression() - correlation_model = ExponentialCorrelation() - K = Kriging(regression_model=regression_model, correlation_model=correlation_model, - optimizer=MinimizeOptimizer('l-bfgs-b'), - optimizations_number=10, correlation_model_parameters=[1, 1], random_state=1) + kernel1 = RBF() + bounds_1 = [[10 ** (-4), 10 ** 3], [10 ** (-3), 10 ** 2], [10 ** (-3), 10 ** 2]] + optimizer1 = MinimizeOptimizer(method='L-BFGS-B', bounds=bounds_1) + gpr = GaussianProcessRegression(kernel=kernel1, hyperparameters=[1, 10 ** (-3), 10 ** (-2)], optimizer=optimizer1, + optimizations_number=50, noise=False, regression_model=LinearRegression(), + random_state=0) # OPTIONS: 'U', 'EFF', 'Weighted-U' learning_function = WeightedUFunction(weighted_u_stop=2) with pytest.raises(NotImplementedError): - a = AdaptiveKriging(distributions=[Normal(loc=0., scale=4.)]*3, runmodel_object=rmodel, surrogate=K, - learning_nsamples=10**3, n_add=1, learning_function=learning_function, + a = AdaptiveKriging(distributions=[Normal(loc=0., scale=4.)] * 3, runmodel_object=rmodel, surrogate=gpr, + learning_nsamples=10 ** 3, n_add=1, learning_function=learning_function, random_state=2, samples=x.samples) def test_akmcs_u_run_from_init(): - from UQpy.surrogates.kriging.regression_models.LinearRegression import LinearRegression - from UQpy.surrogates.kriging.correlation_models.ExponentialCorrelation import ExponentialCorrelation - marginals = [Normal(loc=0., scale=4.), Normal(loc=0., scale=4.)] x = MonteCarloSampling(distributions=marginals, nsamples=20, random_state=1) model = PythonModel(model_script='series.py', model_object_name="series") rmodel = RunModel(model=model) - regression_model = LinearRegression() - correlation_model = ExponentialCorrelation() - K = Kriging(regression_model=regression_model, correlation_model=correlation_model, - optimizer=MinimizeOptimizer('l-bfgs-b'), - optimizations_number=10, correlation_model_parameters=[1, 1], random_state=0) + kernel1 = RBF() + bounds_1 = [[10 ** (-4), 10 ** 3], [10 ** (-3), 10 ** 2], [10 ** (-3), 10 ** 2]] + optimizer1 = MinimizeOptimizer(method='L-BFGS-B', bounds=bounds_1) + gpr = GaussianProcessRegression(kernel=kernel1, hyperparameters=[1, 10 ** (-3), 10 ** (-2)], optimizer=optimizer1, + optimizations_number=100, noise=False, regression_model=LinearRegression(), + random_state=0) # OPTIONS: 'U', 'EFF', 'Weighted-U' learning_function = UFunction(u_stop=2) - a = AdaptiveKriging(distributions=marginals, runmodel_object=rmodel, surrogate=K, - learning_nsamples=10**3, n_add=1, learning_function=learning_function, + a = AdaptiveKriging(distributions=marginals, runmodel_object=rmodel, surrogate=gpr, + learning_nsamples=10 ** 3, n_add=1, learning_function=learning_function, random_state=2, nsamples=25, samples=x.samples) - assert a.samples[23, 0] == -4.141979058326188 - assert a.samples[20, 1] == -1.6476534435429009 + assert a.samples[23, 0] == -3.781937137406927 + assert a.samples[20, 1] == 0.17610325620498946 From 2818b256a9edf049914a2a0943d1f2db52f42af6 Mon Sep 17 00:00:00 2001 From: Dimitris Tsapetis Date: Mon, 21 Nov 2022 14:08:51 -0500 Subject: [PATCH 3/6] Fixes RSS tests --- .../RefinedStratifiedSampling.py | 2 +- .../refinement/baseclass/Refinement.py | 2 +- .../sampling/test_adaptive_kriging.py | 4 +- .../sampling/test_refined_stratified.py | 45 ++++++++----------- 4 files changed, 22 insertions(+), 31 deletions(-) diff --git a/src/UQpy/sampling/stratified_sampling/RefinedStratifiedSampling.py b/src/UQpy/sampling/stratified_sampling/RefinedStratifiedSampling.py index 6a233ea1e..e891651c1 100644 --- a/src/UQpy/sampling/stratified_sampling/RefinedStratifiedSampling.py +++ b/src/UQpy/sampling/stratified_sampling/RefinedStratifiedSampling.py @@ -45,7 +45,7 @@ def __init__( self.random_state = random_state if isinstance(self.random_state, int): - self.random_state = np.random.RandomState(self.random_state) + self.random_state = np.random.default_rng(self.random_state) elif not isinstance(self.random_state, (type(None), np.random.RandomState)): raise TypeError('UQpy: random_state must be None, an int or an np.random.Generator object.') if self.random_state is None: diff --git a/src/UQpy/sampling/stratified_sampling/refinement/baseclass/Refinement.py b/src/UQpy/sampling/stratified_sampling/refinement/baseclass/Refinement.py index 9e3f37f78..0ca80e523 100644 --- a/src/UQpy/sampling/stratified_sampling/refinement/baseclass/Refinement.py +++ b/src/UQpy/sampling/stratified_sampling/refinement/baseclass/Refinement.py @@ -46,7 +46,7 @@ def finalize(self, samples, samples_per_iteration): def identify_bins(strata_metrics, points_to_add, random_state): bins2break = np.array([]) points_left = points_to_add - while (np.where(strata_metrics == strata_metrics.max())[0].shape[0] < points_left): + while np.where(strata_metrics == strata_metrics.max())[0].shape[0] < points_left: bin = np.where(strata_metrics == strata_metrics.max())[0] bins2break = np.hstack([bins2break, bin]) strata_metrics[bin] = 0 diff --git a/tests/unit_tests/sampling/test_adaptive_kriging.py b/tests/unit_tests/sampling/test_adaptive_kriging.py index 73070dba5..9e75f9003 100644 --- a/tests/unit_tests/sampling/test_adaptive_kriging.py +++ b/tests/unit_tests/sampling/test_adaptive_kriging.py @@ -52,7 +52,7 @@ def test_akmcs_u(): random_state=2) a.run(nsamples=25, samples=x.samples) - assert a.samples[23, 0] == -3.781937137406927 + assert a.samples[23, 0] == 4.027342825480197 assert a.samples[20, 1] == 0.17610325620498946 @@ -75,7 +75,7 @@ def test_akmcs_expected_feasibility(): random_state=2) a.run(nsamples=25, samples=x.samples) - assert a.samples[23, 0] == 5.423754197908594 + assert a.samples[23, 0] == 4.553078100499578 assert a.samples[20, 1] == 2.0355505295053384 diff --git a/tests/unit_tests/sampling/test_refined_stratified.py b/tests/unit_tests/sampling/test_refined_stratified.py index 6915add90..7389d7101 100644 --- a/tests/unit_tests/sampling/test_refined_stratified.py +++ b/tests/unit_tests/sampling/test_refined_stratified.py @@ -1,6 +1,7 @@ import pytest from beartype.roar import BeartypeCallHintPepParamException +from UQpy import GaussianProcessRegression, RBF, LinearRegression from UQpy.run_model.model_execution.PythonModel import PythonModel from UQpy.utilities.MinimizeOptimizer import MinimizeOptimizer from UQpy.sampling.stratified_sampling.refinement.GradientEnhancedRefinement import GradientEnhancedRefinement @@ -71,31 +72,22 @@ def test_rect_gerss(): x = TrueStratifiedSampling(distributions=marginals, strata_object=strata, nsamples_per_stratum=1) model = PythonModel(model_script='python_model_function.py', model_object_name="y_func") rmodel = RunModel(model=model) - from UQpy.surrogates.kriging.regression_models import LinearRegression - from UQpy.surrogates.kriging.correlation_models import ExponentialCorrelation - K = Kriging(regression_model=LinearRegression(), correlation_model=ExponentialCorrelation(), optimizations_number=20, random_state=0, - correlation_model_parameters=[1, 1], optimizer=MinimizeOptimizer('l-bfgs-b'), ) - K.fit(samples=x.samples, values=rmodel.qoi_list) + kernel1 = RBF() + bounds_1 = [[10 ** (-4), 10 ** 3], [10 ** (-3), 10 ** 2], [10 ** (-3), 10 ** 2]] + optimizer1 = MinimizeOptimizer(method='L-BFGS-B', bounds=bounds_1) + gpr = GaussianProcessRegression(kernel=kernel1, hyperparameters=[1, 10 ** (-3), 10 ** (-2)], optimizer=optimizer1, + optimizations_number=100, noise=False, regression_model=LinearRegression(), + random_state=0) + # gpr.fit(samples=x.samples, values=rmodel.qoi_list) refinement = GradientEnhancedRefinement(strata=x.strata_object, runmodel_object=rmodel, - surrogate=K, nearest_points_number=4) + surrogate=gpr, nearest_points_number=4) z = RefinedStratifiedSampling(stratified_sampling=x, random_state=2, refinement_algorithm=refinement) z.run(nsamples=6) assert np.allclose(z.samples, np.array([[0.417022, 0.36016225], [1.00011437, 0.15116629], [0.14675589, 0.5461693], [1.18626021, 0.67278036], - [1.51296312, 0.77483124], [0.74237455, 0.66026822]])) - # assert np.allclose(z.samples, np.array([[0.417022, 0.36016225], [1.00011437, 0.15116629], - # [0.14675589, 0.5461693], [1.18626021, 0.67278036], - # [1.59254104, 0.96577043], [1.97386531, 0.24237455]])) - # assert np.allclose(z.samples, np.array([[0.417022, 0.36016225], [1.00011437, 0.15116629], - # [0.14675589, 0.5461693], [1.18626021, 0.67278036], - # [1.59254104, 0.96577043], [1.7176612, 0.2101839]])) - # assert np.allclose(z.samplesU01, np.array([[0.208511, 0.36016225], [0.50005719, 0.15116629], - # [0.07337795, 0.5461693], [0.59313011, 0.67278036], - # [0.79627052, 0.96577043], [0.98693265, 0.24237455]])) - # assert np.allclose(z.samplesU01, np.array([[0.208511, 0.36016225], [0.50005719, 0.15116629], - # [0.07337795, 0.5461693], [0.59313011, 0.67278036], - # [0.79627052, 0.96577043], [0.8588306 , 0.2101839]])) + [1.64924557, 0.90711287], + [0.54595797, 0.30005026]])) def test_vor_rss(): @@ -123,19 +115,18 @@ def test_vor_gerss(): marginals = [Uniform(loc=0., scale=2.), Uniform(loc=0., scale=1.)] strata_vor = VoronoiStrata(seeds_number=4, dimension=2, random_state=10) x_vor = TrueStratifiedSampling(distributions=marginals, strata_object=strata_vor, nsamples_per_stratum=1, ) - from UQpy.surrogates.kriging.regression_models.LinearRegression import LinearRegression - from UQpy.surrogates.kriging.correlation_models.ExponentialCorrelation import ExponentialCorrelation model = PythonModel(model_script='python_model_function.py', model_object_name="y_func") rmodel = RunModel(model=model) - K_ = Kriging(regression_model=LinearRegression(), correlation_model=ExponentialCorrelation(), optimizations_number=20, - optimizer=MinimizeOptimizer('l-bfgs-b'), random_state=0, - correlation_model_parameters=[1, 1]) - - K_.fit(samples=x_vor.samples, values=rmodel.qoi_list) + kernel1 = RBF() + bounds_1 = [[10 ** (-4), 10 ** 3], [10 ** (-3), 10 ** 2], [10 ** (-3), 10 ** 2]] + optimizer1 = MinimizeOptimizer(method='L-BFGS-B', bounds=bounds_1) + gpr = GaussianProcessRegression(kernel=kernel1, hyperparameters=[1, 10 ** (-3), 10 ** (-2)], optimizer=optimizer1, + optimizations_number=100, noise=False, regression_model=LinearRegression(), + random_state=0) z_vor = RefinedStratifiedSampling(stratified_sampling=x_vor, nsamples=6, random_state=x_vor.random_state, refinement_algorithm=GradientEnhancedRefinement(strata=x_vor.strata_object, runmodel_object=rmodel, - surrogate=K_, + surrogate=gpr, nearest_points_number=4)) assert np.allclose(z_vor.samples, np.array([[1.78345908, 0.01640854], [1.46201137, 0.70862104], [0.4021338, 0.05290083], [0.1062376, 0.88958226], From 38c96f502c99785550489d261ac74a3d4bf19ddf Mon Sep 17 00:00:00 2001 From: Dimitris Tsapetis Date: Mon, 21 Nov 2022 15:27:53 -0500 Subject: [PATCH 4/6] Fixes RSS tests --- src/UQpy/sampling/SimplexSampling.py | 2 +- .../sampling/test_adaptive_kriging.py | 8 ++-- .../sampling/test_refined_stratified.py | 40 +++++++++---------- 3 files changed, 25 insertions(+), 25 deletions(-) diff --git a/src/UQpy/sampling/SimplexSampling.py b/src/UQpy/sampling/SimplexSampling.py index 2b7426635..613bba256 100644 --- a/src/UQpy/sampling/SimplexSampling.py +++ b/src/UQpy/sampling/SimplexSampling.py @@ -10,7 +10,7 @@ def __init__( self, nodes: Union[list, Numpy2DFloatArray] = None, nsamples: PositiveInteger = None, - random_state: RandomStateType = None, + random_state: Union[RandomStateType, np.random.Generator] = None, ): """ Generate uniform random samples inside an n-dimensional simplex. diff --git a/tests/unit_tests/sampling/test_adaptive_kriging.py b/tests/unit_tests/sampling/test_adaptive_kriging.py index 9e75f9003..4fcf3b1e6 100644 --- a/tests/unit_tests/sampling/test_adaptive_kriging.py +++ b/tests/unit_tests/sampling/test_adaptive_kriging.py @@ -43,7 +43,7 @@ def test_akmcs_u(): bounds_1 = [[10 ** (-4), 10 ** 3], [10 ** (-3), 10 ** 2], [10 ** (-3), 10 ** 2]] optimizer1 = MinimizeOptimizer(method='L-BFGS-B', bounds=bounds_1) gpr = GaussianProcessRegression(kernel=kernel1, hyperparameters=[1, 10 ** (-3), 10 ** (-2)], optimizer=optimizer1, - optimizations_number=10, noise=False, regression_model=LinearRegression(), + optimizations_number=100, noise=False, regression_model=LinearRegression(), random_state=0) # OPTIONS: 'U', 'EFF', 'Weighted-U' learning_function = UFunction(u_stop=2) @@ -52,7 +52,7 @@ def test_akmcs_u(): random_state=2) a.run(nsamples=25, samples=x.samples) - assert a.samples[23, 0] == 4.027342825480197 + assert a.samples[23, 0] == -3.781937137406927 assert a.samples[20, 1] == 0.17610325620498946 @@ -66,7 +66,7 @@ def test_akmcs_expected_feasibility(): bounds_1 = [[10 ** (-4), 10 ** 3], [10 ** (-3), 10 ** 2], [10 ** (-3), 10 ** 2]] optimizer1 = MinimizeOptimizer(method='L-BFGS-B', bounds=bounds_1) gpr = GaussianProcessRegression(kernel=kernel1, hyperparameters=[1, 10 ** (-3), 10 ** (-2)], optimizer=optimizer1, - optimizations_number=20, noise=False, regression_model=LinearRegression(), + optimizations_number=100, noise=False, regression_model=LinearRegression(), random_state=0) # OPTIONS: 'U', 'EFF', 'Weighted-U' learning_function = ExpectedFeasibility(eff_a=0, eff_epsilon=2, eff_stop=0.001) @@ -75,7 +75,7 @@ def test_akmcs_expected_feasibility(): random_state=2) a.run(nsamples=25, samples=x.samples) - assert a.samples[23, 0] == 4.553078100499578 + assert a.samples[23, 0] == 5.423754197908594 assert a.samples[20, 1] == 2.0355505295053384 diff --git a/tests/unit_tests/sampling/test_refined_stratified.py b/tests/unit_tests/sampling/test_refined_stratified.py index 7389d7101..18b617cc7 100644 --- a/tests/unit_tests/sampling/test_refined_stratified.py +++ b/tests/unit_tests/sampling/test_refined_stratified.py @@ -23,10 +23,10 @@ def test_rss_simple_rectangular(): samples_per_iteration=2, refinement_algorithm=algorithm, random_state=2) - assert y.samples[16, 0] == 0.06614276178462988 - assert y.samples[16, 1] == 0.7836449863362334 - assert y.samples[17, 0] == 0.1891972651582183 - assert y.samples[17, 1] == 0.2961099664117288 + assert y.samples[16, 0] == 0.22677821757428504 + assert y.samples[16, 1] == 0.2729789855337742 + assert y.samples[17, 0] == 0.07501256574570675 + assert y.samples[17, 1] == 0.9321401317029486 def test_rss_simple_voronoi(): @@ -40,10 +40,10 @@ def test_rss_simple_voronoi(): samples_per_iteration=2, refinement_algorithm=algorithm, random_state=2) - assert np.round(y.samples[16, 0], 6) == 0.363793 - assert np.round(y.samples[16, 1], 6) == 0.467625 - assert np.round(y.samples[17, 0], 6) == 0.424586 - assert np.round(y.samples[17, 1], 6) == 0.217301 + assert np.round(y.samples[16, 0], 6) == 0.324738 + assert np.round(y.samples[16, 1], 6) == 0.488029 + assert np.round(y.samples[17, 0], 6) == 0.349367 + assert np.round(y.samples[17, 1], 6) == 0.132426 def test_rect_rss(): @@ -57,10 +57,10 @@ def test_rect_rss(): refinement_algorithm=RandomRefinement(strata=strata)) assert np.allclose(y.samples, np.array([[0.417022, 0.36016225], [1.00011437, 0.15116629], [0.14675589, 0.5461693], [1.18626021, 0.67278036], - [0.77483124, 0.7176612], [1.7101839, 0.66516741]])) + [1.90711287, 0.04595797], [0.80005026, 0.86428026]])) assert np.allclose(np.array(y.samplesU01), np.array([[0.208511, 0.36016225], [0.50005719, 0.15116629], [0.07337795, 0.5461693], [0.59313011, 0.67278036], - [0.38741562, 0.7176612], [0.85509195, 0.66516741]])) + [0.95355644, 0.04595797], [0.40002513, 0.86428026]])) def test_rect_gerss(): @@ -123,17 +123,17 @@ def test_vor_gerss(): gpr = GaussianProcessRegression(kernel=kernel1, hyperparameters=[1, 10 ** (-3), 10 ** (-2)], optimizer=optimizer1, optimizations_number=100, noise=False, regression_model=LinearRegression(), random_state=0) - z_vor = RefinedStratifiedSampling(stratified_sampling=x_vor, nsamples=6, random_state=x_vor.random_state, + z_vor = RefinedStratifiedSampling(stratified_sampling=x_vor, nsamples=6, random_state=0, refinement_algorithm=GradientEnhancedRefinement(strata=x_vor.strata_object, runmodel_object=rmodel, surrogate=gpr, nearest_points_number=4)) assert np.allclose(z_vor.samples, np.array([[1.78345908, 0.01640854], [1.46201137, 0.70862104], [0.4021338, 0.05290083], [0.1062376, 0.88958226], - [0.61246269, 0.47160095], [1.16609055, 0.30832536]])) + [0.66730342, 0.46988084], [1.37411577, 0.39064685]])) assert np.allclose(z_vor.samplesU01, np.array([[0.89172954, 0.01640854], [0.73100569, 0.70862104], [0.2010669, 0.05290083], [0.0531188, 0.88958226], - [0.30623134, 0.47160095], [0.58304527, 0.30832536]])) + [0.33365171, 0.46988084], [0.68705789, 0.39064685]])) def test_rss_random_state(): @@ -155,17 +155,17 @@ def test_rss_runmodel_object(): marginals = [Uniform(loc=0., scale=2.), Uniform(loc=0., scale=1.)] strata = RectangularStrata(strata_number=[2, 2]) x = TrueStratifiedSampling(distributions=marginals, strata_object=strata, nsamples_per_stratum=1, random_state=1) - from UQpy.surrogates.kriging.regression_models import LinearRegression - from UQpy.surrogates.kriging.correlation_models import ExponentialCorrelation - - K = Kriging(regression_model=LinearRegression(), correlation_model=ExponentialCorrelation(), optimizations_number=20, - correlation_model_parameters=[1, 1], optimizer=MinimizeOptimizer('l-bfgs-b'), ) + kernel1 = RBF() + bounds_1 = [[10 ** (-4), 10 ** 3], [10 ** (-3), 10 ** 2], [10 ** (-3), 10 ** 2]] + optimizer1 = MinimizeOptimizer(method='L-BFGS-B', bounds=bounds_1) + gpr = GaussianProcessRegression(kernel=kernel1, hyperparameters=[1, 10 ** (-3), 10 ** (-2)], optimizer=optimizer1, + optimizations_number=100, noise=False, regression_model=LinearRegression(), + random_state=0) model = PythonModel(model_script='python_model_function.py', model_object_name="y_func") rmodel = RunModel(model=model) - K.fit(samples=x.samples, values=rmodel.qoi_list) with pytest.raises(BeartypeCallHintPepParamException): refinement = GradientEnhancedRefinement(strata=x.strata_object, runmodel_object='abc', - surrogate=K) + surrogate=gpr) RefinedStratifiedSampling(stratified_sampling=x, samples_number=6, samples_per_iteration=2, refinement_algorithm=refinement) From 1ed054000f8080d5e76d5abbbda8d944113709d7 Mon Sep 17 00:00:00 2001 From: Dimitris Tsapetis Date: Mon, 21 Nov 2022 16:20:43 -0500 Subject: [PATCH 5/6] Fixes RSS tests --- tests/unit_tests/sampling/test_refined_stratified.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/unit_tests/sampling/test_refined_stratified.py b/tests/unit_tests/sampling/test_refined_stratified.py index 18b617cc7..695bf8383 100644 --- a/tests/unit_tests/sampling/test_refined_stratified.py +++ b/tests/unit_tests/sampling/test_refined_stratified.py @@ -130,7 +130,7 @@ def test_vor_gerss(): nearest_points_number=4)) assert np.allclose(z_vor.samples, np.array([[1.78345908, 0.01640854], [1.46201137, 0.70862104], [0.4021338, 0.05290083], [0.1062376, 0.88958226], - [0.66730342, 0.46988084], [1.37411577, 0.39064685]])) + [0.66730342, 0.46988084], [1.5015904 , 0.97050966]])) assert np.allclose(z_vor.samplesU01, np.array([[0.89172954, 0.01640854], [0.73100569, 0.70862104], [0.2010669, 0.05290083], [0.0531188, 0.88958226], [0.33365171, 0.46988084], [0.68705789, 0.39064685]])) From 32bdd299f9f9f94df7d96d4e1187a8e962549c62 Mon Sep 17 00:00:00 2001 From: Dimitris Tsapetis Date: Mon, 21 Nov 2022 17:19:35 -0500 Subject: [PATCH 6/6] Fixes RSS tests --- tests/unit_tests/sampling/test_refined_stratified.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/unit_tests/sampling/test_refined_stratified.py b/tests/unit_tests/sampling/test_refined_stratified.py index 695bf8383..13afea15f 100644 --- a/tests/unit_tests/sampling/test_refined_stratified.py +++ b/tests/unit_tests/sampling/test_refined_stratified.py @@ -133,7 +133,7 @@ def test_vor_gerss(): [0.66730342, 0.46988084], [1.5015904 , 0.97050966]])) assert np.allclose(z_vor.samplesU01, np.array([[0.89172954, 0.01640854], [0.73100569, 0.70862104], [0.2010669, 0.05290083], [0.0531188, 0.88958226], - [0.33365171, 0.46988084], [0.68705789, 0.39064685]])) + [0.33365171, 0.46988084], [0.7507952 , 0.97050966]])) def test_rss_random_state():