Visualizing Forney Factor Graphs when using RxInfer (Part 3)

In Part 2 we provided a configurable layout scheme suitable for Forney Factor Graphs (FFGs). We also moved to using TikZ rather than GraphViz (please see Part 2 for an overview of the limitations of this package for our purposes). In this part we:

• add decorations for a Constrained Forney Factor Graph (CFFG)
  • currently, constraint decorations are added blindly, rather than making their addition dependent on information from the RxInfer-derived graph
  • in a future part, we will rectify this
• renamed some variables in generate_tikz() to improve consistency and simplicity
• remove the seqlup() function
• remove the str_facs() function
• remove the str_individual_vars() function
• remove the str_vector_vars() function
• add decorations for delta (also known as clamped) factors

One of the distinguishing features of the RxInfer Julia package is that it uses reactive message passing (RMP). This means node updates in the factor graph happens asynchronously. A user usually tends to gloss over these details due to the increased complexity. When node updates happen in the more traditional way, i.e. synchronously (for example as with the PyMDP Python package), users tend to have the option of being aware of which node is currently being udpated. It is quite beneficial to conceptualize the workings of factor graphs in general and the visualization of the factor graph itself is often helpful. The existing visualization code is rather limited and that is the reason this project was undertaken. The ultimate goal is reach a visualization which conforms to the Constrained Forney Factor Graph (CFFG) as defined in

https://biaslab.github.io/pdf/arxiv/2306.08014.pdf

https://biaslab.github.io/pdf/arxiv/2306.02733.pdf

GraphPPL.jl exports the @model macro for model specification allowing the acceptance of a regular Julia function and then builds a Forney Factor Graph (FFG) under the hood. This graph is bipartite where nodes belong to one of two sets: one for factors and one for variables. Note that we use the terms node and link rather than vertex and edge to align with RxInfer terminology.

Node Taxonomy

• nodes
  • facs (FFG/RxInfer factor nodes)
  • vars (FFG/RxInfer variable nodes)
    • parvars (parameter variables, e.g. \(\alpha, a, A_t, B, \gamma_t, \mu\))
    • sysvars (system variables, e.g. \(x_t, s_t, u_t\))

Node Visualization

• a fac node is visualized by a
  • square (style fac) if it is a factor
  • small black square (style dlt) if it is a delta (\(\delta\)) factor
    • delta factors connect to their associated parvars and sysvars
• a var node is visualized by a
  • point (style var) if connected to a single fac
  • equals-square (style eql) if connected to multiple facs
• rather than converting the bipartite graph from RxInfer to a graph that only contains a single kind of node, we keep the bipartite graph and visualize the nodes differently

For each of the examples in this project the following steps happen:

• A @model is created

• The @model is conditioned on some data

  • ususally a few points so that the visualization does not become unwieldy

• A RxInfer model is created (rxi_model)

• A GraphPPL model is obtained (gppl_model)

• A meta graph is obtained (meta_graph)

• The graph is rendered by the existing implementation (GraphPlot.gplot())

• Some vertex labels are inspected for the sake of interest

• Setup global parameters (indicated by a leading underscore character) to influence the behavior of generate_tikz()

  • _san_node_ids is provided for node names that needs to be sanitized
  • _lup_enabled is set to true/false to enable lookup of math names
  • _raw_links_enabled is set to true/false to enable the drawing of raw links for reference
  • _CFFG is set to true/false to enable the drawing of beads to turn a FFG into a CFFG
  • _lup_dict is provided for all variable names needing math equivalents

• The TikZ string is generated (generate_tikz())

  • heading of the FFG
  • GraphPPL model
  • layers for the FFG
  • layer names (if applicable)
    • cntrl (control layer)
    • paramB1 (parameter B1 layer)
    • paramB2 (parameter B2 layer)
    • state1 (state 1 layer)
    • state2 (state 2 layer)
    • state3 (state 3 layer)
    • obser (observation layer)
    • prefr (preference layer, for setpoints)
    • paramA (parameter A layer)

• The TikZ string is printed

• The FFG is saved to a pdf file (show_tikz())

Setup notes

• If additional Ubuntu packages are needed for a Julia package
  • sudo apt-get install package
• If you want to run bash on the devcontainer in a local terminal; kobus@Kobuss-Mac-mini ~ %
  • docker exec -it b20ca0308c9b2698fd7c19fdac794941cd50ecc7a11e900b7e47e5018383e033 /bin/bash

Development cycle notes

• execute & save the generated tikz/pgf code
• copy the generated .tex file from the devcontainer to the local folder; kobus@Kobuss-Mac-mini ~ %
  • use vscode command line (View -> Terminal) to discover the folder name on the devcontainer
  • use the docker ui to copy the id of the container
  • identify the local folder name from where the .tex will be rendered
  • overlay the next command with what you discovered above
  • docker cp b20ca0308c9b2698fd7c19fdac794941cd50ecc7a11e900b7e47e5018383e033:/workspaces/2024-07-05^FFGViz2/myoutput.tex /Users/kobus/2024_LaTeX
• in the /Users/kobus/2024_LaTeX project on the mini, the pdf will be automatically rendered when a change is detected in the copied .tex

versioninfo() ## Julia version

Julia Version 1.10.0
Commit 3120989f39b (2023-12-25 18:01 UTC)
Build Info:
  Official https://julialang.org/ release
Platform Info:
  OS: Linux (x86_64-linux-gnu)
  CPU: 12 × Intel(R) Core(TM) i7-8700B CPU @ 3.20GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-15.0.7 (ORCJIT, skylake)
  Threads: 1 on 12 virtual cores
Environment:
  JULIA_NUM_THREADS = 

import Pkg; Pkg.instantiate(); Pkg.activate()
Pkg.add(Pkg.PackageSpec(;name="RxInfer", version="3.0.0"))
Pkg.add(Pkg.PackageSpec(;name="MetaGraphsNext", version="0.7.0"))
## 
# Pkg.add(Pkg.PackageSpec(;name="Distributions", version="0.25.108"))
# Pkg.add(Pkg.PackageSpec(;name="GraphPlot", version="0.5.2"))
# Pkg.add(Pkg.PackageSpec(;name="Cairo", version="1.0.5"))
# Pkg.add(Pkg.PackageSpec(;name="Graphs", version="1.9.0"))
# Pkg.add(Pkg.PackageSpec(;name="MetaGraphsNext", version="0.7.0"))
# Pkg.add(Pkg.PackageSpec(;name="GraphPPL", version="4.0.2"))
# Pkg.add(Pkg.PackageSpec(;name="GraphViz", version="0.2.0"))
# Pkg.add(Pkg.PackageSpec(;name="Dictionaries", version="0.4.2"))
# Pkg.add(Pkg.PackageSpec(;name="Plots", version="1.40.4"))
# Pkg.add(Pkg.PackageSpec(;name="StableRNGs", version="1.0.2"))
# Pkg.add(Pkg.PackageSpec(;name="StatsPlots", version="0.15.7"))
# Pkg.add(Pkg.PackageSpec(;name="LaTeXStrings", version="1.3.1"))
# Pkg.add(Pkg.PackageSpec(;name="DataFrames", version="1.6.1"))
# Pkg.add(Pkg.PackageSpec(;name="CSV", version="0.10.14"))
# Pkg.add(Pkg.PackageSpec(;name="GLM", version="1.9.0"))

  Activating project at `~/.julia/environments/v1.10`
   Resolving package versions...
  No Changes to `~/.julia/environments/v1.10/Project.toml`
  No Changes to `~/.julia/environments/v1.10/Manifest.toml`
   Resolving package versions...
  No Changes to `~/.julia/environments/v1.10/Project.toml`
  No Changes to `~/.julia/environments/v1.10/Manifest.toml`

## using RxInfer, Distributions, Random, GraphPlot, Cairo, Graphs, MetaGraphsNext, GraphPPL, GraphViz, Dictionaries, Plots, StableRNGs, LinearAlgebra, StatsPlots, LaTeXStrings, DataFrames, CSV, GLM
using RxInfer, MetaGraphsNext

Pkg.status()

Status `~/.julia/environments/v1.10/Project.toml`
  [fa8bd995] MetaGraphsNext v0.7.0
⌃ [86711068] RxInfer v3.0.0
Info Packages marked with ⌃ have new versions available and may be upgradable.

"""
Takes the tikz string given by generate_tikz() and 
writes it to a .tex file to produce a TikZ visualisation. 
"""
function show_tikz(tikz_code_graph::String)
    ## eval(Meta.parse(dot_code_graph)) ## for GraphViz

    ## see problem in ^v1 under TikzPictures section 
    ## tikz_picture = TikzPicture(tikz_code_graph)
    ## return tikz_picture
    file_path = "myoutput.tex" ## write to file rather than show in notebook
    open(file_path, "w") do file
        write(file, tikz_code_graph)
    end
end

show_tikz

##
# # If your PDF contains a TikZ picture (a vector graphic created using LaTeX), 
# # you can still use the method I mentioned earlier to display it in a Julia notebook. 
# # Here’s the method again for your reference:
# import Pkg; Pkg.add("PGFPlotsX")
# using PGFPlotsX

# # Create a struct that represents a PDF file
# struct PDF
#     file::String
# end

# # Define a function to show the PDF as inline SVG
# function Base.show(io::IO, ::MIME"image/svg+xml", pdf::PDF)
#     svg = first(splitext(pdf.file)) * ".svg"
#     PGFPlotsX.convert_pdf_to_svg(pdf.file, svg)
#     write(io, read(svg))
#     try; rm(svg; force=true); catch e; end
#     return nothing
# end

# # Now you can display a PDF file like this:
# display(PDF("myoutput.pdf"))

# # Internal Error: cairo context error: invalid matrix (not invertible)<0a>
# # cairo error: invalid matrix (not invertible)

function sanitized(vertex::String)
    if haskey(_san_node_ids, vertex)
        return _san_node_ids[vertex]
    else
        return vertex
    end
end

function lup(label::String)
    if _lup_enabled
        lup_label = get(_lup_dict, label, replace(label, "_" => "\\_"))
    else
        return replace(label, "_" => "\\_")
    end
end

function show_meta_graph(meta_graph) ## just for info
    str = ""
    for i in MetaGraphsNext.vertices(meta_graph)
        label = MetaGraphsNext.label_for(meta_graph, i)
        str *= "    $(label);\n"
    end
    for edge in MetaGraphsNext.edges(meta_graph)
        source_vertex = MetaGraphsNext.label_for(meta_graph, edge.src)
        dest_vertex = MetaGraphsNext.label_for(meta_graph, edge.dst)
        str *= "    $(source_vertex) -- $(dest_vertex);\n"
    end
    print(str)
end

show_meta_graph (generic function with 1 method)

function spaces(n::Int)
    d = 1/(n + 1)
    return [i*d for i in 1:n]
end
function setup_divn(n)
    divn = Dict()
    for i in 1:n
        divn[i] = spaces(i)
    end
    return divn
end
function radius(constrained::Bool)
    if constrained
        return "\\ar" ## arc radius of a fac node
    else
        return "\\nr" ## node radius of a fac node
    end
end

radius (generic function with 1 method)

_ytop = 18 ## grid y-value at the top of picture
_show_grid = false ## displays a grid with the FFG if true
_node_size = "15mm" ## the length of the sides of factor node boxes (in millimeters)

"15mm"

