Solving Expectation Problems

solve(exprob::ExpectationProblem, expalg::MonteCarlo)

Estimate the expectation defined by exprob with independent Monte Carlo samples.

Arguments

• exprob: Expectation problem whose distribution supplies the samples and whose observable supplies the evaluated quantity.
• expalg: Monte Carlo algorithm. Its trajectories field determines the number of samples or ensemble trajectories.

Returns

An ExpectationSolution whose u field is the sample mean. For a SystemMap or ProcessNoiseSystemMap, the samples are evaluated through an ensemble solve.

Example

using Distributions, SciMLExpectations

prob = ExpectationProblem((x, p) -> x[1]^2, GenericDistribution(Uniform(-1, 1)), nothing)
sol = solve(prob, MonteCarlo(10_000))
sol.u # approximately 1 / 3

solve(exprob::ExpectationProblem, expalg::Koopman;
      maxiters = 1_000_000, batch = nothing, quadalg = HCubatureJL(),
      ireltol = 1e-2, iabstol = 1e-2, kwargs...)

Compute the expectation defined by exprob with deterministic quadrature.

Arguments

• exprob: Expectation problem to integrate. Its distribution must provide finite integration bounds through extrema and a density through pdf.
• expalg: Koopman algorithm that selects the differentiation strategy used by differentiable expectation solves.

Keyword Arguments

• maxiters: Maximum number of integrand evaluations accepted by the quadrature algorithm.
• batch: Number of samples evaluated together for system-map integrands. Leave as nothing for scalar evaluation.
• quadalg: Integrals.jl quadrature algorithm. The default is HCubatureJL().
• ireltol, iabstol: Relative and absolute integration tolerances.
• kwargs...: Additional keyword arguments forwarded to the Integrals.jl solve.

Returns

An ExpectationSolution containing the integral value in u, the integration residual in resid, and the underlying Integrals.jl solution in original.

Example

using Distributions, SciMLExpectations

prob = ExpectationProblem((x, p) -> x[1]^2, GenericDistribution(Uniform(-1, 1)), nothing)
sol = solve(prob, Koopman())
sol.u # approximately 1 / 3

build_integrand(prob::ExpectationProblem, expalg::Koopman, mid, p, batch)

Construct the integral function used by a Koopman expectation solve.

This helper lowers an ExpectationProblem to the IntegralFunction or BatchIntegralFunction consumed by Integrals.jl. It is primarily useful for advanced users who need to inspect or customize the integrand before integration.

Arguments

• prob: Expectation problem to lower.
• expalg: Koopman expectation algorithm.
• mid: Midpoint of the integration domain, used to infer output shape for batched integrands.
• p: Parameters passed to the integrand.
• batch: nothing for scalar integration or an integer batch size.

Returns

An IntegralFunction for unbatched integration or a BatchIntegralFunction for batched system-map integration.

centralmoment(n, exprob::ExpectationProblem, expalg::AbstractExpectationAlgorithm; kwargs...)

Compute the first n central moments of the observable defined in exprob.

Returns a vector [μ₁, μ₂, ..., μₙ] where μₖ is the k-th central moment. Note that the first central moment is always 0 (by definition).

The central moments are computed using the binomial expansion: μₙ = E[(X - μ)ⁿ] = Σₖ₌₀ⁿ C(n,k) (-μ)^(n-k) E[X^k]

where μ = E[X] is the mean.

Arguments

• n: The highest order central moment to compute (n ≥ 1)
• exprob: An ExpectationProblem with a scalar-valued observable
• expalg: The algorithm to use (Koopman() or MonteCarlo(trajectories))
• kwargs...: Additional keyword arguments passed to solve

Returns

A vector of central moments from order 1 through order n.

Examples

using Distributions, SciMLExpectations

gd = GenericDistribution(Uniform(-1, 1))
g(u, p) = u[1]  # Identity function
exprob = ExpectationProblem(g, gd, nothing)

# Compute first 4 central moments
moments = centralmoment(4, exprob, Koopman())
# moments[1] ≈ 0 (1st central moment is always 0)
# moments[2] ≈ 1 / 3 (variance of Uniform(-1, 1))
# moments[3] ≈ 0 (symmetry)
# moments[4] ≈ 1 / 5

centralmoment(n, exprob::ExpectationProblem{F}, expalg::AbstractExpectationAlgorithm; kwargs...) where {F <: Union{SystemMap, ProcessNoiseSystemMap}}

Compute the first n central moments for differential equation systems.

See centralmoment(n, exprob, expalg) for details.

ExpectationSolution

Container returned by solve for an ExpectationProblem.

Fields

• u: Computed expectation value.
• resid: Residual reported by the underlying integration algorithm, or nothing when unavailable.
• original: Original solver result returned by the underlying algorithm, or nothing when unavailable.