{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# 5 - Compare one-step rate distortion against general three-variable case\n",
"This notebook is part of the supplementary material for: \n",
"Genewein T., Leibfried F., Grau-Moya J., Braun D.A. (2015) *Bounded rationality, abstraction and hierarchical decision-making: an information-theoretic optimality principle*, Frontiers in Robotics and AI. \n",
"\n",
"More information on how to run the notebook on the accompanying [github repsitory](https://github.com/tgenewein/BoundedRationalityAbstractionAndHierarchicalDecisionMaking) where you can also find updated versions of the code and notebooks.\n",
"\n",
"This notebook corresponds to Sections 4 and 5 and reproduces Figure 11 of the paper."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Compare any three-variable case against the one-step case\n",
"\n",
"This notebook allows to compare the solutions of any three-variable case (general, parallel, serial) against the solution of the one-step rate-distortion case (Section 2 of the paper). By default, it illustrates how information processing in the parallel hierarchical case can be split up on different information processing pathways (see Section 4.3 of the paper).\n",
"\n",
"Below we show (in the following order)\n",
"* Utility function (i.e. task-setup\n",
"* Solution for the one-step case (Section 2 of the paper)\n",
"* Solution for the general case (corrsponding to Section 5 of the paper)\n",
"* Comparison between both solutions (one-step and general)"
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"WARNING: replacing module RateDistortionDecisionMaking\n"
]
},
{
"data": {
"text/plain": [
"RateDistortionDecisionMaking"
]
},
"execution_count": 15,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"#only run this once\n",
"include(\"RateDistortionDecisionMaking.jl\")"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"using RateDistortionDecisionMaking, Distances, DataFrames, Colors, Gadfly, Distributions, Interact, Reactive\n",
"\n",
"#make the default plot size a bit larger\n",
"set_default_plot_size(15cm, 12cm)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Setup task $U(a,w)$\n",
"\n",
"By default, we re-use the medical system example (see main paper Section 4 or the Jupyter notebook number 4 for more details).\n",
"\n",
"### [Interact] Use any of the examples shown in the other notebooks by changing the utility function in the code-cell below"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"image/svg+xml": [
"\n",
"\n"
],
"text/html": [
"\n",
"\n"
],
"text/plain": [
"Plot(...)"
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"#set up taxonomy example\n",
"#include(\"TaxonomyExample.jl\")\n",
"#w_vec, w_strings, a_vec, a_strings, p_w, U = setuptaxonomy()\n",
"#a_values = a_vec\n",
"#w_values = w_vec\n",
"\n",
"#set up medical system example\n",
"include(\"MedicalSystemExample.jl\")\n",
"w_values, w_strings, a_values, a_strings, p_w, U = setup_medical_example(uniform_w=true)\n",
"num_acts = length(a_values)\n",
"a_vec = collect(1:num_acts)\n",
"num_worldstates = length(w_values)\n",
"w_vec = collect(1:num_worldstates)\n",
"\n",
"\n",
"\n",
"\n",
"#pre-compute utilities, find maxima\n",
"U_pre, Umax = setuputilityarrays(a_values,w_values,U)\n",
"\n",
"#visualize utility\n",
"plt_utility = visualizeMatrix(U_pre, w_values, a_values, w_strings, a_strings, xlabel=\"Disease type w\",\n",
" ylabel=\"Treatment a\", legendlabel=\"U(a,w)\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## [Interact] Set parameters for one-step and general case here\n",
"\n",
"Make sure to re-evaluate all the subsequent cells after changing parameters here. The one-step case corresponds to Section 2 of the paper. The general case is found in Section 5 and by choosing certain settings of the inverse temperatures (see Table 2 in the paper), the parallel, serial or even the one-step case can be recovered."
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"2000"
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"#one-step Blahut-Arimoto\n",
"β = 10 #inverse temperature\n",
"\n",
"#general case\n",
"\n",
"#β1=β3 and β2=∞ is identical to the one-variable case\n",
"#β1 = β #inverse temperature for p(x)→p(x|w)\n",
"#β2 = Inf #inverse temperature for p(a)→p(a|x)\n",
"#β3 = β #inverse temperature for p(a|x)→p(a|x,w)\n",
"\n",
"#parallel case: β2=∞\n",
"β1 = 15 #inverse temperature for p(x)→p(x|w)\n",
"β2 = Inf #17.5 #inverse temperature for p(a)→p(a|x)\n",
"β3 = 14.99 #3.999 #inverse temperature for p(a|x)→p(a|x,w)\n",
"\n",
"#sequential case: β3=0\n",
"#β1 = 1 #inverse temperature for p(o)→p(x|w)\n",
"#β2 = 8 #inverse temperature for p(a)→p(a|x)\n",
"#β3 = 0 #inverse temperature for p(a|x)→p(a|x,w)\n",
"\n",
"\n",
"#cardinality of percept\n",
"num_x = 3#length(a_vec) #cardinality of X is equal to action-space\n",
"\n",
"\n",
"ε = 0.0001 #convergence critetion for BAiterations\n",
"maxiter = 2000 #maximum number of BA iterations"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# One-step Blahut-Arimoto (rate-distortion) - as in Section 2 of the paper"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/html": [
"\n",
"\n"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"\n",
"\n"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"#solve the example with one-step Blahut-Arimoto\n",
"\n",
"#initialize p(a) uniformly\n",
"num_acts = length(a_vec)\n",
"pa_init = ones(num_acts)/num_acts \n",
"\n",
"#Blahut-Arimotot iterations\n",
"p_agw_one, p_a_one, perf_df_one = BAiterations(pa_init, β, U_pre, p_w, ε, maxiter,\n",
" compute_performance=true, performance_as_dataframe=true)\n",
"\n",
"#visualize solution\n",
"plt_marg, plt_cond = visualizeBAsolution(p_a_one, p_agw_one, a_vec, w_vec, a_strings, w_strings,\n",
" wlabel=\"Disease type w\", alabel=\"Treatment a\",\n",
" legendlabel_marginal=\"p(a)\", legendlabel_conditional=\"p*(a|w)\",\n",
" suppress_vis=true);\n",
"\n",
"display(plt_utility)\n",
"display(plt_cond)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Three-variable general case as in Section 5 of the paper"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"x_vec = collect(1:num_x)\n",
"x_strings = map((x)->\"x\"*string(x), x_vec)\n",
"\n",
"num_worldstates = length(w_vec)\n",
"\n",
"\n",
"\n",
"\n",
"#This function performs Blahut-Arimoto iterations for the three-variable general case\n",
"px, pa, pxgw, pagx, pagxw, pagw, performance_df = threevarBAiterations(num_x, β1, β2, β3, U_pre, p_w, ε, maxiter,\n",
" compute_performance=true, performance_per_iteration=true,\n",
" performance_as_dataframe=true,\n",
" init_pogw_sparse = true, init_pogw_uniformly = false,\n",
" init_pagow_uniformly = true);"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Show convergence behavior of iteration of self consistent equations by inspecting information terms, expected utility and the value of the objective\n",
"\n",
"**! Note that in the legends below $O$ is used instead of $X$**"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {
"collapsed": false,
"scrolled": false
},
"outputs": [
{
"data": {
"text/html": [
"\n",
"\n"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"#plot convergence (i.e. evolution across BA iterations)\n",
"plot_convergence_MI, plot_convergence_Util = plot_three_var_BA_convergence(performance_df);"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Solution of general case\n",
"\n",
"**Note that $X$ is labeled here with the general term \"Observation\" but in case of the medical example it indicates the first diagnosis (or model)**"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {
"collapsed": false,
"scrolled": false
},
"outputs": [
{
"data": {
"text/html": [
"\n",
"\n"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"\n",
"\n"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"\n",
"\n"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"\n",
"\n"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"#call the routine that creates plots for all the probability distributions\n",
"plt_px, plt_pa, plt_pxgw, plt_pagx, plt_pagw, dpdown, plt_pagxw_vis = visualize_three_var_BAsolution(px, pa,\n",
" pxgw, pagx, pagxw, pagw,\n",
" x_vec, a_vec, w_vec,\n",
" x_strings, a_strings, w_strings, \n",
" olabel_string=\"x\", alabel_string=\"a\", wlabel_string=\"w\")\n",
"\n",
"#visualize p(x) and p(a)\n",
"display(hstack(plt_px,plt_pa))\n",
"\n",
"#visualize p(x|w)\n",
"display(plt_pxgw)\n",
"#visualize p(a|x)\n",
"display(plt_pagx)\n",
"\n",
"#visualize p(a|w)\n",
"display(plt_pagw)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### $p(a|x,w)$ of general case"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/html": [],
"text/plain": [
"Interact.Options{:Dropdown,ASCIIString}([Reactive.Input{ASCIIString}] h1,\"World state w\",\"h1\",\"h1\",Interact.OptionDict(Any[\"h1\",\"h2\",\"l1\",\"l2\",\"l3\",\"l4\"],Dict{Any,Any}(\"h1\"=>\"h1\",\"h2\"=>\"h2\",\"l2\"=>\"l2\",\"l4\"=>\"l4\",\"l3\"=>\"l3\",\"l1\"=>\"l1\")),Any[],Any[])"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"image/svg+xml": [
"\n",
"\n"
],
"text/html": [
""
],
"text/plain": [
"Plot(...)"
]
},
"metadata": {
"comm_id": "f6a45565-1a15-428d-87ce-b4529251d3ac",
"reactive": true
},
"output_type": "display_data"
}
],
"source": [
"#visualize p(a|x,w)\n",
"display(dpdown) #show the dropdown box to interactively select wk in p(a|x,w=wk)\n",
"display(plt_pagxw_vis)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"If you see the drowpdown-box but not the plot and ``\"Javascript error adding output!\"``, try selecting another $w$ with the dropdwon once, then the plot should appear.\n",
"\n",
"If this does not work (i.e. kernel is busy with this cell but nothing ever happens), try ckecking out the most current version of Interact with: `Pkg.checkout(\"Interact\")` in a Julia console. You can undo this and go back to the latest released version with `Pkg.free(\"Interact\")`."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Performance measures of general case\n",
"\n",
"**! Note that in the legends below $O$ is used instead of $X$**"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"collapsed": false,
"scrolled": false
},
"outputs": [
{
"data": {
"text/html": [
"\n",
"\n"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"\n",
"\n"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"p_MI, p_composed, p_perf=plot_three_var_performancemeasures(performance_df, maximum(U_pre), β1, β2, β3)\n",
"\n",
"display(hstack(p_MI, p_composed))\n",
"display(p_perf)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Comparison between one-step case (left) and general case (right)\n",
"\n",
"Whenever the same quantity is shown twice in any of the plots, the left quantity corresponds to the one-step case (denoted by subindex ``_one`` and the right quantity corresponds to the general case."
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/html": [
"\n",
"\n"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"#compare solution against one-step case\n",
"\n",
"#compare performance measures, EU, J, I(A;W)\n",
"\n",
"util_all = [perf_df_one[end,:E_U], performance_df[end,:E_U],\n",
" perf_df_one[end,:RD_obj], performance_df[end,:Objective_value]]\n",
"\n",
"color_label = [\"E[U]_one\", \"E[U]\", \"J_one\", \"J\"]\n",
"\n",
"colors = BAdiscretecolorscale(4).f(4)\n",
"colscale = Scale.color_discrete_manual(colors[1],colors[1],colors[end],colors[end])\n",
"\n",
"p_perf_comp = plot(x=color_label, y=util_all, color=color_label, Geom.bar(),\n",
" Guide.ylabel(\"[utils]\"), Guide.xlabel(\"\"), Guide.colorkey(\"\"),\n",
" Guide.title(\"One-step solution β=$β vs. \\n multi-step β1=$β1, β2=$β2, β3 = $β3\"),\n",
" colscale, Scale.y_continuous(minvalue=0, maxvalue=maximum(U_pre)),\n",
" BAtheme(key_position = :none))\n",
"\n",
"\n",
"#compare mutual informations I(A;W)\n",
"MI_comp = [perf_df_one[end,:I_aw], performance_df[end,:I_aw]]\n",
"color_label = [\"I(A;W)_one\", \"I(A;W)\"]\n",
"\n",
"colors = Scale.color_discrete_hue().f(4)\n",
"colscale = Scale.color_discrete_manual(colors[end],colors[end])\n",
"\n",
"p_perf_comp_MI = plot(x=color_label, y=MI_comp, color=color_label, Geom.bar(),\n",
" Guide.ylabel(\"[bits]\"), Guide.xlabel(\"\"), Guide.colorkey(\"\"),\n",
" Guide.title(\"I(A;W) for β=$β vs. \\n multi-step β1=$β1, β2=$β2, β3 = $β3\"),\n",
" colscale, Scale.y_continuous(minvalue=0), BAtheme(key_position = :none))\n",
"\n",
"display(hstack(p_perf_comp,p_perf_comp_MI))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Compare overall behavior $p(a|w)$\n",
"\n",
"Top: one-step case\n",
"Bottom: general case"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/html": [
"\n",
"\n"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"\n",
"\n"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"#compare solutions p(a|w)\n",
"plt_cond_one = visualizeBAconditional(p_agw_one, a_vec, w_vec, a_strings, w_strings,\n",
" wlabel=\"w\", alabel=\"a\",legendlabel=\"p*(a|w)_one\")\n",
"\n",
"plt_cond_multi = visualizeBAconditional(pagw, a_vec, w_vec, a_strings, w_strings,\n",
" wlabel=\"w\", alabel=\"a\",legendlabel=\"p*(a|w)\")\n",
"\n",
"display(plt_cond_one)\n",
"display(plt_cond_multi)\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Difference between both final solutions"
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"\"Difference: norm(p*(a|w)_one - p*(a|w))=6.118579364182079e-5\""
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"image/svg+xml": [
"\n",
"\n"
],
"text/html": [
"\n",
"\n"
],
"text/plain": [
"Plot(...)"
]
},
"execution_count": 13,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"display(\"Difference: norm(p*(a|w)_one - p*(a|w))=$(norm(p_agw_one-pagw))\")\n",
"\n",
"plt_cond_diff = visualizeMatrix(p_agw_one - pagw, w_vec, a_vec, w_strings, a_strings,\n",
" xlabel=\"Disease type w\", ylabel=\"Treatment a\",legendlabel=\"p*(a|w)_one - p*(a|w)\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Compose final plot for paper\n",
"\n",
"**! Note that in the legends below $O$ is used instead of $X$**"
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/html": [
"\n",
"\n"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"Compose.Measure{Compose.MeasureNil,Compose.MeasureNil}(150.0,Compose.MeasureNil(),Compose.MeasureNil(),0.0,0.0)"
]
},
"execution_count": 14,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"#for the paper\n",
"plt_final = vstack((hstack(plt_cond_multi, p_perf_comp_MI )), hstack(p_MI, p_perf_comp))\n",
"display(plt_final)\n",
"\n",
"w = 18cm\n",
"h = 15cm\n",
"#draw(SVG(\"Figures/Hierarchical_OneStep_Comparison.svg\", w, h), plt_final) #uncomment to store figure#\n"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Julia 0.4.0",
"language": "julia",
"name": "julia-0.4"
},
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "0.4.0"
}
},
"nbformat": 4,
"nbformat_minor": 0
}