function get_preamble_and_postamble(heading)
    iob = IOBuffer() ## preamble
    write(iob, "\\documentclass{standalone}\n")
    write(iob, "\\usepackage{tikz}\n")
    write(iob, "\\usetikzlibrary{calc}\n")    
    write(iob, "\\usepackage{amsmath}\n")
    ## write(iob, "\\newcommand\\ns{10mm} %% define node (minimum) size\n")
    write(iob, "\\newcommand\\ns{$(_node_size)} %% define node (minimum) size\n")
    write(iob, "\\newcommand\\nr{.5*\\ns} %% define node radius\n")
    write(iob, "\\newcommand\\ar{.4*\\ns} %% define arc radius\n")
    write(iob, "\\newcommand\\br{.19*\\ns} %% define bead radius\n")
    write(iob, "\\DeclareUnicodeCharacter{03B8}{\\theta}\n")
    write(iob, "\\DeclareUnicodeCharacter{03B3}{\\gamma}\n")
    ## write(iob, "\\DeclareUnicodeCharacter{2080}{\\0}\n") ## theta_0
    write(iob, "\\begin{document}\n")
    write(iob, "\\begin{tikzpicture}[\n")
    write(iob, "    scale=1, \n")
    write(iob, "    draw=black, \n")
    write(iob, "    fac/.style={rectangle, draw=black!100, fill=black!0, minimum size=\\ns, inner sep=0pt},\n")
    write(iob, "    dlt/.style={rectangle, draw=black!100, fill=black!100, minimum size=.2*\\ns},\n")
    write(iob, "    var/.style={circle, draw=gray!100, fill=gray!100, minimum size=1mm, inner sep=0pt},\n")
    write(iob, "    bus/.style={rectangle, draw=black!100, fill=black!100, minimum width=30mm, minimum height=.5mm, inner sep=0pt},\n")
    write(iob, "    eql/.style={rectangle, draw=black!100, fill=black!0, minimum size=.5*\\ns},\n")
    ## write(iob, "    text=gray,\n")
    write(iob, "    text=black,\n")
    write(iob, "    ]\n")
    write(iob, "\\node[font=\\Large] at (4,$(_ytop-1)) {\\textbf{$(heading)}};\n")
    if _show_grid write(iob, "\\draw[lightgray] (0,0) grid ($(_ytop), $(_ytop));\n") end
    preamble = String(take!(iob))
    print(preamble)
    println("⋮"); println("⋮"); println("⋮") ## [\]vdots[tab]

    iob = IOBuffer() ## postamble
    ## Adjust bounding box. Must be here, i.e. after the drawing
    write(iob, "\\path[use as bounding box] ([xshift=-1cm, yshift=-1cm]current bounding box.south west) rectangle ([xshift=1cm, yshift=2cm]current bounding box.north east);\n")
    write(iob, "\\end{tikzpicture}\n")
    write(iob, "\\end{document}\n")
    postamble = String(take!(iob))
    print(postamble)

    return preamble, postamble
end

get_preamble_and_postamble (generic function with 1 method)

