Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
227 changes: 114 additions & 113 deletions active_subspaces/domains.py

Large diffs are not rendered by default.

27 changes: 14 additions & 13 deletions active_subspaces/gradients.py
Original file line number Diff line number Diff line change
@@ -1,32 +1,33 @@
"""Utilities for approximating gradients."""
import numpy as np
from numbers import Integral
from utils.misc import process_inputs
from utils.simrunners import SimulationRunner

def local_linear_gradients(X, f, p=None, weights=None):
"""Estimate a collection of gradients from input/output pairs.

Given a set of input/output pairs, choose subsets of neighboring points and
build a local linear model for each subset. The gradients of these local
linear models comprise estimates of sampled gradients.

Parameters
----------
X : ndarray
X : ndarray
M-by-m matrix that contains the m-dimensional inputs
f : ndarray
f : ndarray
M-by-1 matrix that contains scalar outputs
p : int, optional
how many nearest neighbors to use when constructing the local linear
how many nearest neighbors to use when constructing the local linear
model (default 1)
weights : ndarray, optional
M-by-1 matrix that contains the weights for each observation (default
M-by-1 matrix that contains the weights for each observation (default
None)

Returns
-------
df : ndarray
M-by-m matrix that contains estimated partial derivatives approximated
M-by-m matrix that contains estimated partial derivatives approximated
by the local linear models

Notes
Expand All @@ -39,12 +40,12 @@ def local_linear_gradients(X, f, p=None, weights=None):

if p is None:
p = int(np.minimum(np.floor(1.7*m), M))
elif not isinstance(p, int):
elif not isinstance(p, Integral):
raise TypeError('p must be an integer.')

if p < m+1 or p > M:
raise Exception('p must be between m+1 and M')

if weights is None:
weights = np.ones((M, 1)) / M

Expand All @@ -67,18 +68,18 @@ def finite_difference_gradients(X, fun, h=1e-6):

Parameters
----------
X : ndarray
M-by-m matrix that contains the points to estimate the gradients with
X : ndarray
M-by-m matrix that contains the points to estimate the gradients with
finite differences
fun : function
function that returns the simulation's quantity of interest given inputs
h : float, optional
h : float, optional
the finite difference step size (default 1e-6)

Returns
-------
df : ndarray
M-by-m matrix that contains estimated partial derivatives approximated
df : ndarray
M-by-m matrix that contains estimated partial derivatives approximated
by finite differences
"""
X, M, m = process_inputs(X)
Expand Down
95 changes: 47 additions & 48 deletions active_subspaces/integrals.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,7 @@

import numpy as np
import utils.quadrature as gq
from numbers import Integral
from utils.misc import conditional_expectations
from utils.designs import maximin_design
from utils.simrunners import SimulationRunner
Expand All @@ -12,30 +13,30 @@

def integrate(fun, avmap, N, NMC=10):
"""Approximate the integral of a function of m variables.

Parameters
----------
fun : function
an interface to the simulation that returns the quantity of interest
fun : function
an interface to the simulation that returns the quantity of interest
given inputs as an 1-by-m ndarray
avmap : ActiveVariableMap
avmap : ActiveVariableMap
a domains.ActiveVariableMap
N : int
N : int
the number of points in the quadrature rule
NMC : int, optional
the number of points in the Monte Carlo estimates of the conditional
the number of points in the Monte Carlo estimates of the conditional
expectation and conditional variance (default 10)

Returns
-------
mu : float
an estimate of the integral of the function computed against the weight
mu : float
an estimate of the integral of the function computed against the weight
function on the simulation inputs
lb : float
a central-limit-theorem 95% lower confidence from the Monte Carlo part
lb : float
a central-limit-theorem 95% lower confidence from the Monte Carlo part
of the integration
ub : float
a central-limit-theorem 95% upper confidence from the Monte Carlo part
ub : float
a central-limit-theorem 95% upper confidence from the Monte Carlo part
of the integration

See Also
Expand All @@ -51,7 +52,7 @@ def integrate(fun, avmap, N, NMC=10):
if not isinstance(avmap, ActiveVariableMap):
raise TypeError('avmap should be an ActiveVariableMap.')

if not isinstance(N, int):
if not isinstance(N, Integral):
raise TypeError('N should be an integer')

# get the quadrature rule
Expand Down Expand Up @@ -84,16 +85,16 @@ def av_integrate(avfun, avmap, N):

Parameters
----------
avfun : function
avfun : function
a function of the active variables
avmap : ActiveVariableMap
avmap : ActiveVariableMap
a domains.ActiveVariableMap
N : int
N : int
the number of points in the quadrature rule

Returns
-------
mu : float
mu : float
an estimate of the integral

Notes
Expand All @@ -104,7 +105,7 @@ def av_integrate(avfun, avmap, N):
if not isinstance(avmap, ActiveVariableMap):
raise TypeError('avmap should be an ActiveVariableMap.')

if not isinstance(N, int):
if not isinstance(N, Integral):
raise TypeError('N should be an integer.')

Yp, Yw = av_quadrature_rule(avmap, N)
Expand All @@ -117,27 +118,27 @@ def av_integrate(avfun, avmap, N):

def quadrature_rule(avmap, N, NMC=10):
"""Get a quadrature rule on the space of simulation inputs.

