Initial designs for quad (#745)

Author: mmahsereci

docs/source/api/quad/solvers.initial_designs.rst

@@ -0,0 +1,5 @@
+Initial Designs
+---------------
+.. automodapi:: probnum.quad.solvers.initial_designs
+    :no-heading:
+    :headings: "*"

docs/source/api/quad/solvers.rst

@@ -19,3 +19,8 @@ probnum.quad.solvers
     :hidden:
 
     solvers.stopping_criteria
+
+.. toctree::
+    :hidden:
+
+    solvers.initial_designs

src/probnum/quad/_bayesquad.py

@@ -18,23 +18,19 @@
 from probnum.quad.typing import DomainLike, DomainType
 from probnum.randprocs.kernels import Kernel
 from probnum.randvars import Normal
-from probnum.typing import FloatLike, IntLike
+from probnum.typing import IntLike
 
 
 def bayesquad(
     fun: Callable,
-    input_dim: int,
+    input_dim: IntLike,
     kernel: Optional[Kernel] = None,
-    domain: Optional[DomainLike] = None,
     measure: Optional[IntegrationMeasure] = None,
+    domain: Optional[DomainLike] = None,
     policy: Optional[str] = "bmc",
-    scale_estimation: Optional[str] = "mle",
-    max_evals: Optional[IntLike] = None,
-    var_tol: Optional[FloatLike] = None,
-    rel_tol: Optional[FloatLike] = None,
-    batch_size: IntLike = 1,
+    initial_design: Optional[str] = None,
     rng: Optional[np.random.Generator] = None,
-    jitter: FloatLike = 1.0e-8,
+    options: Optional[dict] = None,
 ) -> Tuple[Normal, BQIterInfo]:
     r"""Infer the solution of the uni- or multivariate integral
     :math:`\int_\Omega f(x) d \mu(x)`
@@ -66,44 +62,56 @@ def bayesquad(
         Input dimension of the integration problem.
     kernel
         The kernel used for the GP model. Defaults to the ``ExpQuad`` kernel.
+    measure
+        The integration measure. Defaults to the Lebesgue measure on ``domain``.
     domain
         The integration domain. Contains lower and upper bound as scalar or
         ``np.ndarray``. Obsolete if ``measure`` is given.
-    measure
-        The integration measure. Defaults to the Lebesgue measure on ``domain``.
     policy
         Type of acquisition strategy to use. Defaults to 'bmc'. Options are
 
         ==========================  =======
          Bayesian Monte Carlo [2]_  ``bmc``
-        ==========================  =======
-
-        ==========================  =======
          van Der Corput points      ``vdc``
         ==========================  =======
 
-    scale_estimation
-        Estimation method to use to compute the scale parameter. Defaults to 'mle'.
-        Options are
-
-        ==============================  =======
-         Maximum likelihood estimation  ``mle``
-        ==============================  =======
-
-    max_evals
-        Maximum number of function evaluations.
-    var_tol
-        Tolerance on the variance of the integral.
-    rel_tol
-        Tolerance on consecutive updates of the integral mean.
-    batch_size
-        Number of new observations at each update. Defaults to 1.
+    initial_design
+        The type of initial design to use. If ``None`` is given, no initial design is
+        used. Options are
+
+        ==========================  =========
+         Samples from measure       ``mc``
+         Latin hypercube [3]_       ``latin``
+        ==========================  =========
+
     rng
-        Random number generator. Used by Bayesian Monte Carlo other random sampling
-        policies.
-    jitter
-        Non-negative jitter to numerically stabilise kernel matrix inversion.
-        Defaults to 1e-8.
+        The random number generator used for random methods.
+
+    options
+        A dictionary with the following optional solver settings
+
+            scale_estimation : Optional[str]
+                Estimation method to use to compute the scale parameter. Defaults
+                to 'mle'. Options are
+
+                ==============================  =======
+                 Maximum likelihood estimation  ``mle``
+                ==============================  =======
+
+            max_evals : Optional[IntLike]
+                Maximum number of function evaluations.
+            var_tol : Optional[FloatLike]
+                Tolerance on the variance of the integral.
+            rel_tol : Optional[FloatLike]
+                Tolerance on consecutive updates of the integral mean.
+            jitter : Optional[FloatLike]
+                Non-negative jitter to numerically stabilise kernel matrix
+                inversion. Defaults to 1e-8.
+            batch_size : Optional[IntLike]
+                Number of new observations at each update. Defaults to 1.
+            num_initial_design_nodes : Optional[IntLike]
+                The number of nodes created by the initial design. Defaults to
+                ``input_dim * 5`` if an initial design is given.
 
     Returns
     -------
@@ -119,7 +127,8 @@ def bayesquad(
 
     Warns
     -----
-    When ``domain`` is given but not used.
+    UserWarning
+        When ``domain`` is given but not used.
 
     Notes
     -----
@@ -138,6 +147,8 @@ def bayesquad(
         computation?, *Statistical Science 34.1*, 2019, 1-22, 2019
     .. [2] Rasmussen, C. E., and Z. Ghahramani, Bayesian Monte Carlo, *Advances in
         Neural Information Processing Systems*, 2003, 505-512.
+    .. [3] Mckay et al., A Comparison of Three Methods for Selecting Values of Input
+        Variables in the Analysis of Output from a Computer Code, *Technometrics*, 1979.
 
     Examples
     --------
@@ -162,12 +173,8 @@ def bayesquad(
         measure=measure,
         domain=domain,
         policy=policy,
-        scale_estimation=scale_estimation,
-        max_evals=max_evals,
-        var_tol=var_tol,
-        rel_tol=rel_tol,
-        batch_size=batch_size,
-        jitter=jitter,
+        initial_design=initial_design,
+        options=options,
     )
 
     # Integrate
@@ -182,10 +189,9 @@ def bayesquad_from_data(
     nodes: np.ndarray,
     fun_evals: np.ndarray,
     kernel: Optional[Kernel] = None,
-    domain: Optional[DomainLike] = None,
     measure: Optional[IntegrationMeasure] = None,
-    scale_estimation: Optional[str] = "mle",
-    jitter: FloatLike = 1.0e-8,
+    domain: Optional[DomainLike] = None,
+    options: Optional[dict] = None,
 ) -> Tuple[Normal, BQIterInfo]:
     r"""Infer the value of an integral from a given set of nodes and function
     evaluations.
@@ -199,16 +205,25 @@ def bayesquad_from_data(
         *shape=(n_eval,)* -- Function evaluations at ``nodes``.
     kernel
         The kernel used for the GP model. Defaults to the ``ExpQuad`` kernel.
+    measure
+        The integration measure. Defaults to the Lebesgue measure.
     domain
         The integration domain. Contains lower and upper bound as scalar or
         ``np.ndarray``. Obsolete if ``measure`` is given.
-    measure
-        The integration measure. Defaults to the Lebesgue measure.
-    scale_estimation
-        Estimation method to use to compute the scale parameter. Defaults to 'mle'.
-    jitter
-        Non-negative jitter to numerically stabilise kernel matrix inversion.
-        Defaults to 1e-8.
+    options
+        A dictionary with the following optional solver settings
+
+            scale_estimation : Optional[str]
+                Estimation method to use to compute the scale parameter. Defaults
+                to 'mle'. Options are
+
+                ==============================  =======
+                 Maximum likelihood estimation  ``mle``
+                ==============================  =======
+
+            jitter : Optional[FloatLike]
+                Non-negative jitter to numerically stabilise kernel matrix
+                inversion. Defaults to 1e-8.
 
     Returns
     -------
@@ -224,7 +239,8 @@ def bayesquad_from_data(
 
     Warns
     -----
-    When ``domain`` is given but not used.
+    UserWarning
+        When ``domain`` is given but not used.
 
     See Also
     --------
@@ -256,8 +272,8 @@ def bayesquad_from_data(
         measure=measure,
         domain=domain,
         policy=None,
-        scale_estimation=scale_estimation,
-        jitter=jitter,
+        initial_design=None,
+        options=options,
     )
 
     # Integrate
@@ -274,6 +290,8 @@ def _check_domain_measure_compatibility(
     measure: Optional[IntegrationMeasure],
 ) -> Tuple[int, Optional[DomainType], IntegrationMeasure]:
 
+    input_dim = int(input_dim)
+
     # Neither domain nor measure given
     if domain is None and measure is None:
         raise ValueError(

src/probnum/quad/solvers/__init__.py

@@ -1,6 +1,6 @@
 """Bayesian quadrature methods and their components."""
 
-from . import belief_updates, policies, stopping_criteria
+from . import belief_updates, initial_designs, policies, stopping_criteria
 from ._bayesian_quadrature import BayesianQuadrature
 from ._bq_state import BQIterInfo, BQState