function generate_tikz(;
heading::String,
Model::GraphPPL.Model,
layers::Vector{String}, 

paramB1_iniparvars=nothing,
paramB1_inifac=nothing,
paramB1_parvars=nothing,

paramB2_iniparvars=nothing,
paramB2_inifac=nothing,
paramB2_parvars=nothing,

state1_iniparvars=nothing,
state1_inifac=nothing,
state1_sysvars=nothing,
state1_facs=nothing,

state2_iniparvars=nothing,
state2_inifac=nothing,
state2_sysvars=nothing,
state2_facs=nothing,
state2_parvars=nothing,

state3_iniparvars=nothing,
state3_inifac=nothing,
state3_sysvars=nothing,
state3_facs=nothing,
state3_parvars=nothing,

obser_sysvars=nothing,
obser_facs=nothing,
obser_parvars=nothing,

paramA_iniparvars=nothing,
paramA_inifac=nothing,
paramA_parvars=nothing,
)    
    preamble, postamble = get_preamble_and_postamble(heading)
    tikz = preamble
    divn = setup_divn(15)
    iob = IOBuffer() ## graph
    ytop = _ytop; xleft = 0; ysep = 3; xsep = 2; y = ytop; x = xleft
    write(iob, "%% --------------------------- NODES ---------------------------\n")
    if "paramB1" in layers
        y -= ysep; x = xleft + xsep; y_paramB1 = y
        write(iob, "\\node($(paramB1_inifac))[fac] at ($(x), $(y)) {\$$(lup(paramB1_inifac))\$};\n")
        n = length(paramB1_iniparvars)
        for (i,v) in enumerate(paramB1_iniparvars)
            write(iob, "\\node($(v))[var] at (\$ ($(paramB1_inifac).south west) + (-.5*\\ns, $(divn[n][i])*\\ns) \$) {};\n")
            write(iob, "\\node[above] at ($(v)) {\$$(lup(v))\$};\n") ## label
            write(iob, "\\node($(v))[dlt] at (\$ ($(v).center) + (-.4*\\ns, 0) \$) {};\n")
        end
        x += xsep
        for (i,v) in enumerate(paramB1_parvars)
            write(iob, "\\node($(v))[eql] at ($(x), $(y)) {\$=\$};\n")
            write(iob, "\\node[xshift=.4*\\ns, yshift=.4*\\ns] at ($(v)) {\$$(lup(v))\$};\n") ## label
        end
    end
    if "paramB2" in layers
        y -= ysep; x = xleft + xsep; y_paramB2 = y
        write(iob, "\\node($(paramB2_inifac))[fac] at ($(x), $(y)) {\$$(lup(paramB2_inifac))\$};\n")
        n = length(paramB2_iniparvars)
        for (i,v) in enumerate(paramB2_iniparvars)
            write(iob, "\\node($(v))[var] at (\$ ($(paramB2_inifac).south west) + (-.5*\\ns, $(divn[n][i])*\\ns) \$) {};\n")
            write(iob, "\\node[above] at ($(v)) {\$$(lup(v))\$};\n") ## label
            write(iob, "\\node($(v))[dlt] at (\$ ($(v).center) + (-.4*\\ns, 0) \$) {};\n")
        end
        x += xsep
        for (i,v) in enumerate(paramB2_parvars)
            write(iob, "\\node($(v))[eql] at ($(x), $(y)) {\$=\$};\n")
            write(iob, "\\node[xshift=.4*\\ns, yshift=.4*\\ns] at ($(v)) {\$$(lup(v))\$};\n") ## label
        end
    end
    if "state1" in layers
        y -= ysep; x = xleft + xsep; y_state1 = y
        write(iob, "\\node($(state1_inifac))[fac] at ($(x), $(y)) {\$$(lup(state1_inifac))\$};\n")
        n = length(state1_iniparvars)
        for (i,v) in enumerate(state1_iniparvars)
            write(iob, "\\node($(v))[var] at (\$ ($(state1_inifac).south west) + (-.5*\\ns, $(divn[n][i])*\\ns) \$) {};\n")
            ##- write(iob, "\\node($(v))[var] at ( ($(inistate1_fac).south west) ++ (-.5*\\ns, $(divn[n][i])*\\ns) ) {};\n")
            write(iob, "\\node[above] at ($(v)) {\$$(lup(v))\$};\n") ## label
            write(iob, "\\node($(v))[dlt] at (\$ ($(v).center) + (-.4*\\ns, 0) \$) {};\n")
        end
        x += xsep
        write(iob, "\\node($(state1_sysvars[1]))[var] at ($(x), $(y)) {};\n")
        ## write(iob, "\\node[above] at ($(state10_var)) {\$s_0\$};\n") ## label
        write(iob, "\\node[above] at ($(state1_sysvars[1])) {\$$(lup(state1_sysvars[1]))\$};\n") ## label
        x += xsep
        for (i,v) in enumerate(state1_facs)
            write(iob, "\\node($(v))[fac] at ($(x), $(y)) {\$$(lup(v))\$};\n")
            if _CFFG
                ## write(iob, "\\draw ($(v))++(\\ar,0mm) arc[start angle=0, end angle=180, radius=\\ar];\n") ##TOP
                write(iob, "\\draw ($(v))++(\\ar,0mm) arc[start angle=0, end angle=-180, radius=\\ar];\n") ##BOTTOM
            end
            x += 2*xsep
        end
        x = 4*xsep     
        for (i,v) in enumerate(state1_sysvars[2:end]) 
            if i < length(state1_sysvars) - 1
                write(iob, "\\node($(v))[eql] at ($(x), $(y)) {\$=\$};\n")            
                write(iob, "\\node[xshift=.4*\\ns, yshift=.4*\\ns] at ($(v)) {\$$(lup(v))\$};\n") ## label
            else
                write(iob, "\\node($(v))[var] at ($(x), $(y)) {};\n")
                write(iob, "\\node[above] at ($(v)) {\$$(lup(v))\$};\n") ## label
            end
            x += 2*xsep
        end
    end
    if "state2" in layers
        y -= ysep; x = xleft; y_state2 = y
        x = 4*xsep
        for (i,v) in enumerate(state2_facs)
            write(iob, "\\node($(v))[fac] at ($(x), $(y)) {\$$(lup(v))\$};\n")
            x += 2*xsep
        end
        x = 4*xsep             
        for (i,v) in enumerate(state2_parvars)
            write(iob, "\\node($(v))[var] at (\$ ($(state2_facs[i]).south west) + (-.5*\\ns, $(divn[1][1])*\\ns) \$) {};\n")
            write(iob, "\\node[above] at ($(v)) {\$$(lup(v))\$};\n") ## label
            write(iob, "\\node($(v))[dlt] at (\$ ($(v).center) + (-.4*\\ns, 0) \$) {};\n")
            x += 2*xsep
        end
        x = 4*xsep     
        for (i,v) in enumerate(state2_sysvars[2:end])
            if i < length(state2_sysvars) - 1
                write(iob, "\\node($(v))[eql] at ($(x), $(y)) {\$=\$};\n")            
                write(iob, "\\node[xshift=.4*\\ns, yshift=.4*\\ns] at ($(v)) {\$$(lup(v))\$};\n") ## label
            else
                write(iob, "\\node($(v))[var] at ($(x), $(y)) {};\n")
                write(iob, "\\node[above] at ($(v)) {\$$(lup(v))\$};\n") ## label
            end
            x += 2*xsep
        end
    end
    if "state3" in layers
        y -= ysep; x = xleft; x = 4*xsep; y_state3 = y
        for (i,v) in enumerate(state3_facs)
            write(iob, "\\node($(v))[fac] at ($(x), $(y)) {\$$(lup(v))\$};\n")
            x += 2*xsep
        end
        x = 4*xsep             
        for (i,v) in enumerate(state3_parvars)
            write(iob, "\\node($(v))[var] at ($(x)-1, $(y)) {};\n")
            write(iob, "\\node[left] at ($(v)) {\$$(lup(v))\$};\n") ## label
            x += 2*xsep
        end
        x = 4*xsep     
        for (i,v) in enumerate(state3_sysvars)
            if i < length(state3_sysvars) - 1
                write(iob, "\\node($(v))[eql] at ($(x), $(y)) {\$=\$};\n")
                write(iob, "\\node[xshift=.4*\\ns, yshift=.4*\\ns] at ($(v)) {\$$(lup(v))\$};\n") ## label
            else
                write(iob, "\\node($(v))[var] at ($(x), $(y)) {};\n")
                write(iob, "\\node[above] at ($(v)) {\$$(lup(v))\$};\n") ## label
            end
            x += 2*xsep
        end
    end
    if "obser" in layers 
        y -= ysep; x = 4*xsep; y_obser = y
        for (i,v) in enumerate(obser_facs)
            write(iob, "\\node($(v))[fac] at ($(x), $(y)) {\$$(lup(v))\$};\n")
            x += 2*xsep
        end
        x = 4*xsep             
        for (i,v) in enumerate(obser_parvars)
            write(iob, "\\node($(v))[var] at (\$ ($(obser_facs[i]).south west) + (-.5*\\ns, $(divn[1][1])*\\ns) \$) {};\n")
            write(iob, "\\node[above] at ($(v)) {\$$(lup(v))\$};\n") ## label
            write(iob, "\\node($(v))[dlt] at (\$ ($(v).center) + (-.4*\\ns, 0) \$) {};\n")
            x += 2*xsep
        end
        y -= 2; x = 4*xsep
        for (i,v) in enumerate(obser_sysvars)
            write(iob, "\\node($(v))[var] at ($(x), $(y)) {};\n")
            write(iob, "\\node[above right] at ($(v)) {\$$(lup(v))\$};\n") ## label
            ## write(iob, "\\node($(v))[dlt] at (\$ ($(obser_facs[i]).center) + (0, -2.0*\\ns) \$) {};\n")
            write(iob, "\\node($(v))[dlt] at (\$ ($(v).center) + (0, -.4*\\ns) \$) {};\n")
            x += 2*xsep
        end
    end
    if "paramA" in layers
        y -= ysep; x = xleft + xsep; y_paramA = y
        write(iob, "\\node($(paramA_inifac))[fac] at ($(x), $(y)) {\$$(lup(paramA_inifac))\$};\n")
        n = length(paramA_iniparvars)
        for (i,v) in enumerate(paramA_iniparvars)
            write(iob, "\\node($(v))[var] at (\$ ($(paramA_inifac).south west) + (-.5*\\ns, $(divn[n][i])*\\ns) \$) {};\n")
            write(iob, "\\node[above] at ($(v)) {\$$(lup(v))\$};\n") ## label
            write(iob, "\\node($(v))[dlt] at (\$ ($(v).center) + (-.4*\\ns, 0) \$) {};\n")
        end
        x = 4*xsep
        ## x += xsep
        for (i,v) in enumerate(paramA_parvars)
            ## write(iob, "\\node($(v))[bus, minimum width=30mm, minimum height=.5mm, inner sep=0pt] at (8, $(y)) {};\n")
            write(iob, "\\node($(v))[eql] at ([shift={(-.5*\\ns, -$(y_obser-y_paramA))}]$(obser_facs[i]).west) {\$=\$};\n")
            write(iob, "\\node[xshift=.4*\\ns, yshift=.4*\\ns] at ($(v)) {\$$(lup(v))\$};\n") ## label
        end
    end
    
    write(iob, "%% --------------------------- LINKS ---------------------------\n")
    y = ytop; x = xleft
    if "paramB1" in layers
        n = length(paramB1_iniparvars)
        for (i,v) in enumerate(paramB1_iniparvars)
            write(iob, "\\draw ($(v)) -- ([shift={(0, $(divn[n][i])*\\ns)}]$(paramB1_inifac).south west);\n")
        end; x += xsep
        write(iob, "\\draw ($(paramB1_inifac)) -| ($(paramB1_parvars[1]));\n"); x += 2*xsep
        if "state1" in layers ## connect paramB1 layer to state1 layer
            for (i,v) in enumerate(state1_facs)
                write(iob, "\\draw ($(paramB1_parvars[1])) -| ([shift={(-.2*\\ns, 0)}]$(state1_facs[i]).north east);\n")
            end
        end
    end
    if "paramB2" in layers        
        y -= ysep; x = xleft
        n = length(paramB2_iniparvars)
        for (i,v) in enumerate(paramB2_iniparvars)
            write(iob, "\\draw ($(v)) -- ([shift={(0, $(divn[n][i])*\\ns)}]$(paramB2_inifac).south west);\n")
        end; x += xsep
        write(iob, "\\draw ($(paramB2_inifac)) -| ($(paramB2_parvars[1]));\n"); x += 2*xsep
        if "state1" in layers ## connect paramB2 layer to state1 layer
            for (i,v) in enumerate(state1_facs)  
                write(iob, "\\draw ($(paramB2_parvars[1])) -| ([shift={(.2*\\ns, 0)}]$(state1_facs[i]).north west);\n")
            end
        end            
    end
    if "state1" in layers ## state vars in this layer
        n = length(state1_iniparvars)
        for (i,v) in enumerate(state1_iniparvars)
            write(iob, "\\draw ($(v)) -- ([shift={(0, $(divn[n][i])*\\ns)}]$(state1_inifac).south west);\n")
        end
        write(iob, "\\draw ($(state1_inifac)) -| ($(state1_sysvars[1]));\n")
        for (i,v) in enumerate(state1_facs)
            write(iob, "\\draw ([shift={($(radius(_CFFG)),0)}]$(v).center) -- ($(state1_sysvars[i+1]).west);\n")
        end
        for (i,v) in enumerate(state1_facs)
            write(iob, "\\draw ([shift={(-$(radius(_CFFG)),0)}]$(v).center) -- ($(state1_sysvars[i]).east);\n")
        end
        if "state2" in layers ## connect state1 layer to state2 layer
            for (i,v) in enumerate(state1_sysvars[2:end])
                write(iob, "\\draw ($(v)) -- ($(state2_facs[i]));\n")
            end
        elseif "state3" in layers ## connect state1 layer to state2 layer
            for (i,v) in enumerate(state1_vars)
                write(iob, "\\draw ($(state1_vars[i])) -- ($(trans3_facs[i]));\n")
            end
        else
            for (i,v) in enumerate(state1_sysvars[2:end])
                write(iob, "\\draw ($(v)) -- ($(obser_facs[i]).north);\n")
            end
        end
    end
    if "state2" in layers
        n = length(state2_iniparvars)
        for (i,v) in enumerate(state2_iniparvars)
            write(iob, "\\draw ($(v)) -- ([shift={(0, $(divn[n][i])*\\ns)}]$(inistate2_fac).south west);\n")
        end
        for (i,v) in enumerate(state2_parvars)
            write(iob, "\\draw ($(v)) -- ([shift={(0, $(divn[1][1])*\\ns)}]$(state2_facs[i]).south west);\n")
        end
        if "state3" in layers ## connect state1 layer to state2 layer
            for (i,v) in enumerate(state2_facs)
                write(iob, "\\draw ($(state2_facs[i])) -- ($(state3_facs[i]));\n")
            end
        end
    end
    if "state3" in layers
        if "obser" in layers ## connect state1 layer to state2 layer
            for (i,v) in enumerate(state3_facs)
                write(iob, "\\draw ($(state3_facs[i])) -- ($(obser_facs[i]));\n")
            end
        end
    end
    if "obser" in layers
        for (i,v) in enumerate(obser_parvars)
            write(iob, "\\draw ($(v)) -- ([shift={(0, $(divn[1][1])*\\ns)}]$(obser_facs[i]).south west);\n")
        end
        for (i,v) in enumerate(obser_sysvars)
            write(iob, "\\draw ($(obser_facs[i])) -- ($(obser_sysvars[i]));\n")
        end
    end
    if "paramA" in layers
        n = length(paramA_iniparvars)
        for (i,v) in enumerate(paramA_iniparvars)
            ## write(iob, "\\draw[blue, line width=3pt] ($(v)) -- ([shift={(0, $(divn[length(iniparamA_vars)][i]))}]$(iniparamA_fac).south west);\n")
            write(iob, "\\draw ($(v)) -- ([shift={(0, $(divn[n][i])*\\ns)}]$(paramA_inifac).south west);\n")
        end  
        write(iob, "\\draw ($(paramA_inifac)) -| ($(paramA_parvars[1]));\n")
        for i in obser_facs
            write(iob, "\\draw ($(paramA_parvars[1])) -| ([shift={(-.5*\\ns,0)}]$(i).west) -- ($(i));\n")
        end
    end

    ## DRAW ALL RAW LINKS FOR REFERENCE
    if _raw_links_enabled
        G = Model.graph
        for edge in MetaGraphsNext.edges(G)
            source_vertex = MetaGraphsNext.label_for(G, edge.src)
            dest_vertex = MetaGraphsNext.label_for(G, edge.dst)
            write(iob, "\\draw[dotted, red, line width=1pt] ($(source_vertex)) -- ($(sanitized(string(dest_vertex))));\n")
        end
    end
    ## WRITE EQUAL NODES AGAIN
    if "paramB1" in layers
        for (i,v) in enumerate(paramB1_parvars)
            write(iob, "\\node($(v))[eql] at ($(v)) {\$=\$};\n")
        end
    end
    if "paramB2" in layers
        for (i,v) in enumerate(paramB2_parvars)
            write(iob, "\\node($(v))[eql] at ($(v)) {\$=\$};\n")
        end
    end
    if "paramA" in layers
        for (i,v) in enumerate(paramA_parvars)
            write(iob, "\\node($(v))[eql] at ($(v)) {\$=\$};\n")
        end        
    end
    ## OVERLAY BEADS
    if _CFFG
        if "paramB1" in layers write(iob, "\\draw ($(paramB1_inifac).east) circle(\\br)[fill=white];\n") end
        if "paramB2" in layers write(iob, "\\draw ($(paramB2_inifac).east) circle(\\br)[fill=white];\n") end        
        if "state1" in layers write(iob, "\\draw ($(state1_inifac).east) circle(\\br)[fill=white];\n") end
        if "paramA" in layers write(iob, "\\draw ($(paramA_inifac).east) circle(\\br)[fill=white];\n") end
        for (i,v) in enumerate(obser_facs)
            write(iob, "\\draw ($(v).north) circle(\\br)[fill=white];\n")
            write(iob, "\\draw ($(v).south)++(-1*\\br,-1*\\br) rectangle([shift={(\\br,\\br)}]$(v).south)[fill=white];\n")
        end
    end

    graph = String(take!(iob))
    tikz *= graph
    tikz *= postamble
    return tikz
end

generate_tikz (generic function with 1 method)

Coin-toss Model

In this example, we are going to perform an exact inference for a coin-toss model that can be represented as:

\[\begin{aligned} p(\theta) &= \mathrm{Beta}(\theta|a, b),\\ p(y_i|\theta) &= \mathrm{Bern}(y_i|\theta),\\ \end{aligned}\]

where \(y_i \in \{0, 1\}\) is a binary observation induced by Bernoulli likelihood while \(\theta\) is a Beta prior distribution on the parameter of Bernoulli. We are interested in inferring the posterior distribution of \(\theta\).

The generative model is:

\[ \begin{aligned} p(\boldsymbol{y},\theta) &= p(\boldsymbol{y} \mid \theta) \cdot p(\theta) \\ &= \prod_{i=1}^N{p(y_i \mid \theta)} \cdot p(\theta) \\ &= \prod_{i=1}^N{\mathrm{Bern}(y_i \mid \theta)} \cdot \mathrm{Beta}(\theta \mid a,b) \end{aligned} \]

See Fig 4.1 in 2023_Bagaev_RPPfSBI_Thesis. So, \(N\) \(\delta\) factors, 1 Beta factor, \(N\) Ber factors.

@model function coin_model(y, a, b)
    ## We endow θ parameter of our model with some prior
    θ ~ Beta(a, b)

    ## We assume that outcome of each coin flip is governed by the Bernoulli distribution
    for i in eachindex(y)
        y[i] ~ Bernoulli(θ)
    end
end

