Homec4science
Diffusion RationalROMPy (develop) Browse

RROMPy/rrompy/utilities/poly_fitting/gaussian_process/gaussian_process_interpolator.py
d6451a36626adevelop

Last Change
Raw File

gaussian_process_interpolator.py
View Options

# Copyright (C) 2018-2020 by the RROMPy authors
#
# This file is part of RROMPy.
#
# RROMPy is free software: you can redistribute it and/or modify
# it under the terms of the GNU Lesser General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# RROMPy is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Lesser General Public License for more details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with RROMPy. If not, see <http://www.gnu.org/licenses/>.
#
import warnings
from copy import deepcopy as copy
import numpy as np
from collections.abc import Iterable
from sklearn.exceptions import ConvergenceWarning
from sklearn.gaussian_process import GaussianProcessRegressor as GPR
from rrompy.utilities.base.types import List, ListAny, DictAny, Np1D, paramList
from rrompy.utilities.poly_fitting.interpolator import GenericInterpolator
from rrompy.utilities.poly_fitting.custom_fit import customFit
from .kernel import kernels
from rrompy.utilities.poly_fitting.polynomial import (polyval, polyvander,
polyfitname)
from rrompy.utilities.numerical.degree import mapTotalToTensor
from rrompy.utilities.exception_manager import RROMPyException, RROMPyAssert
from rrompy.parameter import parameterList, checkParameterList
__all__ = ['GaussianProcessInterpolator']
class GaussianProcessInterpolator(GenericInterpolator):
"""Function class with setup by polynomial-biased GPR."""
def __init__(self, other = None):
if other is None: return
self.coeffsGlobal = other.coeffsGlobal
self.GP = other.GP
self.npar = other.npar
self.polybasis = other.polybasis
@property
def shape(self):
if self.coeffsGlobal.ndim > self.npar:
return self.coeffsGlobal.shape[self.npar :]
else:
return 1
@property
def deg(self):
return [x - 1 for x in self.coeffsGlobal.shape[: self.npar]]
def __call__(self, mu:paramList, der : List[int] = None,
scl : Np1D = None, **gpvalPars):
if der is not None and np.max(der) > 0:
raise RROMPyException(("Derivatives of GPR not implemented. "
"Manual access to self.GP is required."))
if isinstance(mu, (parameterList,)):
RROMPyAssert(mu.shape[1], self.npar, "Number of parameters")
muEff = mu.data
else:
muEff = np.reshape(mu, (-1, self.npar))
gpval = self.GP.predict(muEff, **gpvalPars)
if not isinstance(gpval, (tuple,)): gpval = (gpval,)
if len(self.coeffsGlobal) > 0:
pval = polyval(mu, self.coeffsGlobal, self.polybasis, der, scl)
else:
pval = 0.
return (pval + gpval[0],) + gpval[1 :]
def __copy__(self):
return GaussianProcessInterpolator(self)
def __deepcopy__(self, memo):
other = GaussianProcessInterpolator()
(other.coeffsGlobal, other.GP, other.npar, other.polybasis) = copy(
(self.coeffsGlobal, self.GP,
self.npar, self.polybasis), memo)
return other
def setupByInterpolation(self, support:paramList, values:ListAny,
deg:int, polybasis : str = "MONOMIAL_GAUSSIAN",
verbose : bool = True, degreetype : str = "TOTAL",
vanderCoeffs : DictAny = {},
kernelCoeffs : DictAny = {},
fitCoeffs : DictAny = {},
trainCoeffs : DictAny = {}):
support = checkParameterList(support)
RROMPyAssert(len(support), len(values), "Number of support values")
degreetype = degreetype.upper().strip().replace(" ","")
if degreetype not in ["TENSOR", "TOTAL"]:
raise RROMPyException("Degree type not recognized.")
self.npar = support.shape[1]
if "reorder" in vanderCoeffs.keys():
support = support[vanderCoeffs["reorder"]]
self.polybasis, kernel = polybasis.split("_")
if not isinstance(deg, Iterable): deg = [deg] * self.npar
msg = None
if deg[0] is None or deg[0] < 0:
valuesEff = values
self.coeffsGlobal, wellcond = np.empty(0), True
else:
outDim = values.shape[1:]
valuesEff = np.reshape(values, (values.shape[0], -1))
vander = polyvander(support, deg, basis = self.polybasis,
degreetype = degreetype, **vanderCoeffs)
fitOut = customFit(vander, valuesEff, full = True, **fitCoeffs)
P, wellcond = fitOut[0], fitOut[1][1] == vander.shape[1]
if verbose:
degM = deg[0] if self.npar == 1 else deg
msg = ("Fitting {} samples with degree {} through {}... "
"Conditioning of LS system: {:.4e}.").format(
len(support), degM,
polyfitname(self.polybasis),
fitOut[1][2][0] / fitOut[1][2][-1])
valuesEff = (valuesEff - vander.dot(P)).reshape((-1,) + outDim)
if degreetype == "TENSOR":
self.coeffsGlobal = P.reshape(tuple([d + 1 for d in deg])
+ outDim)
else: #if degreetype == "TOTAL":
self.coeffsGlobal = mapTotalToTensor(
tuple([d + 1 for d in deg]) + outDim,
self.npar, P)
kernel = kernels[kernel.upper()](**kernelCoeffs)
with warnings.catch_warnings():
warnings.filterwarnings('ignore', category = ConvergenceWarning)
self.GP = GPR(kernel, **trainCoeffs).fit(support.data, valuesEff)
return wellcond, msg
c4science · Help