mean(expectation::ClosedWilliamsProduct, f, q::ExponentialFamilyDistr…

…ibution{NormalMeanVariance}) (#18)

src/ClosedFormExpectations.jl

@@ -55,6 +55,8 @@ include("Exponential/Exponential.jl")
 include("Normal/expectation.jl")
 include("Normal/williams/normal.jl")
 include("Normal/williams/normal_mean_variance.jl")
+include("Normal/williams/ef_parametrization.jl")
+
 # gamma
 include("Gamma/Gamma.jl")
 

src/Normal/williams/ef_parametrization.jl

@@ -0,0 +1,10 @@
+using StaticArrays
+import Distributions: var, mean, std
+import ExponentialFamily: ExponentialFamilyDistribution, NormalMeanVariance, getnaturalparameters
+
+function mean(expectation::ClosedWilliamsProduct, f, q::ExponentialFamilyDistribution{T}) where {T <: NormalMeanVariance}
+    η = getnaturalparameters(q)
+    jacobian = @SMatrix [-inv(2*η[2]) η[1]/(2*η[2]^2); 0.0 (-1/η[2])^(3/2)/(2*sqrt(2))]
+    normal = Normal(mean(q), std(q))
+    return mean(expectation, f, normal)' * jacobian
+end

test/Normal/williams/ef_parametrization/ef_utils.jl

@@ -0,0 +1,5 @@
+include("../normal_utils.jl")
+
+using ExponentialFamily
+
+score(q::ExponentialFamilyDistribution{NormalMeanVariance}, x) = sufficientstatistics(q, x) .- gradlogpartition(q)

test/Normal/williams/ef_parametrization/williams_tests.jl

@@ -0,0 +1,31 @@
+@testitem "mean(::ClosedWilliamsProduct, p::Logpdf{Normal}, q::ExponentialFamilyDistribution{NormalMeanVariance})" begin
+    include("ef_utils.jl")
+    rng = StableRNG(123)
+    for _ in 1:10
+        μ1, σ1 = rand(rng)*10, rand(rng)*5
+        μ2, σ2 = rand(rng)*10, rand(rng)*5
+        ef = convert(ExponentialFamilyDistribution, Normal(μ1, σ1))
+        central_limit_theorem_test(ClosedWilliamsProduct(), Logpdf(Normal(μ2, σ2)), ef, score)
+    end
+end
+
+@testitem "mean(::ClosedWilliamsProduct, p::Abs, q::ExponentialFamilyDistribution{NormalMeanVariance})" begin
+    include("ef_utils.jl")
+    rng = StableRNG(123)
+    for _ in 1:10
+        μ1, σ1 = rand(rng)*10, rand(rng)*5
+        ef = convert(ExponentialFamilyDistribution, Normal(μ1, σ1))
+        central_limit_theorem_test(ClosedWilliamsProduct(), Abs(), ef, score)
+    end
+end
+
+@testitem "mean(::ClosedWilliamsProduct, p::Logpdf{Laplace}, q::ExponentialFamilyDistribution{NormalMeanVariance})" begin
+    include("ef_utils.jl")
+    rng = StableRNG(123)
+    for _ in 1:10
+        μ1, σ1 = rand(rng)*10, rand(rng)*5
+        loc, θ = rand(rng)*10, rand(rng)*10
+        ef = convert(ExponentialFamilyDistribution, Normal(μ1, σ1))
+        central_limit_theorem_test(ClosedWilliamsProduct(), Logpdf(Laplace(loc, θ)), ef, score)
+    end
+end

test/interface/gauss_tests.jl

@@ -27,13 +27,16 @@ end
 
     @testset "ClosedWilliamsProduct interface" begin
         using Distributions
+        using ExponentialFamily
 
         nmv = NormalMeanVariance(0.0, 1.0)
         nmp = NormalMeanPrecision(0.0, 1.0)
+        ef = convert(ExponentialFamilyDistribution, nmv)
 
         @test mean(ClosedWilliamsProduct(), Logpdf(nmv), Normal(0, 1)) isa AbstractArray
         @test mean(ClosedWilliamsProduct(), Logpdf(nmp), Normal(0, 1)) isa AbstractArray
         @test mean(ClosedWilliamsProduct(), Logpdf(nmv), Normal(0, 1)) ≈ mean(ClosedWilliamsProduct(), Logpdf(nmp), Normal(0, 1))
         @test mean(ClosedWilliamsProduct(), Logpdf(nmv), nmv) ≈ mean(ClosedWilliamsProduct(), Logpdf(nmp), nmv)
+        @test mean(ClosedWilliamsProduct(), Logpdf(nmv), ef) isa AbstractArray
     end
 end