coin_conditioned = coin_model(a = 2.0, b = 7.0) | (y = [ true, false, true ], )
coin_rxi_model = RxInfer.create_model(coin_conditioned)
coin_gppl_model = RxInfer.getmodel(coin_rxi_model)
coin_meta_graph = coin_gppl_model.graph

Meta graph based on a Graphs.SimpleGraphs.SimpleGraph{Int64} with vertex labels of type GraphPPL.NodeLabel, vertex metadata of type GraphPPL.NodeData, edge metadata of type GraphPPL.EdgeLabel, graph metadata given by GraphPPL.Context(0, coin_model, "", nothing, {Bernoulli = 3, Beta = 1}, {}, {(Bernoulli, 3) = Bernoulli_12, (Bernoulli, 1) = Bernoulli_8, (Bernoulli, 2) = Bernoulli_10, (Beta, 1) = Beta_6}, {:constvar_2 = constvar_2_3, :θ = θ_1, :constvar_4 = constvar_4_5}, {:y = ResizableArray{GraphPPL.NodeLabel,1}(GraphPPL.NodeLabel[y_7, y_9, y_11])}, {}, {}, Base.RefValue{Any}(nothing)), and default weight 1.0

for i in MetaGraphsNext.vertices(coin_meta_graph)
    label = MetaGraphsNext.label_for(coin_meta_graph, i)
    println("$i: $label")
end

1: θ_1
2: constvar_2_3
3: constvar_4_5
4: Beta_6
5: y_7
6: Bernoulli_8
7: y_9
8: Bernoulli_10
9: y_11
10: Bernoulli_12

##
## fieldnames(GraphPPL.Context)
coin_gppl_model_context = GraphPPL.getcontext(coin_gppl_model)
coin_gppl_model_context.individual_variables
coin_gppl_model_context.vector_variables
coin_gppl_model_context.factor_nodes

4-element Dictionaries.UnorderedDictionary{GraphPPL.FactorID, GraphPPL.NodeLabel}
 (Bernoulli, 3) │ Bernoulli_12
 (Bernoulli, 1) │ Bernoulli_8
 (Bernoulli, 2) │ Bernoulli_10
      (Beta, 1) │ Beta_6

show_meta_graph(coin_meta_graph)

    θ_1;
    constvar_2_3;
    constvar_4_5;
    Beta_6;
    y_7;
    Bernoulli_8;
    y_9;
    Bernoulli_10;
    y_11;
    Bernoulli_12;
    θ_1 -- Beta_6;
    θ_1 -- Bernoulli_8;
    θ_1 -- Bernoulli_10;
    θ_1 -- Bernoulli_12;
    constvar_2_3 -- Beta_6;
    constvar_4_5 -- Beta_6;
    y_7 -- Bernoulli_8;
    y_9 -- Bernoulli_10;
    y_11 -- Bernoulli_12;

_san_node_ids = Dict("dummy" => "dummy")
_lup_enabled = true
_raw_links_enabled = false
_CFFG = true
_lup_dict = Dict{String, String}(
    "constvar_2_3" => "a",
    "constvar_4_5" => "b",
    "Beta_6" => "\\mathcal{B}eta",
    "Bernoulli_8" => "\\mathcal{B}er",
    "Bernoulli_10" => "\\mathcal{B}er",
    "Bernoulli_12" => "\\mathcal{B}er",
    "θ_1" => "\\theta",
    "y_7" => "y_1",
    "y_9" => "y_2",
    "y_11" => "y_3",
)
coin_tikz = generate_tikz(
    heading = "Coin Toss",
    Model = coin_gppl_model,
    layers = [ ## fac layers, i.e. not var (vars are handled with assoced layer)
        ## "cntrl", 
        ## "paramB1",
        ## "paramB2",
        ## "state1",
        ## "state2",
        ## "state3", 
        "obser", 
        ## "prefr", 
        "paramA"
    ],
    obser_sysvars = ["y_7", "y_9", "y_11"],
    obser_facs = ["Bernoulli_8", "Bernoulli_10", "Bernoulli_12"],
    obser_parvars = [],
    
    paramA_iniparvars = ["constvar_2_3", "constvar_4_5"],
    paramA_inifac = "Beta_6",
    paramA_parvars = ["θ_1"]
)
print(coin_tikz)
show_tikz(coin_tikz) ## write to file rather than show in notebook

\documentclass{standalone}
\usepackage{tikz}
\usetikzlibrary{calc}
\usepackage{amsmath}
\newcommand\ns{15mm} %% define node (minimum) size
\newcommand\nr{.5*\ns} %% define node radius
\newcommand\ar{.4*\ns} %% define arc radius
\newcommand\br{.19*\ns} %% define bead radius
\DeclareUnicodeCharacter{03B8}{\theta}
\DeclareUnicodeCharacter{03B3}{\gamma}
\begin{document}
\begin{tikzpicture}[
    scale=1, 
    draw=black, 
    fac/.style={rectangle, draw=black!100, fill=black!0, minimum size=\ns, inner sep=0pt},
    dlt/.style={rectangle, draw=black!100, fill=black!100, minimum size=.2*\ns},
    var/.style={circle, draw=gray!100, fill=gray!100, minimum size=1mm, inner sep=0pt},
    bus/.style={rectangle, draw=black!100, fill=black!100, minimum width=30mm, minimum height=.5mm, inner sep=0pt},
    eql/.style={rectangle, draw=black!100, fill=black!0, minimum size=.5*\ns},
    text=black,
    ]
\node[font=\Large] at (4,17) {\textbf{Coin Toss}};
⋮
⋮
⋮
\path[use as bounding box] ([xshift=-1cm, yshift=-1cm]current bounding box.south west) rectangle ([xshift=1cm, yshift=2cm]current bounding box.north east);
\end{tikzpicture}
\end{document}
\documentclass{standalone}
\usepackage{tikz}
\usetikzlibrary{calc}
\usepackage{amsmath}
\newcommand\ns{15mm} %% define node (minimum) size
\newcommand\nr{.5*\ns} %% define node radius
\newcommand\ar{.4*\ns} %% define arc radius
\newcommand\br{.19*\ns} %% define bead radius
\DeclareUnicodeCharacter{03B8}{\theta}
\DeclareUnicodeCharacter{03B3}{\gamma}
\begin{document}
\begin{tikzpicture}[
    scale=1, 
    draw=black, 
    fac/.style={rectangle, draw=black!100, fill=black!0, minimum size=\ns, inner sep=0pt},
    dlt/.style={rectangle, draw=black!100, fill=black!100, minimum size=.2*\ns},
    var/.style={circle, draw=gray!100, fill=gray!100, minimum size=1mm, inner sep=0pt},
    bus/.style={rectangle, draw=black!100, fill=black!100, minimum width=30mm, minimum height=.5mm, inner sep=0pt},
    eql/.style={rectangle, draw=black!100, fill=black!0, minimum size=.5*\ns},
    text=black,
    ]
\node[font=\Large] at (4,17) {\textbf{Coin Toss}};
%% --------------------------- NODES ---------------------------
\node(Bernoulli_8)[fac] at (8, 15) {$\mathcal{B}er$};
\node(Bernoulli_10)[fac] at (12, 15) {$\mathcal{B}er$};
\node(Bernoulli_12)[fac] at (16, 15) {$\mathcal{B}er$};
\node(y_7)[var] at (8, 13) {};
\node[above right] at (y_7) {$y_1$};
\node(y_7)[dlt] at ($ (y_7.center) + (0, -.4*\ns) $) {};
\node(y_9)[var] at (12, 13) {};
\node[above right] at (y_9) {$y_2$};
\node(y_9)[dlt] at ($ (y_9.center) + (0, -.4*\ns) $) {};
\node(y_11)[var] at (16, 13) {};
\node[above right] at (y_11) {$y_3$};
\node(y_11)[dlt] at ($ (y_11.center) + (0, -.4*\ns) $) {};
\node(Beta_6)[fac] at (2, 10) {$\mathcal{B}eta$};
\node(constvar_2_3)[var] at ($ (Beta_6.south west) + (-.5*\ns, 0.3333333333333333*\ns) $) {};
\node[above] at (constvar_2_3) {$a$};
\node(constvar_2_3)[dlt] at ($ (constvar_2_3.center) + (-.4*\ns, 0) $) {};
\node(constvar_4_5)[var] at ($ (Beta_6.south west) + (-.5*\ns, 0.6666666666666666*\ns) $) {};
\node[above] at (constvar_4_5) {$b$};
\node(constvar_4_5)[dlt] at ($ (constvar_4_5.center) + (-.4*\ns, 0) $) {};
\node(θ_1)[eql] at ([shift={(-.5*\ns, -5)}]Bernoulli_8.west) {$=$};
\node[xshift=.4*\ns, yshift=.4*\ns] at (θ_1) {$\theta$};
%% --------------------------- LINKS ---------------------------
\draw (Bernoulli_8) -- (y_7);
\draw (Bernoulli_10) -- (y_9);
\draw (Bernoulli_12) -- (y_11);
\draw (constvar_2_3) -- ([shift={(0, 0.3333333333333333*\ns)}]Beta_6.south west);
\draw (constvar_4_5) -- ([shift={(0, 0.6666666666666666*\ns)}]Beta_6.south west);
\draw (Beta_6) -| (θ_1);
\draw (θ_1) -| ([shift={(-.5*\ns,0)}]Bernoulli_8.west) -- (Bernoulli_8);
\draw (θ_1) -| ([shift={(-.5*\ns,0)}]Bernoulli_10.west) -- (Bernoulli_10);
\draw (θ_1) -| ([shift={(-.5*\ns,0)}]Bernoulli_12.west) -- (Bernoulli_12);
\node(θ_1)[eql] at (θ_1) {$=$};
\draw (Beta_6.east) circle(\br)[fill=white];
\draw (Bernoulli_8.north) circle(\br)[fill=white];
\draw (Bernoulli_8.south)++(-1*\br,-1*\br) rectangle([shift={(\br,\br)}]Bernoulli_8.south)[fill=white];
\draw (Bernoulli_10.north) circle(\br)[fill=white];
\draw (Bernoulli_10.south)++(-1*\br,-1*\br) rectangle([shift={(\br,\br)}]Bernoulli_10.south)[fill=white];
\draw (Bernoulli_12.north) circle(\br)[fill=white];
\draw (Bernoulli_12.south)++(-1*\br,-1*\br) rectangle([shift={(\br,\br)}]Bernoulli_12.south)[fill=white];
\path[use as bounding box] ([xshift=-1cm, yshift=-1cm]current bounding box.south west) rectangle ([xshift=1cm, yshift=2cm]current bounding box.north east);
\end{tikzpicture}
\end{document}

3465

Coin Toss

Hidden Markov Model with Control (Part 1)

There are 1-hot encodings on all random variables. Using regular bold font to indicate the 1-hot encodings, and also setting \(t=1\), we have:

\[\begin{aligned} p(\boldsymbol{x}, \boldsymbol{s}) &\propto \underbrace{p(\mathbf s_0)}_{\substack{\text{Initial} \\ \text{state} \\ \text{prior}}} \prod_{k=1}^{T+1} \underbrace{p(\mathbf x_k \mid \mathbf s_k)}_{\substack{\text{Observation} \\ \text{model}}} \, \underbrace{p(\mathbf s_k \mid \mathbf s_{k-1}}_{\substack{\text{Transition} \\ \text{model}}}) \end{aligned}\]

where