Parameters
----------
avmap : ActiveVariableMap
avmap : ActiveVariableMap
a domains.ActiveVariableMap
N : int
N : int
the number of quadrature nodes in the active variables
NMC : int, optional
the number of samples in the simple Monte Carlo over the inactive
NMC : int, optional
the number of samples in the simple Monte Carlo over the inactive
variables (default 10)

Returns
-------
Xp : ndarray
(N*NMC)-by-m matrix containing the quadrature nodes on the simulation
Xp : ndarray
(N*NMC)-by-m matrix containing the quadrature nodes on the simulation
input space
Xw : ndarray
(N*NMC)-by-1 matrix containing the quadrature weights on the simulation
Xw : ndarray
(N*NMC)-by-1 matrix containing the quadrature weights on the simulation
input space
ind : ndarray
array of indices identifies which rows of `Xp` correspond to the same
ind : ndarray
array of indices identifies which rows of `Xp` correspond to the same
fixed value of the active variables

See Also
Expand All @@ -163,10 +164,10 @@ def quadrature_rule(avmap, N, NMC=10):
if not isinstance(avmap, ActiveVariableMap):
raise TypeError('avmap should be an ActiveVariableMap.')

if not isinstance(N, int):
if not isinstance(N, Integral):
raise TypeError('N should be an integer.')

if not isinstance(NMC, int):
if not isinstance(NMC, Integral):
raise TypeError('NMC should be an integer.')

# get quadrature rule on active variables
Expand All @@ -182,14 +183,14 @@ def av_quadrature_rule(avmap, N):

Parameters
----------
avmap : ActiveVariableMap
avmap : ActiveVariableMap
a domains.ActiveVariableMap
N : int
N : int
the number of quadrature nodes in the active variables

Returns
-------
Yp : ndarray
Yp : ndarray
quadrature nodes on the active variables
Yw : ndarray
quadrature weights on the active variables
Expand All @@ -215,25 +216,25 @@ def av_quadrature_rule(avmap, N):

def interval_quadrature_rule(avmap, N, NX=10000):
"""Quadrature rule on a one-dimensional interval.

Quadrature when the dimension of the active subspace is 1 and the
simulation parameter space is bounded.

Parameters
----------
avmap : ActiveVariableMap
avmap : ActiveVariableMap
a domains.ActiveVariableMap
N : int
N : int
the number of quadrature nodes in the active variables
NX : int, optional
NX : int, optional
the number of samples to use to estimate the quadrature weights (default
10000)

Returns
-------
Yp : ndarray
Yp : ndarray
quadrature nodes on the active variables
Yw : ndarray
Yw : ndarray
quadrature weights on the active variables

See Also
Expand Down Expand Up @@ -261,25 +262,25 @@ def interval_quadrature_rule(avmap, N, NX=10000):

def zonotope_quadrature_rule(avmap, N, NX=10000):
"""Quadrature rule on a zonotope.

Quadrature when the dimension of the active subspace is greater than 1 and
the simulation parameter space is bounded.

Parameters
----------
avmap : ActiveVariableMap
avmap : ActiveVariableMap
a domains.ActiveVariableMap
N : int
N : int
the number of quadrature nodes in the active variables
NX : int, optional
NX : int, optional
the number of samples to use to estimate the quadrature weights (default
10000)

Returns
-------
Yp : ndarray
quadrature nodes on the active variables
Yw : ndarray
Yw : ndarray
quadrature weights on the active variables

See Also
Expand Down Expand Up @@ -310,5 +311,3 @@ def zonotope_quadrature_rule(avmap, N, NX=10000):

Yp, Yw = points.reshape((T.nsimplex,n)), weights.reshape((T.nsimplex,1))
return Yp, Yw


Loading