# Problems fitting a simple Multinomial model in RxInfer
**URL:** https://discourse.julialang.org/t/problems-fitting-a-simple-multinomial-model-in-rxinfer/95442
**Category:** Probabilistic Programming
**Tags:** question, package, bayesian-inference, rxinfer
**Created:** [March 2, 2023, 12:50pm UTC](https://discourse.julialang.org/t/problems-fitting-a-simple-multinomial-model-in-rxinfer/95442 "2023-03-02T12:50:27Z")
**Posts on this page:** 6
**Page:** 1
### Author:  [@DoktorMike](https://discourse.julialang.org/u/DoktorMike)
#### Post date: [March 2, 2023, 12:50pm UTC](https://discourse.julialang.org/t/problems-fitting-a-simple-multinomial-model-in-rxinfer/95442/1 "2023-03-02T12:50:27Z")
Hello everyone,
I’m experimenting with the awesome package `RxInfer` but apparently I’m having some difficulties in understanding some basic principles.
I tried to adapt one of the first simple examples of `RxInfer` to build a small model for Net promotor score (NPS). Basically each data point is a three dimensional vector where each dimension is the count of Detractors, Neutral and Promoters. Here is the very simple model. I can run it through until it’s actually time to call infer. Then I get a message that I need to define the inference rules. Any hints on how to do that?
```julia
using RxInfer
using Random
ns = rand(Binomial(300, 0.3), 100)
trueθ = [0.3, 0.5, 0.2] # detractors, neutral, promoters
dataset = float.([rand(Multinomial(ns[i], trueθ)) for i in eachindex(ns)])
#dataset = reduce(hcat, dataset)
@model function npsmodel(n)
y = datavar(Vector{Float64}, n)
# We endow θ parameter of our model with a conjugate prior
θ ~ Dirichlet([10, 10, 10])
# We assume that outcome of each trial
# is governed by the Multinomial distribution
for i in 1:n
y[i] ~ Multinomial(90, θ)
end
end
result = inference(
model=npsmodel(size(dataset, 2)),
data=(y=dataset,)
)
```
The stacktrace:
```julia
ERROR: `Multinomial` is not available as a node in the inference engine. Used in `y ~ Multinomial(...)` expression.
Use `@node` macro to add a custom factor node corresponding to `Multinomial`. See `@node` macro for additional documentation and examples.
Stacktrace:
[1] error(s::String)
@ Base ./error.jl:35
[2] make_node(::Type, ::FactorNodeCreationOptions{Tuple{Tuple{Int64}, Tuple{Int64}, Tuple{Int64}}, Nothing, Nothing}, ::DataVariable{PointMass{Vector{Float64}}, Rocket.RecentSubjectInstance{Message{PointMass{Vector{Float64}}}, Subject{Message{PointMass{Vector{Float64}}}, AsapScheduler, AsapScheduler}}}, ::ConstVariable{PointMass{Int64}, SingleObservable{Message{PointMass{Int64}, Nothing}, AsapScheduler}}, ::RandomVariable)
@ ReactiveMP ~/.julia/packages/ReactiveMP/vWHNY/src/node.jl:980
[3] make_node(::FactorGraphModel, ::FactorNodeCreationOptions{Nothing, Nothing, Nothing}, ::Type, ::DataVariable{PointMass{Vector{Float64}}, Rocket.RecentSubjectInstance{Message{PointMass{Vector{Float64}}}, Subject{Message{PointMass{Vector{Float64}}}, AsapScheduler, AsapScheduler}}}, ::ConstVariable{PointMass{Int64}, SingleObservable{Message{PointMass{Int64}, Nothing}, AsapScheduler}}, ::RandomVariable)
@ RxInfer ~/.julia/packages/RxInfer/NiAqM/src/model.jl:340
[4] macro expansion
@ ./REPL[163]:12 [inlined]
[5] macro expansion
@ ~/.julia/packages/GraphPPL/n5QGe/src/model.jl:445 [inlined]
[6] var"##npsmodel#355"(model#352::FactorGraphModel, n::Int64)
@ Main ~/.julia/packages/RxInfer/NiAqM/src/graphppl.jl:33
[7] ModelGenerator
@ ~/.julia/packages/RxInfer/NiAqM/src/model.jl:248 [inlined]
[8] #create_model#98
@ ~/.julia/packages/RxInfer/NiAqM/src/model.jl:265 [inlined]
[9] inference(; model::RxInfer.ModelGenerator{var"###npsmodel#355", Tuple{Int64}, NamedTuple{(), Tuple{}}}, data::NamedTuple{(:y,), Tuple{Vector{Vector{Float64}}}}, initmarginals::Nothing, initmessages::Nothing, constraints::Nothing, meta::Nothing, options::Nothing, returnvars::Nothing, iterations::Nothing, free_energy::Bool, free_energy_diagnostics::Tuple{BetheFreeEnergyCheckNaNs, BetheFreeEnergyCheckInfs}, showprogress::Bool, callbacks::Nothing, addons::Nothing, postprocess::DefaultPostprocess, warn::Bool)
@ RxInfer ~/.julia/packages/RxInfer/NiAqM/src/inference.jl:499
[10] top-level scope
@ REPL[164]:1
```
Update: I realized that I need to make my own @node and corresponding rules. That seems like a lot of stuff that I have no idea how to do. Is there a quick way to achieve this?
---
### Author:  [@cscherrer](https://discourse.julialang.org/u/cscherrer)
#### Post date: [March 7, 2023, 2:57pm UTC](https://discourse.julialang.org/t/problems-fitting-a-simple-multinomial-model-in-rxinfer/95442/2 "2023-03-07T14:57:10Z")
I don’t know how closely he watches this, so for RxInfer questions it might be good to tag @bvdmitri (unless he prefers otherwise 🙂 )
---
### Author:  [@albertpod](https://discourse.julialang.org/u/albertpod)
#### Post date: [March 7, 2023, 3:55pm UTC](https://discourse.julialang.org/t/problems-fitting-a-simple-multinomial-model-in-rxinfer/95442/3 "2023-03-07T15:55:19Z")
> [@DoktorMike](#):
>
> ```julia
> ERROR: `Multinomial` is not available as a node in the inference engine. Used in `y ~ Multinomial(...)` expression.
> Use `@node` macro to add a custom factor node corresponding to `Multin
>
> ```
Hi @DoktorMike!
Thank you for trying out `RxInfer.jl`. Unfortunately, the Multinomial node is currently not available in `ReactiveMP.jl`, which is the inference engine used in `RxInfer.jl`.
I can imagine that it’s frustrating that you are not able to run inference in such a simple model. However, there are analytical rules available for this model, so it is possible to implement it in principle. You could open an issue on `ReactiveMP.jl` requesting this feature.
Alternatively, you could implement the node and corresponding rules yourself. However, I cannot provide you with a link to the documentation as we are currently working on an example of _how to derive rules_.
At the moment, the tutorials can be found in [this thesis](https://research.tue.nl/files/126071163/20190612_Laar.pdf) (see Appendix A), but I understand that they may be difficult to follow. I will open an issue on `RxInfer.jl` to include instructions on how to implement a node in the documentation.
Technically, what needs to be done is to define the node:
```julia
@node Multinomial Stochastic [out, n, k]
```
and then specify the update rule for k interface (θ in your model):
```julia
@rule Multinomial(:k, Marginalisation) (q_out::PointMass, q_n::PointMass, ) = begin
return Dirichlet(probvec(q_out) .+ one(eltype(probvec(q_out))))
end
```
Given all that, your snippet can look as follows:
```julia
using RxInfer
using Random
ns = rand(Binomial(300, 0.3), 100)
trueθ = [0.3, 0.5, 0.2] # detractors, neutral, promoters
dataset = float.([rand(Multinomial(ns[i], trueθ)) for i in eachindex(ns)])
@node Multinomial Stochastic [out, n, k]
@rule Multinomial(:k, Marginalisation) (q_out::PointMass, q_n::PointMass, ) = begin
return Dirichlet(probvec(q_out) .+ one(eltype(probvec(q_out))))
end
@model function npsmodel(n)
y = datavar(Vector{Float64}, n)
# We endow θ parameter of our model with a conjugate prior
θ ~ Dirichlet([10, 10, 10])
# We assume that outcome of each trial
# is governed by the Multinomial distribution
for i in 1:n
y[i] ~ Multinomial(90, θ)
end
end
result = inference(
model=npsmodel(size(dataset, 1)),
data=(y=dataset,)
)
@show mean(result.posteriors[:θ])
mean(result.posteriors[:θ]) = [0.29851203608373783, 0.5084326126391137, 0.1930553512771486]
```
P.S. There might be a way of hacking this model with RxInfer approximations, but I doubt that you want to do it for such a simple case.
---
### Author:  [@DoktorMike](https://discourse.julialang.org/u/DoktorMike)
#### Post date: [March 7, 2023, 4:46pm UTC](https://discourse.julialang.org/t/problems-fitting-a-simple-multinomial-model-in-rxinfer/95442/4 "2023-03-07T16:46:22Z")
Great, thank you for the pointers. I will read up on the methodology and see if I can contribute to the rules. 🙏🏻☺
---
### Author:  [@albertpod](https://discourse.julialang.org/u/albertpod)
#### Post date: [March 7, 2023, 6:40pm UTC](https://discourse.julialang.org/t/problems-fitting-a-simple-multinomial-model-in-rxinfer/95442/5 "2023-03-07T18:40:11Z")
You are welcome. @bvdmitri pointed out at this small [example](https://biaslab.github.io/RxInfer.jl/stable/examples/Advanced%20Tutorial/#Custom-messages-computation-rules) on nodes and rules that will be extended in the near future.
UPD:
The rule I’ve provided appears to be correct.
---
### Author:  [@DoktorMike](https://discourse.julialang.org/u/DoktorMike)
#### Post date: [March 26, 2023, 2:33pm UTC](https://discourse.julialang.org/t/problems-fitting-a-simple-multinomial-model-in-rxinfer/95442/6 "2023-03-26T14:33:10Z")
That perfectly solves it! Thank you @albertpod 🙏🏻