\[\begin{aligned} p(\mathbf s_{0}) &= \mathcal{Cat}(\mathbf s_{0} ∣ \mathbf d) \\ p(\mathbf x_k ∣ \mathbf s_k) &= \mathcal{Cat}(\mathbf x_k ∣ \mathbf A \mathbf s_k) \\ p(\mathbf s_k ∣ \mathbf s_{k-1}) &= \mathcal{Cat}(\mathbf s_k ∣ \mathbf B \mathbf s_{k-1}) \\ \end{aligned}\]

@model function hidden_markov_model(x)
    B ~ MatrixDirichlet(ones(3, 3))
    A ~ MatrixDirichlet([10.0 1.0 1.0; 
                         1.0 10.0 1.0; 
                         1.0 1.0 10.0 ])
    d = fill(1.0/3.0, 3)
    s₀ ~ Categorical(d)
    
    sₖ₋₁ = s₀
    for k in eachindex(x)
        s[k] ~ Transition(sₖ₋₁, B)
        x[k] ~ Transition(s[k], A)
        sₖ₋₁ = s[k]
    end
end

hmm_conditioned = hidden_markov_model() | (x = [[1.0, 0.0, 0.0], [0.0, 0.0, 1.0]],)
hmm_rxi_model = RxInfer.create_model(hmm_conditioned)
hmm_gppl_model = RxInfer.getmodel(hmm_rxi_model)
hmm_meta_graph = hmm_gppl_model.graph

Meta graph based on a Graphs.SimpleGraphs.SimpleGraph{Int64} with vertex labels of type GraphPPL.NodeLabel, vertex metadata of type GraphPPL.NodeData, edge metadata of type GraphPPL.EdgeLabel, graph metadata given by GraphPPL.Context(0, hidden_markov_model, "", nothing, {Categorical{P} where P<:Real = 1, Transition = 4, MatrixDirichlet = 2}, {}, {(Transition, 4) = Transition_20, (MatrixDirichlet, 2) = MatrixDirichlet_8, (Transition, 3) = Transition_18, (MatrixDirichlet, 1) = MatrixDirichlet_4, (Transition, 2) = Transition_16, (Categorical{P} where P<:Real, 1) = Categorical{P} where P<:Real_12, (Transition, 1) = Transition_14}, {:A = A_5, :s₀ = s₀_9, :B = B_1, :constvar_6 = constvar_6_7, :constvar_10 = constvar_10_11, :constvar_2 = constvar_2_3}, {:s = ResizableArray{GraphPPL.NodeLabel,1}(GraphPPL.NodeLabel[s_13, s_17]), :x = ResizableArray{GraphPPL.NodeLabel,1}(GraphPPL.NodeLabel[x_15, x_19])}, {}, {}, Base.RefValue{Any}(nothing)), and default weight 1.0

for i in MetaGraphsNext.vertices(hmm_meta_graph)
    label = MetaGraphsNext.label_for(hmm_meta_graph, i)
    println("$i: $label")
end

1: B_1
2: constvar_2_3
3: MatrixDirichlet_4
4: A_5
5: constvar_6_7
6: MatrixDirichlet_8
7: s₀_9
8: constvar_10_11
9: Categorical{P} where P<:Real_12
10: s_13
11: Transition_14
12: x_15
13: Transition_16
14: s_17
15: Transition_18
16: x_19
17: Transition_20

##
## fieldnames(GraphPPL.Context)
hmm_gppl_model_context = GraphPPL.getcontext(hmm_gppl_model)
hmm_gppl_model_context.individual_variables
hmm_gppl_model_context.vector_variables
hmm_gppl_model_context.factor_nodes

7-element Dictionaries.UnorderedDictionary{GraphPPL.FactorID, GraphPPL.NodeLabel}
                   (Transition, 4) │ Transition_20
              (MatrixDirichlet, 2) │ MatrixDirichlet_8
                   (Transition, 3) │ Transition_18
              (MatrixDirichlet, 1) │ MatrixDirichlet_4
                   (Transition, 2) │ Transition_16
 (Categorical{P} where P<:Real, 1) │ Categorical{P} where P<:Real_12
                   (Transition, 1) │ Transition_14

show_meta_graph(hmm_meta_graph)

    B_1;
    constvar_2_3;
    MatrixDirichlet_4;
    A_5;
    constvar_6_7;
    MatrixDirichlet_8;
    s₀_9;
    constvar_10_11;
    Categorical{P} where P<:Real_12;
    s_13;
    Transition_14;
    x_15;
    Transition_16;
    s_17;
    Transition_18;
    x_19;
    Transition_20;
    B_1 -- MatrixDirichlet_4;
    B_1 -- Transition_14;
    B_1 -- Transition_18;
    constvar_2_3 -- MatrixDirichlet_4;
    A_5 -- MatrixDirichlet_8;
    A_5 -- Transition_16;
    A_5 -- Transition_20;
    constvar_6_7 -- MatrixDirichlet_8;
    s₀_9 -- Categorical{P} where P<:Real_12;
    s₀_9 -- Transition_14;
    constvar_10_11 -- Categorical{P} where P<:Real_12;
    s_13 -- Transition_14;
    s_13 -- Transition_16;
    s_13 -- Transition_18;
    x_15 -- Transition_16;
    s_17 -- Transition_18;
    s_17 -- Transition_20;
    x_19 -- Transition_20;

## _san_node_ids = Dict("dummy" => "dummy")
_san_node_ids = Dict("Categorical{P} where P<:Real_12" => "Categorical_12")
_lup_enabled = true
_raw_links_enabled = false
_CFFG = true
_lup_dict = Dict{String, String}(
    "MatrixDirichlet_4" => "\\mathcal{M}Dir",
    "B_1" => "\\mathbf{B}",
    "constvar_10_11" => "\\mathbf{d}",
    "Categorical_12" => "\\mathcal{C}at",
    "Transition_14" => "\\mathcal{T}",
    "Transition_18" => "\\mathcal{T}",
    "Transition_16" => "\\mathcal{T}",
    "Transition_20" => "\\mathcal{T}",
    "MatrixDirichlet_8" => "\\mathcal{M}Dir",
    "A_5" => "\\mathbf{A}",
    "constvar_2_3" => "\\mathbf{B_0}",
    "constvar_6_7" => "\\mathbf{A_0}",
    ## "s₀_9" => "s_0", ## s₀ must always be in _lup_dict
    "s₀_9" => "\\mathbf{s}_0",
    "s_13" => "\\mathbf{s}_1",
    "s_17" => "\\mathbf{s}_2",
    "x_15" => "\\mathbf{x}_1",
    "x_19" => "\\mathbf{x}_2",
)
hmm_tikz = generate_tikz(
    heading = "Hidden Markov Model",
    Model = hmm_gppl_model,
    layers = [ ## fac layers, i.e. not var (vars are handled with assoced layer)
        ## "cntrl", 
        "paramB1",
        ## "paramB2",
        "state1",
        ## "state2",
        ## "state3", 
        "obser", 
        ## "prefr", 
        "paramA"
    ],
    paramB1_iniparvars = ["constvar_2_3"],
    paramB1_inifac = "MatrixDirichlet_4",
    paramB1_parvars = ["B_1"],

    state1_iniparvars = ["constvar_10_11"],
    state1_inifac = sanitized("Categorical{P} where P<:Real_12"),
    state1_sysvars = ["s₀_9", "s_13", "s_17"],
    state1_facs = ["Transition_14", "Transition_18"],

    obser_sysvars = ["x_15", "x_19"],
    obser_facs = ["Transition_16", "Transition_20"],
    obser_parvars = [],

    paramA_iniparvars = ["constvar_6_7"],
    paramA_inifac = "MatrixDirichlet_8",
    paramA_parvars = ["A_5"]
)
print(hmm_tikz)
show_tikz(hmm_tikz) ## write to file rather than show in notebook

\documentclass{standalone}
\usepackage{tikz}
\usetikzlibrary{calc}
\usepackage{amsmath}
\newcommand\ns{15mm} %% define node (minimum) size
\newcommand\nr{.5*\ns} %% define node radius
\newcommand\ar{.4*\ns} %% define arc radius
\newcommand\br{.19*\ns} %% define bead radius
\DeclareUnicodeCharacter{03B8}{\theta}
\DeclareUnicodeCharacter{03B3}{\gamma}
\begin{document}
\begin{tikzpicture}[
    scale=1, 
    draw=black, 
    fac/.style={rectangle, draw=black!100, fill=black!0, minimum size=\ns, inner sep=0pt},
    dlt/.style={rectangle, draw=black!100, fill=black!100, minimum size=.2*\ns},
    var/.style={circle, draw=gray!100, fill=gray!100, minimum size=1mm, inner sep=0pt},
    bus/.style={rectangle, draw=black!100, fill=black!100, minimum width=30mm, minimum height=.5mm, inner sep=0pt},
    eql/.style={rectangle, draw=black!100, fill=black!0, minimum size=.5*\ns},
    text=black,
    ]
\node[font=\Large] at (4,17) {\textbf{Hidden Markov Model}};
⋮
⋮
⋮
\path[use as bounding box] ([xshift=-1cm, yshift=-1cm]current bounding box.south west) rectangle ([xshift=1cm, yshift=2cm]current bounding box.north east);
\end{tikzpicture}
\end{document}
\documentclass{standalone}
\usepackage{tikz}
\usetikzlibrary{calc}
\usepackage{amsmath}
\newcommand\ns{15mm} %% define node (minimum) size
\newcommand\nr{.5*\ns} %% define node radius
\newcommand\ar{.4*\ns} %% define arc radius
\newcommand\br{.19*\ns} %% define bead radius
\DeclareUnicodeCharacter{03B8}{\theta}
\DeclareUnicodeCharacter{03B3}{\gamma}
\begin{document}
\begin{tikzpicture}[
    scale=1, 
    draw=black, 
    fac/.style={rectangle, draw=black!100, fill=black!0, minimum size=\ns, inner sep=0pt},
    dlt/.style={rectangle, draw=black!100, fill=black!100, minimum size=.2*\ns},
    var/.style={circle, draw=gray!100, fill=gray!100, minimum size=1mm, inner sep=0pt},
    bus/.style={rectangle, draw=black!100, fill=black!100, minimum width=30mm, minimum height=.5mm, inner sep=0pt},
    eql/.style={rectangle, draw=black!100, fill=black!0, minimum size=.5*\ns},
    text=black,
    ]
