# Transducers for sufficient statistics

**URL:** https://discourse.julialang.org/t/transducers-for-sufficient-statistics/27425
**Category:** Probabilistic Programming
**Tags:** question
**Created:** [August 11, 2019, 4:49pm UTC](https://discourse.julialang.org/t/transducers-for-sufficient-statistics/27425 "2019-08-11T16:49:09Z")
**Posts on this page:** 6
**Page:** 1

<div class="post-metadata">

### Author: ![cscherrer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/cscherrer/32/7631_2.png) [@cscherrer](https://discourse.julialang.org/u/cscherrer)
#### Post date: [August 11, 2019, 4:49pm UTC](https://discourse.julialang.org/t/transducers-for-sufficient-statistics/27425/1 "2019-08-11T16:49:09Z")

</div>

For independent observations from some distribution, we end up doing a lot of

```julia
sum(logpdf.(myDist, x))

```

But in some cases, the algebra simplified a lot, and we can just track sufficient statistics. A common example of this is normal distributions, where we only need `(Σx,Σx²)`. This kind of thing generalizes quite a bit - there’s a whole collection of [exponential families](https://en.wikipedia.org/wiki/Exponential_family) that allow this kind of optimization.

I’ve been looking into this a little bit, current progress is [here](https://github.com/cscherrer/ExponentialFamilies.jl/blob/master/src/types.jl). But I’m a bit stuck on

```julia
function iid(d::NatExpFamDist{P,X} ) where {P,X}
    logh(xs :: AbstractArray{<:X}) = sum(d.fam.logh.(xs))
    t(xs :: AbstractArray{<:X}) = reduce((a,b) -> a .+ b,d.fam.t.(xs))
    a(η::P) = d.fam.a(η)

    fam = NaturalExponentialFamily{P,AbstractArray{X}}(logh, t, a)
    NatExpFamDist(fam,d.η)
end

```

I already need to do some weird `reduce` tricks, and I also traverse the data twice (a sum and a reduce), which is not great.

This seems like a natural place to use @tkf’s Transducers, but I’m not familiar with it enough yet to see how to go about it. Any ideas?

Even with these issues and without being too careful yet about the typing, the performance is pretty good. So I think there could be a big advantage once it’s cleaned up 🙂

![ExponentialFamilies](https://global.discourse-cdn.com/julialang/original/3X/6/3/634a4f3377ffb060a37d3d5933ed807ced377a15.png)

---

<div class="post-metadata">

### Author: ![cscherrer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/cscherrer/32/7631_2.png) [@cscherrer](https://discourse.julialang.org/u/cscherrer)
#### Post date: [August 11, 2019, 6:34pm UTC](https://discourse.julialang.org/t/transducers-for-sufficient-statistics/27425/2 "2019-08-11T18:34:18Z")

</div>

Making some progress, think it will be something like this:

```julia
julia> f(a,b) = (a[1]+b[1], a[2]+b[2])
f (generic function with 1 method)

julia> mapfoldl(Zip(Map(identity),Map(x -> x^2)), Completing(f), 1:10, init=(0.0,0.0))
(55.0, 385.0)

```

---

<div class="post-metadata">

### Author: ![tkf](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tkf/32/17635_2.png) [@tkf](https://discourse.julialang.org/u/tkf)
#### Post date: [August 11, 2019, 6:45pm UTC](https://discourse.julialang.org/t/transducers-for-sufficient-statistics/27425/3 "2019-08-11T18:45:32Z")

</div>

It’s great that you find Transducers.jl useful! But I wonder if OnelineStats.jl already does the job out-of-the-box.

> **[GitHub - joshday/OnlineStats.jl: ⚡ Single-pass algorithms for statistics](https://github.com/joshday/OnlineStats.jl)**
>
> ⚡ Single-pass algorithms for statistics. Contribute to joshday/OnlineStats.jl development by creating an account on GitHub.

You can also [convert OnlineStat to a transducer](https://tkf.github.io/Transducers.jl/dev/manual/#Transducers.Transducer-Tuple%7BOnlineStatsBase.OnlineStat%7D). This could be useful if you want to compose pre-/post-processings that are easier to do in transducers.

---

<div class="post-metadata">

### Author: ![tkf](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tkf/32/17635_2.png) [@tkf](https://discourse.julialang.org/u/tkf)
#### Post date: [August 12, 2019, 1:30am UTC](https://discourse.julialang.org/t/transducers-for-sufficient-statistics/27425/4 "2019-08-12T01:30:09Z")

</div>

> [@cscherrer](#):
>
> ```julia
> julia> f(a,b) = (a[1]+b[1], a[2]+b[2])
> f (generic function with 1 method)
> 
> julia> mapfoldl(Zip(Map(identity),Map(x -> x^2)), Completing(f), 1:10, init=(0.0,0.0))
> (55.0, 385.0)
> 
> ```

FYI it can also be written as

```julia
foldl(right, Zip(Scan(+), Map(x -> x^2) |> Scan(+)), 1:10)

```

---

<div class="post-metadata">

### Author: ![cscherrer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/cscherrer/32/7631_2.png) [@cscherrer](https://discourse.julialang.org/u/cscherrer)
#### Post date: [August 12, 2019, 2:19am UTC](https://discourse.julialang.org/t/transducers-for-sufficient-statistics/27425/5 "2019-08-12T02:19:39Z")

</div>

> [@tkf](#):
>
> I wonder if OnelineStats.jl already does the job out-of-the-box.

For the Normal case yes, but I think I need finer control of the computation, and composability. But OnlineStats looks great, I think it will be useful for other things anyway 🙂

> [@tkf](#):
>
> FYI it can also be written as
> 
> ```julia
> foldl(right, Zip(Scan(+), Map(x -> x^2) |> Scan(+)), 1:10)
> 
> ```

I don’t understand this yet. Is one of them more easily made parallel? I had wondered about rewriting it using `mapreduce` if that could help performance

---

<div class="post-metadata">

### Author: ![tkf](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tkf/32/17635_2.png) [@tkf](https://discourse.julialang.org/u/tkf)
#### Post date: [August 12, 2019, 3:19am UTC](https://discourse.julialang.org/t/transducers-for-sufficient-statistics/27425/6 "2019-08-12T03:19:14Z")

</div>

> [@cscherrer](#):
>
> Is one of them more easily made parallel?

Right, `Scan` does not support parallelism (yet). Supporting `Scan` in parallel `reduce` is rather tricky ([https://www.youtube.com/watch?v=JCvT9Rnhyvk](https://www.youtube.com/watch?v=JCvT9Rnhyvk)) but maybe special-casing when the bottom reducing function is `right` is a good solution. It’s relatively easy to implement and can be done without allocations.