\node[font=\Large] at (4,17) {\textbf{Hidden Markov Model}};
%% --------------------------- NODES ---------------------------
\node(MatrixDirichlet_4)[fac] at (2, 15) {$\mathcal{M}Dir$};
\node(constvar_2_3)[var] at ($ (MatrixDirichlet_4.south west) + (-.5*\ns, 0.5*\ns) $) {};
\node[above] at (constvar_2_3) {$\mathbf{B_0}$};
\node(constvar_2_3)[dlt] at ($ (constvar_2_3.center) + (-.4*\ns, 0) $) {};
\node(B_1)[eql] at (4, 15) {$=$};
\node[xshift=.4*\ns, yshift=.4*\ns] at (B_1) {$\mathbf{B}$};
\node(Categorical_12)[fac] at (2, 12) {$\mathcal{C}at$};
\node(constvar_10_11)[var] at ($ (Categorical_12.south west) + (-.5*\ns, 0.5*\ns) $) {};
\node[above] at (constvar_10_11) {$\mathbf{d}$};
\node(constvar_10_11)[dlt] at ($ (constvar_10_11.center) + (-.4*\ns, 0) $) {};
\node(s₀_9)[var] at (4, 12) {};
\node[above] at (s₀_9) {$\mathbf{s}_0$};
\node(Transition_14)[fac] at (6, 12) {$\mathcal{T}$};
\draw (Transition_14)++(\ar,0mm) arc[start angle=0, end angle=-180, radius=\ar];
\node(Transition_18)[fac] at (10, 12) {$\mathcal{T}$};
\draw (Transition_18)++(\ar,0mm) arc[start angle=0, end angle=-180, radius=\ar];
\node(s_13)[eql] at (8, 12) {$=$};
\node[xshift=.4*\ns, yshift=.4*\ns] at (s_13) {$\mathbf{s}_1$};
\node(s_17)[var] at (12, 12) {};
\node[above] at (s_17) {$\mathbf{s}_2$};
\node(Transition_16)[fac] at (8, 9) {$\mathcal{T}$};
\node(Transition_20)[fac] at (12, 9) {$\mathcal{T}$};
\node(x_15)[var] at (8, 7) {};
\node[above right] at (x_15) {$\mathbf{x}_1$};
\node(x_15)[dlt] at ($ (x_15.center) + (0, -.4*\ns) $) {};
\node(x_19)[var] at (12, 7) {};
\node[above right] at (x_19) {$\mathbf{x}_2$};
\node(x_19)[dlt] at ($ (x_19.center) + (0, -.4*\ns) $) {};
\node(MatrixDirichlet_8)[fac] at (2, 4) {$\mathcal{M}Dir$};
\node(constvar_6_7)[var] at ($ (MatrixDirichlet_8.south west) + (-.5*\ns, 0.5*\ns) $) {};
\node[above] at (constvar_6_7) {$\mathbf{A_0}$};
\node(constvar_6_7)[dlt] at ($ (constvar_6_7.center) + (-.4*\ns, 0) $) {};
\node(A_5)[eql] at ([shift={(-.5*\ns, -5)}]Transition_16.west) {$=$};
\node[xshift=.4*\ns, yshift=.4*\ns] at (A_5) {$\mathbf{A}$};
%% --------------------------- LINKS ---------------------------
\draw (constvar_2_3) -- ([shift={(0, 0.5*\ns)}]MatrixDirichlet_4.south west);
\draw (MatrixDirichlet_4) -| (B_1);
\draw (B_1) -| ([shift={(-.2*\ns, 0)}]Transition_14.north east);
\draw (B_1) -| ([shift={(-.2*\ns, 0)}]Transition_18.north east);
\draw (constvar_10_11) -- ([shift={(0, 0.5*\ns)}]Categorical_12.south west);
\draw (Categorical_12) -| (s₀_9);
\draw ([shift={(\ar,0)}]Transition_14.center) -- (s_13.west);
\draw ([shift={(\ar,0)}]Transition_18.center) -- (s_17.west);
\draw ([shift={(-\ar,0)}]Transition_14.center) -- (s₀_9.east);
\draw ([shift={(-\ar,0)}]Transition_18.center) -- (s_13.east);
\draw (s_13) -- (Transition_16.north);
\draw (s_17) -- (Transition_20.north);
\draw (Transition_16) -- (x_15);
\draw (Transition_20) -- (x_19);
\draw (constvar_6_7) -- ([shift={(0, 0.5*\ns)}]MatrixDirichlet_8.south west);
\draw (MatrixDirichlet_8) -| (A_5);
\draw (A_5) -| ([shift={(-.5*\ns,0)}]Transition_16.west) -- (Transition_16);
\draw (A_5) -| ([shift={(-.5*\ns,0)}]Transition_20.west) -- (Transition_20);
\node(B_1)[eql] at (B_1) {$=$};
\node(A_5)[eql] at (A_5) {$=$};
\draw (MatrixDirichlet_4.east) circle(\br)[fill=white];
\draw (Categorical_12.east) circle(\br)[fill=white];
\draw (MatrixDirichlet_8.east) circle(\br)[fill=white];
\draw (Transition_16.north) circle(\br)[fill=white];
\draw (Transition_16.south)++(-1*\br,-1*\br) rectangle([shift={(\br,\br)}]Transition_16.south)[fill=white];
\draw (Transition_20.north) circle(\br)[fill=white];
\draw (Transition_20.south)++(-1*\br,-1*\br) rectangle([shift={(\br,\br)}]Transition_20.south)[fill=white];
\path[use as bounding box] ([xshift=-1cm, yshift=-1cm]current bounding box.south west) rectangle ([xshift=1cm, yshift=2cm]current bounding box.north east);
\end{tikzpicture}
\end{document}

4802

Hidden Markov Model

Sales Forecasting with Time-Varying Autoregressive Models (Part 4)

Our model for the Variational Bayesian Inference of a latent autoregressive model can be represented mathematically as follows:

\[\begin{aligned} p(\gamma) &= \mathcal{Gam}(\gamma|a, b)\\ p(\boldsymbol{\theta}) &= \mathcal{N}(\boldsymbol{\theta}|\mathbf{m}_{\theta}, \mathbf{V}_{\theta})\\ p(\mathbf{s}_t \mid \mathbf{s}_{t-1}, \boldsymbol{\theta}) &= \mathcal{N}\left(\mathbf{s}_t \mid \mathbf{B}\left(\boldsymbol{\theta}\right) \mathbf{s}_{t-1},\, \mathbf{V}(\gamma)\right)\\ p(x_t \mid \mathbf{s}_t) &= \mathcal{N}(x_t \mid \mathbf{s}_t, \tau^{-1}) \end{aligned}\]

where

• \(M^{AR}\) is the order of the autoregression process

• \(\boldsymbol{\theta}=(\theta_1, \theta_2, ..., \theta_{M^{AR}})^{\top}\) is a vector of transition coefficients

• \(\mathbf{s}_{t-1} = (s_{t-1}, s_{t-2}, ..., s_{t-M^{AR}})^\top\) is the previous state

• \(\mathbf{s}_t = (s_{t}, s_{t-1}, ..., s_{t-M^{AR}+1})^\top\) is the current state

• \(\gamma\) is the transition precision

• \(\tau\) is the observation precision

• \(\mathbf{v}_u = (1,0,...,0)^\top\) is an \(M^{AR}\)-dimensional unit vector selecting the first component of the vector that follows

• \(\mathbf{V}(\gamma) = (1/ \gamma)\mathbf{v}_u\mathbf{v}_u^\top\) is the transition covariance

• \(\mathbf{B}(\boldsymbol{\theta}) = \begin{bmatrix} \boldsymbol{\qquad\theta}^\top \\ \mathbf{I}_{M^{AR}-1} & \mathbf{0} \end{bmatrix}\) is the transition matrix

The following diagram represents an autoregressive process:

Autoregressive process

@model function lar_model(
    x, ##. data/observations 
    𝚃ᴬᴿ, ##. Uni/Multi variate 
    Mᴬᴿ, ##. AR order
    vᵤ, ##. unit vector 
    τ) ##. observation precision     
        ## Priors
        γ  ~ Gamma(α = 1.0, β = 1.0) ##. for transition precision    
        if 𝚃ᴬᴿ === Multivariate
            θ  ~ MvNormal(μ = zeros(Mᴬᴿ), Λ = diageye(Mᴬᴿ)) ##.kw μ,Λ only work inside macro
            s₀ ~ MvNormal(μ = zeros(Mᴬᴿ), Λ = diageye(Mᴬᴿ)) ##.kw μ,Λ only work inside macro
        else ## Univariate
            θ  ~ Normal(μ = 0.0, γ = 1.0)
            s₀ ~ Normal(μ = 0.0, γ = 1.0)
        end
        sₜ₋₁ = s₀
        for t in eachindex(x)
            s[t] ~ AR(sₜ₋₁, θ, γ) #.Eq (2b)
            if 𝚃ᴬᴿ === Multivariate
                x[t] ~ Normal(μ = dot(vᵤ, s[t]), γ = τ) #.Eq (2c)
            else
                x[t] ~ Normal(μ = vᵤ*s[t], γ = τ) #.Eq (2c)
            end
            sₜ₋₁ = s[t]
        end
    end

𝚃ᴬᴿ = Univariate
m = 1
τ̃ = 0.001 ## assumed observation precision
lar_conditioned = lar_model(
    𝚃ᴬᴿ=𝚃ᴬᴿ, 
    Mᴬᴿ=m, 
    vᵤ=ReactiveMP.ar_unit(𝚃ᴬᴿ, m), 
    τ=τ̃
) | (x = [266.0, 145.0, 183.0],)

lar_model(𝚃ᴬᴿ = Univariate, Mᴬᴿ = 1, vᵤ = 1.0, τ = 0.001) conditioned on: 
  x = [266.0, 145.0, 183.0]

lar_rxi_model = RxInfer.create_model(lar_conditioned)
lar_gppl_model = RxInfer.getmodel(lar_rxi_model)
lar_meta_graph = lar_gppl_model.graph

Meta graph based on a Graphs.SimpleGraphs.SimpleGraph{Int64} with vertex labels of type GraphPPL.NodeLabel, vertex metadata of type GraphPPL.NodeData, edge metadata of type GraphPPL.EdgeLabel, graph metadata given by GraphPPL.Context(0, lar_model, "", nothing, {AR = 3, NormalMeanPrecision = 5, * = 3, GammaShapeRate = 1}, {}, {(NormalMeanPrecision, 4) = NormalMeanPrecision_38, (AR, 1) = AR_20, (*, 2) = *_34, (GammaShapeRate, 1) = GammaShapeRate_6, (NormalMeanPrecision, 2) = NormalMeanPrecision_18, (NormalMeanPrecision, 1) = NormalMeanPrecision_12, (AR, 3) = AR_40, (*, 3) = *_44, (NormalMeanPrecision, 5) = NormalMeanPrecision_48, (NormalMeanPrecision, 3) = NormalMeanPrecision_28, (*, 1) = *_24, (AR, 2) = AR_30}, {:γ = γ_1, :constvar_42 = constvar_42_43, :constvar_14 = constvar_14_15, :constvar_16 = constvar_16_17, :anonymous_var_graphppl = anonymous_var_graphppl_41, :constvar_22 = constvar_22_23, :constvar_26 = constvar_26_27, :s₀ = s₀_13, :constvar_32 = constvar_32_33, :constvar_36 = constvar_36_37, :constvar_46 = constvar_46_47, :constvar_8 = constvar_8_9, :constvar_2 = constvar_2_3, :θ = θ_7, :constvar_4 = constvar_4_5, :constvar_10 = constvar_10_11}, {:s = ResizableArray{GraphPPL.NodeLabel,1}(GraphPPL.NodeLabel[s_19, s_29, s_39]), :x = ResizableArray{GraphPPL.NodeLabel,1}(GraphPPL.NodeLabel[x_25, x_35, x_45])}, {}, {}, Base.RefValue{Any}(nothing)), and default weight 1.0

for i in MetaGraphsNext.vertices(lar_meta_graph)
    label = MetaGraphsNext.label_for(lar_meta_graph, i)
    println("$i: $label")
end

1: γ_1
2: constvar_2_3
3: constvar_4_5
4: GammaShapeRate_6
5: θ_7
6: constvar_8_9
7: constvar_10_11
8: NormalMeanPrecision_12
9: s₀_13
10: constvar_14_15
11: constvar_16_17
12: NormalMeanPrecision_18
13: s_19
14: AR_20
15: anonymous_var_graphppl_21
16: constvar_22_23
17: *_24
18: x_25
19: constvar_26_27
20: NormalMeanPrecision_28
21: s_29
22: AR_30
23: anonymous_var_graphppl_31
24: constvar_32_33
25: *_34
26: x_35
27: constvar_36_37
28: NormalMeanPrecision_38
29: s_39
30: AR_40
31: anonymous_var_graphppl_41
32: constvar_42_43
33: *_44
34: x_45
35: constvar_46_47
36: NormalMeanPrecision_48

##
## fieldnames(GraphPPL.Context)
lar_gppl_model_context = GraphPPL.getcontext(lar_gppl_model)
lar_gppl_model_context.individual_variables
lar_gppl_model_context.vector_variables
lar_gppl_model_context.factor_nodes

12-element Dictionaries.UnorderedDictionary{GraphPPL.FactorID, GraphPPL.NodeLabel}
 (NormalMeanPrecision, 4) │ NormalMeanPrecision_38
                  (AR, 1) │ AR_20
                   (*, 2) │ *_34
      (GammaShapeRate, 1) │ GammaShapeRate_6
 (NormalMeanPrecision, 2) │ NormalMeanPrecision_18
 (NormalMeanPrecision, 1) │ NormalMeanPrecision_12
                  (AR, 3) │ AR_40
                   (*, 3) │ *_44
 (NormalMeanPrecision, 5) │ NormalMeanPrecision_48
 (NormalMeanPrecision, 3) │ NormalMeanPrecision_28
                   (*, 1) │ *_24
                  (AR, 2) │ AR_30

show_meta_graph(lar_meta_graph)

    γ_1;
    constvar_2_3;
    constvar_4_5;
    GammaShapeRate_6;
    θ_7;
    constvar_8_9;
    constvar_10_11;
    NormalMeanPrecision_12;
    s₀_13;
    constvar_14_15;
    constvar_16_17;
    NormalMeanPrecision_18;
    s_19;
    AR_20;
    anonymous_var_graphppl_21;
    constvar_22_23;
    *_24;
    x_25;
    constvar_26_27;
    NormalMeanPrecision_28;
    s_29;
    AR_30;
    anonymous_var_graphppl_31;
    constvar_32_33;
    *_34;
    x_35;
    constvar_36_37;
    NormalMeanPrecision_38;
    s_39;
    AR_40;
    anonymous_var_graphppl_41;
    constvar_42_43;
    *_44;
    x_45;
    constvar_46_47;
    NormalMeanPrecision_48;
    γ_1 -- GammaShapeRate_6;
    γ_1 -- AR_20;
    γ_1 -- AR_30;
    γ_1 -- AR_40;
    constvar_2_3 -- GammaShapeRate_6;
    constvar_4_5 -- GammaShapeRate_6;
    θ_7 -- NormalMeanPrecision_12;
    θ_7 -- AR_20;
    θ_7 -- AR_30;
    θ_7 -- AR_40;
    constvar_8_9 -- NormalMeanPrecision_12;
    constvar_10_11 -- NormalMeanPrecision_12;
    s₀_13 -- NormalMeanPrecision_18;
    s₀_13 -- AR_20;
    constvar_14_15 -- NormalMeanPrecision_18;
    constvar_16_17 -- NormalMeanPrecision_18;
    s_19 -- AR_20;
    s_19 -- *_24;
    s_19 -- AR_30;
    anonymous_var_graphppl_21 -- *_24;
    anonymous_var_graphppl_21 -- NormalMeanPrecision_28;
    constvar_22_23 -- *_24;
    x_25 -- NormalMeanPrecision_28;
    constvar_26_27 -- NormalMeanPrecision_28;
    s_29 -- AR_30;
    s_29 -- *_34;
    s_29 -- AR_40;
    anonymous_var_graphppl_31 -- *_34;
    anonymous_var_graphppl_31 -- NormalMeanPrecision_38;
    constvar_32_33 -- *_34;
    x_35 -- NormalMeanPrecision_38;
    constvar_36_37 -- NormalMeanPrecision_38;
    s_39 -- AR_40;
    s_39 -- *_44;
    anonymous_var_graphppl_41 -- *_44;
    anonymous_var_graphppl_41 -- NormalMeanPrecision_48;
    constvar_42_43 -- *_44;
    x_45 -- NormalMeanPrecision_48;
    constvar_46_47 -- NormalMeanPrecision_48;

_san_node_ids = Dict("dummy" => "dummy")
_lup_enabled = true
_raw_links_enabled = false
_CFFG = true
_lup_dict = Dict{String, String}(
    "θ_7" => "\\theta",
    "constvar_8_9" => "\\mu",
    "constvar_10_11" => "\\gamma",
    ## "NormalMeanPrecision_12" => "\\mathcal{N}\\_12",
    "NormalMeanPrecision_12" => "\\mathcal{N}",

    "constvar_2_3" => "\\alpha",
    "constvar_4_5" => "\\beta",
    "GammaShapeRate_6" => "\\mathcal{G}am",
    "γ_1" => "\\gamma",

    "constvar_14_15" => "\\mu",
    "constvar_16_17" => "\\gamma",
    "NormalMeanPrecision_18" => "\\mathcal{N}",
    ## "s₀_13" => "s_0", ## s₀ must always be in _lup_dict
    "s₀_13" => "\\mathbf{s}_0",
    "s_19" => "\\mathbf{s}_1",
    "s_29" => "\\mathbf{s}_2",
    "s_39" => "\\mathbf{s}_3",

    "AR_20" => "AR",
    "AR_30" => "AR",
    "AR_40" => "AR",

    "constvar_22_23" => "\\mathbf{B}",
    "constvar_32_33" => "\\mathbf{B}",
    "constvar_42_43" => "\\mathbf{B}",

    "*_24" => "*",
    "*_34" => "*",
    "*_44" => "*",
    
    "anonymous_var_graphppl_21" => "\\mathbf{v}_u",
    "anonymous_var_graphppl_31" => "\\mathbf{v}_u",
    "anonymous_var_graphppl_41" => "\\mathbf{v}_u",

    "constvar_26_27" => "\\tau",
    "constvar_36_37" => "\\tau",
    "constvar_46_47" => "\\tau",

    "NormalMeanPrecision_28" => "\\mathcal{N}",
    "NormalMeanPrecision_38" => "\\mathcal{N}",
    "NormalMeanPrecision_48" => "\\mathcal{N}",

    "x_25" => "x_1",
    "x_35" => "x_2",
    "x_45" => "x_3",
)
lar_tikz = generate_tikz(
    heading = "Time-Varying Autoregressive Model (Univariate)",
    Model = lar_gppl_model,
    layers = [ ## fac layers, i.e. not var (vars are handled with assoced layer)
        ## "cntrl", 
        "paramB1",
        "paramB2",
        "state1", 
        "state2",       
        "state3", 
        "obser", 
        ## "prefr", 
        ## "paramA"
    ],
    paramB1_iniparvars = ["constvar_8_9", "constvar_10_11"],
    paramB1_inifac = "NormalMeanPrecision_12",   
    paramB1_parvars = ["θ_7"],

    paramB2_iniparvars = ["constvar_2_3", "constvar_4_5"],
    paramB2_inifac = "GammaShapeRate_6",         
    paramB2_parvars =  ["γ_1"],

    state1_iniparvars = ["constvar_14_15", "constvar_16_17"],
    state1_inifac = "NormalMeanPrecision_18",
    state1_sysvars = ["s₀_13", "s_19", "s_29", "s_39"],
    state1_facs = ["AR_20", "AR_30", "AR_40"],

    state2_iniparvars = [], 
    state2_inifac = "",
    state2_sysvars = [],      
    state2_facs = ["*_24", "*_34", "*_44", ],         
    state2_parvars = ["constvar_22_23", "constvar_32_33", "constvar_42_43"],

    state3_iniparvars = [],
    state3_inifac = "",    
    state3_sysvars = [],        
    state3_facs = ["anonymous_var_graphppl_21", "anonymous_var_graphppl_31", "anonymous_var_graphppl_41"],       
    state3_parvars = [], 
        
    obser_sysvars = ["x_25", "x_35", "x_45"],
    obser_facs = ["NormalMeanPrecision_28", "NormalMeanPrecision_38", "NormalMeanPrecision_48"],    
    obser_parvars = ["constvar_26_27", "constvar_36_37", "constvar_46_47"],
)
print(lar_tikz)
show_tikz(lar_tikz) ## write to file rather than show in notebook

\documentclass{standalone}
\usepackage{tikz}
\usetikzlibrary{calc}
\usepackage{amsmath}
\newcommand\ns{15mm} %% define node (minimum) size
\newcommand\nr{.5*\ns} %% define node radius
\newcommand\ar{.4*\ns} %% define arc radius
\newcommand\br{.19*\ns} %% define bead radius
\DeclareUnicodeCharacter{03B8}{\theta}
\DeclareUnicodeCharacter{03B3}{\gamma}
\begin{document}
\begin{tikzpicture}[
    scale=1, 
    draw=black, 
    fac/.style={rectangle, draw=black!100, fill=black!0, minimum size=\ns, inner sep=0pt},
    dlt/.style={rectangle, draw=black!100, fill=black!100, minimum size=.2*\ns},
    var/.style={circle, draw=gray!100, fill=gray!100, minimum size=1mm, inner sep=0pt},
    bus/.style={rectangle, draw=black!100, fill=black!100, minimum width=30mm, minimum height=.5mm, inner sep=0pt},
    eql/.style={rectangle, draw=black!100, fill=black!0, minimum size=.5*\ns},
    text=black,
    ]
\node[font=\Large] at (4,17) {\textbf{Time-Varying Autoregressive Model (Univariate)}};
⋮
⋮
⋮
\path[use as bounding box] ([xshift=-1cm, yshift=-1cm]current bounding box.south west) rectangle ([xshift=1cm, yshift=2cm]current bounding box.north east);
\end{tikzpicture}
\end{document}
\documentclass{standalone}
\usepackage{tikz}
\usetikzlibrary{calc}
\usepackage{amsmath}
\newcommand\ns{15mm} %% define node (minimum) size
\newcommand\nr{.5*\ns} %% define node radius
\newcommand\ar{.4*\ns} %% define arc radius
\newcommand\br{.19*\ns} %% define bead radius
\DeclareUnicodeCharacter{03B8}{\theta}
\DeclareUnicodeCharacter{03B3}{\gamma}
\begin{document}
\begin{tikzpicture}[
    scale=1, 
    draw=black, 
    fac/.style={rectangle, draw=black!100, fill=black!0, minimum size=\ns, inner sep=0pt},
    dlt/.style={rectangle, draw=black!100, fill=black!100, minimum size=.2*\ns},
    var/.style={circle, draw=gray!100, fill=gray!100, minimum size=1mm, inner sep=0pt},
    bus/.style={rectangle, draw=black!100, fill=black!100, minimum width=30mm, minimum height=.5mm, inner sep=0pt},
    eql/.style={rectangle, draw=black!100, fill=black!0, minimum size=.5*\ns},
    text=black,
    ]
\node[font=\Large] at (4,17) {\textbf{Time-Varying Autoregressive Model (Univariate)}};
%% --------------------------- NODES ---------------------------
\node(NormalMeanPrecision_12)[fac] at (2, 15) {$\mathcal{N}$};
\node(constvar_8_9)[var] at ($ (NormalMeanPrecision_12.south west) + (-.5*\ns, 0.3333333333333333*\ns) $) {};
\node[above] at (constvar_8_9) {$\mu$};
\node(constvar_8_9)[dlt] at ($ (constvar_8_9.center) + (-.4*\ns, 0) $) {};
\node(constvar_10_11)[var] at ($ (NormalMeanPrecision_12.south west) + (-.5*\ns, 0.6666666666666666*\ns) $) {};
\node[above] at (constvar_10_11) {$\gamma$};
\node(constvar_10_11)[dlt] at ($ (constvar_10_11.center) + (-.4*\ns, 0) $) {};
\node(θ_7)[eql] at (4, 15) {$=$};
\node[xshift=.4*\ns, yshift=.4*\ns] at (θ_7) {$\theta$};
\node(GammaShapeRate_6)[fac] at (2, 12) {$\mathcal{G}am$};
\node(constvar_2_3)[var] at ($ (GammaShapeRate_6.south west) + (-.5*\ns, 0.3333333333333333*\ns) $) {};
\node[above] at (constvar_2_3) {$\alpha$};
\node(constvar_2_3)[dlt] at ($ (constvar_2_3.center) + (-.4*\ns, 0) $) {};
\node(constvar_4_5)[var] at ($ (GammaShapeRate_6.south west) + (-.5*\ns, 0.6666666666666666*\ns) $) {};
\node[above] at (constvar_4_5) {$\beta$};
\node(constvar_4_5)[dlt] at ($ (constvar_4_5.center) + (-.4*\ns, 0) $) {};
\node(γ_1)[eql] at (4, 12) {$=$};
\node[xshift=.4*\ns, yshift=.4*\ns] at (γ_1) {$\gamma$};
\node(NormalMeanPrecision_18)[fac] at (2, 9) {$\mathcal{N}$};
\node(constvar_14_15)[var] at ($ (NormalMeanPrecision_18.south west) + (-.5*\ns, 0.3333333333333333*\ns) $) {};
\node[above] at (constvar_14_15) {$\mu$};
\node(constvar_14_15)[dlt] at ($ (constvar_14_15.center) + (-.4*\ns, 0) $) {};
\node(constvar_16_17)[var] at ($ (NormalMeanPrecision_18.south west) + (-.5*\ns, 0.6666666666666666*\ns) $) {};
\node[above] at (constvar_16_17) {$\gamma$};
\node(constvar_16_17)[dlt] at ($ (constvar_16_17.center) + (-.4*\ns, 0) $) {};
\node(s₀_13)[var] at (4, 9) {};
\node[above] at (s₀_13) {$\mathbf{s}_0$};
\node(AR_20)[fac] at (6, 9) {$AR$};
\draw (AR_20)++(\ar,0mm) arc[start angle=0, end angle=-180, radius=\ar];
\node(AR_30)[fac] at (10, 9) {$AR$};
\draw (AR_30)++(\ar,0mm) arc[start angle=0, end angle=-180, radius=\ar];
\node(AR_40)[fac] at (14, 9) {$AR$};
\draw (AR_40)++(\ar,0mm) arc[start angle=0, end angle=-180, radius=\ar];
\node(s_19)[eql] at (8, 9) {$=$};
\node[xshift=.4*\ns, yshift=.4*\ns] at (s_19) {$\mathbf{s}_1$};
\node(s_29)[eql] at (12, 9) {$=$};
\node[xshift=.4*\ns, yshift=.4*\ns] at (s_29) {$\mathbf{s}_2$};
\node(s_39)[var] at (16, 9) {};
\node[above] at (s_39) {$\mathbf{s}_3$};
\node(*_24)[fac] at (8, 6) {$*$};
\node(*_34)[fac] at (12, 6) {$*$};
\node(*_44)[fac] at (16, 6) {$*$};
\node(constvar_22_23)[var] at ($ (*_24.south west) + (-.5*\ns, 0.5*\ns) $) {};
\node[above] at (constvar_22_23) {$\mathbf{B}$};
\node(constvar_22_23)[dlt] at ($ (constvar_22_23.center) + (-.4*\ns, 0) $) {};
\node(constvar_32_33)[var] at ($ (*_34.south west) + (-.5*\ns, 0.5*\ns) $) {};
\node[above] at (constvar_32_33) {$\mathbf{B}$};
\node(constvar_32_33)[dlt] at ($ (constvar_32_33.center) + (-.4*\ns, 0) $) {};
\node(constvar_42_43)[var] at ($ (*_44.south west) + (-.5*\ns, 0.5*\ns) $) {};
\node[above] at (constvar_42_43) {$\mathbf{B}$};
\node(constvar_42_43)[dlt] at ($ (constvar_42_43.center) + (-.4*\ns, 0) $) {};
\node(anonymous_var_graphppl_21)[fac] at (8, 3) {$\mathbf{v}_u$};
\node(anonymous_var_graphppl_31)[fac] at (12, 3) {$\mathbf{v}_u$};
\node(anonymous_var_graphppl_41)[fac] at (16, 3) {$\mathbf{v}_u$};
\node(NormalMeanPrecision_28)[fac] at (8, 0) {$\mathcal{N}$};
\node(NormalMeanPrecision_38)[fac] at (12, 0) {$\mathcal{N}$};
\node(NormalMeanPrecision_48)[fac] at (16, 0) {$\mathcal{N}$};
\node(constvar_26_27)[var] at ($ (NormalMeanPrecision_28.south west) + (-.5*\ns, 0.5*\ns) $) {};
\node[above] at (constvar_26_27) {$\tau$};
\node(constvar_26_27)[dlt] at ($ (constvar_26_27.center) + (-.4*\ns, 0) $) {};
\node(constvar_36_37)[var] at ($ (NormalMeanPrecision_38.south west) + (-.5*\ns, 0.5*\ns) $) {};
\node[above] at (constvar_36_37) {$\tau$};
\node(constvar_36_37)[dlt] at ($ (constvar_36_37.center) + (-.4*\ns, 0) $) {};
\node(constvar_46_47)[var] at ($ (NormalMeanPrecision_48.south west) + (-.5*\ns, 0.5*\ns) $) {};
\node[above] at (constvar_46_47) {$\tau$};
\node(constvar_46_47)[dlt] at ($ (constvar_46_47.center) + (-.4*\ns, 0) $) {};
\node(x_25)[var] at (8, -2) {};
\node[above right] at (x_25) {$x_1$};
\node(x_25)[dlt] at ($ (x_25.center) + (0, -.4*\ns) $) {};
\node(x_35)[var] at (12, -2) {};
\node[above right] at (x_35) {$x_2$};
\node(x_35)[dlt] at ($ (x_35.center) + (0, -.4*\ns) $) {};
\node(x_45)[var] at (16, -2) {};
\node[above right] at (x_45) {$x_3$};
\node(x_45)[dlt] at ($ (x_45.center) + (0, -.4*\ns) $) {};
%% --------------------------- LINKS ---------------------------
\draw (constvar_8_9) -- ([shift={(0, 0.3333333333333333*\ns)}]NormalMeanPrecision_12.south west);
\draw (constvar_10_11) -- ([shift={(0, 0.6666666666666666*\ns)}]NormalMeanPrecision_12.south west);
\draw (NormalMeanPrecision_12) -| (θ_7);
\draw (θ_7) -| ([shift={(-.2*\ns, 0)}]AR_20.north east);
\draw (θ_7) -| ([shift={(-.2*\ns, 0)}]AR_30.north east);
\draw (θ_7) -| ([shift={(-.2*\ns, 0)}]AR_40.north east);
\draw (constvar_2_3) -- ([shift={(0, 0.3333333333333333*\ns)}]GammaShapeRate_6.south west);
\draw (constvar_4_5) -- ([shift={(0, 0.6666666666666666*\ns)}]GammaShapeRate_6.south west);
\draw (GammaShapeRate_6) -| (γ_1);
\draw (γ_1) -| ([shift={(.2*\ns, 0)}]AR_20.north west);
\draw (γ_1) -| ([shift={(.2*\ns, 0)}]AR_30.north west);
\draw (γ_1) -| ([shift={(.2*\ns, 0)}]AR_40.north west);
\draw (constvar_14_15) -- ([shift={(0, 0.3333333333333333*\ns)}]NormalMeanPrecision_18.south west);
\draw (constvar_16_17) -- ([shift={(0, 0.6666666666666666*\ns)}]NormalMeanPrecision_18.south west);
\draw (NormalMeanPrecision_18) -| (s₀_13);
\draw ([shift={(\ar,0)}]AR_20.center) -- (s_19.west);
\draw ([shift={(\ar,0)}]AR_30.center) -- (s_29.west);
\draw ([shift={(\ar,0)}]AR_40.center) -- (s_39.west);
\draw ([shift={(-\ar,0)}]AR_20.center) -- (s₀_13.east);
\draw ([shift={(-\ar,0)}]AR_30.center) -- (s_19.east);
\draw ([shift={(-\ar,0)}]AR_40.center) -- (s_29.east);
\draw (s_19) -- (*_24);
\draw (s_29) -- (*_34);
\draw (s_39) -- (*_44);
\draw (constvar_22_23) -- ([shift={(0, 0.5*\ns)}]*_24.south west);
\draw (constvar_32_33) -- ([shift={(0, 0.5*\ns)}]*_34.south west);
\draw (constvar_42_43) -- ([shift={(0, 0.5*\ns)}]*_44.south west);
\draw (*_24) -- (anonymous_var_graphppl_21);
\draw (*_34) -- (anonymous_var_graphppl_31);
\draw (*_44) -- (anonymous_var_graphppl_41);
\draw (anonymous_var_graphppl_21) -- (NormalMeanPrecision_28);
\draw (anonymous_var_graphppl_31) -- (NormalMeanPrecision_38);
\draw (anonymous_var_graphppl_41) -- (NormalMeanPrecision_48);
\draw (constvar_26_27) -- ([shift={(0, 0.5*\ns)}]NormalMeanPrecision_28.south west);
\draw (constvar_36_37) -- ([shift={(0, 0.5*\ns)}]NormalMeanPrecision_38.south west);
\draw (constvar_46_47) -- ([shift={(0, 0.5*\ns)}]NormalMeanPrecision_48.south west);
\draw (NormalMeanPrecision_28) -- (x_25);
\draw (NormalMeanPrecision_38) -- (x_35);
\draw (NormalMeanPrecision_48) -- (x_45);
\node(θ_7)[eql] at (θ_7) {$=$};
\node(γ_1)[eql] at (γ_1) {$=$};
\draw (NormalMeanPrecision_12.east) circle(\br)[fill=white];
\draw (GammaShapeRate_6.east) circle(\br)[fill=white];
\draw (NormalMeanPrecision_18.east) circle(\br)[fill=white];
\draw (NormalMeanPrecision_28.north) circle(\br)[fill=white];
\draw (NormalMeanPrecision_28.south)++(-1*\br,-1*\br) rectangle([shift={(\br,\br)}]NormalMeanPrecision_28.south)[fill=white];
\draw (NormalMeanPrecision_38.north) circle(\br)[fill=white];
\draw (NormalMeanPrecision_38.south)++(-1*\br,-1*\br) rectangle([shift={(\br,\br)}]NormalMeanPrecision_38.south)[fill=white];
\draw (NormalMeanPrecision_48.north) circle(\br)[fill=white];
\draw (NormalMeanPrecision_48.south)++(-1*\br,-1*\br) rectangle([shift={(\br,\br)}]NormalMeanPrecision_48.south)[fill=white];
\path[use as bounding box] ([xshift=-1cm, yshift=-1cm]current bounding box.south west) rectangle ([xshift=1cm, yshift=2cm]current bounding box.north east);
\end{tikzpicture}
\end{document}

9051

Time-Varying Autoregressive Model (Univariate)